An intelligent extraction method and system for irrigation area of an irrigation district based on satellite remote sensing images

Abstract

The application provides a kind of based on satellite remote sensing image irrigation area intelligent extraction method and system of irrigation area, it is related to agricultural remote sensing and irrigation monitoring technical field, the method includes: obtaining mulched drip irrigation cotton field multi-temporal satellite image in spring irrigation period, and it is registered with irrigation area engineering BIM planning model in space, obtain the original image set for fine design orientation;According to the original image set for fine design orientation, extract the sequence of change of time of pixel thermal infrared temperature, obtain the pixel temperature time series for supporting irrigation area planning design;According to the pixel temperature time series, calculate the radiation temperature difference value of adjacent transit time in the afternoon, obtain the per-hour cooling amplitude;According to the per-hour cooling amplitude, calculate the cooling amplitude difference value of adjacent time, obtain the cooling rate step amount.The application realizes the intelligent extraction of irrigation area irrigation area by capturing the spatial turning position of cooling amplitude from slow drop to sudden drop along the vertical direction of drip irrigation belt.

Description

Technical Field

[0001] This invention relates to the field of agricultural remote sensing and irrigation monitoring technology, and in particular to a method and system for intelligent extraction of irrigation area in irrigation districts based on satellite remote sensing images. Background Technology

[0002] Accurate determination of the irrigated area is a fundamental element for irrigation district planning and design, optimal water resource allocation, and modernization. Traditional methods for surveying irrigated areas largely rely on manual field surveys or basic optical remote sensing interpretation. However, these existing methods have the following shortcomings in practical engineering applications:

[0003] Manual field surveys are time-consuming and costly, making it difficult to achieve large-scale, high-frequency dynamic monitoring. Existing conventional remote sensing extraction methods often remain at the level of pure surveying theory, only acquiring two-dimensional static information about the land surface. This pure surveying method lacks deep integration with the three-dimensional spatial layout and water conservancy design of irrigation area projects, and the data obtained is difficult to directly support the refined planning and design of modern irrigation areas.

[0004] Existing remote sensing methods rely heavily on the spectral characteristics of vegetation cover when extracting irrigated area, neglecting key physical processes such as actual infiltration, lateral diffusion, and heat dissipation of irrigation water in the soil. This results in extraction results that fail to reflect the true effective irrigated area.

[0005] Data obtained through traditional methods cannot be directly input into and integrated into the BIM system for irrigation district projects. This results in a lack of accurate, dynamic, and integrated underlying data support during the planning and design phase. Consequently, designers lack a basis for optimizing drip irrigation tape layout and adjusting irrigation regimes, hindering the modernization and refined design process of irrigation district projects. Summary of the Invention

[0006] This invention provides a method and system for intelligent extraction of irrigation area in irrigation districts based on satellite remote sensing imagery. By capturing the spatial transition point where the temperature drop changes from a slow to a sudden drop along the direction perpendicular to the drip irrigation belt, the boundary of the irrigated wet area is identified by the spatial second-order change of the temperature drop rate step. This overcomes the shortcomings of fixed temperature thresholds being unable to adapt to differences in soil texture and uniform decay models being unable to identify abrupt changes in decay rate, thus realizing intelligent extraction of irrigation area in irrigation districts.

[0007] To solve the above-mentioned technical problems, the technical solution of the present invention is as follows:

[0008] Firstly, a method for intelligent extraction of irrigation area in irrigation districts based on satellite remote sensing imagery, the method comprising:

[0009] Step 1: Acquire multi-temporal satellite images of cotton fields under mulching and drip irrigation during the spring irrigation season, and spatially register them with the BIM planning model of the irrigation area project to obtain the original image set for refined design; Based on the original image set for refined design, extract the pixel thermal infrared temperature change sequence over time to obtain the pixel temperature time series used to support the irrigation area planning and design.

[0010] Step 2: Based on the pixel temperature time series, calculate the radiation temperature difference between adjacent transit times in the afternoon to obtain the hourly temperature drop.

[0011] Step 3: Based on the hourly temperature drop, calculate the difference in temperature drop between adjacent time intervals to obtain the step rate of temperature drop.

[0012] Step 4: Based on the step size of the cooling rate and the vertical distance between the pixel and the drip irrigation tape, construct a transverse spatial attenuation profile perpendicular to the direction of the drip irrigation tape to obtain the profile attenuation sequence.

[0013] Step 5: Based on the profile attenuation sequence, determine the spatial coordinates of the abrupt change in the sequence slope to obtain the geometric boundary line of water lateral expansion; based on the geometric boundary line of water lateral expansion, extract the set of continuous high-step pixels inside the boundary line to obtain the active pixels of water evaporation.

[0014] Step 6: Based on the active pixels of water evaporation, gradually add the cooling rate step along the drip irrigation tape to obtain the cumulative amount of heat dissipation in the strip.

[0015] Step 7: Calculate the irrigation water infiltration probability of the pixels based on the cumulative heat dissipation of the strip and the short-wave infrared reflectivity to obtain the irrigation status; based on the irrigation status, calculate the actual surface coverage area corresponding to the total number of irrigated pixels to obtain the irrigated area of ​​the irrigation district.

[0016] Secondly, the intelligent extraction system for irrigation area based on satellite remote sensing imagery includes:

[0017] The acquisition module is used to acquire multi-temporal satellite images of cotton fields under mulching and drip irrigation during the spring irrigation season, and spatially register them with the BIM planning model of the irrigation area project to obtain the original image set for refined design; based on the original image set for refined design, the pixel thermal infrared temperature change sequence over time is extracted to obtain the pixel temperature time series used to support the irrigation area planning and design.

[0018] The calculation module is used to calculate the difference in radiation temperature between adjacent afternoon transit times based on the pixel temperature time series, and obtain the hourly temperature drop amplitude; based on the hourly temperature drop amplitude, it calculates the difference in temperature drop amplitude between adjacent times, and obtains the step amount of the temperature drop rate.

[0019] The module is used to construct a transverse spatial attenuation profile perpendicular to the direction of the drip irrigation tape based on the step amount of cooling rate and the vertical distance of the pixel from the drip irrigation tape, and obtain the profile attenuation sequence.

[0020] The interception module is used to determine the spatial coordinates of the abrupt change in the slope of the sequence based on the profile attenuation sequence, and obtain the geometric boundary line of water lateral expansion; based on the geometric boundary line of water lateral expansion, it intercepts the set of continuous high-step pixels inside the boundary line to obtain the active pixels of water evaporation.

[0021] The accumulation module is used to accumulate the cooling rate step by step along the drip irrigation tape based on the active water evaporation pixels to obtain the cumulative heat dissipation of the strip.

[0022] The extraction module is used to calculate the irrigation water infiltration probability of pixels based on the cumulative amount of heat dissipation in the strip and the short-wave infrared reflectivity, so as to obtain the irrigation status; based on the irrigation status, the actual surface coverage area corresponding to the total number of irrigated pixels is calculated to obtain the irrigated area of ​​the irrigation district.

[0023] The above-described solution of the present invention has at least the following beneficial effects:

[0024] By calculating the cooling rate step value using thermal infrared time-series data, the evaporative cooling effect of irrigation is quantified into a thermal signal parameter, overcoming the problem that optical remote sensing cannot penetrate the film to perceive the wet body, and achieving accurate characterization of the intensity of irrigation water evaporative cooling. By constructing a spatial attenuation profile of the cooling rate step value along the transverse direction of the drip irrigation strip, the lateral expansion boundary is located based on the slope change and the inner high step value pixels are extracted, overcoming the problem that fixed threshold cannot adapt to the difference in wet distance, and achieving adaptive high-precision positioning of the drip irrigation wet body boundary. By accumulating the cooling rate step value along the direction of the drip irrigation strip to generate the cumulative heat dissipation, and combining it with short-wave infrared reflectivity to calculate the water penetration probability to screen irrigation pixels and calculate the area, the problem of signal loss due to weather fluctuations or shading of a single thermal signal is overcome, and dual-source mutual verification of thermal signal and water absorption characteristics and accurate extraction of irrigation area are achieved. Attached Figure Description

[0025] Figure 1 This is a flowchart illustrating an intelligent extraction method for irrigation area based on satellite remote sensing imagery provided by an embodiment of the present invention.

[0026] Figure 2 This is a schematic diagram of an intelligent irrigation area extraction system based on satellite remote sensing imagery provided by an embodiment of the present invention.

[0027] Figure 3 This is a flowchart illustrating the process of calculating the irrigation water infiltration probability of a pixel based on the cumulative amount of heat dissipation in a strip and the short-wave infrared reflectivity, thereby obtaining the irrigation status; and calculating the actual surface coverage area corresponding to the total number of irrigated pixels based on the irrigation status, thereby obtaining the irrigated area of ​​the irrigation district. Detailed Implementation

[0028] Exemplary embodiments of the present disclosure will now be described in more detail with reference to the accompanying drawings. While exemplary embodiments of the present disclosure are shown in the drawings, it should be understood that the present disclosure may be implemented in various forms and should not be limited to the embodiments set forth herein. Rather, these embodiments are provided so that this disclosure will be thorough and complete, and will fully convey the scope of the disclosure to those skilled in the art.

[0029] like Figure 1 As shown, an embodiment of the present invention proposes a method for intelligent extraction of irrigation area in irrigation districts based on satellite remote sensing imagery. The method includes the following steps:

[0030] Step 1: Acquire multi-temporal satellite images of cotton fields under mulching and drip irrigation during the spring irrigation season, and spatially register them with the BIM planning model of the irrigation area project to obtain the original image set for refined design; Based on the original image set for refined design, extract the pixel thermal infrared temperature change sequence over time to obtain the pixel temperature time series used to support the irrigation area planning and design.

[0031] Step 2: Based on the pixel temperature time series, calculate the radiation temperature difference between adjacent transit times in the afternoon to obtain the hourly temperature drop.

[0032] Step 3: Based on the hourly temperature drop, calculate the difference in temperature drop between adjacent time intervals to obtain the step rate of temperature drop.

[0033] Step 4: Based on the step size of the cooling rate and the vertical distance between the pixel and the drip irrigation tape, construct a transverse spatial attenuation profile perpendicular to the direction of the drip irrigation tape to obtain the profile attenuation sequence.

[0034] Step 5: Based on the profile attenuation sequence, determine the spatial coordinates of the abrupt change in the sequence slope to obtain the geometric boundary line of water lateral expansion; based on the geometric boundary line of water lateral expansion, extract the set of continuous high-step pixels inside the boundary line to obtain the active pixels of water evaporation.

[0035] Step 6: Based on the active pixels of water evaporation, gradually add the cooling rate step along the drip irrigation tape to obtain the cumulative amount of heat dissipation in the strip.

[0036] Step 7: Calculate the irrigation water infiltration probability of the pixels based on the cumulative heat dissipation of the strip and the short-wave infrared reflectivity to obtain the irrigation status; based on the irrigation status, calculate the actual surface coverage area corresponding to the total number of irrigated pixels to obtain the irrigated area of ​​the irrigation district.

[0037] In this embodiment of the invention, the cooling rate step is calculated using thermal infrared time-series data, quantifying the irrigation evaporative cooling effect into a thermal signal parameter. This overcomes the problem that optical remote sensing cannot penetrate the film to perceive the humid body, achieving accurate characterization of the intensity of irrigation water evaporative cooling. By constructing a spatial attenuation profile of the cooling rate step along the transverse direction of the drip irrigation strip, the lateral expansion boundary is located based on the slope change and high-step pixels on the inner side are extracted. This overcomes the problem that a fixed threshold cannot adapt to differences in humidification distance, achieving adaptive high-precision positioning of the drip irrigation humid body boundary. By accumulating the cooling rate step along the drip irrigation strip direction to generate the cumulative heat dissipation, and combining it with short-wave infrared reflectivity to calculate the water penetration probability to screen irrigation pixels and calculate the area, this overcomes the problem that a single thermal signal is affected by weather fluctuations or shading, resulting in signal loss. This achieves dual-source mutual verification of thermal signals and water absorption characteristics and accurate extraction of irrigation area.

[0038] In a preferred embodiment of the present invention, step 1 above involves acquiring multi-temporal satellite images of cotton fields under mulching and drip irrigation during the spring irrigation season, and spatially registering these images with the BIM planning model of the irrigation district project to obtain an original image set for refined design; based on the original image set for refined design, extracting the time sequence of pixel thermal infrared temperature to obtain the pixel temperature time series used to support the irrigation district planning and design, specifically including:

[0039] Step 1.1: Based on the irrigation district planning and design requirements, determine the geographical scope of the target irrigation district using the vector boundary file provided by the irrigation district management unit. The geographical scope is represented by the longitude and latitude coordinates of the four corners of the circumscribed rectangle, in degrees. Select multi-temporal multispectral satellite remote sensing images that cover the above geographical scope and have a cloud cover rate of less than 10% for each scene. The image source is the Gaofen series satellite, Landsat, or Sentinel-2 satellite platform. The spatial resolution of the image is not less than 10 meters to ensure that a single irrigated plot (usually 0.1 to 1 hectare) occupies enough pixels on the image for reliable identification. The image bands must include at least the red, green, and blue bands of visible light, the near-infrared band, and the shortwave infrared band to meet the needs of subsequent vegetation index and water body index calculations.

[0040] The image acquisition period covers the main irrigation season of the irrigation area. Taking the Hetao Irrigation Area as an example, the main irrigation season is from April to September each year. The time interval between two adjacent images is set according to the irrigation cycle of the irrigation area, with a typical value of ten to sixteen days, corresponding to the revisit cycle of the Sentinel-2 dual-satellite network. The images of each period are arranged in chronological order of imaging time throughout the irrigation season.

[0041] Step 1.2: The recorded values ​​of each pixel in the original satellite remote sensing image are dimensionless digital quantization values. Radiometric calibration converts the values ​​of each band into physically meaningful apparent radiance, in watts per square meter per steradian per micrometer. The original digital quantization value of each pixel in each band is multiplied by the radiometric calibration gain coefficient of that band, and then added to the radiometric calibration bias coefficient of that band to obtain the apparent radiance of that pixel.

[0042] Step 1.3, Atmospheric correction aims to eliminate the influence of atmospheric effects such as atmospheric molecular scattering, aerosol absorption, and water vapor absorption on the surface reflectance signal, converting apparent radiance into true surface reflectance. An atmospheric correction method based on a radiative transfer model is employed. Required input parameters include the solar zenith angle and solar azimuth angle at the time of image imaging, both read from image metadata; atmospheric aerosol optical thickness, obtained from contemporaneous aerosol product data or retrieved from the image itself using the dark target method; and atmospheric water vapor content, obtained from image metadata or contemporaneous meteorological reanalysis data. After atmospheric correction, the value of each pixel in each band represents the surface reflectance, dimensionless, ranging from zero to one. Surface reflectance represents the proportion of incident solar radiation reflected at a specific location on the Earth's surface in a given band; this value depends only on the spectral characteristics of the surface object itself and is unaffected by atmospheric conditions or imaging geometry.

[0043] Step 1.4: Using the evenly distributed ground control points within the irrigation area as a benchmark, the typical number of ground control points is twenty to thirty. The ground control points are permanent features within the irrigation area that are fixed in location and clearly identifiable on the image, including canal intersections, road junctions, and pump station building corners. Their precise geographic coordinates are obtained by on-site measurement using differential positioning equipment or extracted from existing irrigation area engineering survey results.

[0044] Based on the correspondence between image coordinates and their precise geographic coordinates at ground control points, a polynomial correction model is used to establish a transformation function from image coordinates to geographic coordinates, eliminating geometric distortions caused by factors such as sensor attitude, terrain undulation, and Earth curvature. After correction, all pixels in the image have a unified geographic coordinate reference.

[0045] Using the image from the middle of the irrigation period as the reference image, images from other periods are registered with the reference image. The registration method is automatic registration based on normalized mutual matching. Fifty to one hundred pairs of corresponding image points are automatically extracted between the reference image and the image to be registered. The registration accuracy is controlled within half a pixel. All registered images satisfy the one-to-one spatial correspondence at the pixel level, that is, the same row and column coordinates on images from different periods correspond to the same physical location on the surface of the irrigation area.

[0046] Step 1.5: Based on the vector boundary file provided by the irrigation district management unit, perform boundary cropping on all registered images one by one, retaining only the pixels inside the vector boundary and removing the pixels outside the boundary. After cropping, the number of pixel rows and columns in each image is reduced from the original full-scene image size to the size of the smallest bounding rectangle area covering the irrigation district boundary.

[0047] The cropped images from each period are arranged chronologically to obtain the basic dataset of remote sensing images for the irrigation district. This dataset not only serves as the input basis for subsequent calculations of multi-temporal vegetation indices and water body indices, but also as the base map data source for the three-dimensional information planning model of the irrigation district. After spatial registration, it can be overlaid onto the topographic layer of the three-dimensional model of the irrigation district, providing highly timely and unified spatial benchmark surface spectral information support for canal system layout optimization, plot flatness analysis, and spatial planning of irrigation facilities. Compared with the traditional manual household-by-household survey method, which requires several months to half a year to complete a survey of the irrigated area of ​​the entire district, this method can complete the image update of the entire district within 10 to 16 days using satellite remote sensing imagery, increasing the data update frequency by more than 10 times. This provides a timeliness guarantee for dynamic decision-making on the optimal allocation of water resources in the irrigation district that cannot be achieved by traditional methods.

[0048] In a preferred embodiment of the present invention, step 2 above, which calculates the radiation temperature difference between adjacent afternoon transit times based on the pixel temperature time series to obtain the hourly temperature drop, specifically includes:

[0049] Step 2.1: For the surface temperature variation sequence of each pixel in the pixel temperature time series, extract the temperature observation values ​​between 12:00 noon and 5:00 pm local time when the satellite passes over the time series, and obtain the afternoon temperature time series subset; traverse all afternoon temperature observation records of each pixel in the afternoon temperature time series subset within each observation date, and pair every two adjacent passing observations according to the order of passing times into a passing time pair, and record the earlier time in each passing time pair as the preceding passing time and the later time as the following passing time, to obtain the daily adjacent passing time pair sequence of each pixel.

[0050] Step 2.2: For each pair of adjacent transit time pairs in the sequence, extract the surface temperature value corresponding to the preceding transit time and the surface temperature value corresponding to the following transit time. Calculate the surface temperature value of the following transit time minus the surface temperature value of the preceding transit time, and use the difference as the afternoon cooling amount of that pixel within the time interval corresponding to that transit time pair. Arrange the afternoon cooling amounts corresponding to each transit time pair for the same pixel in chronological order across all valid observation dates to form the daily afternoon cooling amount sequence for that pixel. Use the set of the daily afternoon cooling amount sequences for all valid pixels as the hourly cooling amplitude.

[0051] In a preferred embodiment of the present invention, step 3 above, which calculates the difference in cooling rate between adjacent time intervals based on the hourly cooling rate, to obtain the cooling rate step, specifically includes:

[0052] Step 3.1: For the daily afternoon cooling amount sequence of each pixel in the hourly cooling range, pair the afternoon cooling amount of each observation date with the afternoon cooling amount of the previous date. In each pair, the afternoon cooling amount of the previous date is recorded as the previous day's cooling amount, and the afternoon cooling amount of the current date is recorded as the current day's cooling amount, thus obtaining the adjacent date cooling amount pair sequence of each pixel.

[0053] Step 3.2: For each pair of adjacent date cooling amount pairs in the sequence, calculate the difference between the current day's cooling amount and the previous day's cooling amount, and use the difference as the cooling rate change of the pixel between the adjacent dates; arrange the cooling rate change of the same pixel between all adjacent date pairs in chronological order to form the cooling rate change sequence of the pixel; and take the set of the cooling rate change sequences of all valid pixels as the cooling rate step.

[0054] In a preferred embodiment of the present invention, step 4 above, based on the cooling rate step and the vertical distance between the pixel and the drip irrigation tape, constructs a transverse spatial attenuation profile perpendicular to the direction of the drip irrigation tape, and obtains a profile attenuation sequence, including:

[0055] Step 4.1: Based on the cooling rate step, extract the vertical projection distance of each pixel center point in the direction perpendicular to the drip irrigation tape, to obtain the lateral spatial distance sequence, including:

[0056] Step 4.1.1: Based on the step change in cooling rate, extract the center points of continuous pixels whose values ​​reach local valleys to obtain the spatial center line of the drip irrigation tape. Specifically, this includes:

[0057] In cotton fields irrigated by drip irrigation under plastic film, the soil directly above the drip irrigation strip experiences the strongest evaporative cooling effect due to continuous watering. The surface temperature is significantly lower in the afternoon than the dry areas on both sides far from the drip irrigation strip. This temperature difference is represented by a low-temperature valley line extending along the direction of the drip irrigation strip in the cooling rate step data. The cooling rate step value at the location of the drip irrigation strip is lower than that of its adjacent pixels on the left and right. By performing local valley detection, connectivity tracking, and skeleton extraction on the cooling rate step data, the spatial location of the drip irrigation strip can be obtained from the thermal infrared image itself without relying on field measurement data.

[0058] Step 4.1.1.1: Extract discrete pixels with local minima whose values ​​are less than their neighboring pixels to obtain a candidate valley pixel set, specifically including:

[0059] For each valid pixel in the cooling rate step data, obtain its eight direct neighboring pixels along the row and column directions in the two-dimensional pixel plane, forming a spatial eight-neighbor pixel group. Traverse each neighboring pixel in the spatial eight-neighbor pixel group: if the neighboring pixel is a cloud-covered invalid pixel marked in step 1.4, skip it and do not include it in the comparison; if it is a valid pixel, read its cooling rate step value. If all neighboring pixels of a central pixel are invalid pixels, the central pixel does not participate in the local valley value determination and is directly excluded. Compare the cooling rate step value of the central pixel with the cooling rate step values ​​of all valid neighboring pixels one by one:

[0060] If the value of the central pixel is less than the value of every valid neighboring pixel, then the central pixel is determined to be a local minimum discrete pixel and added to the candidate valley value pixel set; otherwise, it is not added. After completing the traversal, the candidate valley value pixel set is obtained.

[0061] Step 4.1.1.2: Perform connectivity tracing based on the spatial eight-neighbor relationship to obtain candidate connected components for the drip irrigation tape, specifically including:

[0062] For the candidate valley value cell set, initialize an empty set of visited markers and a list of connected components. Select an unmarked cell from the candidate valley value cell set as the seed cell, mark it as visited and add it to the current connected component. Check each of the seed cell's eight neighboring cells. If a neighboring cell simultaneously satisfies:

[0063] If a cell belongs to the candidate valley value set and has not yet been marked as visited, it is marked as visited and added to the current connected group. This newly added cell is then used as the new seed cell to continue checking its eight neighbors. The recursion continues until no new cells meeting the criteria are added to the boundary of the current connected group. This connected group is then added as a connected component to the connected component list. The next unmarked cell is selected, and connectivity tracing is performed repeatedly until all cells in the candidate valley value set are marked as visited. The resulting list of connected components constitutes the candidate connected components of the drip irrigation tape.

[0064] Step 4.1.1.3: Remove isolated noise connected components with an area smaller than the background noise area threshold to obtain the candidate skeleton pixel set for the drip irrigation tape, specifically including:

[0065] For each connected component in the candidate connected components of the drip irrigation tape, the total number of candidate valley value pixels contained therein is counted as the pixel area of ​​the connected component. A background noise area threshold is set, and the value of the background noise area threshold ranges from 10 pixel area to 50 pixel area: when the spatial resolution of the satellite image is 30 meters, the value is 10 to 20 pixel area; when the spatial resolution is 10 meters, the value is 30 to 50 pixel area. Connected components with a pixel area less than the threshold are identified as isolated noise connected components and removed. Connected components with a pixel area not less than the threshold are retained. The set of candidate valley value pixels in all retained connected components constitutes the candidate skeleton pixel set of the drip irrigation tape.

[0066] Step 4.1.1.4: Extract the central skeleton line of the connected components to obtain the spatial center line of the drip irrigation tape, specifically including:

[0067] For each connected component in the candidate skeleton cell set of drip irrigation tape, perform the following operations:

[0068] The sum of row coordinates and column coordinates of all candidate valley value pixels within the connected component is calculated. The ratio of the sum of row coordinates to the total number of pixels is taken as the average row coordinate, and the ratio of the sum of column coordinates to the total number of pixels is taken as the average column coordinate. This average value determines the overall geometric center position of the connected component. The offset of each candidate valley value pixel within the connected component relative to the above average value is calculated. For each candidate valley value pixel, the difference between its actual row coordinate and the average row coordinate is taken as the row offset of that pixel, and the difference between its actual column coordinate and the average column coordinate is taken as the column offset of that pixel. The row offset and column offset together constitute the two-dimensional deviation vector of that pixel. Based on the two-dimensional deviation vectors of all candidate valley value pixels, the dispersion of the connected component along the row and column directions is calculated. The squares of all row offsets are summed, and the ratio of the sum to the difference between the total number of pixels and one is taken as the row variance. Square the offsets of all columns and sum them up. The ratio of this sum to the difference between the total number of offsets and one is taken as the column direction variance. Multiply the row direction offset of each cell by its column direction offset and sum them up. The ratio of this sum to the difference between the total number of offsets and one is taken as the row and column direction covariance.

[0069] The row direction variance, column direction variance, and row-column direction covariance are arranged into a 2x2 matrix. The elements in the first row and first column of this matrix are the row direction variance values, the elements in the first row and second column and the first row and second column are the row-column direction covariance values, and the elements in the second row and second column are the column direction variance values. This matrix describes the overall spatial distribution of candidate valley pixels within the connected component on the two-dimensional pixel plane. The dominant extension direction of the spatial distribution of the connected component is extracted from this matrix. Two characteristic directions of this matrix are calculated, where the larger characteristic direction corresponds to the direction with the greatest dispersion of the candidate valley pixel projection, i.e., the longitudinal direction of the drip irrigation tape. The direction of this largest characteristic direction is denoted as the direction direction vector, and the row direction component and column direction component of this vector represent the unit extension component of the drip irrigation tape in the row direction and column direction, respectively.

[0070] Along the stated direction, with a step size of a single pixel width, the process gradually advances from one end of the connected domain to the other. At each step size position, a cross-sectional line perpendicular to the direction is drawn, and all candidate valley value pixels whose spatial positions intersect with the cross-sectional line within the connected domain are collected. The average row coordinates and average column coordinates of these candidate valley value pixels are calculated, and the two average values ​​are rounded to the nearest integer row and column number according to the rounding rules. This integer coordinate is taken as the center pixel of the cross-section at that cross-sectional position.

[0071] Connecting the center pixels of each cross-section along the direction sequentially yields the central skeleton line of the connected domain. Merging all the central skeleton lines of the connected domains forms the spatial centerline of the drip irrigation tape. This central skeleton line completely describes the actual spatial orientation of each drip irrigation tape in the satellite imagery, providing a spatial reference benchmark for subsequent calculations of the vertical distance of each pixel from the drip irrigation tape and for constructing lateral spatial attenuation profiles. Simultaneously, it can be directly imported into the BIM planning model of the irrigation area project as the spatial positioning baseline for the drip irrigation tape laying path.

[0072] Step 4.1.2: Based on the spatial centerline of the drip irrigation tape, extract the extension direction between adjacent center points, construct a lateral reference axis perpendicular to this extension direction, and obtain the lateral projection reference line, specifically including:

[0073] For each pair of adjacent center points on the centerline of the drip irrigation tape, process them one by one according to the extension order of the centerline. Take the first center point in each pair as the starting point and the second center point as the ending point. Calculate the difference between the row coordinates of the ending point and the row coordinates of the starting point, as the row direction extension component of the pair of adjacent center points; calculate the difference between the column coordinates of the ending point and the column coordinates of the starting point, as the column direction extension component of the pair of adjacent center points. These row and column direction extension components together constitute the local extension direction of the pair of adjacent center points.

[0074] Rotating the aforementioned local extension direction 90 degrees clockwise within the two-dimensional pixel plane yields the lateral reference axis corresponding to the pair of adjacent center points. The row direction extension component is negatively evaluated and used as the column direction component of the lateral reference axis, while the column direction extension component is directly used as the row direction component of the lateral reference axis. The direction of the lateral reference axis is always perpendicular to the local orientation of the drip irrigation tape at the pair of center points, pointing laterally to one side of the drip irrigation tape.

[0075] Check the consistency of the local extension direction between adjacent pairs of center points one by one. Calculate the angle between the local extension direction of the previous pair of center points and the local extension direction of the next pair of center points. If the angle exceeds a preset direction change threshold (ranging from 15 to 30 degrees), it indicates that there is a sharp turn in direction at the junction of the drip irrigation tape spatial centerline, which may be caused by random fluctuations or noise in local valley pixels. Smoothing is required. Taking the center point of the junction as the center position, take two center points forward and two center points backward along the centerline, for a total of five center points. Sum the row coordinates of these five center points and take the average as the corrected row coordinates; sum the column coordinates of these five center points and take the average as the corrected column coordinates. Replace the original coordinates of the junction center point with the corrected coordinates, and then recalculate the local extension direction and lateral reference axis of the two pairs of adjacent center points before and after the junction position based on the corrected coordinates.

[0076] After completing the above processing for all adjacent center point pairs along the center line of the drip irrigation tape, the lateral reference axes corresponding to each pair of center points are converged to form a transverse projection reference line.

[0077] Step 4.1.3: Based on the spatial coordinates of the lateral projection reference line and the center point of each pixel, calculate the vertical drop from the center point of each pixel along the lateral reference axis to the spatial center line of the drip irrigation tape, obtaining the lateral spatial distance sequence, specifically including:

[0078] To determine the nearest center point of a pixel on the centerline of the drip irrigation tape, iterate through all center points on the centerline and calculate the Euclidean distance between the pixel's center point and each other. Select the center point with the smallest distance from each of these distances as the nearest center point of the pixel, and determine the lateral reference axis for the location of this nearest center point. In the lateral projection reference line, find the lateral reference axis vector corresponding to the pair of adjacent center points of the nearest center point. Extract the row and column components of this lateral reference axis. If the nearest center point is the endpoint center point (i.e., the beginning or end) of the drip irrigation tape's centerline, select the lateral reference axis corresponding to that endpoint and its only adjacent center point.

[0079] Calculate the offset of the pixel's center point relative to the nearest center point. The difference between the row coordinates of the pixel and the row coordinates of the nearest center point is used as the row offset, and the difference between the column coordinates of the pixel and the column coordinates of the nearest center point is used as the column offset.

[0080] Calculate the vertical projection distance of the pixel along the lateral reference axis to the center line of the drip irrigation tape. Multiply the row direction offset by the row direction component of the lateral reference axis, and add the column direction offset multiplied by the column direction component of the lateral reference axis to obtain the projection length of the offset vector along the lateral reference axis. Take the absolute value of this projection length, and then normalize it by comparing it with the magnitude of the lateral reference axis. Add the squares of the row direction component and the squares of the column direction component of the lateral reference axis, and take the square root. The normalized vertical projection distance is the lateral vertical distance of the pixel from the drip irrigation tape.

[0081] The vertical distance of a pixel from the centerline of the drip irrigation strip is the distance along a lateral direction perpendicular to the local orientation of the strip. A smaller vertical distance indicates that the pixel is closer to the drip irrigation strip and the irrigation water is more abundant; a larger vertical distance indicates that the pixel is farther from the drip irrigation strip and the lateral infiltration effect of irrigation water is weaker. In mulched drip-irrigated cotton fields, the step change in the cooling rate of the soil moist body on both sides of the drip irrigation strip exhibits a nonlinear decay characteristic as the vertical distance increases: a slow initial drop followed by a sudden drop at the edge of the moist body. Therefore, this vertical distance provides a key spatial variable for subsequently constructing a lateral spatial attenuation profile and identifying the geometric boundary of lateral water spread. The vertical distance values ​​corresponding to all effective pixels within the target irrigation area are aggregated to form a lateral spatial distance sequence.

[0082] Step 4.2: Based on the horizontal spatial distance sequence, arrange the pixels in ascending order of vertical projection distance to obtain an ordered distance pixel set, specifically including:

[0083] From the horizontal spatial distance sequence, extract the vertical distance value and its row and column coordinates corresponding to each cell. Using the vertical distance value as the primary sorting key, the row coordinates as the first secondary sorting key, and the column coordinates as the second secondary sorting key, sort all cells in ascending order of vertical distance value. After sorting, the cell closest to the drip irrigation tape is placed at the beginning, and the cell farthest is placed at the end. The resulting sequence is the ordered distance cell set.

[0084] Step 4.3: Based on the ordered distance pixel set, extract the cooling rate step corresponding to each pixel to obtain the profile attenuation sequence, specifically including:

[0085] For each pixel in the ordered distance pixel set, find the corresponding cooling rate change sequence in the cooling rate step data based on its row and column coordinates. Extract the cooling rate change value corresponding to the most recent valid observation date in this sequence as the representative step value for that pixel.

[0086] If a pixel has no valid value in the cooling rate change sequence (e.g., no record for that date due to cloud cover or missing data), the median of the cooling rate change for that pixel across all valid observation dates is used as the representative value of the step change. The cooling rate change values ​​for that pixel on each valid observation date are arranged in ascending order. If the total number of valid observation dates is odd, the value in the middle position after sorting is taken as the median; if the total number of valid observation dates is even, the average of the two values ​​in the middle position after sorting is taken as the median. This process is repeated for all ordered distance pixels, and the representative value of the step change for each pixel is arranged in order of the pixel's position number in the ordered distance pixel set, forming a one-dimensional sequence of cooling rate step changes distributed from near to far along the vertical distance direction.

[0087] In this sequence, the value at the beginning corresponds to the step value of the pixel closest to the drip irrigation strip, and the value at the end corresponds to the step value of the pixel farthest from the drip irrigation strip. This sequence fully depicts the spatial variation of the cooling rate step value as it decreases with increasing distance from the drip irrigation strip outwards. Inside the irrigated wet zone (near the drip irrigation strip), the step value is higher and decays slowly; at the boundary zone transitioning from the wet zone to the dry bare soil zone, the step value changes from a slow drop to a sudden drop; after entering the dry bare soil zone, the step value drops to a low value and tends to stabilize.

[0088] The above one-dimensional sequence is used as a profile decay sequence and output to step 5.

[0089] In a preferred embodiment of the present invention, step 5 above, based on the profile attenuation sequence, determines the spatial coordinates of the abrupt change in the sequence slope to obtain the geometric boundary line of water lateral expansion; based on the geometric boundary line of water lateral expansion, a set of continuous high-step pixels inside the boundary line is extracted to obtain pixels active in water evaporation, including:

[0090] Step 5.1: Based on the profile attenuation sequence, calculate the difference in cooling rate step between adjacent ordered pixels to obtain the profile spatial change sequence; based on the profile spatial change sequence, calculate the difference in spatial change between adjacent profiles to obtain the profile spatial change acceleration sequence, specifically including:

[0091] For each pair of adjacent positions in the profile attenuation sequence, the difference between the representative value of the step amount of the later pixel and the representative value of the step amount of the earlier pixel is calculated as the profile spatial variation between the pair of adjacent positions. This profile spatial variation reflects the magnitude of the cooling rate step between two adjacent ordered pixel positions along the direction perpendicular to the drip irrigation belt. A negative difference indicates that the step amount decreases in that segment away from the drip irrigation belt, while a positive difference indicates that the step amount increases in that segment. The profile spatial variation values ​​corresponding to all adjacent positions are arranged sequentially in the direction away from the drip irrigation belt to form a profile spatial variation sequence.

[0092] For each pair of adjacent positions in the sequence of spatial changes in the profile, the difference between the spatial change that follows and the spatial change that precedes is calculated, serving as a measure of acceleration of spatial change between the adjacent positions. This acceleration measure describes the degree of acceleration in the rate of drop of the cooling rate step in the lateral space: a positive acceleration measure indicates that the downward trend of the step is decelerating in this segment, i.e., the attenuation profile is becoming gentler, corresponding to the irrigated and moist area near the drip irrigation belt, where the soil moisture content is sufficient, the cooling rate step remains at a high level and only decreases slowly with increasing distance; a negative acceleration measure indicates that the downward trend of the step is accelerating in this segment, i.e., the attenuation profile is becoming steeper, corresponding to the physical transition zone from the irrigated and moist area to the dry bare soil area, where the soil moisture content drops sharply, and the cooling rate step falls rapidly with increasing distance.

[0093] From a numerical perspective, the acceleration measure of spatial change in the profile is equivalent to the second-order change of the profile decay sequence at three consecutive adjacent positions. That is, the representative value of the step amount at the farthest position is subtracted by twice the representative value of the step amount at the middle position, and the representative value of the step amount at the nearest position is added. The rate of change of the first-order decay rate of the step amount along the spatial distance direction is extracted, and the second-order change is used to sensitively capture the inflection point of the decay profile from a slow descent to a sharp descent.

[0094] All acceleration measures of spatial profile change are arranged sequentially in the direction away from the drip irrigation tape to form a sequence of spatial profile change accelerations. The positive value segments of this sequence correspond to the interior of the irrigated wet zone, and the negative value segments correspond to the transition zone from the wet zone to the dry zone. The zero-crossing position where the positive value turns into the negative value is the candidate position for determining the geometric boundary of water lateral expansion in subsequent steps.

[0095] Step 5.2: Based on the spatial change acceleration sequence of the profile, extract the pixels corresponding to the zero-crossing positions where the values ​​change from positive to negative, and obtain the spatial coordinates of the abrupt changes in the sequence slope. Specifically, this includes:

[0096] For the spatial acceleration sequence of the profile, scan position by position starting from the beginning of the sequence. At each position, check whether the acceleration metric at that position and the immediately following position satisfy the sign change from positive to negative. That is, the acceleration metric at the current position is positive and the acceleration metric at the immediately following position is negative. This sign reversal means that between this pair of positions, the decay behavior of the cooling rate step change has undergone a qualitative change from deceleration to acceleration. The former corresponds to the characteristic of sufficient moisture in the irrigated wet area and the cooling effect slowly decaying with distance, while the latter corresponds to the characteristic of a sudden decrease in moisture and a sharp decay of the cooling effect when the wet area transitions to the dry bare soil area.

[0097] If multiple position pairs in the acceleration sequence satisfy the aforementioned sign reversal condition (positive to negative), the optimal abrupt change position that best reflects the true physical boundary is selected according to the following priority rule: among all position pairs satisfying the condition, the position pair with the largest absolute value of the positive acceleration metric is selected first, while also requiring that the negative acceleration metric of this position pair has the largest absolute value among all position pairs satisfying the condition. The profile segment corresponding to the position pair satisfying these two conditions experiences a dramatic transition from maximum deceleration to maximum acceleration, and this transition process most strongly reflects the physical boundary position of the transition from the drip-irrigated wet zone to the dry bare soil zone.

[0098] If there are no position pairs in the acceleration sequence that satisfy the above sign reversal condition, it indicates that there is no obvious attenuation mode transition at this profile location. There are two possible reasons:

[0099] First, the entire profile is located in a dry, bare soil area, and the cooling rate step change remains consistently low without significant fluctuations. Second, the entire profile is located within an irrigated, moist area, and the cooling rate step change gradually and uniformly decreases from near the drip irrigation tape to far away from it, without any abrupt drops. In either case, it is determined that there is no geometric boundary for lateral water expansion at the location of this profile, and this profile is not involved in the subsequent determination of the boundary.

[0100] For a selected pair of positions that meet the zero-crossing condition, the abrupt change position needs to be further determined between the two profile positions crossed by the pair. The absolute value of the positive acceleration metric corresponding to the current position is compared with the absolute value of the negative acceleration metric corresponding to the next position: if the absolute value of the positive acceleration metric is less than or equal to the absolute value of the negative acceleration metric, it indicates that the acceleration zero-crossing point is numerically closer to the side of the profile position that comes first in the pair, and the pixel with the ordered distance pixel cluster that comes first is taken as the abrupt change position pixel; if the absolute value of the positive acceleration metric is greater than the absolute value of the negative acceleration metric, it indicates that the acceleration zero-crossing point is numerically closer to the side of the profile position that comes second in the pair, and the pixel with the ordered distance pixel cluster that comes second is taken as the abrupt change position pixel.

[0101] The zero-crossing detection is performed independently at each longitudinal profile position along the direction of the drip irrigation tape on both sides. If the abrupt change pixels identified at adjacent longitudinal profile positions on the same side of the drip irrigation tape show a slight misalignment in the row coordinates, i.e., the difference between the row coordinates of the abrupt change pixels identified at two adjacent profile positions exceeds the width of a single pixel, then the row coordinates of the abrupt change pixels on that side are uniformly smoothed and corrected. The row coordinate values ​​of all successfully identified valid abrupt change pixels on that side are collected, and the median of these row coordinates is calculated as the reference row coordinate. The row coordinates of each abrupt change pixel are compared with the reference row coordinates. If the deviation of the row coordinates of a certain abrupt change pixel from the reference row coordinates exceeds a preset tolerance range (the tolerance range is 1 to 3 pixel widths), then the row coordinates of that abrupt change pixel are corrected to the reference row coordinate value.

[0102] After the above corrections are completed, the row and column coordinates of each aberration pixel in the original two-dimensional pixel plane together constitute the spatial coordinates of the slope aberration position of that side sequence. These spatial coordinates define the endpoint of the lateral diffusion of water in the drip irrigation wetting body in the direction perpendicular to the drip irrigation tape, providing spatial anchor points for connecting each aberration pixel sequentially along the drip irrigation tape in step 5.3 to delineate the geometric boundary of water lateral diffusion. These anchor points, together with the extracted spatial centerline of the drip irrigation tape, jointly define the complete spatial range of the irrigated wetting area on both sides of each drip irrigation tape.

[0103] Step 5.3: Based on the spatial coordinates of the abrupt change in sequence slope, determine the geometric boundary line of water lateral expansion on both sides of the drip irrigation tape to obtain the lateral expansion boundary range; based on the lateral expansion boundary range, extract continuous pixels within the range where the cooling rate step is greater than zero to obtain the inner candidate pixel set, specifically including:

[0104] For each longitudinal section along the direction of the drip irrigation tape, the spatial coordinates of the abrupt change in slope are determined. These coordinate points are then connected sequentially along the direction of the drip irrigation tape to form the geometric boundary line of water lateral expansion on that side. The boundary line on the left side of the drip irrigation tape and the boundary line on the right side of the drip irrigation tape together enclose a strip-shaped closed area extending longitudinally along the drip irrigation tape on the two-dimensional pixel plane. This strip-shaped area is used as the lateral expansion boundary range.

[0105] For each valid pixel within the lateral expansion boundary, its representative value of the cooling rate step is extracted. To eliminate temperature benchmark differences caused by different soil textures, an adaptive background threshold is used instead of a fixed zero threshold: the representative values ​​of the cooling rate step of all valid pixels within the lateral expansion boundary are sorted from smallest to largest, and the value at the 20th percentile after sorting is taken as the adaptive background threshold. Pixels within the lateral expansion boundary whose representative value of the cooling rate step is greater than the adaptive background threshold are selected to obtain the inner candidate pixel set.

[0106] Step 5.4: Based on the inner candidate pixel set, remove edge pixels affected by lateral heat conduction from the damaged edge of the plastic film to obtain pixels with active water evaporation, specifically including:

[0107] Within the inner candidate pixel set, there may be pseudo-irrigated pixels originating from the edge of a damaged mulch film. At the damaged area, the exposed soil directly exchanges water vapor with the atmosphere due to the lack of mulch film obstruction. The evaporative cooling effect exhibits a positive value in the step rate of cooling, similar to that of the irrigated wet area. However, the heat transfer direction at the damaged edge is disordered and divergent, while in the true irrigated wet area, the heat transfer direction decreases orderly from the drip irrigation tape outwards. The two can be distinguished by analyzing the consistency of the spatial variation direction of the step rate of cooling within the local neighborhood of each candidate pixel.

[0108] Step 5.4.1: Calculate the spatial drop direction of the cooling rate step in the two-dimensional neighborhood of each pixel to obtain the heat flow direction in the pixel space, specifically including:

[0109] For each target pixel in the inner candidate pixel set, determine its spatial neighborhood on the two-dimensional pixel plane consisting of four directly adjacent pixels, namely the left adjacent pixel, the right adjacent pixel, the upper adjacent pixel, and the lower adjacent pixel. Check whether each of these four adjacent pixels belongs to the inner candidate pixel set. If it does not belong, skip it and do not include it in the calculation. If it does belong, retain it.

[0110] The lateral spatial drop is calculated, reflecting the spatial gradient of the cooling rate step along the pixel row direction (i.e., perpendicular to the column direction). If the target pixel has both left and right adjacent pixels, the lateral drop is calculated by subtracting the representative value of the cooling rate step of the left adjacent pixel from the representative value of the cooling rate step of the right adjacent pixel. If the target pixel has only right adjacent pixels and no left adjacent pixels, the difference between the target pixel's own representative value and the representative value of the right adjacent pixel is multiplied by two to obtain an approximate value of the lateral spatial drop. If the target pixel has only left adjacent pixels and no right adjacent pixels, the difference between the target pixel's own representative value and the representative value of the left adjacent pixel (i.e., the target pixel value minus the left pixel value) is multiplied by two. If there are no effective adjacent pixels on either the left or right, the lateral spatial drop is set to zero.

[0111] The vertical spatial drop is calculated, reflecting the spatial gradient of the cooling rate step along the column direction (i.e., perpendicular to the row direction). If a target pixel has both an upper and lower adjacent pixel, the vertical drop is calculated by subtracting the representative value of the cooling rate step of the upper adjacent pixel from the representative value of the cooling rate step of the lower adjacent pixel. The approximation rules for unilateral or bilateral missing pixels are the same as for the lateral pixels.

[0112] The lateral spatial drop is taken as the row direction component of the two-dimensional spatial pointing vector, and the vertical spatial drop is taken as the column direction component. Together, they constitute the spatial heat flow direction vector of the target pixel. The direction of this vector is the direction in which the cooling rate step at the target pixel increases the most rapidly in the two-dimensional spatial neighborhood. In the irrigated wet area, the cooling rate step gradually decreases from the highest value directly below the drip irrigation strip to the outside. Therefore, the heat flow direction of each pixel approximately points away from the drip irrigation strip, that is, the direction of heat transfer from the high-temperature area near the drip irrigation strip to the low-temperature area on the outside. If both the row direction component and the column direction component are zero, it indicates that there is no directional step spatial change in the local neighborhood of the pixel.

[0113] Step 5.4.2 calculates the standard deviation of the angle between the spatial heat flow direction vectors of each pixel in the local neighborhood to obtain the heat flow direction dispersion, specifically including:

[0114] For each target pixel with an effective spatial heat flow direction in the inner candidate pixel set, a demarcation area of ​​size centered on that pixel is defined. A local square window of one pixel. The rules for determining the value are as follows:

[0115] When the spatial resolution of satellite imagery is 30 meters When the spatial resolution of satellite imagery is 10 meters , The number must be odd to ensure the target pixel is centered within the window. Collect the heat flow direction vectors of all pixels within the window that belong to the inner candidate pixel set and have a valid spatial heat flow direction. Let the total number of collected vectors be... Let there be vectors, denoted as . .like If the directional information within the window is insufficient to assess consistency, then the heat flow direction dispersion of the target pixel is... Set it directly to zero.

[0116] and Calculate the included angle pair by pair Output The dimension is radians, and the calculation is performed. The angle between each pair of vectors (common mean of (pairs) The standard deviation of the included angle is calculated as the heat flow direction dispersion of the target pixel. :

[0117] ;

[0118] The heat flow direction dispersion of the target pixel is defined as follows: the cooling rate step of each pixel in the irrigated wet area decreases uniformly from the drip irrigation tape outwards; the heat flow direction of each pixel approximately points in the same direction away from the drip irrigation tape, with concentrated angles. Smaller; at the edge of the damaged plastic film, heat conduction is disorganized, and the heat flow direction of each pixel is chaotic and dispersed, with discrete angles. Relatively large.

[0119] Step 5.4.3: Based on the heat flow direction dispersion, remove pixels from the inner candidate pixel set whose heat flow direction dispersion is greater than a preset dispersion threshold to obtain pixels with active water evaporation, specifically including:

[0120] For each of the inner candidate cell sets that has been calculated The target pixel will With preset dispersion threshold Compare them. The value ranges from 30 to 60 degrees. When the cotton mulch film is relatively intact and has few breaks during the spring irrigation period, a smaller value (30 to 40 degrees) is used to more strictly exclude a small number of damaged artifacts. When the mulch film is widely aged and damaged, a larger value (45 to 60 degrees) is used to retain sufficient pixels in the irrigated area. Comparison rules:

[0121] like If the heat flow direction of the pixel is determined to be disordered, the pixel is marked as a false pixel at the edge of the damaged mulch film and removed from the inner candidate pixel set; if If the heat flow direction of the target pixel is consistent with the interior of the irrigated and moistened area, it will be retained; if the target pixel's heat flow direction is consistent with the interior of the irrigated and moistened area, it will be retained. In step 5.4.2, because If set to zero, then The retention criteria are automatically met. After traversing all pixels in the inner candidate pixel set and completing the above judgment and elimination, the retained pixels constitute the set of pixels with active water evaporation. Each pixel in this set simultaneously meets three conditions: it is located inside the geometric boundary of water lateral expansion, the representative value of the cooling rate step is higher than the adaptive background threshold, and the spatial heat flow direction is consistent with the heat diffusion direction inside the irrigated wet area. The pixels with active water evaporation are used as the final output of step 5 and passed to step 6 to accumulate the cooling rate step segment by segment along the drip irrigation tape.

[0122] In a preferred embodiment of the present invention, step 6 above, which involves progressively accumulating the cooling rate step along the drip irrigation tape based on the active water evaporation pixels to obtain the cumulative heat dissipation in the strip, includes:

[0123] Step 6.1: Based on the active pixels of water evaporation, extract the longitudinal projection coordinates of each active pixel along the spatial centerline of the drip irrigation tape to obtain the longitudinal position sequence of the active pixels, specifically including:

[0124] Determine the nearest center point of the active pixel on the center line of the drip irrigation tape space. Traverse all center points on the center line and calculate the Euclidean straight-line distance between the active pixel and each center point to determine the local extension direction of the nearest center point. Find the local extension direction corresponding to the pair of adjacent center points on the center line of the drip irrigation tape space where the nearest center point is located. This direction is calculated in step 4.1.2. The row direction component of the local extension direction is the row coordinate of the endpoint of the pair of adjacent center points minus the row coordinate of the starting point, and the column direction component is the column coordinate of the endpoint minus the column coordinate of the starting point. If the nearest center point is the first endpoint of the center line of the drip irrigation tape space, then take the direction corresponding to the first endpoint and its unique next center point as the local extension direction; if it is the last endpoint, then take the direction corresponding to the last endpoint and its unique previous center point.

[0125] Calculate the offset of the active pixel relative to the nearest center point. The difference between the row coordinates of the active pixel and the row coordinates of the nearest center point is taken as the row offset, and the difference between the column coordinates of the active pixel and the column coordinates of the nearest center point is taken as the column offset. Calculate the projection length of this offset along the local extension direction, which is taken as the longitudinal projection displacement of the active pixel at the nearest center point. Multiply the row offset by the row component of the local extension direction, and add the column offset multiplied by the column component of the local extension direction to obtain the dot product of the offset vector and the extension direction. Ratio this dot product to the magnitude of the local extension direction vector yields the longitudinal projection displacement, in pixels. The physical meaning of this longitudinal projection displacement is:

[0126] The degree of longitudinal offset of an active pixel relative to its nearest center point along the direction of the drip irrigation tape is indicated by a positive displacement, which means the pixel is downstream of the nearest center point, i.e., in front of the drip irrigation tape; and a negative displacement, which means the pixel is upstream of the nearest center point, i.e., behind the drip irrigation tape.

[0127] The total mileage of the active pixel along the centerline of the drip irrigation tape from the starting point to the current longitudinal position is determined. The centerline of the drip irrigation tape is formed by connecting several center points sequentially along the direction. Starting from the first center point at the starting point, the straight-line distances between adjacent center points are accumulated segment by segment along the centerline. The cumulative mileage from the starting point to the nearest center point is equal to the sum of the straight-line distances between all adjacent center points from the first center point to the nearest center point. The total longitudinal projected coordinates of the active pixel along the drip irrigation tape direction are calculated. The cumulative mileage from the starting point to the nearest center point is added to the longitudinal projected displacement of the active pixel at the nearest center point. The result is the total longitudinal projected coordinates of the active pixel along the drip irrigation tape direction. The above longitudinal projected coordinate calculation is performed on each active pixel in the set of active pixels for water evaporation. The longitudinal projected coordinate values ​​corresponding to all active pixels are collected to form a sequence of active pixel longitudinal positions. This sequence completely describes the spatial distribution of all active pixels for water evaporation along the drip irrigation tape direction.

[0128] Step 6.2: Based on the longitudinal position sequence of active pixels, the drip irrigation tape is divided into continuously accumulating segments according to a fixed spatial step size, resulting in a set of segmented intervals, specifically including:

[0129] The active pixel longitudinal position sequence provides the one-dimensional longitudinal coordinates of each active pixel along the drip irrigation belt. To group discrete pixels by spatial location for segment-by-segment aggregation, it is necessary to divide the drip irrigation belt direction into sub-intervals with a fixed step size. The fixed spatial step size is physically defined as the length of adjacent cumulative segments, and its value ranges from the ground resolution of a single pixel to 10 times the ground resolution of the pixel. When the satellite image spatial resolution is 30 meters, the step size ranges from 30 meters to 150 meters; when the spatial resolution is 10 meters, the step size ranges from 10 meters to 50 meters, i.e., 1 to 5 pixel widths. A smaller step size results in finer segments and higher accuracy in calculating the irrigated area along the drip irrigation belt direction, but the increased number of segments leads to a greater computational burden. A larger step size results in coarser segments and a corresponding decrease in calculation accuracy, but better robustness to local noise and individual pixel anomalies. In practical irrigation projects, the selection of step length needs to take into account the degree of fragmentation of the irrigation area, the uniformity of drip irrigation tape laying, and the spatial accuracy requirements of the BIM planning model for irrigation area data.

[0130] When determining the segmentation scheme, first find the minimum and maximum values ​​of all longitudinal projected coordinates from the active pixel longitudinal position sequence, and use them as the starting and ending points of the drip irrigation tape along the longitudinal axis, respectively. The difference between the starting and ending points is the total effective length of the drip irrigation tape carrying the effective irrigation signal upwards.

[0131] Following the drip irrigation tape, starting from the beginning, the tape advances segment by segment at a predetermined fixed step length, dividing the total effective length into several continuous intervals. The number of intervals is equal to the ratio of the total effective length to the step length. If this ratio is not an integer, it is rounded up to ensure that the last interval covers the end. The first interval starts from the beginning and covers the length of one step length; the second interval follows the first interval and covers the length of another step length; and so on. At the last interval, if the tape extends beyond the end after advancing by the step length, the upper boundary of the interval is truncated to the end position to ensure that the last interval does not exceed the effective range. All intervals are then combined to form a strip-shaped set of segmented intervals.

[0132] Step 6.3: Based on the strip-shaped segmented interval set, filter the active pixels for water evaporation that fall within each segment to obtain a subset of active pixels for each segment, specifically including:

[0133] For each segment interval in the strip segment interval set, traverse each active pixel in the active pixel vertical position sequence and check whether the vertical projection coordinates of the active pixel fall within the vertical coordinate range of the current segment interval. The determination rule is as follows:

[0134] The longitudinal projection coordinates of active pixels must be greater than or equal to the lower bound of the segment interval and less than the upper bound of the segment interval. For the last segment interval, since the upper bound has been truncated to the end position according to the rule in step 6.2, the judgment condition is adjusted to: the longitudinal projection coordinates of active pixels are greater than or equal to the lower bound of the segment interval and less than or equal to the end position.

[0135] Active pixels that meet the above-mentioned longitudinal coordinate range determination conditions are assigned to the active pixel subset corresponding to the current segment interval. After traversal, each segment interval corresponds to an active pixel subset, which contains all active pixels for water evaporation whose spatial location falls within the range of that segment along the drip irrigation tape. If the longitudinal projection coordinates of no active pixels in a certain segment interval fall within its range, then the active pixel subset corresponding to that segment interval is an empty set. The appearance of an empty set usually means that no effective active pixels for water evaporation were detected at that drip irrigation tape segment due to insufficient irrigation water penetration, intact and undamaged mulch film, or cloud cover during satellite transit.

[0136] Step 6.4: Based on the active pixel subsets of each segment, cumulatively add the cooling rate step of all active water evaporation pixels in the corresponding segment to obtain the cumulative heat dissipation of the band, specifically including:

[0137] For each active pixel subset, process them segment by segment. For each non-empty active pixel subset, count the total number of active pixels in the subset, and sum the representative values ​​of the cooling rate step of all active pixels in the subset. The sum is the cumulative heat dissipation of the corresponding segment. If the active pixel subset of a certain segment is empty, the cumulative heat dissipation of the segment is directly recorded as zero, indicating that there is no heat dissipation signal contributed by the evaporative cooling effect of irrigation water at that drip irrigation strip segment.

[0138] Within a certain number of drip irrigation strip segments, this is the sum of the cooling rate steps of all pixels identified as irrigated wet areas. This value comprehensively reflects the overall intensity of the afternoon cooling effect caused by irrigation at that segment. The larger the cumulative heat dissipation in the strip, the more sufficient the irrigation water, the more obvious the evaporative cooling effect, and the better the irrigation effect at that segment; the smaller or zero the cumulative heat dissipation in the strip, the less water penetration or the absence of effective irrigation at that segment.

[0139] The cumulative heat dissipation of all drip irrigation tape segments is arranged sequentially according to the segmentation order along the drip irrigation tape direction to form the cumulative heat dissipation of the strip.

[0140] like Figure 3 As shown, in a preferred embodiment of the present invention, step 7 above calculates the pixel irrigation water infiltration probability based on the cumulative amount of strip heat dissipation and shortwave infrared reflectivity to obtain the irrigation status; based on the irrigation status, the actual surface coverage area corresponding to the total number of irrigated pixels is calculated to obtain the irrigated area of ​​the irrigation district, including:

[0141] Step 7.1: Based on the cumulative amount of heat dissipation in the band, extract the thermal response intensity level corresponding to each pixel to obtain the heat dissipation feature sequence; based on the short-wave infrared reflectance, extract the spectral absorption depth corresponding to each pixel to obtain the moisture absorption feature sequence, specifically including:

[0142] The cumulative heat dissipation of the strip is obtained by aggregating the drip irrigation strip segments. All active water evaporation pixels in the same segment share the same cumulative heat dissipation value, which reflects the total irrigation cooling intensity along the direction of the drip irrigation strip in that segment. This value is then assigned back to each active pixel in that segment to form the heat response intensity level of each pixel.

[0143] For each drip irrigation tape segment, if the segment has been marked as a segment without irrigation signal (i.e., the cumulative heat dissipation of the strip is zero and the active pixel subset is empty), then the segment is skipped, and its pixels do not participate in subsequent heat dissipation feature extraction. The segment is directly determined as having no irrigation, and there is no need to calculate the irrigation water infiltration probability. If the segment has a valid irrigation signal (i.e., the cumulative heat dissipation of the strip is greater than zero or the active pixel subset is not empty), then the thermal response intensity level of each active water evaporation pixel within the segment is assigned the value of the cumulative heat dissipation of the strip for that segment. The physical meaning of the thermal response intensity level is:

[0144] The change in the total cooling rate of the drip irrigation strip segment where the pixel is located due to the evaporative cooling effect of irrigation. The higher the thermal response intensity level, the more abundant the irrigation water and the stronger the irrigation cooling signal at the location of the pixel; the lower the thermal response intensity level, the weaker the irrigation cooling signal at the location of the pixel.

[0145] For each pixel with active water evaporation, the raw reflectance value in the shortwave infrared band is obtained from the images output in step 1.5 based on its row and column coordinates. The raw reflectance is dimensionless and typically ranges from 0 to 1. Its physical meaning is the proportion of solar radiation reflected by the Earth's surface in the shortwave infrared band.

[0146] Dry bare soil has high short-wave infrared reflectance, while moist soil has low reflectance due to the absorption of short-wave infrared radiation by moisture. To convert short-wave infrared reflectance into spectral absorption depth, which is positively correlated with soil moisture content (i.e., the more moisture, the greater the absorption depth), the original reflectance value is subtracted from 1. The difference is the spectral absorption depth of that pixel. This spectral absorption depth is dimensionless and ranges from zero to one: a spectral absorption depth close to one indicates high soil moisture content, where a large amount of short-wave infrared radiation is absorbed and reflectance is low; a spectral absorption depth close to zero indicates dry soil, where short-wave infrared reflectance is high and absorption is low. For areas with intact plastic film coverage, the high reflectivity of the film itself causes the spectral absorption depth to approach zero, which is consistent with physical reality, as the moisture state of the soil beneath the film cannot be directly perceived by the short-wave infrared band. However, the spectral absorption depth of exposed soil at locations with damaged film and soil in uncovered areas accurately reflects soil moisture content.

[0147] For each active pixel in the set of active pixels for water evaporation, its thermal response intensity level and spectral absorption depth are calculated according to the above rules. The thermal response intensity levels of all active pixels are arranged in sequence according to the pixel number to form a thermal dissipation characteristic sequence. The spectral absorption depths of all active pixels are arranged in the same sequence to form a water absorption characteristic sequence.

[0148] Step 7.2: Based on the heat dissipation characteristic sequence and the water absorption characteristic sequence, the heat response intensity level and the spectral absorption depth are scaled proportionally to a unified numerical range. The product of the scaled values ​​within the same pixel is calculated as the spatial matching degree to obtain the pixel irrigation water infiltration probability. Specifically, this includes:

[0149] The thermal response intensity levels of all active pixels in the heat dissipation feature sequence are normalized and scaled to uniformly map them to the range of 0 to 1. The minimum and maximum values ​​are found from all thermal response intensity levels. For each pixel, the difference between its thermal response intensity level and the minimum value is calculated, and then compared with the difference between the maximum and minimum values. The result is the scaled thermal response intensity value. If the thermal response intensity levels of all active pixels in the entire irrigation area are exactly the same (maximum value equals minimum value), it indicates that there is no difference in the heat signal among the pixels. In this case, the scaled thermal response intensity values ​​of all pixels are uniformly assigned to 0.5, degenerating into equal-weight matching. The scaled thermal response intensity values ​​are dimensionless. A value of 1 corresponds to the pixel containing the drip irrigation strip segment with the strongest cooling heat dissipation in the entire irrigation area, i.e., the location of the most sufficient irrigation; a value of zero corresponds to the pixel containing the weakest cooling heat dissipation but still a positive value in the entire irrigation area, i.e., the location at the edge of irrigation.

[0150] The same normalization scaling operation is performed on the spectral absorption depth values ​​of all active pixels in the water absorption feature sequence, mapping them to a value range of zero to one. The scaled spectral absorption depth values ​​are dimensionless, with a value of one corresponding to the pixel with the strongest spectral absorption and the highest soil moisture content in the entire irrigation area; a value of zero corresponds to the pixel with the weakest spectral absorption, the driest soil, or complete mulch coverage in the entire irrigation area. If the spectral absorption depth values ​​are exactly the same throughout the entire irrigation area, it is also degraded to an equal-weighted matching of 0.5.

[0151] For each pixel with active water evaporation, its scaled thermal response intensity value is multiplied by its scaled spectral absorption depth value. The product is used as the spatial matching degree of that pixel. This product is dimensionless and ranges from zero to one. The product represents the degree of consistency between the cooling heat dissipation evidence from the thermal infrared channel and the soil moisture absorption evidence from the short-wave infrared channel at the spatial location of that pixel. This product is directly used as the irrigation water infiltration probability of that pixel. The automatic discrimination effect of the product coupling mechanism on three typical cases is as follows:

[0152] When the thermal response intensity value is high and the spectral absorption depth value is also high, the product result is large. The thermal infrared channel indicates that the evaporative cooling effect is strong at this location, and the short-wave infrared channel indicates that the soil moisture content is high at this location. The evidence from the two independent physical channels consistently points to the irrigated wet area. The irrigation water penetration probability of this pixel is high, and it is correctly identified as an irrigated pixel.

[0153] When the thermal response intensity value is high but the spectral absorption depth value is low, the product result is small. The thermal infrared channel detects a cooling signal, but the shortwave infrared channel does not detect the corresponding moisture absorption characteristics. This situation often occurs in areas with complete mulch film coverage. The high reflectivity of the mulch film reduces the dynamic range of shortwave infrared reflectivity, causing the spectral absorption depth to approach zero. Even if the soil under the film is actually moist, the shortwave infrared channel cannot detect it. The evidence from the two channels is inconsistent. Product coupling automatically suppresses the irrigation water penetration probability of this pixel to a low level, avoiding false irrigation misjudgments caused by relying solely on the thermal infrared channel.

[0154] When the thermal response intensity value is low but the spectral absorption depth value is high, the product result is also small. The shortwave infrared channel detects soil moisture absorption characteristics, but the thermal infrared channel does not detect the corresponding cooling heat dissipation effect. This situation often occurs in temporarily moist areas after rainfall. Rainwater increases the surface soil moisture content and decreases the shortwave infrared reflectance. However, the evaporative cooling effect of rainfall itself and the continuous evaporative cooling effect of irrigation water show different cooling rate steps. The thermal infrared channel cannot identify it as an irrigation signal. The evidence from the two channels is inconsistent. The product coupling also automatically suppresses the irrigation water infiltration probability to a low level, avoiding false irrigation misjudgments caused by relying solely on the shortwave infrared channel.

[0155] Step 7.3: Based on the pixel irrigation water infiltration probability, filter pixels with probability values ​​greater than the infiltration probability threshold to obtain the set of pixels in the irrigation state, specifically including:

[0156] The irrigation water infiltration probability reflects the confidence level that each pixel with active water evaporation belongs to a truly irrigated wet zone. To convert continuous confidence levels into a binary irrigation status determination—that is, whether the pixel is irrigated or non-irrigated—an infiltration probability threshold is set as the determination threshold. The infiltration probability threshold ranges from 0.2 to 0.5, and its selection follows these principles:

[0157] When the physical data quality of both the thermal infrared and shortwave infrared channels is high (i.e., the image is unobstructed by clouds, the mulch film is intact, and the atmospheric correction effect is good), a lower threshold value is adopted, ranging from 0.2 to 0.3, to reduce the false positive rate of irrigation pixels and ensure the integrity of the irrigation coverage area. When the data quality is average (i.e., there are thin cloud remnants, significant mulch film damage, or large atmospheric correction residuals in some areas), a higher threshold value is adopted, ranging from 0.4 to 0.5, to reduce the false positive rate of spurious irrigation pixels and ensure the accuracy of irrigation area extraction. In actual irrigation area engineering applications, this threshold value can be adaptively adjusted by the irrigation area management unit based on the satellite image quality assessment results of the current season and the historical irrigation area verification data.

[0158] For each pixel with active water evaporation, its irrigation water infiltration probability is compared with the infiltration probability threshold. If the irrigation water infiltration probability is greater than the infiltration probability threshold, the pixel is determined to be in an irrigated state and is included in the set of irrigated pixels. If the irrigation water infiltration probability does not exceed the infiltration probability threshold, the pixel is determined not to be in an irrigated state and is not included in the set of irrigated pixels.

[0159] After traversing all pixels with active water evaporation and completing the above screening, a set of irrigated state pixels is obtained. Each pixel in this set satisfies all four progressive judgment conditions: the pixel is located inside the geometric boundary of the water lateral expansion, i.e., within the spatial boundary of the drip irrigation wet body; the cooling rate step value of the pixel is positive, indicating the presence of an afternoon cooling effect caused by irrigation evaporation cooling; the spatial heat flow direction of the pixel coincides with the consistent heat dissipation direction within the irrigated wet area, excluding pseudo-irrigated pixels at the edge of mulch film damage; the dual-channel product coupling probability of the pixel, i.e., the matching degree between thermal infrared cooling evidence and short-wave infrared water absorption evidence, exceeds the infiltration probability critical value, passing the dual-source cross-validation of thermal and water signals. After the above four conditions are progressively screened, the pixels in the set of irrigated state pixels are considered as irrigated pixels in the target irrigation area that are actually affected by the evaporation cooling effect of irrigation water in this irrigation season.

[0160] Step 7.4: Based on the set of irrigated state pixels, sum the number of irrigated state pixels within the set to obtain the total number of target irrigated pixels; based on the total number of target irrigated pixels, calculate the actual land cover area corresponding to the total number of target irrigated pixels to obtain the irrigated area of ​​the irrigation district, specifically including:

[0161] The set of irrigation status pixels contains all pixels that are determined to be in an irrigation state. The total number of pixels in this set is the total number of effective irrigation pixels in the target irrigation area that are determined to be affected by the evaporation and cooling effect of irrigation water.

[0162] Based on the total number of irrigation pixels mentioned above, the irrigated area of ​​the irrigation district is calculated by multiplying the total number of irrigation pixels by the ground area corresponding to a single pixel. The ground area corresponding to a single pixel is determined by the spatial resolution of the satellite imagery used in this study. If the spatial resolution is 30 meters, the ground area corresponding to a single pixel is 900 square meters; if the spatial resolution is 10 meters, the ground area corresponding to a single pixel is 100 square meters. This product is the final irrigated area extracted value of the target irrigation district for the cotton field under mulching and drip irrigation during this spring irrigation season.

[0163] The extracted irrigation area value is used as the final output of the entire methodology. This output is exported as a vector boundary file or a geotagged raster format, which can be directly imported into the spatial data layer of the irrigation district's BIM planning model. In the BIM platform, the extracted irrigation area is overlaid and merged with the existing canal system orientation layer, plot boundary layer, and topographic elevation layer in the model. Designers can then compare the differences between the planned irrigation area and the actual remotely sensed irrigation area on a plot-by-plot basis. Based on this, they can adjust the drip irrigation tape spacing for the next irrigation season, optimize the canal valve grouping scheme, and correct irrigation system parameters, thus achieving a data loop from satellite remote sensing monitoring to refined irrigation district planning and design.

[0164] In this embodiment of the invention, the step value of cooling rate is used instead of the absolute temperature threshold as the parameter for identifying irrigated wet areas, overcoming the defect that a fixed temperature threshold cannot adapt to the spatial heterogeneity of soil texture. By constructing a transverse spatial attenuation profile perpendicular to the direction of the drip irrigation belt and calculating its spatial second difference, the zero-crossing point where the step value of cooling rate changes from a slow drop to a sudden drop is automatically captured as the geometric boundary line of water lateral expansion, overcoming the defect that the uniform attenuation model cannot identify the spatial abrupt change location of the attenuation rate, resulting in blurred boundaries. By detecting local valleys and extracting the skeleton, the spatial center line of the drip irrigation belt is obtained from the thermal infrared image itself, realizing pure data-driven extraction of the drip irrigation belt location and reducing the dependence on field measured data. By using the spatial heat flow direction dispersion mechanism to remove false pixels at the edge of the damaged mulch film, the mutual verification of dual-channel evidence and the automatic suppression of false irrigation judgment are realized. The step value of cooling rate is accumulated in segments along the direction of the drip irrigation belt at a fixed spatial step size, preserving the spatial non-uniformity information of irrigation intensity along the direction of the drip irrigation belt, realizing integrated automated processing from multi-temporal satellite image input to irrigation area output.

[0165] like Figure 2 As shown, embodiments of the present invention also provide an intelligent system for extracting irrigation area in irrigation districts based on satellite remote sensing imagery, comprising:

[0166] The acquisition module is used to acquire multi-temporal satellite images of cotton fields under mulching and drip irrigation during the spring irrigation season, and spatially register them with the BIM planning model of the irrigation area project to obtain the original image set for refined design; based on the original image set for refined design, the pixel thermal infrared temperature change sequence over time is extracted to obtain the pixel temperature time series used to support the irrigation area planning and design.

[0167] The calculation module is used to calculate the difference in radiation temperature between adjacent afternoon transit times based on the pixel temperature time series, and obtain the hourly temperature drop amplitude; based on the hourly temperature drop amplitude, it calculates the difference in temperature drop amplitude between adjacent times, and obtains the step amount of the temperature drop rate.

[0168] The module is used to construct a transverse spatial attenuation profile perpendicular to the direction of the drip irrigation tape based on the step amount of cooling rate and the vertical distance of the pixel from the drip irrigation tape, and obtain the profile attenuation sequence.

[0169] The interception module is used to determine the spatial coordinates of the abrupt change in the slope of the sequence based on the profile attenuation sequence, and obtain the geometric boundary line of water lateral expansion; based on the geometric boundary line of water lateral expansion, it intercepts the set of continuous high-step pixels inside the boundary line to obtain the active pixels of water evaporation.

[0170] The accumulation module is used to accumulate the cooling rate step by step along the drip irrigation tape based on the active water evaporation pixels to obtain the cumulative heat dissipation of the strip.

[0171] The extraction module is used to calculate the irrigation water infiltration probability of pixels based on the cumulative amount of heat dissipation in the strip and the short-wave infrared reflectivity, so as to obtain the irrigation status; based on the irrigation status, the actual surface coverage area corresponding to the total number of irrigated pixels is calculated to obtain the irrigated area of ​​the irrigation district.

[0172] It should be noted that this system is a system corresponding to the above method. All implementation methods in the above method embodiments are applicable to this embodiment and can achieve the same technical effect.

[0173] Embodiments of the present invention also provide a computing device, including: a processor and a memory storing a computer program, wherein the computer program, when executed by the processor, performs the method described above. All implementations in the above method embodiments are applicable to this embodiment and can achieve the same technical effects.

[0174] Embodiments of the present invention also provide a computer-readable storage medium storing instructions that, when executed on a computer, cause the computer to perform the method described above. All implementations in the above method embodiments are applicable to this embodiment and can achieve the same technical effects.

[0175] The above description represents the preferred embodiments of the present invention. It should be noted that those skilled in the art can make various improvements and modifications without departing from the principles of the present invention, and these improvements and modifications should also be considered within the scope of protection of the present invention.

Claims

1. A method for intelligent extraction of irrigation area in irrigation districts based on satellite remote sensing imagery, characterized in that, The method includes: Step 1: Acquire multi-temporal satellite images of cotton fields under mulching and drip irrigation during the spring irrigation season, and spatially register them with the BIM planning model of the irrigation area project to obtain the original image set for refined design; Based on the original image set for refined design, extract the pixel thermal infrared temperature change sequence over time to obtain the pixel temperature time series used to support the irrigation area planning and design. Step 2: Based on the pixel temperature time series, calculate the radiation temperature difference between adjacent transit times in the afternoon to obtain the hourly temperature drop. Step 3: Based on the hourly temperature drop, calculate the difference in temperature drop between adjacent time intervals to obtain the step rate of temperature drop. Step 4: Based on the step size of the cooling rate and the vertical distance between the pixel and the drip irrigation tape, construct a transverse spatial attenuation profile perpendicular to the direction of the drip irrigation tape to obtain the profile attenuation sequence. Step 5: Based on the profile attenuation sequence, determine the spatial coordinates of the abrupt change in the sequence slope to obtain the geometric boundary line of water lateral expansion; based on the geometric boundary line of water lateral expansion, extract the set of continuous high-step pixels inside the boundary line to obtain the active pixels of water evaporation. Step 6: Based on the active pixels of water evaporation, gradually add the cooling rate step along the drip irrigation tape to obtain the cumulative amount of heat dissipation in the strip. Step 7: Calculate the irrigation water infiltration probability of the pixels based on the cumulative heat dissipation of the strip and the short-wave infrared reflectivity to obtain the irrigation status; based on the irrigation status, calculate the actual surface coverage area corresponding to the total number of irrigated pixels to obtain the irrigated area of ​​the irrigation district.

2. The intelligent extraction method for irrigation area of ​​irrigation districts based on satellite remote sensing imagery according to claim 1, characterized in that, Based on the cooling rate step and the vertical distance of the pixel from the drip irrigation tape, a lateral spatial attenuation profile perpendicular to the direction of the drip irrigation tape is constructed, resulting in a profile attenuation sequence, including: Based on the step of the cooling rate, the vertical projection distance of the center point of each pixel in the direction perpendicular to the drip irrigation tape is extracted to obtain the horizontal spatial distance sequence. Based on the horizontal spatial distance sequence, the pixels are arranged from near to far according to the vertical projection distance from small to large, resulting in an ordered distance pixel set. Based on the ordered distance pixel set, the cooling rate step value corresponding to each pixel is extracted to obtain the profile decay sequence.

3. The intelligent extraction method for irrigation area of ​​irrigation district based on satellite remote sensing imagery according to claim 2, characterized in that, Based on the cooling rate step, the vertical projection distance of each pixel center point in the direction perpendicular to the drip irrigation tape is extracted to obtain the lateral spatial distance sequence, including: Based on the step rate of cooling, the center points of continuous pixels that reach local valley values ​​are extracted to obtain the spatial center line of the drip irrigation tape; Based on the center line of the drip irrigation tape, the extension direction between adjacent center points is extracted, and a lateral reference axis perpendicular to this extension direction is constructed to obtain the lateral projection reference line. Based on the spatial coordinates of the lateral projection reference line and the center point of each pixel, the vertical drop from the center point of each pixel along the lateral reference axis to the center line of the drip irrigation tape is calculated to obtain the lateral spatial distance sequence.

4. The intelligent extraction method for irrigation area of ​​irrigation district based on satellite remote sensing imagery according to claim 3, characterized in that, Based on the step change in cooling rate, the center points of continuous pixels that reach local valley values ​​are extracted to obtain the spatial center line of the drip irrigation tape, including: Based on the step rate of cooling, local minimum discrete pixels with values ​​smaller than their surrounding neighboring pixels are extracted to obtain a set of candidate valley value pixels. Based on the candidate valley value pixel set, connectivity tracing is performed according to the spatial eight-neighbor relationship to obtain the candidate connected domain of the drip irrigation tape; Based on the candidate connected components of the drip irrigation tape, isolated noise connected components with an area smaller than the background noise area threshold are removed to obtain the candidate skeleton pixel set of the drip irrigation tape. Based on the candidate skeleton pixel set of the drip irrigation tape, the center skeleton line of the connected domain is extracted to obtain the spatial center line of the drip irrigation tape.

5. The intelligent extraction method for irrigation area of ​​irrigation district based on satellite remote sensing imagery according to claim 4, characterized in that, Based on the profile attenuation sequence, the spatial coordinates of the abrupt change in the sequence slope are determined, thus obtaining the geometric boundary line for water lateral expansion. Based on this boundary line, a set of continuous high-step intensity pixels inside the boundary line is extracted to obtain the pixels with active water evaporation, including: Based on the profile attenuation sequence, the difference in cooling rate step between adjacent ordered pixels is calculated to obtain the profile spatial change sequence; based on the profile spatial change sequence, the difference in spatial change between adjacent profiles is calculated to obtain the profile spatial change acceleration sequence. Based on the spatial change acceleration sequence of the profile, the pixels corresponding to the zero-crossing positions where the values ​​change from positive to negative are extracted to obtain the spatial coordinates of the positions where the sequence slope changes abruptly. Based on the spatial coordinates of the abrupt change in the sequence slope, the geometric boundary line of water lateral expansion on both sides of the drip irrigation tape is determined, and the lateral expansion boundary range is obtained. Based on the lateral expansion boundary range, continuous pixels with a cooling rate step greater than zero within the range are extracted to obtain the inner candidate pixel set. Based on the inner candidate pixel set, edge pixels affected by lateral heat conduction interference from the damaged edge of the mulch film are removed to obtain pixels with active water evaporation.

6. The intelligent extraction method for irrigation area based on satellite remote sensing imagery according to claim 5, characterized in that, Based on the inner candidate pixel set, edge pixels affected by lateral heat conduction interference from the damaged edge of the plastic film are removed, resulting in pixels with active water evaporation, including: Based on the inner candidate pixel set, calculate the spatial drop direction of the cooling rate step in the two-dimensional spatial neighborhood of each pixel to obtain the heat flow direction in the pixel space. Based on the heat flow direction in the pixel space, the standard deviation of the angle between the heat flow direction vectors of each pixel space in the local neighborhood is calculated to obtain the heat flow direction dispersion. Based on the heat flow direction dispersion, pixels with heat flow direction dispersion greater than a preset dispersion threshold are removed from the inner candidate pixel set to obtain pixels with active water evaporation.

7. The intelligent extraction method for irrigation area of ​​irrigation district based on satellite remote sensing imagery according to claim 6, characterized in that, Based on the inner candidate pixel set, the spatial drop direction of the cooling rate step in the two-dimensional spatial neighborhood of each pixel is calculated to obtain the spatial heat flow direction of the pixel, including: Based on the inner candidate pixel set, extract the direct neighboring pixels of each target pixel in the horizontal and vertical directions to obtain the spatial neighborhood pixel group; based on the spatial neighborhood pixel group, calculate the difference in cooling rate step between the target pixel and the horizontally adjacent pixels to obtain the horizontal spatial drop. Based on the spatial neighborhood pixel group, calculate the difference in the cooling rate step between the target pixel and its longitudinally adjacent pixels to obtain the longitudinal spatial drop. Based on the horizontal and vertical spatial drop, a two-dimensional spatial pointing vector is constructed to obtain the direction of heat flow in the pixel space.

8. The intelligent extraction method for irrigation area of ​​irrigation district based on satellite remote sensing imagery according to claim 7, characterized in that, Based on the pixels with active water evaporation, the cooling rate step change is accumulated segment by segment along the drip irrigation strip to obtain the cumulative strip heat dissipation, including: Based on the active pixels of water evaporation, the longitudinal projection coordinates of each active pixel along the spatial center line of the drip irrigation tape are extracted to obtain the longitudinal position sequence of the active pixels; Based on the longitudinal position sequence of active pixels, the drip irrigation tape is divided into continuous cumulative segments according to a fixed spatial step size, resulting in a strip segmented interval set. Based on the strip segment interval set, the active water evaporation pixels falling within each segment are selected to obtain the active pixel subset of each segment; Based on the active pixel subsets of each segment, the cooling rate step of all active water evaporation pixels in the corresponding segment is accumulated segment by segment to obtain the cumulative amount of heat dissipation in the strip.

9. The intelligent extraction method for irrigation area of ​​irrigation district based on satellite remote sensing imagery according to claim 8, characterized in that, The irrigation status is obtained by calculating the pixel irrigation water infiltration probability based on the cumulative amount of heat dissipation in the strip and the short-wave infrared reflectivity. Based on the irrigation status, the actual surface coverage area corresponding to the total number of irrigated pixels is calculated to obtain the irrigated area of ​​the irrigation district, including: Based on the cumulative amount of heat dissipation in the band, the thermal response intensity level corresponding to each pixel is extracted to obtain the heat dissipation feature sequence; based on the short-wave infrared reflectance, the spectral absorption depth corresponding to each pixel is extracted to obtain the moisture absorption feature sequence. Based on the heat dissipation characteristic sequence and the water absorption characteristic sequence, the heat response intensity level and the spectral absorption depth are scaled proportionally to a unified numerical range. The product of the scaled values ​​within the same pixel is calculated as the spatial matching degree to obtain the pixel irrigation water penetration probability. Based on the pixel irrigation water infiltration probability, pixels with probability values ​​greater than the infiltration probability threshold are selected to obtain the set of pixels in the irrigation state; Based on the set of irrigated state pixels, the number of irrigated state pixels in the set is accumulated to obtain the total number of target irrigated pixels; based on the total number of target irrigated pixels, the actual surface coverage area corresponding to the total number of target irrigated pixels is calculated to obtain the irrigated area of ​​the irrigation district.

10. An intelligent system for extracting irrigation area of ​​irrigation districts based on satellite remote sensing imagery, wherein the system implements the method as described in any one of claims 1 to 9, characterized in that, include: The acquisition module is used to acquire multi-temporal satellite images of cotton fields under mulching and drip irrigation during the spring irrigation season, and spatially register them with the BIM planning model of the irrigation area project to obtain the original image set for refined design; based on the original image set for refined design, the pixel thermal infrared temperature change sequence over time is extracted to obtain the pixel temperature time series used to support the irrigation area planning and design. The calculation module is used to calculate the difference in radiation temperature between adjacent afternoon transit times based on the pixel temperature time series, and obtain the hourly temperature drop amplitude; based on the hourly temperature drop amplitude, it calculates the difference in temperature drop amplitude between adjacent times, and obtains the step amount of the temperature drop rate. The module is used to construct a transverse spatial attenuation profile perpendicular to the direction of the drip irrigation tape based on the step amount of cooling rate and the vertical distance of the pixel from the drip irrigation tape, and obtain the profile attenuation sequence. The interception module is used to determine the spatial coordinates of the abrupt change in the slope of the sequence based on the profile attenuation sequence, and obtain the geometric boundary line of water lateral expansion; based on the geometric boundary line of water lateral expansion, it intercepts the set of continuous high-step pixels inside the boundary line to obtain the active pixels of water evaporation. The accumulation module is used to accumulate the cooling rate step by step along the drip irrigation tape based on the active water evaporation pixels to obtain the cumulative heat dissipation of the strip. The extraction module is used to calculate the irrigation water infiltration probability of a pixel based on the cumulative amount of heat dissipation in the strip and the short-wave infrared reflectivity, so as to obtain the irrigation status. Based on the irrigation status, the actual surface coverage area corresponding to the total number of irrigated pixels is calculated to obtain the irrigated area of ​​the irrigation district.

Images