Agricultural film residue pollution monitoring method based on satellite-ground remote sensing data

Abstract

The invention discloses an agricultural film residual pollution monitoring method based on satellite-ground remote sensing data, and relates to the technical field of remote sensing monitoring and image processing, in particular to an agricultural film residual pollution monitoring method. The objective of the invention is to solve the problem that the existing satellite and unmanned aerial vehicle monitoring system cannot evolve automatically during operation and is always difficult to find an optimal balance point between the coverage range and the monitoring precision, so that the monitoring accuracy of agricultural film residue pollution is low. The method comprises the following steps: 1, collecting ground feature samples to construct a spectrum library, and determining a characteristic wave band for distinguishing agricultural film residues; 2, fusing satellite images and regional statistical data, and carrying out macroscopic inversion to delimit a high-risk key monitoring region; 3, planning a route according to the key area, and controlling the unmanned aerial vehicle to obtain hyperspectral image data; 4, analyzing the hyperspectral data by taking a laboratory spectrum library as a benchmark, and calculating the residual abundance of the agricultural film; and 5, performing scale conversion on the microcosmic data to serve as a true value, and feeding back and correcting satellite macroscopic model parameters so as to realize dynamic optimization.

Description

Technical Field

[0001] This invention relates to the field of remote sensing monitoring and image processing technology, specifically to a method for monitoring agricultural film residue pollution. Background Technology

[0002] my country is a major agricultural producer, and plastic film mulching technology has been widely used in agricultural production due to its significant effects in warming and conserving soil moisture, suppressing weed growth, and increasing crop yields. This technology has greatly improved crop yields and quality, ensuring national food security and the supply of important agricultural products. However, with the continuous increase in the amount of plastic film used and the extension of its service life, the problem of agricultural non-point source pollution caused by plastic film residues has become increasingly prominent. Waste plastic film remaining in the soil, due to its main components such as polyethylene and other high molecular compounds, is difficult to degrade naturally and can remain in the soil for decades or even centuries. These residues not only damage the soil's aggregate structure and block soil pores, affecting soil aeration and permeability, but also hinder the growth and development of crop roots, making it difficult for crops to absorb water and nutrients, ultimately leading to reduced crop yields. Therefore, establishing a scientific, efficient, widely covered, and reliable monitoring system for plastic film residues is the prerequisite and foundation for implementing comprehensive management of plastic film pollution and protecting the agricultural ecological environment.

[0003] Current agricultural film residue monitoring technologies mainly rely on manual ground surveys, satellite remote sensing, and low-altitude drone monitoring. However, in actual agricultural environmental supervision and pollution control, single monitoring methods often fail to meet increasingly complex management needs, exhibiting significant scale gaps and efficiency bottlenecks. While manual ground surveys provide accurate data considered true values, their limitations are most pronounced. This method relies heavily on manpower to go deep into the fields for digging, sampling, washing, and weighing, resulting in extremely low efficiency. For monitoring needs at the county, city, or even provincial levels, manual surveys can only extrapolate the overall situation from a very small number of sample points, leading to serious problems of extrapolating the whole from local sampling points, a lack of spatial continuity in the data, and an inability to conduct high-frequency dynamic monitoring. While existing satellite remote sensing technology possesses the capability for macroscopic, large-scale observation and is suitable for regional surveys, it faces the challenge of insufficient resolution in specific applications of agricultural film monitoring. Agricultural film residues are often fragmented and dispersed, frequently buried in the soil or obscured by crop straw. Existing civilian multispectral satellites have limited spatial resolution, making it difficult to directly identify the fine textures of agricultural film. They are also susceptible to interference from different soil types, farming methods, and surface moisture. This results in satellite inversion often providing only a vague risk level, failing to offer precise location and quantitative residual data, and thus hindering direct guidance for specific recycling and remediation efforts.

[0004] While low-altitude drone monitoring technology can acquire high-resolution images and clearly identify agricultural film residue on the ground, its operation suffers from severe fragmentation and high costs. Drones have short flight times and small coverage areas per flight. Using drones for comprehensive area scanning would incur enormous time, equipment wear and tear, and manpower costs, making it economically unfeasible. Furthermore, current drone monitoring is often isolated; the high-precision data acquired is typically used only to generate local orthophotos and then archived, failing to effectively translate into knowledge or experience to improve satellite monitoring capabilities. This lack of data synergy leads to a waste of data resources. More critically, current agricultural film monitoring lacks a complete collaborative mechanism. Satellite and drone data are often acquired at different times by different departments, lacking unified scheduling and verification standards. When satellites detect suspected contamination areas, rapid ground verification is often lacking; drones collect high-precision data but cannot feed it back to satellite models for correction. This fragmented approach prevents the monitoring system from evolving during operation, making it difficult to find the optimal balance between coverage and monitoring accuracy.

[0005] Therefore, there is an urgent need to research a collaborative monitoring method that can organically combine the area survey capabilities of satellites with the point-based precision survey capabilities of drones. This method should not only be able to automatically select key areas to guide drone operations, but also be able to feed back the high-precision monitoring results from drones to satellite models, forming a dynamic optimization system of monitoring, verification, and correction, thereby achieving efficient, accurate, and dynamic supervision of regional agricultural film residue pollution. Summary of the Invention

[0006] The purpose of this invention is to address the problem that existing satellite and UAV monitoring systems cannot evolve on their own during operation and always struggle to find the optimal balance between coverage and monitoring accuracy, resulting in low accuracy in monitoring agricultural film residue pollution. Therefore, this invention proposes a method for monitoring agricultural film residue pollution based on satellite-ground remote sensing data.

[0007] The specific process of a method for monitoring agricultural film residue pollution based on satellite-ground remote sensing data is as follows:

[0008] Step 1: Manually collect samples of agricultural film residue, soil surface, crop leaves and stems, and weeds from the target monitoring area; obtain the raw spectral data of each sample; construct a ground cover spectral database based on the raw spectral data; and select key feature bands that distinguish agricultural film from background ground cover from all the pure raw spectral data contained in the constructed ground cover spectral database.

[0009] Step 2: Acquire satellite remote sensing images of the target monitoring area;

[0010] Obtain agricultural statistics for the target monitoring area from the agricultural sector;

[0011] Spectral features are extracted from the acquired satellite remote sensing images of the target monitoring area to obtain the spectral feature values ​​of the satellite remote sensing images;

[0012] The background risk coefficient is calculated based on agricultural statistical data of the target monitoring area obtained from the agricultural sector.

[0013] Spectral features based on satellite remote sensing imagery and background risk coefficient Calculate the residual pollution index ;

[0014] Set pollution warning thresholds ;

[0015] Extract all The contiguous pixels are designated as low-risk pollution areas;

[0016] Extract all The contiguous pixels are designated as medium-risk pollution areas;

[0017] Extract all The contiguous pixels are designated as high-risk pollution areas;

[0018] Obtain the geographic coordinate vector file of high-risk pollution areas;

[0019] Step 3: Use the convex decomposition algorithm to divide the high-risk contamination area in the vector file obtained in Step 2 into multiple convex polygon sub-regions, and record the geographic coordinates of the vertices of each convex polygon sub-region;

[0020] Based on residual pollution index Determine the UAV's flight altitude by combining each convex polygonal sub-region. ;

[0021] Based on the drone's flight altitude By combining the pixel size, camera focal length, and field of view of the shortwave infrared hyperspectral camera carried by the UAV, the ground resolution and track spacing of the flight path operation are calculated.

[0022] For each convex polygonal sub-region, a bow-shaped reciprocating scanning path is generated based on the track spacing;

[0023] Control the drone to perform automated spectral acquisition operations according to its flight altitude and scanning path;

[0024] Step 4: Perform radiometric calibration on the raw hyperspectral image acquired by the UAV to obtain the wavelength of the target pixels in the raw hyperspectral image. Actual reflectivity at ;

[0025] Extracting endmember matrices from raw hyperspectral image data collected by drones The basic components, namely endmembers, are used to determine the types of endmember substances using spectral libraries;

[0026] Based on actual reflectance and end-member matrix A mixed pixel decomposition equation was established; the mixed pixel decomposition equation was solved using the fully constrained least squares method to obtain the abundance values ​​of various material endmembers in each pixel;

[0027] Based on the abundance values ​​of endmembers of various substances in each pixel, endmember abundance data for the substance category "agricultural film" was extracted. ;

[0028] Based on the end-member abundance data of the extracted end-member material category "agricultural film" Generate an abundance distribution map of agricultural film;

[0029] Step 5: Convert the agricultural film abundance distribution map generated in Step 4 into UAV pixel ground values ​​with the same resolution as satellites. ;

[0030] The residual pollution index obtained in step two As a predicted value, the true value of the drone pixels As the objective truth value, construct the optimization objective function. ;

[0031] Using the objective function For weight vector Adaptive adjustments are made to obtain the optimal weight vector. ;

[0032] The calculated optimal weight vector substitution formula replace ;

[0033] The calculated optimal weight vector substitution formula replace and ;

[0034] The target area is recalculated to generate a true-value corrected distribution map of agricultural film residue.

[0035] The beneficial effects of this invention are as follows:

[0036] The purpose of this invention is to provide a method for monitoring agricultural film residue pollution based on satellite-ground remote sensing data. This method constructs a satellite-UAV collaborative operation process, using satellites for macroscopic guidance and UAVs for microscopic verification. It also combines laboratory spectral library physical benchmarks and multi-scale data feedback mechanisms to solve the contradiction between efficiency and accuracy in large-scale monitoring and achieve dynamic optimization of the monitoring system.

[0037] 1. This invention achieves an optimal balance between full coverage and local precision. It breaks through the limitations of traditional, singular monitoring methods by utilizing the macroscopic perspective of satellite data to quickly pinpoint key areas, followed by targeted, detailed verification using drones. This tiered monitoring model avoids the high costs of blanket drone operations across the entire area while also solving the problem of insufficient satellite remote sensing resolution, significantly improving the overall cost-effectiveness and execution efficiency of monitoring tasks.

[0038] 2. A quantitative inversion system with physical benchmarks was constructed. By building a localized ground object spectral library in the laboratory and selecting short-wave infrared characteristic bands, this invention provides a solid physical basis for hyperspectral inversion of UAVs. This allows monitoring results to go beyond the color or texture features of images and delve into the chemical bond absorption characteristics of substances, thereby effectively distinguishing agricultural film from similarly colored saline-alkali land, white pollution, and other interfering objects, significantly improving the accuracy of identification and anti-interference capabilities.

[0039] 3. A dynamic optimization mechanism for satellite-UAV collaboration was established. This invention goes beyond simple multi-source data acquisition; more importantly, it establishes a feedback loop between data sources. The high-precision data acquired by the UAV is not only used as monitoring results but also as ground truth labels to back-train and correct the large-scale inversion model of the satellite. As monitoring missions continue, the richer the sample data accumulated by the system, the more accurate the inversion parameters of the satellite model become, thereby realizing the self-learning and dynamic evolution of the monitoring system.

[0040] 4. Enhanced ability to detect concealed pollution sources. Macroscopic initial screening combining multi-source information can fully utilize prior knowledge such as historical mulching years and crop types to identify areas with extremely high potential residual risk, even if the surface mulch is not clearly exposed. Combined with short-wave infrared hyperspectral technology from drones, even severely damaged or partially buried mulch can be effectively detected through spectral characteristics, thus significantly reducing the false negative rate.

[0041] 5. It provides scientific decision support for agricultural environmental management. The monitoring results generated by this invention include both macroscopic distribution trends and microscopic fine distribution in key areas. This provides agricultural management departments with multi-dimensional and high-precision quantitative data support for formulating agricultural film recycling policies, evaluating the effectiveness of remediation, and conducting law enforcement supervision, resulting in significant social and ecological benefits. Attached Figure Description

[0042] Figure 1 This is a flowchart illustrating the implementation of the present invention. Detailed Implementation

[0043] Specific Implementation Method 1: The specific process of this implementation method for monitoring agricultural film residue pollution based on satellite-to-ground remote sensing data is as follows:

[0044] Step 1: Manually collect samples of agricultural film residue, soil surface, crop leaves and stems, and weeds from the target monitoring area; obtain the raw spectral data of each sample; construct a ground cover spectral database based on the raw spectral data; and select the key feature bands most suitable for distinguishing agricultural film from background ground cover from all the pure raw spectral data contained in the constructed ground cover spectral database.

[0045] Step 2: Acquire satellite remote sensing images of the target monitoring area. The satellite remote sensing images include, but are not limited to, Sentinel-2, Landsat series, Gaofen series, or other multispectral and hyperspectral satellite images with shortwave infrared bands.

[0046] Obtain agricultural statistics for the target monitoring area from the agricultural sector;

[0047] Spectral features are extracted from the acquired satellite remote sensing images of the target monitoring area to obtain the spectral feature values ​​of the satellite remote sensing images;

[0048] The background risk coefficient is calculated based on agricultural statistical data of the target monitoring area obtained from the agricultural sector.

[0049] Spectral features based on satellite remote sensing imagery and background risk coefficient Calculate the residual pollution index ;

[0050] Set pollution warning thresholds ;

[0051] Extract all The contiguous pixels are designated as low-risk suspected contamination areas;

[0052] Extract all The contiguous pixels are designated as medium-risk suspected contamination areas;

[0053] Extract all The contiguous pixels are designated as high-risk suspected contamination areas;

[0054] Obtain precise geographic coordinates (latitude and longitude) vector files of high-risk suspected contamination areas (satellite images have built-in geographic coordinates) as target objects for subsequent UAV collaborative operations, thereby achieving precise fixed-point monitoring;

[0055] Step 3: Use the convex decomposition algorithm to divide the high-risk suspected contamination area in the vector file obtained in Step 2 into multiple convex polygon sub-regions, and record the geographical coordinates of the vertices of each convex polygon sub-region;

[0056] Based on residual pollution index Determine the UAV's flight altitude by combining each convex polygonal sub-region. ;

[0057] Based on the drone's flight altitude By combining the pixel size, camera focal length, and field of view of the shortwave infrared hyperspectral camera carried by the UAV, the ground resolution (GSD) and track spacing of the flight path operation are calculated.

[0058] For each convex polygonal sub-region, a bow-shaped reciprocating scanning path is generated based on the track spacing;

[0059] Control the drone to perform automated spectral acquisition operations according to the path;

[0060] Step 4: Perform radiometric calibration on the raw hyperspectral image acquired by the UAV to obtain the wavelength of the target pixels in the raw hyperspectral image. Actual reflectivity at ;

[0061] To construct accurate unmixing equations, the endmember matrix is ​​first extracted from the raw hyperspectral image data collected by the UAV. The basic components, namely endmembers, are used to determine the types of endmember substances using spectral libraries;

[0062] Based on actual reflectance and end-member matrix A mixed pixel decomposition equation was established; the fully constrained least squares (FCLS) method was used to solve the mixed pixel decomposition equation to obtain the abundance values ​​of various endmembers in each pixel.

[0063] Based on the abundance values ​​of endmembers of various substances in each pixel, endmember abundance data for the substance category "agricultural film" was extracted. Ignoring soil and vegetation abundance;

[0064] Based on the end-member abundance data of the extracted end-member material category "agricultural film" Generate sub-pixel level agricultural film abundance distribution maps;

[0065] Step 5: Convert the agricultural film abundance distribution map generated in Step 4 into UAV pixel ground values ​​with the same resolution as satellites. ;

[0066] The residual pollution index obtained in step two As a predicted value, the true value of the drone pixels As the objective truth value, construct the optimization objective function. ;

[0067] Using the objective function For weight vector Adaptive adjustments are made to obtain the optimal weight vector. ;

[0068] The calculated optimal weight vector substitution formula replace ;

[0069] The calculated optimal weight vector substitution formula replace and ;

[0070] The target area is recalculated to generate a high-precision map of agricultural film residue distribution that has been corrected for true values. This map combines the macroscopic view of the satellite with the microscopic precision of the UAV, providing the final monitoring results.

[0071] Specific Implementation Method Two: This implementation method differs from Specific Implementation Method One in that, in step one, samples of agricultural film residue, soil surface, crop leaves and stems, and weeds are manually collected from the target monitoring area; the original spectral data of each sample are obtained; a ground cover spectral database is constructed based on the original spectral data; and the key feature bands most suitable for distinguishing agricultural film from background ground cover are selected from all the pure original spectral data contained in the constructed ground cover spectral database.

[0072] The specific process is as follows:

[0073] Step 11: Manually collect samples of agricultural film residue, soil surface, crop leaves and stems, and weeds from the target monitoring area; the specific process is as follows:

[0074] Manually collect agricultural film residue samples of different materials (polyethylene PE, polypropylene PP, etc.) and different aging levels (new film, 3-month-old film, 6-month-old film and completely broken agricultural film) in the target monitoring area;

[0075] Soil surface samples of different textures (sand, loam, clay) and different moisture contents (air-dried, field holding capacity) were manually collected from the target monitoring area.

[0076] Leaf and stem samples of major crops (such as cotton, corn, and wheat) at different growth stages (seedling stage, maturity stage, and withering stage) in the target monitoring area were collected manually.

[0077] Samples of common field weeds at different growth stages were manually collected from the target monitoring area;

[0078] Steps 1 and 2: Obtain the original spectral data of each sample based on the samples obtained in Step 1; the specific process is as follows:

[0079] Step 121: In a darkroom environment in the laboratory, under a halogen light source, a hyperspectral imager covering the short-wave infrared band (900nm-2500nm) is used to perform spectral measurements on agricultural film residue samples of different materials and aging degrees collected manually in the target monitoring area, so as to obtain the original spectral data of agricultural film residue samples of different materials and aging degrees in the target monitoring area.

[0080] Step 122: In a darkroom environment in the laboratory, under a halogen light source, a hyperspectral imager covering the short-wave infrared band (900nm-2500nm) is used to perform spectral measurements on soil surface samples with different textures and water contents collected manually in the target monitoring area, so as to obtain the original spectral data of soil surface samples with different textures and water contents in the target monitoring area.

[0081] Steps 1, 2, and 3: In a darkroom environment in the laboratory, under a halogen light source, a hyperspectral imager covering the short-wave infrared band (900nm-2500nm) is used to perform spectral measurements on leaf and stem samples of major crops at different growth stages in the manually collected target monitoring area, so as to obtain the original spectral data of leaf and stem samples of major crops at different growth stages in the target monitoring area.

[0082] Steps 1, 2, and 4: In a darkroom environment in the laboratory, under a halogen light source, a hyperspectral imager covering the short-wave infrared band (900nm-2500nm) is used to perform spectral measurements on common field weed samples at different growth stages in the manually collected target monitoring area to obtain the original spectral data of common field weed samples at different growth stages in the target monitoring area.

[0083] Step 13: Construct a ground cover spectral database based on the original spectral data from Step 1 and Step 2; the specific process is as follows:

[0084] Step 131: The raw spectral data of agricultural film residue samples of different materials and aging degrees in the target monitoring area are sequentially subjected to noise reduction, smoothing and resampling preprocessing (outlier removal), and the preprocessed pure raw spectral data are entered into the ground cover spectral database.

[0085] Resampling addresses the issue that the wavelength center and sampling interval of the acquired spectral data points may deviate slightly, for example, by 1.1 nm in some areas and 0.9 nm in others. Furthermore, the sampling points may differ from those obtained from data acquired by different devices, or sometimes the data may contain dead pixels. To facilitate subsequent data processing, the data can be reprocessed, similar to mathematical interpolation, to ensure that the sampling points are consistent.

[0086] Step 1, 3, 2: Perform noise reduction, smoothing, and resampling preprocessing (outlier removal) on the raw spectral data of soil surface samples with different textures and moisture contents in the target monitoring area in sequence, and enter the clean raw spectral data after preprocessing into the ground cover spectral database.

[0087] Step 133: Perform noise reduction, smoothing and resampling preprocessing (outlier removal) on the raw spectral data of leaf and stem samples of major crops at different growth stages in the target monitoring area, and enter the clean raw spectral data after preprocessing into the ground cover spectral database.

[0088] Steps 1, 3, and 4: The raw spectral data of common field weeds at different growth stages in the target monitoring area are sequentially subjected to noise reduction, smoothing, and resampling preprocessing (outlier removal). The preprocessed clean raw spectral data are then entered into the ground cover spectral database.

[0089] Step 135: Based on Step 131, Step 132, Step 133, and Step 134, complete the construction of a high-precision ground object spectral database;

[0090] Step 14: Select the key feature bands most suitable for distinguishing agricultural film from background land cover from all the pure raw spectral data contained in the land cover spectral database constructed in Step 13.

[0091] The other steps and parameters are the same as in Specific Implementation Method 1.

[0092] Specific Implementation Method Three: This implementation method differs from Specific Implementation Method One or Two in that, in step one-four, the key feature bands most suitable for distinguishing agricultural film from background land cover are selected from all the pure raw spectral data contained in the land cover spectral database constructed in step one-three; the specific process is as follows:

[0093] Step 141: In order to eliminate the drift of the background baseline and highlight the absorption characteristics of the spectrum, the reflectance curves of all pure original spectral data contained in the ground cover spectral database constructed in Step 13 are subjected to envelope removal processing to obtain all pure original spectral data after envelope removal processing.

[0094] The formula for calculating envelope removal is as follows:

[0095]

[0096] In the formula, Wavelength representing pure, raw spectral data The envelope removal value at the location;

[0097] Wavelengths of pure, raw spectral data Reflectance at that location;

[0098] Wavelengths of pure, raw spectral data The value on the convex circumscribed envelope of the spectral curve;

[0099] Step 142: Based on all the pure original spectral data obtained in Step 141 after envelope removal, extract the characteristic absorption peaks from all the pure original spectral curves after envelope removal to explore the characteristic absorption patterns of agricultural film materials. The specific process is as follows:

[0100] The first-order differential spectrum of the original spectral data of agricultural film residue samples was calculated using the central difference method. The formula is as follows:

[0101]

[0102] in, Wavelength of the original spectral data for agricultural film residue samples The envelope removal value at the location; Wavelength of the original spectral data for agricultural film residue samples The envelope removal value at the location; wavelength place, wavelength Place;

[0103] Identify the zero point position in the first-order differential spectrum; find the inflection point where the slope of the original spectrum curve of the agricultural film residue sample changes from negative to positive, corresponding to the zero point of the first-order differential. The inflection point is the minimum point, which is the potential absorption peak center.

[0104] To eliminate interference from minute fluctuations caused by spectral noise, a pre-defined significance discrimination condition is introduced to screen potential absorption peak center wavelengths.

[0105] The significance criteria include:

[0106] If the wavelength of the potential absorption peak center is lower than the preset intensity threshold (indicating that the absorption is deep enough), the potential absorption peak center is retained.

[0107] Analyze the rate of change of the derivative in the neighborhood of the potential absorption peak center (1-3 sampling intervals on both sides of the potential absorption peak center). If the rate of change of the derivative on both sides of zero is greater than the preset rate of change threshold (indicating that the trough is sharp enough and not smooth noise), then the potential absorption peak center is retained.

[0108] By performing the above double verification, false feature zeros that do not meet the significance condition are eliminated, thereby accurately determining the characteristic absorption peak position of the agricultural film.

[0109] The center wavelength corresponding to the center position of the potential absorption peak retained after screening is used as the candidate position of the key characteristic band for the spectral separability assessment in step one, four, and three.

[0110] Step 143: Calculate the spectral separability index (SDI); obtain key characteristic bands based on the SDI; the specific process is as follows:

[0111] Step 1431: To quantitatively assess the ability of different wavelength bands to identify agricultural film, the spectral separability index (SDI) is introduced as an evaluation indicator. The SDI measures the statistical difference between the agricultural film spectrum and the background spectra of soil, vegetation, etc., at a specific wavelength. The SDI is calculated using the following formula:

[0112]

[0113] In the formula,

[0114] Indicates wavelength The separability index at the location;

[0115] This indicates the use of the wavelength acquired in step one. The average value of the envelope removal value of the agricultural film sample set (the data obtained from point 141 are divided into two types according to the material type: agricultural film and background, and the average value is calculated for each).

[0116] This indicates the use of the wavelength acquired in step one. The average value of the envelope removal value of the background sample set (the data obtained from point 141 are divided into two types according to the material type: agricultural film and background, and the average value of each is calculated separately);

[0117] Indicates wavelength The standard deviation of the envelope removal values ​​in the agricultural film sample set at that location;

[0118] Indicates wavelength The standard deviation of the envelope removal values ​​of the background sample set at that location;

[0119] wavelength ;

[0120] Step 1432: Traverse the entire spectral range and use the SDI formula to calculate the spectral separability of the agricultural film residue sample set and the background ground cover sample set (including soil, crops, weeds, etc.) in each spectral band.

[0121] From the candidate locations of the key feature bands determined in steps one, four, and two, select band intervals with SDI values ​​greater than a preset threshold; if a candidate location has obvious absorption by agricultural film but a low SDI value (indicating that background ground features also have similar absorption in this band), it is removed.

[0122] Finally, P band intervals with high discriminative power were selected as the key feature band set.

[0123] This band set not only represents the specific spectral characteristics of agricultural film materials, but also minimizes interference from background objects, serving as a direct basis for subsequent UAV sensor channel configuration and data acquisition.

[0124] Other steps and parameters are the same as in specific implementation method one or two.

[0125] Specific Implementation Method Four: This implementation method differs from one of the specific implementation methods one to three in that, in step two, satellite remote sensing images of the target monitoring area are acquired. The satellite remote sensing images include, but are not limited to, Sentinel-2, Landsat series, Gaofen series, or other multispectral and multi-spectral satellite images with shortwave infrared bands.

[0126] Obtain agricultural statistics for the target monitoring area from the agricultural sector;

[0127] Spectral features are extracted from the acquired satellite remote sensing images of the target monitoring area to obtain the spectral feature values ​​of the satellite remote sensing images;

[0128] The background risk coefficient is calculated based on agricultural statistical data of the target monitoring area obtained from the agricultural sector.

[0129] Spectral features based on satellite remote sensing imagery and background risk coefficient Calculate the residual pollution index ;

[0130] Set pollution warning thresholds ;

[0131] Extract all The contiguous pixels are designated as low-risk suspected contamination areas;

[0132] Extract all The contiguous pixels are designated as medium-risk suspected contamination areas;

[0133] Extract all The contiguous pixels are designated as high-risk suspected contamination areas;

[0134] Obtain precise geographic coordinates (latitude and longitude) vector files of high-risk suspected contamination areas (satellite images have built-in geographic coordinates) as target objects for subsequent UAV collaborative operations, thereby achieving precise fixed-point monitoring;

[0135] The specific process is as follows:

[0136] Step 21: Acquire satellite remote sensing images of the target monitoring area;

[0137] Obtain agricultural statistics for the target monitoring area from the agricultural sector;

[0138] Agricultural statistics for the target monitoring area include crop types (cotton, corn, etc.), the historical years of film covering in each plot (e.g., 5 years, 10 years), the amount of agricultural film used in each plot in the current year, and the method of agricultural film recycling (manual collection or mechanical recycling).

[0139] Multi-temporal satellite remote sensing images of the target monitoring area were acquired through a remote sensing data service platform. Image data sources were selected from multispectral satellites with wide swath coverage, such as Sentinel-2 or Landsat-8 / 9, with a focus on images taken before spring planting or after autumn harvest when the ground was bare to minimize interference from vegetation cover on the surface spectrum. Simultaneously, agricultural statistics for the area were obtained from the agricultural sector, including crop planting structure (cotton, corn, etc.), historical mulching duration (e.g., 5 years, 10 years), current year's mulch usage, and recycling methods (manual collection or mechanical recycling).

[0140] Based on this, a risk feature set of agricultural film residue is constructed, and macroscopic inversion and residual pollution index (RPI) calculation are performed. This is a process of spatially fusing physical spectral characteristics with statistical prior knowledge;

[0141] Step 22: Extract spectral features from the satellite remote sensing images of the target monitoring area obtained in Step 21 to obtain the spectral feature values ​​of the satellite remote sensing images;

[0142] The specific process is as follows:

[0143] Step 221: Perform radiometric calibration and atmospheric correction preprocessing on the satellite remote sensing images of the target monitoring area obtained in Step 21 to obtain the preprocessed satellite remote sensing images of the target monitoring area;

[0144] Step Two Two Two:

[0145] From the bands of the preprocessed satellite remote sensing image of the target monitoring area, select one band from the key feature bands obtained in step one as the sensitive band, with a reflectance of [missing information]. And a band outside the key feature band set interval is used as a reference band, with a reflectance of Based on reflectivity and reflectivity Constructing spectral eigenvalues The calculation formula is as follows:

[0146]

[0147] This index can preliminarily extract the bright areas of suspected plastic residue in ground cover by calculating the difference in reflectance under different absorption intensities.

[0148] Step Two and Three: Calculate the background risk coefficient based on the agricultural statistical data of the target monitoring area obtained from the agricultural sector in Step Two and One; the specific process is as follows:

[0149] Agricultural statistics for the target monitoring area include crop types (cotton, corn, etc.), the historical years of film covering in each plot (e.g., 5 years, 10 years), the amount of agricultural film used in each plot in the current year, and the method of agricultural film recycling (manual collection or mechanical recycling).

[0150] Background risk coefficient as follows:

[0151]

[0152] In the formula, For the first The weights of statistical variables (such as years of mulching, crop type coefficient);

[0153] The statistical variables are four: crop type (cotton, corn, etc.), the historical years of mulching in each plot (e.g., 5 years, 10 years), the amount of agricultural film used in each plot in the current year, and the method of agricultural film recycling (manual collection or mechanical recycling).

[0154] For the first The normalized values ​​of statistical variables (such as years of mulching, crop type coefficient);

[0155] For example, for cotton fields that have been continuously covered with plastic film for more than 10 years, their The value will be assigned a higher weight;

[0156] Step 24: Spectral feature values ​​of the satellite remote sensing image obtained in Step 22 and the background risk coefficient obtained in steps two and three Calculate the residual pollution index ;

[0157] To eliminate the vulnerability of single-satellite spectral inversion to interference from various factors, this method weights and fuses spectral features with background risk.

[0158] The specific process is as follows:

[0159]

[0160] In the formula, For spectral feature weights, and The background risk weights are adaptively adjusted based on the actual conditions of different regions.

[0161] Step 25: Set pollution warning thresholds ;

[0162] Extract all The contiguous pixels are designated as low-risk suspected contamination areas;

[0163] Extract all The contiguous pixels are designated as medium-risk suspected contamination areas;

[0164] Extract all The contiguous pixels are designated as high-risk suspected contamination areas;

[0165] Obtain precise geographic coordinates (latitude and longitude) vector files of high-risk suspected contamination areas to serve as target objects for subsequent UAV collaborative operations, thereby achieving precise point-to-point monitoring.

[0166] The other steps and parameters are the same as those in one of the specific implementation methods one to three.

[0167] Specific Implementation Method 5: This implementation method differs from one of the specific implementation methods 1 to 4 in that, in step 3, a convex decomposition algorithm is used to divide the high-risk suspected contamination area in the vector file obtained in step 2 into multiple convex polygon sub-regions, and the geographical coordinates of the vertices of each convex polygon sub-region are recorded.

[0168] Based on residual pollution index The drone's flight altitude is determined by each convex polygonal sub-region. ;

[0169] Based on the drone's flight altitude By combining the pixel size, camera focal length, and field of view of the shortwave infrared hyperspectral camera carried by the UAV, the ground resolution (GSD) and track spacing of the flight path operation are calculated.

[0170] For each convex polygonal sub-region, a bow-shaped reciprocating scanning path is generated based on the track spacing;

[0171] Control the drone to perform automated spectral acquisition operations according to parameters such as path.

[0172] The specific process is as follows:

[0173] Step 31: Use the convex decomposition algorithm to divide the high-risk suspected contamination areas in the vector file obtained in Step 2 into multiple convex polygon sub-regions, and record the geographic coordinates of the vertices of each convex polygon sub-region; the specific process is as follows:

[0174] Since suspected contaminated areas identified by satellite (RPI > T) often present as irregular polygons, which are inconvenient for UAV flight path planning, convex decomposition algorithms, such as vertex-based greedy convex decomposition algorithms, are used to process them.

[0175] The high-risk suspected contamination area obtained in step two is an irregular polygon. Identify all concave vertices in the irregular polygon and emit rays along the angle bisector or diagonal at the concave vertices to cut the irregular polygon into multiple sub-regions.

[0176] Determine whether all sub-regions are convex polygonal sub-regions;

[0177] If so, record the geographic coordinates of the vertices of each convex polygon sub-region;

[0178] If not, for the non-convex polygon sub-region, emit rays along the angle bisector or diagonal direction at the concave vertex to cut the non-convex polygon sub-region into sub-regions until all sub-regions are convex polygon sub-regions, and record the geographic coordinates of the vertices of each convex polygon sub-region.

[0179] The concave vertex is a vertex with an interior angle greater than 180 degrees;

[0180] Step 3.2: Residual pollution index obtained in Step 2 Determine the UAV's flight altitude by combining each convex polygonal sub-region. ;

[0181] The specific process is as follows:

[0182] Based on the residual pollution index obtained in step two Determine the contaminated area for each convex polygon sub-region;

[0183] If the residual pollution index of the convex polygon sub-region The corresponding convex polygon sub-region is a low-risk suspected contamination area;

[0184] If the residual pollution index of the convex polygon sub-region Then the corresponding convex polygon sub-region is a medium-risk suspected contamination area;

[0185] If the residual pollution index of the convex polygon sub-region Then the corresponding convex polygon sub-region is a high-risk suspected contamination area;

[0186] For high-risk suspected contamination areas, where the residual film is severely fragmented and densely distributed, extremely high-resolution imagery is required; the drone's flight altitude should be adjusted accordingly. The range is set as follows: (Low-altitude range with civil aviation height restrictions);

[0187] For areas suspected of being medium-risk contamination, the residual film distribution is determined to be relatively uniform or the patches are relatively large; to balance efficiency, the drone's flight altitude is set at a certain level. The range is set as follows: (Higher altitude ranges limited by civilian use) to increase the ground coverage width of a single image and improve operational efficiency;

[0188] For low-risk suspected contamination areas, the drone flight altitude will be adjusted. The range is set as follows: ;

[0189] Step 33: Based on the drone flight altitude determined in Step 32 By combining the pixel size, camera focal length, and field of view of the shortwave infrared hyperspectral camera carried by the UAV, the ground resolution (GSD) and track spacing of the flight path operation are calculated. For each convex polygon sub-region, a bow-shaped reciprocating scanning path is generated based on the track spacing.

[0190] Steps 3 and 4: Control the drone to perform automated spectral acquisition operations according to the flight altitude and scanning path.

[0191] The other steps and parameters are the same as those in one of the specific implementation methods one to four.

[0192] Specific Implementation Method Six: This implementation method differs from Specific Implementation Methods One to Five in that the UAV flight altitude determined in step Three-Three is based on the flight altitude of the UAV determined in step Three-Two. By combining the pixel size, focal length, and field of view of the shortwave infrared hyperspectral camera carried by the UAV, the ground resolution (GSD) and track spacing of the flight path are calculated, and the UAV is controlled to perform automated spectral acquisition operations according to key parameters; the specific process is as follows:

[0193] Step 331: Calculate the ground resolution ; indicates as:

[0194]

[0195] in, For pixel size, The focal length of the camera;

[0196] Step 332: Calculate the ground cover width of a single image ; indicates as:

[0197]

[0198] in, This refers to the camera's field of view.

[0199] Step 333: Based on the ground coverage width of a single image Calculate track spacing ; indicates as:

[0200]

[0201] in, To meet the preset overlap requirements (e.g., 60%–80%) of the puzzle.

[0202] The other steps and parameters are the same as those in one of the specific implementation methods one to five.

[0203] Specific Implementation Method Seven: This implementation method differs from Specific Implementation Methods One through Six in that, in steps three and four, for each convex polygonal sub-region, a bow-shaped reciprocating scanning path is generated based on the track spacing; the specific process is as follows:

[0204] The longest side direction of each convex polygon sub-region is selected as the reference direction of the main flight path. Translation sampling is performed along the perpendicular direction of the reference direction at intervals D to generate a series of parallel flight paths.

[0205] Buffer waypoints (i.e., auxiliary turning points located outside the boundary of the work area) are automatically added at the end of the extension line of each parallel flight line scanning segment to compensate for the attitude oscillation of the UAV during the turning process and ensure that the camera is in a stable orthogonal state when entering the effective work area.

[0206] Combining the geographic coordinates of the vertices of each convex polygon sub-region determined in step 31, the flight altitude H determined in step 32, the track spacing D calculated in step 33, and the preset heading overlap rate, the planar coordinates of the endpoints of each parallel flight segment and the buffer waypoint (the geographic coordinates in the vector file of the precise geographic coordinates (latitude and longitude) of the high-risk suspected pollution area obtained in step 2) are spatially mapped and superimposed with the flight altitude H to obtain the three-dimensional coordinate information containing the endpoints of each parallel flight segment and the buffer waypoint.

[0207] The flight speed and camera trigger command are determined based on the coverage rate of the ground coverage area of ​​a single image by the camera sampling frequency and the preset forward overlap rate.

[0208] The UAV is controlled to perform automated spectral acquisition operations based on the three-dimensional coordinate information of the endpoints and buffer waypoints of each parallel flight segment, flight speed, and camera trigger commands.

[0209] The other steps and parameters are the same as those in one of the specific implementation methods one to six.

[0210] Specific Implementation Method Eight: This implementation method differs from Specific Implementation Methods One to Seven in that, in step four, radiometric calibration is performed on the original hyperspectral image acquired by the UAV to obtain the wavelength of the target pixels in the original hyperspectral image. Actual reflectivity at ;

[0211] To construct accurate unmixing equations, the endmember matrix is ​​first extracted from the raw hyperspectral image data collected by the UAV. The basic components, namely endmembers, are used to determine the types of endmember substances using spectral libraries;

[0212] Based on actual reflectance and end-member matrix A mixed pixel decomposition equation was established; the fully constrained least squares (FCLS) method was used to solve the mixed pixel decomposition equation to obtain the abundance values ​​of various endmembers in each pixel.

[0213] Based on the abundance values ​​of various substances in each pixel, endmember abundance data for the endmember substance category "agricultural film" was extracted. Ignoring soil and vegetation abundance;

[0214] Based on the end-member abundance data of the extracted end-member material category "agricultural film" Generate sub-pixel level agricultural film abundance distribution maps;

[0215] The specific process is as follows:

[0216] Step 41: Perform radiometric calibration on the raw hyperspectral image acquired by the UAV to obtain the wavelength of the target pixels in the raw hyperspectral image. Actual reflectivity at ;

[0217] The purpose of radiometric calibration of raw hyperspectral images acquired by UAVs is to convert the dimensionless digital quantization (DN) values ​​recorded by the sensors into physically meaningful surface reflectance data, which is the basis for all subsequent spectral analyses.

[0218] The calculation formula is:

[0219]

[0220] In the formula,

[0221] For target pixels in the original hyperspectral image at wavelength The actual reflectivity at that location The observations are directly used as input to the linear spectral mixture model (LSMM).

[0222] For target pixels in the original hyperspectral image at wavelength The original observations at the location;

[0223] These are the observations from a standard whiteboard;

[0224] This represents the dark current noise value.

[0225] Laboratory calibration of reflectivity for a standard whiteboard;

[0226] Step 42: The radiometrically calibrated images undergo preprocessing such as atmospheric correction and geometric correction to obtain preprocessed full-band images;

[0227] Subsequently, based on the key feature band set determined in step one, the corresponding N feature channel data are extracted from the preprocessed full-band image by wavelength position alignment to obtain the dimension-reduced image data.

[0228] This step involves removing redundant bands from the full-band image that are irrelevant to the target features, thereby reducing computational load while enhancing the agricultural film feature signal, thus obtaining dimension-reduced image data.

[0229] Step 43: To construct an accurate unmixing equation, the endmember matrix is ​​first extracted from the raw hyperspectral image data collected by the UAV. The basic components, namely endmembers, are used to determine the types of endmember substances using spectral libraries;

[0230] Step 44: Based on the actual reflectivity in Step 41 and the endmember matrix in step four three A mixed pixel decomposition equation was established; the fully constrained least squares (FCLS) method was used to solve the mixed pixel decomposition equation to obtain the abundance values ​​of various endmembers in each pixel.

[0231] The mixed pixel decomposition equation is expressed as:

[0232]

[0233] The mixed pixel decomposition equation satisfies: (The sum of the proportions of each component is 1). (Non-negative percentage);

[0234] In the formula,

[0235] Endmember matrix The Middle Abundance of endmembers in pixels (i.e., area ratio); Endmember matrix Total number of mid-range elements; Endmember matrix The Middle Seed endogenous element; For residuals;

[0236] Steps 4 and 5: Based on the abundance values ​​of various substances in each pixel obtained in Step 44, extract the endmember abundance data of the endmember substance category "agricultural film" determined in Step 43. Ignoring soil and vegetation abundance;

[0237] Based on the extracted endmember abundance data of "agricultural film" Generate sub-pixel level agricultural film abundance distribution maps;

[0238] Abundance data of each endmember in the agricultural film abundance distribution map It is used directly as the basis for the ground truth calculation in step five.

[0239] The other steps and parameters are the same as those in any of the specific implementation methods one to seven.

[0240] Specific Implementation Method Nine: This implementation method differs from Specific Implementation Methods One through Eight in that, in step four-three, in order to construct an accurate unmixing equation, the endmember matrix is ​​first extracted from the original hyperspectral image data collected by the UAV. The basic components, namely endmembers, are determined using a spectral library to classify endmember substances; the specific process is as follows:

[0241] Step 431: Extract the basic components constituting the mixed pixels, i.e., candidate endmembers, from the raw hyperspectral image data collected by the UAV; the specific process is as follows:

[0242] Using endmember extraction algorithms, such as the Alternating Volume Maximization (AVmax) algorithm, we find the minimum volume simplex that encloses all pixels in the dimensionality-reduced image data obtained in step four-two, and extract all the minimum volume simplexes. There are vertices, each vertex being a candidate endmember vector. These vectors represent several major pure substances present in the image;

[0243] Step 432: Calculate the vector of each candidate endmember extracted in Step 431. Reference spectra in the ground cover spectral database constructed in step one spectral angles between ;

[0244] Take the spectral angle Minimum reference spectrum As the matching object, the candidate endmember vector The attribute is marked as The corresponding material category;

[0245] spectral angle Represented as:

[0246]

[0247] In the formula, It is the spectral angle;

[0248] For target pixels in the original hyperspectral image at wavelength Each endmember vector at;

[0249] Wavelengths in the ground cover spectral database constructed in step one Reference spectrum at that location;

[0250] Step 433: Obtain the material category corresponding to all candidate endmember vectors;

[0251] Construct an end-member matrix containing different components of agricultural film, soil, and vegetation. , used for subsequent pixel decomposition;

[0252] As with step one, the same substance can have multiple spectra, so it's necessary to identify the substance category of all candidate endmembers. Step four then requires integrating the abundance data of all candidate endmembers belonging to the agricultural film endmember.

[0253] The other steps and parameters are the same as those in one of the specific implementation methods one to eight.

[0254] Specific Implementation Method Ten: This implementation method differs from Specific Implementation Methods One to Nine in that, in step five, the agricultural film abundance distribution map generated in step four is converted into UAV pixel ground values ​​consistent with satellite resolution. ;

[0255] The residual pollution index obtained in step two As a predicted value, the true value of the drone pixels As the objective truth value, construct the optimization objective function. ;

[0256] Using the objective function For weight vector Adaptive adjustments are made to obtain the optimal weight vector. ;

[0257] The calculated optimal weight vector substitution formula replace ;

[0258] The calculated optimal weight vector substitution formula replace and ;

[0259] The target area is recalculated to generate a high-precision distribution map of agricultural film residue after true value correction; this map combines the macroscopic view of the satellite with the microscopic precision of the UAV, which is the final monitoring result;

[0260] The specific process is as follows:

[0261] Step 51: Due to the significant scale difference between satellite pixels and UAV pixels, the agricultural film abundance distribution map generated in Step 4 needs to be converted into UAV pixel ground truth values ​​with the same resolution as the satellite. ;

[0262] The specific process is as follows:

[0263] 1) Spatial alignment: Using a Geographic Information System (GIS), the agricultural film abundance distribution map generated in step four is registered with the geographic coordinates in the vector file of the high-risk suspected pollution area obtained in step two, to obtain the registered satellite pixels and UAV pixels;

[0264] 2) Pixel aggregation: Identify each registered satellite pixel. All K drone pixels included within the spatial coverage area;

[0265] 3) Truth value calculation: The area-weighted average method is used to calculate the true value of the UAV pixel corresponding to each registered satellite pixel. The calculation formula is:

[0266]

[0267] In the formula, For the first Agricultural film abundance value (area percentage) of each drone pixel;

[0268] For each satellite pixel The total number of drone pixels included within the spatial coverage area;

[0269] This yields true values ​​that match the scale of satellite observation.

[0270] Step 52: Use the initial satellite inversion results obtained in Step 2 As a predicted value, the true value of the drone pixels As the objective truth value, construct the optimization objective function. This is used to quantify the deviation of the satellite model under the current parameters; it is expressed as:

[0271]

[0272] In the formula, The total number of satellite pixels in the target monitoring area; Weights, including spectral feature weights Background risk weights and the weights of each statistical factor ; For the current weight vector The residual pollution index obtained below; The true value of the drone pixels obtained in step 51;

[0273] Step 53: Using the objective function For weight vector Adaptive adjustments are made to obtain the optimal weight vector; the specific process is as follows:

[0274] 1) Gradient calculation: The gradient descent algorithm is used to calculate the objective function. Regarding the weight vector The partial derivatives of each component are used to determine the gradient direction in which the error decreases the most.

[0275] 2) Parameter update: Update the weight vector according to the direction of the fastest descent of the determined error. This makes the predicted value Gradually approaching the true value ;

[0276] 3) Convergence criterion: When the objective function... Stop the calculation and output the optimal weight vector when the weight no longer decreases significantly or when the preset number of iterations is reached. ;

[0277] Step 54:

[0278] The calculated optimal weight vector substitution formula replace ;

[0279] The calculated optimal weight vector substitution formula replace and ;

[0280] Step 55: Take the steps from step 54... Return to step 52, and repeat steps 53 and 54 until... convergence;

[0281] Each pixel corresponds to a convergent point. ;

[0282] Based on all converged values ​​corresponding to all pixels in the target monitoring region Obtain a high-precision distribution map of agricultural film residue after truth value correction.

[0283] This image combines the macroscopic view from the satellite with the microscopic precision of the drone, representing the final monitoring results.

[0284] The other steps and parameters are the same as those in any of the specific implementation methods one to nine.

[0285] This invention may have other embodiments. Without departing from the spirit and essence of this invention, those skilled in the art can make various corresponding changes and modifications according to this invention, but these corresponding changes and modifications should all fall within the protection scope of the appended claims.

Claims

1. A method for monitoring agricultural film residue pollution based on satellite-to-ground remote sensing data, characterized in that: The specific process of the method is as follows: Step 1: Manually collect samples of agricultural film residue, soil surface, crop leaves and stems, and weeds from the target monitoring area; obtain the raw spectral data of each sample; construct a ground cover spectral database based on the raw spectral data; and select key feature bands that distinguish agricultural film from background ground cover from all the pure raw spectral data contained in the constructed ground cover spectral database. Step 2: Acquire satellite remote sensing images of the target monitoring area; Obtain agricultural statistics for the target monitoring area from the agricultural sector; Spectral features are extracted from the acquired satellite remote sensing images of the target monitoring area to obtain the spectral feature values ​​of the satellite remote sensing images; The background risk coefficient is calculated based on agricultural statistical data of the target monitoring area obtained from the agricultural sector. Spectral features based on satellite remote sensing imagery and background risk coefficient Calculate the residual pollution index ; Set pollution warning thresholds ; Extract all The contiguous pixels are designated as low-risk pollution areas; Extract all The contiguous pixels are designated as medium-risk pollution areas; Extract all The contiguous pixels are designated as high-risk pollution areas; Obtain the geographic coordinate vector file of high-risk pollution areas; Step 3: Use the convex decomposition algorithm to divide the high-risk contamination area in the vector file obtained in Step 2 into multiple convex polygon sub-regions, and record the geographic coordinates of the vertices of each convex polygon sub-region; Based on residual pollution index Determine the UAV's flight altitude by combining each convex polygonal sub-region. ; Based on the drone's flight altitude By combining the pixel size, camera focal length, and field of view of the shortwave infrared hyperspectral camera carried by the UAV, the ground resolution and track spacing of the flight path operation are calculated. For each convex polygonal sub-region, a bow-shaped reciprocating scanning path is generated based on the track spacing; Control the drone to perform automated spectral acquisition operations according to its flight altitude and scanning path; Step 4: Perform radiometric calibration on the raw hyperspectral image acquired by the UAV to obtain the wavelength of the target pixels in the raw hyperspectral image. Actual reflectivity at ; Extracting endmember matrices from raw hyperspectral image data collected by drones The basic components, namely endmembers, are used to determine the types of endmember substances using spectral libraries; Based on actual reflectance and end-member matrix Establish the hybrid pixel decomposition equation; The fully constrained least squares method was used to solve the mixed pixel decomposition equation to obtain the abundance values ​​of various endmembers in each pixel. Based on the abundance values ​​of various endmembers in each pixel, endmember abundance data for the endmember material category "agricultural film" were extracted. ; Based on the end-member abundance data of the extracted end-member material category "agricultural film" Generate an abundance distribution map of agricultural film; Step 5: Convert the agricultural film abundance distribution map generated in Step 4 into UAV pixel ground values ​​with the same resolution as satellites. ; The residual pollution index obtained in step two As a predicted value, the true value of the drone pixels As the objective truth value, construct the optimization objective function. ; Using the objective function For weight vector Adaptive adjustments are made to obtain the optimal weight vector. ; The calculated optimal weight vector substitution formula replace ; The calculated optimal weight vector substitution formula replace and ; The target area is recalculated to generate a true-value corrected distribution map of agricultural film residue.

2. The method for monitoring agricultural film residue pollution based on satellite-to-ground remote sensing data according to claim 1, characterized in that: In step one, samples of agricultural film residue, soil surface, crop leaves and stems, and weeds are manually collected from the target monitoring area; raw spectral data of each sample are obtained; a ground cover spectral database is constructed based on the raw spectral data; and key feature bands that distinguish agricultural film from background ground cover are selected from all the pure raw spectral data contained in the constructed ground cover spectral database. The specific process is as follows: Step 11: Manually collect samples of agricultural film residue, soil surface, crop leaves and stems, and weeds from the target monitoring area; the specific process is as follows: Manually collect agricultural film residue samples of different materials and aging levels from the target monitoring area; Soil surface samples with different textures and moisture contents were manually collected from the target monitoring area. Manually collect crop leaf and stem samples at different growth stages in the target monitoring area; Weed samples at different growth stages were manually collected from the target monitoring area; Steps 1 and 2: Obtain the original spectral data of each sample based on the samples obtained in Step 1; the specific process is as follows: Step 121: In a darkroom environment in the laboratory, under a halogen light source, a hyperspectral imager covering the short-wave infrared band is used to perform spectral measurements on agricultural film residue samples of different materials and aging degrees collected manually in the target monitoring area, so as to obtain the original spectral data of agricultural film residue samples of different materials and aging degrees in the target monitoring area. Step 122: In a darkroom environment in the laboratory, under a halogen light source, a hyperspectral imager covering the short-wave infrared band is used to perform spectral measurements on soil surface samples with different textures and water contents collected manually in the target monitoring area, so as to obtain the original spectral data of soil surface samples with different textures and water contents in the target monitoring area. Steps 1, 2, and 3: In a darkroom environment in the laboratory, under a halogen light source, a hyperspectral imager covering the short-wave infrared band is used to perform spectral measurements on crop leaf and stem samples at different growth stages in the manually collected target monitoring area to obtain the original spectral data of crop leaf and stem samples at different growth stages in the target monitoring area. Steps 1, 2, and 4: In a darkroom environment in the laboratory, a hyperspectral imager covering the short-wave infrared band is used to perform spectral measurements on weed samples at different growth stages in the manually collected target monitoring area under a halogen light source to obtain the original spectral data of weed samples at different growth stages in the target monitoring area. Step 13: Construct a ground cover spectral database based on the original spectral data from Step 1 and Step 2; the specific process is as follows: Step 131: The raw spectral data of agricultural film residue samples of different materials and aging degrees in the target monitoring area are sequentially subjected to noise reduction, smoothing and resampling preprocessing, and the preprocessed pure raw spectral data are entered into the ground cover spectral database. Step 1, 3, 2: The raw spectral data of soil surface samples with different textures and moisture contents in the target monitoring area are sequentially subjected to noise reduction, smoothing and resampling preprocessing, and the preprocessed pure raw spectral data are entered into the ground cover spectral database. Step 133: Perform noise reduction, smoothing and resampling preprocessing on the raw spectral data of crop leaves and stems at different growth stages in the target monitoring area, and enter the clean raw spectral data after preprocessing into the ground cover spectral database. Steps 1, 3, and 4: The raw spectral data of weed samples at different growth stages in the target monitoring area are sequentially subjected to noise reduction, smoothing, and resampling preprocessing. The preprocessed pure raw spectral data are then entered into the ground cover spectral database. Step 135: Based on Step 131, Step 132, Step 133, and Step 134, complete the construction of the ground cover spectral database; Step 14: Select key feature bands from all the pure raw spectral data contained in the land cover spectral database constructed in Step 13 to distinguish agricultural film from background land cover.

3. The method for monitoring agricultural film residue pollution based on satellite-to-ground remote sensing data according to claim 2, characterized in that: In step one four, key feature bands that distinguish agricultural film from background land cover are selected from all the pure raw spectral data contained in the land cover spectral database constructed in step one three. The specific process is as follows: Step 141: Perform envelope removal processing on the reflectance curves of all pure raw spectral data contained in the ground cover spectral database constructed in Step 13 to obtain all pure raw spectral data after envelope removal processing. The formula for calculating envelope removal is as follows: In the formula, Wavelength representing pure, raw spectral data The envelope removal value at the location; Wavelengths of pure, raw spectral data Reflectance at that location; Wavelengths of pure, raw spectral data The value on the convex circumscribed envelope of the spectral curve; Step 142: Based on all the pure original spectral data after envelope removal processing obtained in Step 141, extract the characteristic absorption peak positions from all the pure original spectral curves after envelope removal processing; the specific process is as follows: The first-order differential spectrum of the original spectral data of agricultural film residue samples was calculated using the central difference method. The formula is as follows: in, Wavelength of the original spectral data for agricultural film residue samples The envelope removal value at the location; Wavelength of the original spectral data for agricultural film residue samples The envelope removal value at the location; wavelength place, wavelength Place; Identify the zero point position in the first-order differential spectrum; find the inflection point where the slope of the original spectrum curve of the agricultural film residue sample changes from negative to positive, corresponding to the zero point of the first-order differential. The inflection point is the minimum point, which is the potential absorption peak center. Pre-defined significance criteria are introduced to screen potential absorption peak center wavelengths; The significance criteria include: If the wavelength of the potential absorption peak center is lower than the preset intensity threshold, the potential absorption peak center is retained. Analyze the rate of change of the derivative in the neighborhood of the potential absorption peak center. If the rate of change of the derivative on both sides of zero is greater than the preset rate of change threshold, the potential absorption peak center is retained. The center wavelength corresponding to the center position of the potential absorption peak retained after screening is used as the candidate position of the key characteristic band. Step 143: Calculate the spectral separability index (SDI); obtain key characteristic bands based on the SDI; the specific process is as follows: Step 1431: Calculate the spectral separability index (SDI); the calculation formula is as follows: In the formula, Indicates wavelength The separability index at the location; This indicates the use of the wavelength acquired in step one. The average value of the envelope removal value of the agricultural film sample set at that location; This indicates the use of the wavelength acquired in step one. The average value of the envelope removal values ​​of the background sample set at that location; Indicates wavelength The standard deviation of the envelope removal values ​​in the agricultural film sample set at that location; Indicates wavelength The standard deviation of the envelope removal values ​​of the background sample set at that location; wavelength ; Step 1432: Traverse the entire spectral range and use the SDI formula to calculate the spectral separability of the agricultural film residue sample set and the background ground cover sample set in each spectral band; From the candidate locations of the key feature bands determined in steps one, four, and two, select the band intervals with SDI values ​​greater than a preset threshold. Finally, P band intervals with high discriminative power were selected as the key feature band set.

4. The method for monitoring agricultural film residue pollution based on satellite-to-ground remote sensing data according to claim 3, characterized in that: In step two, satellite remote sensing images of the target monitoring area are acquired. Obtain agricultural statistics for the target monitoring area from the agricultural sector; Spectral features are extracted from the acquired satellite remote sensing images of the target monitoring area to obtain the spectral feature values ​​of the satellite remote sensing images; The background risk coefficient is calculated based on agricultural statistical data of the target monitoring area obtained from the agricultural sector. Spectral features based on satellite remote sensing imagery and background risk coefficient Calculate the residual pollution index ; Set pollution warning thresholds ; Extract all The contiguous pixels are designated as low-risk pollution areas; Extract all The contiguous pixels are designated as medium-risk pollution areas; Extract all The contiguous pixels are designated as high-risk pollution areas; Obtain the geographic coordinate vector file of high-risk pollution areas; The specific process is as follows: Step 21: Acquire satellite remote sensing images of the target monitoring area; Obtain agricultural statistics for the target monitoring area from the agricultural sector; Agricultural statistics for the target monitoring area include crop types, the historical number of years of plastic film mulching in each plot, the amount of plastic film used in each plot in the current year, and the method of plastic film recycling. Step 22: Extract spectral features from the satellite remote sensing images of the target monitoring area obtained in Step 21 to obtain the spectral feature values ​​of the satellite remote sensing images; The specific process is as follows: Step 221: Perform radiometric calibration and atmospheric correction preprocessing on the satellite remote sensing images of the target monitoring area obtained in Step 21 to obtain the preprocessed satellite remote sensing images of the target monitoring area; Step Two Two Two: From the bands of the preprocessed satellite remote sensing image of the target monitoring area, select one band from the key feature bands obtained in step one as the sensitive band, with a reflectance of [missing information]. And a band outside the key feature band set interval is used as a reference band, with a reflectance of Based on reflectivity and reflectivity Constructing spectral eigenvalues The calculation formula is as follows: Steps 2 and 3: Calculate the background risk coefficient based on the agricultural statistics data of the target monitoring area obtained from the agricultural sector in Step 21; The specific process is as follows: Agricultural statistics for the target monitoring area include crop types, the historical number of years of plastic film mulching in each plot, the amount of plastic film used in each plot in the current year, and the method of plastic film recycling. Background risk coefficient as follows: In the formula, For the first The weights of the statistical variables; The statistical variables are four: crop type, historical years of mulching in each plot, amount of agricultural film used in each plot in the current year, and method of agricultural film recycling. For the first Normalized values ​​of several statistical variables; Step 24: Spectral feature values ​​of the satellite remote sensing image obtained in Step 22 and the background risk coefficient obtained in steps two and three Calculate the residual pollution index ; The specific process is as follows: In the formula, For spectral feature weights, and Background risk weights; Step 25: Set pollution warning thresholds ; Extract all The contiguous pixels are designated as low-risk pollution areas; Extract all The contiguous pixels are designated as medium-risk pollution areas; Extract all The contiguous pixels are designated as high-risk pollution areas; Obtain the geographic coordinate vector file of high-risk pollution areas.

5. The method for monitoring agricultural film residue pollution based on satellite-to-ground remote sensing data according to claim 4, characterized in that: In step three, a convex decomposition algorithm is used to divide the high-risk contamination area in the vector file obtained in step two into multiple convex polygon sub-regions, and the geographical coordinates of the vertices of each convex polygon sub-region are recorded. Based on residual pollution index The drone's flight altitude is determined by each convex polygonal sub-region. ; Based on the drone's flight altitude By combining the pixel size, camera focal length, and field of view of the shortwave infrared hyperspectral camera carried by the UAV, the ground resolution and track spacing of the flight path operation are calculated. For each convex polygonal sub-region, a bow-shaped reciprocating scanning path is generated based on the track spacing; The specific process is as follows: Step 31: Use the convex decomposition algorithm to divide the high-risk contamination area in the vector file obtained in Step 2 into multiple convex polygon sub-regions, and record the geographic coordinates of the vertices of each convex polygon sub-region; the specific process is as follows: The high-risk contamination area obtained in step two is an irregular polygon. Identify all concave vertices in the irregular polygon and emit rays along the angle bisector or diagonal at the concave vertices to cut the irregular polygon into multiple sub-regions. Determine whether all sub-regions are convex polygonal sub-regions; If so, record the geographic coordinates of the vertices of each convex polygon sub-region; If not, for the non-convex polygon sub-region, emit rays along the angle bisector or diagonal direction at the concave vertex to cut the non-convex polygon sub-region into sub-regions until all sub-regions are convex polygon sub-regions, and record the geographic coordinates of the vertices of each convex polygon sub-region. The concave vertex is a vertex with an interior angle greater than 180 degrees; Step 3.2: Residual pollution index obtained in Step 2 Determine the UAV's flight altitude by combining each convex polygonal sub-region. The specific process is as follows: Based on the residual pollution index obtained in step two Determine the contaminated area for each convex polygon sub-region; If the residual pollution index of the convex polygon sub-region The corresponding convex polygon sub-region is a low-risk pollution area; If the residual pollution index of the convex polygon sub-region Then the corresponding convex polygon sub-region is a medium-risk pollution area; If the residual pollution index of the convex polygon sub-region Then the corresponding convex polygon sub-region is a high-risk pollution area; For high-risk pollution areas, the drone flight altitude will be adjusted. The range is set as follows: ; For medium-risk contamination areas, the drone flight altitude will be adjusted. The range is set as follows: ; For low-risk pollution areas, the drone flight altitude will be adjusted. The range is set as follows: ; Step 33: Based on the drone flight altitude determined in Step 32 By combining the pixel size, camera focal length, and field of view of the shortwave infrared hyperspectral camera carried by the UAV, the ground resolution and track spacing of the flight path operation are calculated. For each convex polygonal sub-region, a bow-shaped reciprocating scanning path is generated based on the track spacing; Steps 3 and 4: Control the drone to perform automated spectral acquisition operations according to the flight altitude and scanning path.

6. The method for monitoring agricultural film residue pollution based on satellite-to-ground remote sensing data according to claim 5, characterized in that: In step 33, the drone flight altitude is determined based on step 32. Based on the pixel size, focal length, and field of view of the shortwave infrared hyperspectral camera carried by the UAV, the ground resolution and track spacing of the flight path operation are calculated; the specific process is as follows: Step 331: Calculate the ground resolution ; indicates as: in, For pixel size, The focal length of the camera; Step 332: Calculate the ground cover width of a single image ; indicates as: in, This refers to the camera's field of view. Step 333: Based on the ground coverage width of a single image Calculate track spacing ; indicates as: in, To meet the preset overlap requirements of the jigsaw puzzle.

7. The method for monitoring agricultural film residue pollution based on satellite-to-ground remote sensing data according to claim 6, characterized in that: In steps three and four, for each convex polygonal sub-region, a bow-shaped reciprocating scanning path is generated based on the track spacing; the specific process is as follows: The longest side direction of each convex polygon sub-region is selected as the reference direction of the main flight path. Translation sampling is performed along the perpendicular direction of the reference direction at intervals D to generate a series of parallel flight paths. Automatically add buffer waypoints at the end of the extension lines of each parallel route scan segment; Combining the geographic coordinates of the vertices of each convex polygon sub-region determined in step 31, the flight altitude H determined in step 32, the track spacing D calculated in step 33, and the preset heading overlap rate, the planar coordinates of the endpoints of each parallel segment and the buffer waypoint are spatially mapped and superimposed with the flight altitude H to obtain the three-dimensional coordinate information containing the endpoints of each parallel segment and the buffer waypoint. The flight speed and camera trigger command are determined based on the coverage rate of the ground coverage area of ​​a single image by the camera sampling frequency and the preset forward overlap rate. The UAV is controlled to perform automated spectral acquisition operations based on the three-dimensional coordinate information of the endpoints and buffer waypoints of each parallel flight segment, flight speed, and camera trigger commands.

8. The method for monitoring agricultural film residue pollution based on satellite-to-ground remote sensing data according to claim 7, characterized in that: In step four, radiometric calibration is performed on the raw hyperspectral image acquired by the UAV to obtain the wavelength of the target pixels in the raw hyperspectral image. Actual reflectivity at ; Extracting endmember matrices from raw hyperspectral image data collected by drones The basic components, namely endmembers, are used to determine the types of endmember substances using spectral libraries; Based on actual reflectance and end-member matrix A mixed pixel decomposition equation was established; the mixed pixel decomposition equation was solved using the fully constrained least squares method to obtain the abundance values ​​of various material endmembers in each pixel; Based on the abundance values ​​of various substances in each pixel, endmember abundance data for the endmember substance category "agricultural film" was extracted. ; Based on the end-member abundance data of the extracted end-member material category "agricultural film" Generate an abundance distribution map of agricultural film; The specific process is as follows: Step 41: Perform radiometric calibration on the raw hyperspectral image acquired by the UAV to obtain the wavelength of the target pixels in the raw hyperspectral image. Actual reflectivity at ; The calculation formula is: In the formula, For target pixels in the original hyperspectral image at wavelength The actual reflectivity at that location; For target pixels in the original hyperspectral image at wavelength The original observations at the location; These are the observations from a standard whiteboard; This represents the dark current noise value. Laboratory calibration of reflectivity for a standard whiteboard; Step 42: The radiometrically calibrated images undergo atmospheric correction and geometric correction preprocessing to obtain preprocessed full-band images; Subsequently, based on the key feature band set determined in step one, the corresponding N feature channel data are extracted from the preprocessed full-band image by wavelength position alignment to obtain the dimension-reduced image data. Step 43: Extract the endmember matrix from the raw hyperspectral image data collected by the UAV. The basic components, namely endmembers, are used to determine the types of endmember substances using spectral libraries; Step 44: Based on the actual reflectivity in Step 41 and the endmember matrix in step four three A mixed pixel decomposition equation was established; the mixed pixel decomposition equation was solved using the fully constrained least squares method to obtain the abundance values ​​of various material endmembers in each pixel; The mixed pixel decomposition equation is expressed as: The mixed pixel decomposition equation satisfies: , ; In the formula, Endmember matrix The Middle Abundance of endmembers in pixels; Endmember matrix Total number of mid-range elements; Endmember matrix The Middle Seed endogenous element; For residuals; Step 4-5: Based on the abundance values ​​of various substances in each pixel obtained in Step 4-4, extract the endmember abundance data of "agricultural film" as the endmember substance category determined in Step 4-3. ; Based on the extracted endmember abundance data of "agricultural film" Generate an abundance distribution map of agricultural film.

9. A method for monitoring agricultural film residue pollution based on satellite-to-ground remote sensing data according to claim 8, characterized in that: In step four-three, endmember matrices are extracted from the raw hyperspectral image data collected by the UAV. The basic components, namely endmembers, are determined using a spectral library to classify endmember substances; the specific process is as follows: Step 431: Extract the basic components constituting the mixed pixels, i.e., candidate endmembers, from the raw hyperspectral image data collected by the UAV; the specific process is as follows: Using the endmember extraction algorithm, we find the minimum volume simplex that encloses all pixels in the dimensionality-reduced image data obtained in step four-two, and extract all the minimum volume simplexes. There are vertices, each vertex being a candidate endmember vector. ; Step 432: Calculate the vector of each candidate endmember extracted in Step 431. Reference spectra in the ground cover spectral database constructed in step one spectral angles between ; Take the spectral angle Minimum reference spectrum As the matching object, the candidate endmember vector The attribute is marked as The corresponding material category; spectral angle Represented as: In the formula, It is the spectral angle; For target pixels in the original hyperspectral image at wavelength Each endmember vector at; Wavelengths in the ground cover spectral database constructed in step one Reference spectrum at that location; Step 433: Obtain the material category corresponding to all candidate endmember vectors; Construct an end-member matrix containing different components of agricultural film, soil, and vegetation. .

10. The method for monitoring agricultural film residue pollution based on satellite-to-ground remote sensing data according to claim 9, characterized in that: In step five, the agricultural film abundance distribution map generated in step four is converted into UAV pixel ground values ​​with satellite resolution. ; The residual pollution index obtained in step two As a predicted value, the true value of the drone pixels As the objective truth value, construct the optimization objective function. ; Using the objective function For weight vector Adaptive adjustments are made to obtain the optimal weight vector. ; The calculated optimal weight vector substitution formula replace ; The calculated optimal weight vector substitution formula replace and ; The target area is recalculated to generate a true-value corrected distribution map of agricultural film residue. The specific process is as follows: Step 51: Convert the agricultural film abundance distribution map generated in Step 4 into UAV pixel ground values ​​with the same resolution as satellites. ; The specific process is as follows: 1) Spatial alignment: Using a geographic information system, the agricultural film abundance distribution map generated in step four is registered with the geographic coordinates in the geographic coordinate vector file of the high-risk pollution area obtained in step two to obtain the registered satellite pixels and UAV pixels; 2) Pixel aggregation: Identify each registered satellite pixel. All K drone pixels included within the spatial coverage area; 3) Truth value calculation: The area-weighted average method is used to calculate the true value of the UAV pixel corresponding to each registered satellite pixel. The calculation formula is: In the formula, For the first The abundance value of agricultural film per drone pixel; For each satellite pixel The total number of drone pixels included within the spatial coverage area; Step 52: Use the initial satellite inversion results obtained in Step 2 As a predicted value, the true value of the drone pixels As the objective truth value, construct the optimization objective function. ; indicates as: In the formula, The total number of satellite pixels in the target monitoring area; Weights, including spectral feature weights Background risk weights and the weights of each statistical factor ; For the current weight vector The residual pollution index obtained below; The true value of the drone pixels obtained in step 51; Step 53: Using the objective function For weight vector Adaptive adjustments are made to obtain the optimal weight vector; the specific process is as follows: 1) Gradient calculation: The gradient descent algorithm is used to calculate the objective function. Regarding the weight vector The partial derivatives of each component are used to determine the gradient direction in which the error decreases the most. 2) Parameter update: Update the weight vector according to the direction of the fastest descent of the determined error. This makes the predicted value Gradually approaching the true value ; 3) Convergence criterion: When the objective function... Stop the calculation and output the optimal weight vector when the weight no longer decreases significantly or when the preset number of iterations is reached. ; Step 54: The calculated optimal weight vector substitution formula replace ; The calculated optimal weight vector substitution formula replace and ; Step 55: Take the steps from step 54... Return to step 52, and repeat steps 53 and 54 until... convergence; Each pixel corresponds to a convergent point. ; Based on all converged values ​​corresponding to all pixels in the target monitoring region Obtain a distribution map of agricultural film residues after truth value correction.

Images