An image calibration-based farmland remote sensing monitoring method and system

Abstract

This invention relates to the field of farmland remote sensing monitoring and image processing technology, and discloses a farmland remote sensing monitoring method and system based on image calibration. The method includes acquiring original image sequences and historical image sequences; performing edge extraction and spectral difference calculation on the original image sequences to obtain preliminary labeled images; using texture parameters and Kriging interpolation algorithms to reconstruct pixels in occluded areas of the preliminary labeled images to obtain reconstructed regions; performing spectral verification on the reconstructed regions to obtain de-occluded restored images, and performing illumination consistency calibration on the de-occluded restored images to obtain illumination-balanced images; calculating the normalized vegetation index based on the illumination-balanced images and performing spatial clustering to obtain anomalous distribution images. This method can achieve high-fidelity restoration and accurate health status monitoring of images under complex weather conditions, effectively solving the assessment distortion problems caused by uneven illumination and cloud cover.

Description

Technical Field

[0001] This invention relates to the field of farmland remote sensing monitoring and image processing technology, and in particular to a farmland remote sensing monitoring method and system based on image calibration. Background Technology

[0002] Currently, farmland remote sensing monitoring technology, as a core means of precision management in modern agriculture, obtains vegetation growth status and environmental parameters through multispectral Earth observation, and has become an important pillar for ensuring food security and achieving sustainable agricultural development.

[0003] In existing technologies, large-area farmland images are typically acquired using drones or satellite platforms equipped with multispectral sensors. These images are then combined with pre-defined algorithms to perform image recognition, image classification, or vegetation index inversion to achieve real-time assessment of crop health. However, when monitoring farmland in remote mountainous areas, complex weather conditions (such as cloud and rain) often result in extensive cloud cover and shadow interference in the images, leading to incomplete raw data. Furthermore, the uneven distribution of natural light due to the undulating terrain further distorts the spectral characteristics of the images, making direct comparison of monitoring results from different areas difficult and severely impacting the reliability of vegetation health index calculations. Although existing technologies attempt to address this through simple interpolation or filtering, the lack of precise segmentation of cloud edges and in-depth calibration for illumination consistency makes it highly susceptible to misjudgments in extreme environments, failing to meet the decision-making needs of precision agriculture.

[0004] Existing technologies suffer from low accuracy in monitoring the health of vegetation in remote mountainous farmland. Summary of the Invention

[0005] This invention provides a method and system for remote sensing monitoring of farmland based on image calibration, in order to solve the problem of low accuracy in monitoring the health of vegetation in remote mountain farmland in the prior art.

[0006] In a first aspect, to address the aforementioned technical problems, the present invention provides a method for remote sensing monitoring of farmland based on image calibration, comprising: Obtain the original image sequence, and perform edge extraction and spectral difference calculation on the original image sequence to obtain a preliminary labeled image; Extract the cloud-occluded area and the surrounding unoccluded area from the preliminary labeled image. Calculate the texture parameters of the surrounding unoccluded area using the gray-level co-occurrence matrix algorithm and select reference samples. Based on the reference samples, use the Kriging interpolation algorithm to perform pixel reconstruction on the cloud-occluded area to obtain the reconstructed area. Local spectral smoothing calculations are performed on the reconstructed region, and consistency verification is performed in conjunction with a preset spectral deviation threshold to obtain the de-occluded restored image; Spatial smoothing filtering is performed on the de-occluded and restored image to obtain an optimized clear image; The global brightness deviation is obtained by performing histogram statistical calculation on the optimized clear image, and the brightness-adjusted image is obtained by adjusting the band gain based on the global brightness deviation. Dynamic pixel value adjustment is performed on pixels in the brightness adjustment image that exceed a preset brightness threshold to obtain an illumination equalization image; The near-infrared and red band data of the illumination-balanced image are extracted to calculate the normalized vegetation index. Based on the normalized vegetation index and a preset anomaly threshold, spatial clustering is performed to obtain anomaly distribution images.

[0007] Secondly, the present invention provides a farmland remote sensing monitoring system based on image calibration, comprising: The edge segmentation module is used to acquire the original image sequence and perform edge extraction and spectral difference calculation on the original image sequence to obtain a preliminary labeled image. The texture reconstruction module is used to extract the cloud-occluded area and the surrounding non-occluded area in the preliminary labeled image, calculate the texture parameters of the surrounding non-occluded area using the gray-level co-occurrence matrix algorithm and select reference samples, and perform pixel reconstruction on the cloud-occluded area using the kriging interpolation algorithm based on the reference samples to obtain the reconstructed area. The spectral verification module is used to perform local spectral smoothing calculations on the reconstructed region and perform consistency verification in combination with a preset spectral deviation threshold to obtain the de-occluded restored image. The smoothing and noise reduction module is used to perform spatial smoothing filtering calculations on the de-occlusion restored image to obtain an optimized clear image. The band gain module is used to perform histogram statistical calculations on the optimized clear image to obtain the global brightness deviation, and adjust the band gain of the optimized clear image according to the global brightness deviation to obtain a brightness-adjusted image. The pixel equalization module is used to dynamically adjust the pixel values ​​of pixels in the brightness adjustment image that exceed a preset brightness threshold to obtain a light-balanced image. The clustering and localization module is used to extract near-infrared band data and red band data from the illumination-balanced image to calculate the normalized vegetation index, and to perform spatial clustering calculation based on the normalized vegetation index and a preset anomaly threshold to obtain anomaly distribution images.

[0008] Compared with the prior art, the present invention has the following beneficial effects: (1) This invention utilizes the gray-level co-occurrence matrix to extract texture features and combines kriging interpolation and spectral consistency verification to reconstruct pixels in occluded areas, thus overcoming the edge blurring and spectral deviation problems caused by traditional interpolation algorithms that rely solely on the spatial neighborhood mean. This mechanism achieves high-fidelity restoration of occluded areas in both spatial texture and physical spectral properties, ensuring that the restored image can accurately participate in subsequent vegetation index calculations and significantly improving the usability of remote sensing data under cloudy and rainy weather conditions.

[0009] (2) By performing global brightness deviation statistics and combining them with adaptive pixel value adjustment technology, this invention can effectively eliminate local illumination unevenness caused by mountainous terrain undulations, cloud shadows, and sensor angle fluctuations. This dynamic calibration mechanism weakens the interference of non-vegetation factors on pixel reflectivity from the source, avoids false vegetation decline misjudgments caused by complex mountainous environments, and provides a physically consistent data foundation for cross-regional and cross-time period farmland health monitoring.

[0010] (3) This invention achieves a technological leap from single-point-of-time image recognition to long-cycle state evolution analysis by integrating DBSCAN spatial clustering and ARIMA time series analysis models and introducing external operational data such as meteorological and fertilization data for correction. The scheme constructs anomaly judgment logic based on spatiotemporal multidimensional joint verification, which can automatically eliminate false alarms caused by temporary interference and accurately locate vegetation growth anomalies caused by insufficient water or nutrients. The resulting structured decision report has higher scientific value and engineering application value. Attached Figure Description

[0011] Figure 1 This is a schematic diagram of a farmland remote sensing monitoring method based on image calibration provided in the first embodiment of the present invention; Figure 2 This is a schematic diagram of the structure of a farmland remote sensing monitoring system based on image calibration provided in the second embodiment of the present invention. Detailed Implementation

[0012] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0013] Reference Figure 1 The first embodiment of the present invention provides a method for remote sensing monitoring of farmland based on image calibration, comprising the following steps: S11, Obtain the original image sequence, and perform edge extraction and spectral difference calculation on the original image sequence to obtain a preliminary labeled image; S12, extract the cloud-occluded area and the surrounding unoccluded area in the preliminary labeled image, calculate the texture parameters of the surrounding unoccluded area using the gray-level co-occurrence matrix algorithm and select reference samples, and use the kriging interpolation algorithm to perform pixel reconstruction on the cloud-occluded area to obtain the reconstructed area based on the reference samples. S13, perform local spectral smoothing calculation on the reconstructed area, and perform consistency verification in combination with a preset spectral deviation threshold to obtain the de-occluded restored image; S14, Spatial smoothing filtering calculation is performed on the de-occlusion restored image to obtain an optimized clear image; S15, perform histogram statistical calculation on the optimized clear image to obtain the global brightness deviation, and adjust the band gain of the optimized clear image according to the global brightness deviation to obtain the brightness adjusted image; S16, dynamically adjust the pixel values ​​of pixels in the brightness adjustment image that exceed the preset brightness threshold to obtain a balanced illumination image. S17, extract the near-infrared band data and red band data of the illumination-balanced image to calculate the normalized vegetation index, and perform spatial clustering calculation based on the normalized vegetation index and the preset anomaly threshold to obtain the anomaly distribution image.

[0014] In step S11, the original image sequence is acquired, and edge extraction and spectral difference calculation are performed on the original image sequence to obtain a preliminary labeled image, including: The Canny edge detection algorithm is used to extract gradient features from the original image sequence to obtain a binary edge map. Extract the blue light band reflectance and near-infrared band reflectance of the original image sequence, and calculate the difference between the two to obtain the reflectance difference value; Pixels with reflectance differences greater than a preset difference threshold are marked as occluded pixels. The preliminary marked image is obtained by segmenting and classifying the occluded pixels and the edge binary map using a region growing algorithm.

[0015] In one implementation, a sequence of original images covering multiple bands and channels is acquired from a drone equipped with a multispectral sensor. This original image sequence is then extracted and input into a processing unit. This embodiment utilizes the Canny edge detection algorithm to extract gradient features from the original image sequence. Specifically, a Gaussian smoothing filter is used to perform a weighted average calculation on the neighborhood of each pixel in the original image sequence to eliminate random high-frequency noise. A first-order finite difference operator is used to calculate the partial derivatives in the horizontal and vertical directions respectively. The square root of the sum of the squares of the two partial derivatives is calculated to obtain the gradient magnitude. The arctangent function of the ratio of the two partial derivatives is calculated to obtain the gradient direction. Non-maximum suppression processing is performed on the gradient magnitude along the gradient direction to retain local gradient extrema.

[0016] It should be noted that, when performing dual-threshold edge connectivity, this embodiment extracts a set of gradient magnitudes from historical images under similar weather conditions, calculates the cumulative probability distribution histogram of this set, and extracts the gradient magnitude corresponding to the 95th percentile as the high threshold parameter. The high threshold parameter is then multiplied by one-third to obtain the low threshold parameter. In this embodiment, pixels with gradient magnitudes greater than the high threshold parameter are directly identified as strong edge pixels, while pixels with gradient magnitudes between the low and high threshold parameters are identified as weak edge pixels. The eight-neighborhood of each weak edge pixel is traced; if at least one strong edge pixel exists within its eight-neighborhood, the weak edge pixel is retained; otherwise, it is discarded. Through these operations, a matrix containing only connected edge pixels is output, generating a binary edge map.

[0017] In one implementation, the reflectance of the blue light band with a center wavelength of 450 nm and the reflectance of the near-infrared band with a center wavelength of 850 nm are extracted from the original image sequence. In this embodiment, the absolute difference between the near-infrared band reflectance and the blue light band reflectance is calculated pixel by pixel to form a numerical matrix with the same resolution as the original image, thus obtaining the reflectance difference.

[0018] It is worth noting that regarding the method for determining the preset difference threshold, this embodiment collects historical manually labeled and verified cloud-occluded pixel sample sets and unoccluded bare ground pixel sample sets, and calculates the reflectance difference of all pixels in the two sample sets respectively. Using the Otsu's method, the maximum inter-class variance algorithm is used to traverse all segmentation boundaries within the possible floating-point range, calculating the intra-class variance and inter-class variance of the two classes of pixel difference sets divided by each boundary. The segmentation boundary value that maximizes the inter-class variance is extracted and objectively labeled as the preset difference threshold. This embodiment determines pixel-by-pixel whether the reflectance difference is greater than the preset difference threshold. If it is greater, a specific identifier is assigned to it in the mask matrix, marking it as an occluded pixel.

[0019] For example, assuming that the segmentation value that maximizes the inter-class variance within a certain traversal interval is 0.3, this embodiment strictly sets it as a preset difference threshold. For a pixel whose reflectance difference is 0.45, the system determines that 0.45 is greater than 0.3, and thus the pixel is marked as an occluded pixel.

[0020] In one implementation, the preliminary labeled image is obtained by segmenting and classifying the occluded pixels and the edge binary map using a region growing algorithm. In this embodiment, the intersection of the set of strong edge pixels in the edge binary map and the set of occluded pixels is performed, and the resulting overlapping pixels are used as the initial seed point set for the region growing algorithm. This embodiment iterates through the eight spatial neighbor pixels of each seed point in the initial seed point set, extracts the reflectance difference corresponding to the neighbor pixels, and calculates the absolute distance between the reflectance difference of the neighbor pixels and the reflectance difference of the seed point. If the absolute distance is less than a preset similarity threshold, the neighbor pixel is classified into the occluded region category of the current seed point, and the neighbor pixel is updated as a new seed point and added to the search queue.

[0021] It is worth noting that, regarding the preset similarity limit, this embodiment extracts the reflectance difference distribution of the aforementioned cloud-occluded pixel sample set, calculates its standard deviation, and directly uses this standard deviation as the preset similarity limit to constrain the physical range of growth. This embodiment iteratively executes the above-mentioned neighborhood search and classification process until the search queue is empty and no new pixels meet the classification conditions, at which point growth stops, generating a binarized segmentation mask for the entire spatial domain. This embodiment performs channel overlay mapping between this binarized segmentation mask and the original image sequence, forcibly colorizing pixels classified as occluded areas (e.g., uniformly assigning the maximum value to the red channel and clearing other channels), keeping the pixel values ​​of non-occluded areas unchanged, and reconstructing and outputting the final preliminary labeled image.

[0022] In step S12, the cloud-occluded region and the surrounding unoccluded region in the preliminary labeled image are extracted. The texture parameters of the surrounding unoccluded region are calculated using the gray-level co-occurrence matrix algorithm, and reference samples are selected. Based on the reference samples, the cloud-occluded region is reconstructed using the kriging interpolation algorithm to obtain the reconstructed region, including: The contrast and correlation of a preset pixel window within the surrounding unoccluded area are calculated using the gray-level co-occurrence matrix algorithm to obtain texture feature values; Based on the texture feature values, the Euclidean distance algorithm is used to select non-occluded pixels with textures similar to the pixel to be reconstructed as reference samples; wherein, the texture features of the pixel to be reconstructed are derived based on edge continuity. Based on the spatial coordinates and spectral values ​​of the reference sample, and combined with a preset distance attenuation coefficient, the Kriging interpolation algorithm is used to perform a weighted average calculation on the pixels to be reconstructed within the cloud-covered area to obtain the reconstructed area.

[0023] In one implementation, texture modeling is performed on the edges of occluded regions and their neighborhoods identified in the preliminary labeled image. Using the Gray-Level Co-occurrence Matrix (GLCM) algorithm, a preset pixel window (e.g., a 5×5 pixel window) is slid across the non-occluded region. In this embodiment, the step size is set to 1 pixel, and the gray-level distribution of pixel pairs is statistically analyzed in four directions: 0 degrees, 45 degrees, 90 degrees, and 135 degrees. Texture feature values ​​such as contrast, correlation, and entropy are calculated to quantitatively characterize the texture patterns of different crop plots or mountainous terrains.

[0024] It should be noted that, regarding the method for determining the preset pixel window size, this embodiment extracts the spatial resolution parameters of the original image and the average geometric feature size of the target farmland plot. The ratio of the average geometric feature size to the spatial resolution is used as a benchmark reference value. Through multi-scale window experiments, the coefficient of variation of texture features under different windows is calculated. The minimum window size when the coefficient of variation tends to stabilize is extracted and strictly calibrated as the preset pixel window.

[0025] It is worth noting that this embodiment uses the Euclidean distance algorithm to select reference samples based on the texture feature values. This embodiment constructs texture feature vectors and calculates the Euclidean distance between the predicted feature vector of the pixel to be reconstructed in the occluded area (derived based on edge continuity) and the feature vector of the surrounding unoccluded pixels. By determining whether the Euclidean distance is less than a preset similarity threshold (e.g., 0.2), qualified unoccluded pixels are extracted and determined as the reference samples. Regarding the preset similarity threshold, this embodiment collects the texture feature distribution of historical typical land features (such as terraces, woodlands, and bare land), plots classification ROC curves under different similarity values, and extracts the critical value when the Youden exponent is maximized as the preset similarity threshold.

[0026] It is worth noting that the texture features of the pixel to be reconstructed are derived based on edge continuity. This means extracting the pixels adjacent to the unobstructed area on the boundary of the cloud-occluded area, and directly assigning the texture feature value of the boundary pixel to the pixel to be reconstructed with the nearest Euclidean distance. If the pixel to be reconstructed is equidistant from multiple boundary pixels, the arithmetic mean of the texture feature values ​​of these boundary pixels is taken as the predicted texture feature vector of the pixel to be reconstructed.

[0027] In one implementation, Kriging interpolation based on spatial autocorrelation is performed using the spatial coordinates and spectral values ​​of the reference sample. This embodiment extracts the reflectance values ​​of the reference sample in the visible and near-infrared bands to establish a semi-variogram model. Combined with a preset distance attenuation coefficient (e.g., 0.3), weights are assigned according to the spatial geometric distance between the reference sample and the pixel to be reconstructed, and a weighted average of the spectral values ​​of the pixel to be reconstructed is calculated.

[0028] The preset distance attenuation coefficient is obtained by collecting multiple sets of historical, cloud-free, and similar farmland images, and then performing simulated occlusion reconstruction on known pixels under different attenuation coefficient values. The root mean square error between the reconstructed pixel value and the original real pixel value is calculated. The attenuation coefficient value that minimizes the root mean square error is selected and calibrated as the preset distance attenuation coefficient. For scenarios without historical images, the median value of 0.3 within the empirical range is used as the initial default value.

[0029] For example, assuming a preset distance attenuation coefficient of 0.3, there are three reference samples within a search range of 10 pixels, with distances of 2, 5, and 8 pixels from the pixel to be reconstructed, respectively. In this embodiment, the weight allocation of each sample is calculated based on the reciprocal of the distance combined with the attenuation coefficient. The spectral values ​​of the three samples are then synthesized to obtain the pixel values ​​of the reconstructed region.

[0030] In step S13, local spectral smoothing calculation is performed on the reconstructed region, and consistency verification is performed in conjunction with a preset spectral deviation threshold to obtain the de-occluded restored image, including: The spectral values ​​of the reconstructed region and the reference sample in the visible and near-infrared bands are extracted, and the absolute value of the difference between the spectral values ​​is calculated to obtain the spectral deviation. Determine whether the spectral deviation is greater than the preset spectral deviation threshold; If the spectral deviation is greater than the preset spectral deviation threshold, the corresponding pixel is numerically reconstructed using a local spectral smoothing algorithm; if it is not greater, the original pixel values ​​of the reconstructed region are maintained, thereby obtaining the de-occluded restored image.

[0031] In one implementation, the spectral information of each pixel to be verified is extracted from the aforementioned reconstructed region. Specifically, the spectral reflectance values ​​of each pixel within the reconstructed region are extracted in the visible light band (e.g., the green light band with a center wavelength of 550 nm) and the near-infrared band (e.g., the band with a center wavelength of 800 nm). Simultaneously, this embodiment obtains the average spectral reflectance of the reference sample in the same band as described in the preceding steps. This embodiment subtracts the average reflectance of the reference sample in the corresponding band from the reflectance value of the pixel to be verified in each band, and takes the absolute value of the difference. This calculated absolute value of the difference is determined as the spectral deviation.

[0032] It should be noted that the method for determining the preset spectral deviation threshold (e.g., 0.05) in this embodiment adopts a statistical analysis method based on sensor noise calibration. This embodiment acquires calibration whiteboard test data of the multispectral sensor under an unobstructed standard environment, extracting the electronic noise level and spectral response non-uniformity data across the entire band. By statistically analyzing the spectral fluctuations of homogeneous areas in historically acquired unobstructed farmland images, the standard deviation sequence of reflectance variation is calculated. The mean of the standard deviation sequence is multiplied by a scaling factor of three to serve as the critical value for evaluating whether spectral distortion has occurred in the reconstructed pixels, and this value is objectively calibrated as the preset spectral deviation threshold.

[0033] It is worth noting that this embodiment determines whether the spectral deviation of each reconstructed pixel is greater than the preset spectral deviation threshold. If the determination result is greater than the threshold, it is determined that the pixel has introduced spectral deviation during the reconstruction process, and a numerical reconstruction operation is performed. Specifically, the pixel is recalibrated using a local spectral smoothing algorithm. In this embodiment, a 3×3 local pixel window is established with the pixel as the center, and the spectral reflectance values ​​of the other 8 neighboring pixels within the window (excluding the center pixel) are extracted. This embodiment calculates the arithmetic mean of the reflectance of these 8 neighboring pixels and assigns this arithmetic mean to the center pixel as the reconstructed pixel value.

[0034] The local spectral smoothing algorithm establishes a square neighborhood window of odd size centered on the pixel to be reconstructed; calculates the sum of squares of the differences in spectral reflectance between the pixel and each neighboring pixel in all bands; assigns weights according to the magnitude of the absolute value of the difference, with the smaller the difference, the greater the weight; and uses the assigned weights to perform a weighted average of the spectral values ​​of the neighboring pixels, using the result as the reconstructed value of the pixel, thereby smoothing spectral noise while preserving the original spatial texture structure.

[0035] If the determination result is not greater than, then the spectral features of the pixel are determined to meet the spatial continuity requirement. In this embodiment, the original pixel value of the pixel in the reconstructed area is directly retained. This embodiment traverses all reconstructed pixels and performs the above verification and reconstruction process to generate a de-occluded restored image with continuous global spectral features, and stores it in GeoTIFF format to ensure that the image resolution is consistent with the original data.

[0036] For example, assume a preset spectral deviation threshold of 0.05. If the reflectance of a reconstructed pixel in the near-infrared band is 0.42, while the mean of the surrounding reference samples is 0.35, the calculated spectral deviation is 0.07. Since 0.07 is greater than 0.05, this embodiment extracts the reflectance of other pixels within a 3×3 window surrounding this pixel (e.g., with a mean of 0.36), reconstructs the pixel value to 0.36, and eliminates spectral abrupt changes.

[0037] In step S14, spatial smoothing filtering is performed on the de-occluded restored image to obtain an optimized clear image, including: The edge intensity distribution is obtained by performing gradient calculation on the deoccluded and restored image using the Sobel operator. An adaptive weighting function is constructed using the edge intensity distribution, and a Gaussian filter is used to perform edge-preserving smoothing processing on the de-occluded restored image to obtain the optimized clear image.

[0038] In one implementation, a two-dimensional pixel matrix of the de-occluded and restored image is extracted. The Sobel operator is used to perform two-dimensional discrete convolution on the two-dimensional pixel matrix. Specifically, a predefined third-order horizontal edge detection convolution kernel and a third-order vertical edge detection convolution kernel are used to perform sliding convolution calculations along the row and column directions of the two-dimensional pixel matrix, respectively, to obtain the horizontal and vertical gradient components corresponding to each pixel. In this embodiment, the sum of the squares of the horizontal and vertical gradient components is calculated, and the square root of the sum is taken to obtain the gradient magnitude. The gradient magnitudes of all pixels are arranged according to their original spatial coordinates to generate the edge intensity distribution.

[0039] It should be noted that this embodiment utilizes the edge intensity distribution to construct an adaptive weight function. The gradient magnitude of each pixel in the edge intensity distribution is extracted. The square of the gradient magnitude is negative, divided by a preset attenuation control parameter, and the resulting quotient is used as the exponent of the natural logarithm. The calculated output weight coefficient is between zero and one. Regarding the objective method for determining the preset attenuation control parameter, this embodiment calculates the arithmetic mean of the gradient magnitudes across the entire image, multiplies the square of this arithmetic mean by a constant two, and strictly labels the result as the preset attenuation control parameter. This operator logic ensures that the weight coefficient output by edge regions with larger gradient magnitudes approaches zero, and the weight coefficient output by flat regions with smaller gradient magnitudes approaches one, thus constituting the adaptive weight function.

[0040] It is worth noting that this embodiment combines a Gaussian filter to perform edge-preserving smoothing processing on the de-occluded restored image, resulting in the optimized clear image. This embodiment constructs a two-dimensional Gaussian filter. Regarding the standard deviation parameter of the Gaussian filter, this embodiment extracts the background thermal noise variance data of historical similar sensors under dark field conditions, takes the square root of this variance data, and directly sets it as the standard deviation parameter. According to the three-times-standard-deviation principle, three times the standard deviation parameter is calculated, and the nearest odd number is taken upwards, which is then set as the spatial kernel size of the Gaussian filter. Within the determined spatial kernel sliding window, this embodiment calculates the spatial distance weights from the center point to each neighboring point based on a Gaussian function. Subsequently, this embodiment performs a pixel-by-pixel product operation on the spatial distance weights and the weight coefficients corresponding to the aforementioned calculated adaptive weight function to obtain a joint filtering weight matrix. Finally, the joint filtering weight matrix is ​​used to perform a weighted average calculation on the pixel values ​​of the de-occluded restored image within the sliding window, and the output result is updated with the new value of the center pixel. After traversing the entire image to complete the filtering, an optimized clear image with edge-preserving characteristics is output.

[0041] In step S15, a histogram statistical calculation is performed on the optimized clear image to obtain a global brightness deviation. Based on the global brightness deviation, the band gain of the optimized clear image is adjusted to obtain a brightness-adjusted image, including: The mean brightness of the optimized clear image in different image blocks is calculated using a histogram statistical algorithm, and the mean deviation value is calculated based on the mean brightness of each image block to obtain the global brightness deviation. The correction gain for each band is determined based on the global brightness deviation. The correction gain is then used to perform linear mapping calculations on the pixel values ​​of each band in the optimized clear image to obtain the brightness-adjusted image.

[0042] In one implementation, the spatial resolution parameter of the optimized clear image and the average physical scale parameter of the target monitored farmland plot are extracted. The average physical scale parameter is divided by the spatial resolution parameter to obtain the grid pixel ratio. In this embodiment, the grid pixel ratio is used as the side length parameter to divide the optimized clear image into multiple non-overlapping image blocks. For each image block, this embodiment uses a histogram statistical algorithm to calculate the comprehensive brightness probability density distribution across the entire band. That is, each gray level between zero and the upper limit of sensor sensitivity is used as the independent variable, and the frequency of pixel occurrence corresponding to each gray level is counted and divided by the total number of pixels in the block to obtain the probability distribution. Subsequently, the sum of the products of each gray level value and its corresponding probability (i.e., the mathematical expectation) is calculated, and the calculation result is strictly output as the brightness mean of the image block. This embodiment further calculates the arithmetic mean of the brightness mean of all image blocks, and calculates the square of the difference between the brightness mean of each image block and the arithmetic mean; the square root of the arithmetic mean of the squares of the differences corresponding to all image blocks is taken to obtain the standard deviation value. In this embodiment, the standard deviation value is directly used as the mean deviation value and determined as the global brightness deviation.

[0043] It should be noted that, regarding the step of determining the correction gain for each band based on the global brightness deviation, this embodiment extracts a preset target reference brightness value. Regarding the objective method for determining the target reference brightness value, this embodiment retrieves historical reference images of the same type of multispectral sensor acquired under cloud-free conditions, and calculates the arithmetic mean of the global brightness of these reference images in various bands, such as the near-infrared and red bands. These arithmetic mean values ​​are then objectively calibrated as the target reference brightness values ​​corresponding to each band. In this embodiment, the target reference brightness value for the corresponding band is used as the dividend, and the currently calculated average brightness value for the corresponding band within the image block is used as the divisor to perform a division operation, obtaining the initial gain coefficient.

[0044] It is worth noting that, in order to prevent excessive initial gain due to abnormal shadows and thus spectral saturation, this embodiment extracts the global brightness deviation calculated in the previous steps to construct a smooth attenuation mechanism. This embodiment divides the preset control constant by the sum of the preset control constant and the global brightness deviation to calculate an attenuation weight with a value between zero and one. When the global brightness deviation is extremely large, this attenuation weight approaches zero, causing the final correction gain to conservatively converge to one, thereby effectively avoiding excessive amplification of local noise into pseudo-spectral signals. Subsequently, the gain difference between the initial gain coefficient and one is calculated, and this gain difference is multiplied by the attenuation weight. The product result is then added to one to calculate the final correction gain for each band. Regarding the preset control constant, this embodiment extracts the dynamic range value of the corresponding band sensor, multiplies this dynamic range value by a 10% scaling factor, and objectively calibrates the calculation result as the preset control constant. This algorithm ensures that when the global illumination distortion is extremely large, the attenuation weight approaches zero, so that the final correction gain conservatively converges to a value of one, thus strictly achieving the unity of structured fidelity and overexposure prevention.

[0045] In one implementation, the correction gain is used to perform linear mapping calculations on the pixel values ​​of each band in the optimized clear image. Specifically, in this embodiment, image blocks are used as units. The original pixel value of each independent pixel in each band is extracted, and this value is multiplied pixel-by-pixel with the correction gain determined for that block in the corresponding band. After completing the multiplication mapping for all image blocks, this embodiment performs bilinear interpolation smoothing on pixels crossing block boundaries, merging and outputting the brightness-adjusted image with a natural physical transition and consistent band proportions.

[0046] For example, suppose the initial gain coefficient of a certain block is calculated to be 5.0 (i.e., extremely dark areas need to be magnified 5 times), while the global brightness deviation is extremely large. The attenuation weight calculated by the system according to the attenuation mechanism formula is 0.2. The system first calculates the gain difference (5.0 minus 1 equals 4.0), multiplies it by the attenuation weight 0.2 to get 0.8, and finally adds the value 1 to obtain the final corrected gain of 1.8. At the same time, the system extracts the specific pixel with the original pixel value of 90 in the near-infrared band within the image block, multiplies it by the corrected gain of 1.8, and performs a linear mapping calculation to obtain the updated pixel value of 162. Through this objective calculation process, the system successfully limits the originally dangerous 5x magnification to a safe 1.8x, effectively avoiding the risk of amplifying local noise into pseudo-spectral signals, and compiles the processed data to generate a brightness-adjusted image.

[0047] In step S16, the pixels in the brightness adjustment image that exceed a preset brightness threshold are dynamically adjusted to obtain a balanced illumination image, including: Extract the dark pixels with brightness values ​​below a preset low brightness threshold and the bright pixels with brightness values ​​above a preset high brightness threshold from the brightness adjustment image. The illumination-equalized image is obtained by performing numerical remapping calculations on the extracted overly dark pixels and overly bright pixels according to preset contrast limit parameters.

[0048] In one implementation, the brightness adjustment image output from the preceding step is extracted, and the spectral brightness values ​​of all pixels in the image are traversed. In this embodiment, the brightness values ​​of the brightness adjustment image are compared pixel by pixel, and a set of pixels with brightness values ​​lower than a preset low-brightness threshold is extracted and identified as the overly dark pixels; a set of pixels with brightness values ​​higher than a preset high-brightness threshold is extracted and identified as the overly bright pixels.

[0049] It should be noted that, regarding the objective method for determining the preset low-brightness threshold and the preset high-brightness threshold, this embodiment extracts a high-quality benchmark farmland image sample set with no historical cloud cover and no shadow interference, statistically analyzes the probability density distribution of the brightness values ​​of all pixels in the sample set, and calculates and generates a cumulative distribution function curve. This embodiment extracts the brightness value corresponding to the fifth percentile of this cumulative distribution function curve and strictly defines it as the preset low-brightness threshold; simultaneously, it extracts the brightness value corresponding to the ninety-fifth percentile of the curve and strictly defines it as the preset high-brightness threshold. This determination method eliminates the arbitrariness of subjective settings and strictly defines the abnormal brightness range that needs adjustment.

[0050] It is worth noting that this embodiment employs a performance-driven grid search method to determine the preset contrast limiting parameter. This embodiment acquires validation set images containing local illumination extremes, sets candidate intervals and fixed iteration step sizes for the contrast limiting parameter, and processes the validation set images using adaptive histogram equalization core logic. This embodiment calculates the image information entropy index of the output image under each candidate parameter, extracts the parameter value that maximizes the image information entropy index, and objectively determines it as the preset contrast limiting parameter.

[0051] In one implementation, the extracted overly dark pixels and overly bright pixels are numerically remapped according to the preset contrast limit parameter. Specifically, the calculation logic is as follows: In this embodiment, the spatial domain of the brightness adjustment image is divided into multiple continuous and non-overlapping local window matrices, and local brightness frequency histograms are calculated for each local window matrix. Pixels whose frequency in the local brightness frequency histogram exceeds the preset contrast limit parameter are truncated, and the total number of truncated pixels is evenly distributed across all brightness levels of the local brightness frequency histogram to calculate a restricted histogram. Discrete integration is performed on the restricted histogram to obtain a local cumulative distribution function, which is then used as a nonlinear mapping operator to transform and replace the original brightness values ​​of the overly dark pixels and overly bright pixels located within the corresponding local window. After completing the remapping of all local window matrices, this embodiment uses a bilinear interpolation algorithm to perform weighted smoothing calculations on the pixel values ​​at the boundaries of adjacent local windows to eliminate block stitching artifacts and output the resulting illumination-equalized image.

[0052] For example, assume that the preset low-brightness threshold is 50 and the preset high-brightness threshold is 200, determined by the cumulative distribution of historical samples. For a pixel with a brightness value of 40 in the brightness adjustment image, the system identifies and extracts it as an overly dark pixel and locates its 8x8 local window. The system truncates and uniformly redistributes the histogram of this window based on a preset contrast limit parameter (value 2.0) objectively determined by grid search. The system uses a recalculated local cumulative distribution function to perform a non-linear remapping of the overly dark pixel's value of 40, obtaining an updated value of 105. The system performs the above value remapping calculation on all identified overly dark and overly bright pixels, ultimately merging them to generate a uniformly bright, evenly distributed illumination image.

[0053] In step S17, near-infrared and red band data of the illumination-equalized image are extracted to calculate the normalized vegetation index (NDI). Based on the NDI and a preset anomaly threshold, spatial clustering is performed to obtain anomaly distribution images, including: Extract the near-infrared band data and the red band data from the illumination-balanced image, and calculate the difference between the two and the quotient of their sum to obtain the normalized vegetation index. A vegetation status distribution map is obtained by performing gridding processing on the normalized vegetation index using a spatial interpolation algorithm. Extract the spatial coordinates of pixels in the vegetation status distribution map whose values ​​are lower than the preset anomaly threshold; The abnormal distribution image is obtained by using a density-based spatial clustering algorithm to classify and merge the spatial coordinates of pixels in the region into neighborhoods.

[0054] After obtaining the anomalous distribution image, the process also includes: Acquire historical image sequences and extract historical vegetation index data of the abnormal spatial range corresponding to the abnormal distribution images in the historical image sequences. The numerical rate of change is calculated by fitting time series trends. Meteorological environment data and fertilizer application data are acquired, and the normalized vegetation index is compensated and corrected by using a linear regression algorithm in combination with the meteorological environment data and the fertilizer application data to obtain the corrected state image. The anomalous spatial range in the anomalous distribution image is further determined and eliminated using the corrected state image to obtain the final monitoring result.

[0055] In one implementation, near-infrared band data with wavelengths ranging from 760 to 900 nanometers and red band data with wavelengths ranging from 650 to 680 nanometers are extracted from the illumination-equalized image. In this embodiment, the reflectance values ​​corresponding to the near-infrared band data and the red band data are extracted pixel by pixel. The difference between the near-infrared band reflectance value and the red band reflectance value is calculated, and their sum is calculated. The difference is divided by the sum, and a floating-point matrix with values ​​ranging from -1 to +1 is output, which is determined as the normalized vegetation index.

[0056] It should be noted that this embodiment uses a spatial interpolation algorithm to perform gridding processing on the normalized vegetation index. Specifically, ordinary kriging interpolation is used to extract discrete normalized vegetation index sampling points as input, and a semi-variogram model is constructed based on spatial autocorrelation. This embodiment sets a fixed spatial resolution of 5 meters by 5 meters to construct a two-dimensional physical grid. The obtained kriging weight matrix is ​​used to perform weighted summation prediction on the index values ​​of unknown coordinate points within the grid, generating a globally continuous two-dimensional matrix with uniform spatial resolution, which is then output as the vegetation state distribution map.

[0057] It is worth noting that, regarding the objective calibration method for the preset anomaly threshold, this embodiment extracts the nominal normalized vegetation index set for the same period in the past three years without disasters in the target monitoring area, calculates the arithmetic mean and standard deviation of the index set; subtracts twice the standard deviation from the arithmetic mean, and strictly sets the resulting difference value as the preset anomaly threshold. This embodiment traverses each grid cell in the vegetation status distribution map, determining whether the corresponding value is lower than the preset anomaly threshold; if the determination result is lower, the two-dimensional latitude and longitude coordinates or planar projection coordinates corresponding to the grid cell are extracted and compiled to form a set of pixel spatial coordinates of the region.

[0058] In one implementation, a density-based spatial clustering algorithm (DBSCAN) is used to perform neighborhood classification and merging calculations on the spatial coordinates of the region's pixels. Regarding the neighborhood radius parameter (Eps) upon which this algorithm relies, this embodiment calculates the Euclidean distance (K is fixed at 5) from all coordinate points to their Kth nearest neighbor. The distance values ​​are then arranged in descending order to create a K-distance map. The distance value corresponding to the inflection point of maximum curvature of the curve (e.g., 10 meters) is extracted and objectively set as the neighborhood radius parameter. The minimum number of contained points parameter (MinPts) is directly set to 5. This embodiment performs a density reachability traversal on the set of spatial coordinates of the region's pixels using the above parameters, dividing spatially densely connected anomalous coordinate points into the same anomalous cluster and removing isolated noise points that fail to cluster. This embodiment then overlays the spatial boundary information of the retained anomalous clusters onto the original projected coordinate system to generate the anomalous distribution image.

[0059] It should be noted that, regarding the extracted anomalous spatial range, this embodiment acquires the historical image sequence of the past six months with a temporal resolution of thirty days. Within the same physical spatial boundary, monthly normalized vegetation index (NDI) values ​​are extracted to constitute the time series of the historical vegetation index data. This embodiment uses an autoregressive integral moving average (ARIMA) model to perform time series trend fitting on the time series. Regarding the determination of model parameters p, d, and q, this embodiment performs a stationarity test on the sequence to determine the difference order d (usually 1), and calculates the Akaike Information Criterion (AIC) value by traversing possible parameter combinations, extracting the parameter combination that minimizes the AIC value as the optimal model. The historical vegetation index data is input into this optimal ARIMA model, and the first derivative value of its fitted curve is extracted and output as the rate of change of the value.

[0060] It is worth noting that in this embodiment, the average rainfall and average temperature values ​​matching the monitoring period are used to construct the meteorological environment data, and the nitrogen fertilizer application rate values ​​for the plot are retrieved from the business system to construct the fertilizer application data. This embodiment uses a multiple linear regression algorithm for compensation and correction. Regarding the specific training and fitting process of this multiple linear regression algorithm, this embodiment first collects historical monitoring datasets of the target monitoring area over the past three years, which have been manually verified to be free from disease interference; from the historical monitoring dataset, historical normalized vegetation index, historical meteorological environment data of the same period, and historical fertilizer application data are extracted as joint input feature vectors, and the standard vegetation index corresponding to the physical space under ideal growth conditions is calibrated as the regression target variable.

[0061] The standard vegetation index refers to the normalized vegetation index benchmark value obtained by collecting and calculating data through drones or satellite remote sensing under ideal healthy conditions, where the same geographical area, crop type, and growing season are met, and the area is confirmed by on-site manual verification to be free from pests and diseases, with sufficient water and nutrient supply, and unaffected by extreme weather. In this embodiment, multiple images that meet the above conditions are selected from the historical database, the normalized vegetation index of the corresponding plots in each image is extracted, and their arithmetic mean is calculated. This mean is used as the standard vegetation index of the plot.

[0062] This embodiment constructs a multiple linear regression equation and uses ordinary least squares (OLS) to fit the joint input feature vector to the regression target variable. Specifically, it constructs an objective function with the sum of squared residuals between the predicted values ​​and the true regression target variable as the independent variable. By taking the partial derivatives of each regression coefficient and setting them to zero, it solves the normal equation system to calculate the feature weight coefficients that minimize the sum of squared residuals. This embodiment further uses ten-fold cross-validation to evaluate the regression equation. When the coefficient of determination is greater than a preset fitting threshold (e.g., 0.80), the feature weight coefficients are fixed, and the trained regression equation is output. If the coefficient of determination is not greater than the preset fitting threshold, the model parameter self-adjustment logic is triggered, and the model is refitted by increasing the historical sample sampling step size or introducing a nonlinear polynomial term until the coefficient of determination reaches the threshold. In this embodiment, the current meteorological environment data and the fertilizer application data are input into the trained regression equation, and a linear weighted sum is performed with the corresponding feature weight coefficients to calculate a numerical compensation term. This numerical compensation term is added to the currently measured normalized vegetation index to obtain a true health characterization quantity that is not affected by short-term weather or fertilizer application delays, and it is updated as the corrected state image.

[0063] It is worth noting that, regarding the objective method for determining the preset fitting threshold, this embodiment extracts the determination coefficient data of multiple similar regression models that have been successfully deployed in historical business systems and verified to effectively eliminate false alarms caused by environmental interference, forming a validation set; calculates the arithmetic mean and standard deviation of all determination coefficients in the validation set; subtracts one time the standard deviation from the arithmetic mean, and strictly sets the calculated difference value as the preset fitting threshold. This embodiment determines whether the currently calculated determination coefficient is greater than the preset fitting threshold; if so, it determines that the model's goodness of fit meets the standard, fixes the weight coefficients of each feature at this time, and outputs the regression equation after training is completed.

[0064] For example, suppose that the initial normalized vegetation index (NRI) mean of an anomaly cluster after clustering is 0.28, which is lower than the preset anomaly threshold of 0.35. The system extracts data from the same period and finds that the rainfall in this area over the past week was extremely low, and the fertilization amount was only 60% of the normal level. Using a linear regression model, it is calculated that this extreme weather and fertilization deviation will cause a non-pathological drop of -0.08 in the NRI. The system performs compensation correction, adding 0.08 to 0.28, resulting in a corrected NRI value of 0.36, generating a corrected state image. Since the corrected value of 0.36 is no longer lower than the preset anomaly threshold of 0.35, the system determines in the secondary judgment step that this area is not a true example of vegetation disease decline, and thus removes this area from the anomaly spatial range, finally outputting the final monitoring result filtered out for false positives caused by temporary environmental fluctuations.

[0065] In summary, this invention utilizes edge detection and spectral difference analysis to identify occluded areas in images and generate preliminary labeled images. It then combines gray-level co-occurrence matrix to extract texture features and employs Kriging interpolation and spectral consistency verification to perform pixel reconstruction, generating de-occluded restored images. Furthermore, it employs edge intensity-guided adaptive smoothing filtering and band gain mapping to perform global illumination consistency calibration, generating illumination-balanced images. Finally, it calculates the normalized vegetation index to perform spatial clustering and localization, and integrates multi-temporal trend fitting and multi-source heterogeneous data compensation to perform state correction. This invention achieves high-fidelity physical restoration and accurate quantitative assessment of the health status of remote sensing images of mountain farmland under complex weather conditions, effectively solving the technical challenges of distortion and low accuracy in remote sensing monitoring and assessment caused by occlusion interference and uneven local illumination.

[0066] Reference Figure 2 The second embodiment of the present invention provides a farmland remote sensing monitoring system based on image calibration, comprising: The edge segmentation module is used to acquire the original image sequence and perform edge extraction and spectral difference calculation on the original image sequence to obtain a preliminary labeled image. The texture reconstruction module is used to extract the cloud-occluded area and the surrounding non-occluded area in the preliminary labeled image, calculate the texture parameters of the surrounding non-occluded area using the gray-level co-occurrence matrix algorithm and select reference samples, and perform pixel reconstruction on the cloud-occluded area using the kriging interpolation algorithm based on the reference samples to obtain the reconstructed area. The spectral verification module is used to perform local spectral smoothing calculations on the reconstructed region and perform consistency verification in combination with a preset spectral deviation threshold to obtain the de-occluded restored image. The smoothing and noise reduction module is used to perform spatial smoothing filtering calculations on the de-occlusion restored image to obtain an optimized clear image. The band gain module is used to perform histogram statistical calculations on the optimized clear image to obtain the global brightness deviation, and adjust the band gain of the optimized clear image according to the global brightness deviation to obtain a brightness-adjusted image. The pixel equalization module is used to dynamically adjust the pixel values ​​of pixels in the brightness adjustment image that exceed a preset brightness threshold to obtain a light-balanced image. The clustering and localization module is used to extract near-infrared band data and red band data from the illumination-balanced image to calculate the normalized vegetation index, and to perform spatial clustering calculation based on the normalized vegetation index and a preset anomaly threshold to obtain anomaly distribution images.

[0067] It should be noted that the image calibration-based farmland remote sensing monitoring system provided in this embodiment of the invention is used to execute all the process steps of the image calibration-based farmland remote sensing monitoring method in the above embodiment. The working principles and beneficial effects of the two are one-to-one, so they will not be described again.

[0068] It should be noted that the system embodiments described above are merely illustrative. The units described as separate components may or may not be physically separate, and the components shown as units may or may not be physical units; that is, they may be located in one place or distributed across multiple network units. Some or all of the modules can be selected to achieve the purpose of this embodiment according to actual needs. Furthermore, in the accompanying drawings of the system embodiments provided by this invention, the connection relationships between modules indicate that they have communication connections, which can be specifically implemented as one or more communication buses or signal lines. Those skilled in the art can understand and implement this without any creative effort.

[0069] The specific embodiments described above further illustrate the purpose, technical solution, and beneficial effects of the present invention. It should be understood that the above descriptions are merely specific embodiments of the present invention and are not intended to limit the scope of protection of the present invention. In particular, it should be noted that any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of the present invention should be included within the scope of protection of the present invention for those skilled in the art.

Claims

1. An image calibration-based farmland remote sensing monitoring method, characterized in that, include: Obtain the original image sequence, and perform edge extraction and spectral difference calculation on the original image sequence to obtain a preliminary labeled image; Extract the cloud-occluded area and the surrounding unoccluded area from the preliminary labeled image. Calculate the texture parameters of the surrounding unoccluded area using the gray-level co-occurrence matrix algorithm and select reference samples. Based on the reference samples, use the Kriging interpolation algorithm to perform pixel reconstruction on the cloud-occluded area to obtain the reconstructed area. Local spectral smoothing calculations are performed on the reconstructed region, and consistency verification is performed in conjunction with a preset spectral deviation threshold to obtain the de-occluded restored image; Spatial smoothing filtering is performed on the de-occluded and restored image to obtain an optimized clear image; The global brightness deviation is obtained by performing histogram statistical calculation on the optimized clear image, and the brightness-adjusted image is obtained by adjusting the band gain based on the global brightness deviation. Dynamic pixel value adjustment is performed on pixels in the brightness adjustment image that exceed a preset brightness threshold to obtain an illumination equalization image; The near-infrared and red band data of the illumination-balanced image are extracted to calculate the normalized vegetation index. Based on the normalized vegetation index and a preset anomaly threshold, spatial clustering is performed to obtain anomaly distribution images. 2.The image calibration based farmland remote sensing monitoring method according to claim 1, characterized in that, The process of acquiring the original image sequence and performing edge extraction and spectral difference calculation on the original image sequence to obtain a preliminary labeled image includes: The Canny edge detection algorithm is used to extract gradient features from the original image sequence to obtain a binary edge map. Extract the blue light band reflectance and near-infrared band reflectance of the original image sequence, and calculate the difference between the two to obtain the reflectance difference value; Pixels with reflectance differences greater than a preset difference threshold are marked as occluded pixels. The preliminary marked image is obtained by segmenting and classifying the occluded pixels and the edge binary map using a region growing algorithm. 3.The image calibration based farmland remote sensing monitoring method according to claim 1, characterized in that, The process involves calculating the texture parameters of the surrounding unoccluded area using the gray-level co-occurrence matrix algorithm and selecting reference samples. Based on these reference samples, the kriging interpolation algorithm is used to reconstruct the cloud-occluded area pixel by pixel to obtain the reconstructed area, including: The contrast and correlation of a preset pixel window within the surrounding unoccluded area are calculated using the gray-level co-occurrence matrix algorithm to obtain texture feature values; Based on the texture feature values, the Euclidean distance algorithm is used to select non-occluded pixels with textures similar to the pixel to be reconstructed as reference samples; wherein, the texture features of the pixel to be reconstructed are derived based on edge continuity. Based on the spatial coordinates and spectral values ​​of the reference sample, and combined with a preset distance attenuation coefficient, the Kriging interpolation algorithm is used to perform a weighted average calculation on the pixels to be reconstructed within the cloud-covered area to obtain the reconstructed area.

4. The farmland remote sensing monitoring method based on image calibration according to claim 1, characterized in that, The step of performing local spectral smoothing calculations on the reconstructed region and combining this with a preset spectral deviation threshold for consistency verification to obtain the de-occluded restored image includes: The spectral values ​​of the reconstructed region and the reference sample in the visible and near-infrared bands are extracted, and the absolute value of the difference between the spectral values ​​is calculated to obtain the spectral deviation. Determine whether the spectral deviation is greater than the preset spectral deviation threshold; If the spectral deviation is greater than the preset spectral deviation threshold, the corresponding pixel is numerically reconstructed using a local spectral smoothing algorithm; if it is not greater, the original pixel values ​​of the reconstructed region are maintained, thereby obtaining the de-occluded restored image.

5. The farmland remote sensing monitoring method based on image calibration according to claim 1, characterized in that, The step of performing spatial smoothing filtering calculations on the de-occluded restored image to obtain an optimized clear image includes: The edge intensity distribution is obtained by performing gradient calculation on the deoccluded and restored image using the Sobel operator. An adaptive weighting function is constructed using the edge intensity distribution, and a Gaussian filter is used to perform edge-preserving smoothing processing on the de-occluded restored image to obtain the optimized clear image.

6. The farmland remote sensing monitoring method based on image calibration according to claim 1, characterized in that, The step of performing histogram statistical calculations on the optimized clear image to obtain a global brightness deviation, and then adjusting the band gain of the optimized clear image based on the global brightness deviation to obtain a brightness-adjusted image, includes: The mean brightness of the optimized clear image in different image blocks is calculated using a histogram statistical algorithm, and the mean deviation value is calculated based on the mean brightness of each image block to obtain the global brightness deviation. The correction gain for each band is determined based on the global brightness deviation. The correction gain is then used to perform linear mapping calculations on the pixel values ​​of each band in the optimized clear image to obtain the brightness-adjusted image.

7. The farmland remote sensing monitoring method based on image calibration according to claim 1, characterized in that, The step of dynamically adjusting the pixel values ​​of pixels in the brightness-adjusted image that exceed a preset brightness threshold to obtain a balanced illumination image includes: Extract the dark pixels with brightness values ​​below a preset low brightness threshold and the bright pixels with brightness values ​​above a preset high brightness threshold from the brightness adjustment image. The illumination-equalized image is obtained by performing numerical remapping calculations on the extracted overly dark pixels and overly bright pixels according to preset contrast limit parameters.

8. The farmland remote sensing monitoring method based on image calibration according to claim 1, characterized in that, The process of extracting near-infrared and red band data from the illumination-equalized image to calculate the normalized vegetation index (NDI), and then performing spatial clustering calculations based on the NDI and a preset anomaly threshold to obtain anomaly distribution images, includes: Extract the near-infrared band data and the red band data from the illumination-balanced image, and calculate the difference between the two and the quotient of their sum to obtain the normalized vegetation index. A vegetation status distribution map is obtained by performing gridding processing on the normalized vegetation index using a spatial interpolation algorithm. Extract the spatial coordinates of pixels in the vegetation status distribution map whose values ​​are lower than the preset anomaly threshold; The abnormal distribution image is obtained by using a density-based spatial clustering algorithm to classify and merge the spatial coordinates of pixels in the region into neighborhoods.

9. The farmland remote sensing monitoring method based on image calibration according to claim 1, characterized in that, After obtaining the anomalous distribution image, the following is also included: Acquire historical image sequences and extract historical vegetation index data of the abnormal spatial range corresponding to the abnormal distribution images in the historical image sequences. The numerical rate of change is calculated by fitting time series trends. Meteorological environment data and fertilizer application data are acquired, and the normalized vegetation index is compensated and corrected by using a linear regression algorithm in combination with the meteorological environment data and the fertilizer application data to obtain the corrected state image. The anomalous spatial range in the anomalous distribution image is further determined and eliminated using the corrected state image to obtain the final monitoring result.

10. A farmland remote sensing monitoring system based on image calibration, characterized in that, include: The edge segmentation module is used to acquire the original image sequence and perform edge extraction and spectral difference calculation on the original image sequence to obtain a preliminary labeled image. The texture reconstruction module is used to extract the cloud-occluded area and the surrounding non-occluded area in the preliminary labeled image, calculate the texture parameters of the surrounding non-occluded area using the gray-level co-occurrence matrix algorithm and select reference samples, and perform pixel reconstruction on the cloud-occluded area using the kriging interpolation algorithm based on the reference samples to obtain the reconstructed area. The spectral verification module is used to perform local spectral smoothing calculations on the reconstructed region and perform consistency verification in combination with a preset spectral deviation threshold to obtain the de-occluded restored image. The smoothing and noise reduction module is used to perform spatial smoothing filtering calculations on the de-occlusion restored image to obtain an optimized clear image. The band gain module is used to perform histogram statistical calculations on the optimized clear image to obtain the global brightness deviation, and adjust the band gain of the optimized clear image according to the global brightness deviation to obtain a brightness-adjusted image. The pixel equalization module is used to dynamically adjust the pixel values ​​of pixels in the brightness adjustment image that exceed a preset brightness threshold to obtain a light-balanced image. The clustering and localization module is used to extract near-infrared band data and red band data from the illumination-balanced image to calculate the normalized vegetation index, and to perform spatial clustering calculation based on the normalized vegetation index and a preset anomaly threshold to obtain anomaly distribution images.

Images