Remote sensing estimation method and system for separating groundwater fluxes of infiltration and lateral inflow

Abstract

The present application belongs to the technical field of cross-basin groundwater circulation flux estimation, and particularly relates to a remote sensing estimation method and system for separating groundwater flux of infiltration and lateral inflow. The method comprises the following steps: obtaining a multi-source remote sensing data set of a study area; correcting evapotranspiration data of different land use types in the study area; dividing the study area into equal-area units, identifying vertical dominant units and lateral dominant units, and then determining a cross-basin groundwater exchange zone; using the vertical dominant units to solve local vertical infiltration and form a background field by interpolation, and establishing a joint equation set for the lateral dominant units with the background field as the initial value and groundwater dynamic characteristic parameters as additional constraints to solve the external lateral inflow, and then obtain the cross-basin groundwater circulation flux. The present application identifies the vertical dominant units and the lateral dominant units, and realizes quantitative estimation of the cross-basin groundwater circulation flux.

Description

Technical Field

[0001] This invention belongs to the field of cross-basin groundwater circulation flux estimation technology, and particularly relates to a remote sensing estimation method and system for groundwater flux that separates infiltration and lateral inflow. Background Technology

[0002] Inter-basin groundwater circulation flux refers to the lateral exchange of groundwater between natural watershed boundaries. It is of great significance for regional water resource total volume accounting, environmental impact assessment of inter-basin water transfer projects, and sustainable development and utilization of groundwater. However, groundwater flow is concealed, watershed boundary conditions are complex, and traditional monitoring methods are costly and difficult to deploy over large-scale areas, making it difficult to quantitatively estimate inter-basin groundwater circulation flux.

[0003] To address the aforementioned issues, existing technologies have attempted to estimate using remote sensing data: First, remote sensing precipitation data with higher correlation to ground stations are selected. Then, taking the entire watershed as a unit, a globally unified evapotranspiration correction coefficient is calculated based on the principle of water balance to correct the remote sensing evapotranspiration products. Finally, the multi-year average net recharge is calculated by simplifying the water balance equation, and areas with positive net recharge values ​​are considered recharge areas, while areas with negative values ​​are considered discharge areas.

[0004] However, this method still has certain limitations when applied to estimating inter-basin groundwater circulation flux. This method simplifies the complex vertical and lateral coupling system into a vertically independent system, treating the net recharge as the residual amount of local vertical precipitation infiltration. When inter-basin lateral groundwater flow exists in actual hydrogeological conditions, this net recharge is essentially a mixed signal of local vertical infiltration and external lateral inflow, making it impossible to effectively separate the two and thus difficult to directly obtain inter-basin groundwater circulation flux. Summary of the Invention

[0005] To overcome the shortcomings of the prior art, this invention provides a remote sensing estimation method and system for groundwater flux that separates infiltration and lateral inflow. The study area is divided into multiple units, and the response characteristics of land water storage changes and precipitation in each unit are calculated. Based on this, vertically dominant units and laterally dominant units are identified, achieving effective separation of local vertical infiltration and external lateral inflow. Then, based on the vertically dominant units, the local vertical infiltration is solved to form a spatial background field. For the laterally dominant units, a joint equation system is constructed with the background field as the initial value and the response characteristics as constraints to invert the external lateral inflow, thereby realizing quantitative estimation of cross-basin groundwater circulation flux.

[0006] To achieve the above objectives, one or more embodiments of the present invention provide the following technical solutions:

[0007] The first aspect of this invention provides a remote sensing estimation method for groundwater flux that separates infiltration and lateral inflow.

[0008] A remote sensing method for estimating groundwater flux by separating infiltration and lateral inflow includes the following steps:

[0009] Acquire a multi-source remote sensing dataset for the study area, which includes land water storage change data, remote sensing precipitation data, remote sensing evapotranspiration data, digital elevation model data, and land use type data;

[0010] Based on land use type data, the study area was divided into vegetated area and bare soil area. Based on the principle of water balance, the evapotranspiration correction coefficients of vegetated area and bare soil area were calculated respectively. The evapotranspiration correction coefficients were applied to the remote sensing evapotranspiration data of the corresponding area to obtain the corrected spatial distribution data of actual evapotranspiration.

[0011] The study area was divided into units of equal area. For each unit, the time series of changes in terrestrial water storage and the time series of precipitation were extracted, and the correlation coefficient and phase lag time of the two were calculated.

[0012] Based on the correlation coefficient and phase lag time of each unit, the corresponding units are divided into vertically dominant units or laterally dominant units. The boundaries of surface sub-basins are extracted using digital elevation models, and the strip-shaped area composed of laterally dominant units belonging to different surface sub-basins is identified as the cross-basin groundwater exchange zone.

[0013] The local vertical infiltration rate is solved using vertically dominant units and interpolated to form a background field. Combined with the corrected spatial distribution data of actual evapotranspiration, a joint equation system is established for laterally dominant units with the background field as the initial value and the correlation coefficient and phase lag time as constraints. The external lateral inflow rate is then obtained by solving the equation system.

[0014] The external lateral inflows located within the inter-basin groundwater exchange zone and crossing the surface watershed boundary are statistically summarized to obtain the inter-basin groundwater circulation flux.

[0015] A second aspect of the present invention provides a remote sensing estimation system for groundwater flux that separates infiltration and lateral inflow.

[0016] A remote sensing estimation system for groundwater flux that separates infiltration and lateral inflow includes:

[0017] The data acquisition module is configured to acquire a multi-source remote sensing dataset of the study area, which includes land water storage change data, remote sensing precipitation data, remote sensing evapotranspiration data, digital elevation model data, and land use type data.

[0018] The zoning correction module is configured to: divide the study area into vegetated areas and bare soil areas according to land use type data; calculate the evapotranspiration correction coefficients for vegetated areas and bare soil areas respectively based on the water balance principle; apply the evapotranspiration correction coefficients to the remote sensing evapotranspiration data of the corresponding areas; and obtain the corrected spatial distribution data of actual evapotranspiration.

[0019] The unit division module is configured to: divide the study area into units of equal area, extract the time series of changes in terrestrial water storage and precipitation for each unit, and calculate the correlation coefficient and phase lag time between the two.

[0020] The type identification module is configured to: classify the corresponding unit into vertically dominant units or laterally dominant units based on the correlation coefficient and phase lag time of each unit; extract the surface sub-basin boundary using the digital elevation model; and identify the strip-shaped area composed of laterally dominant units belonging to different surface sub-basins as the cross-basin groundwater exchange zone.

[0021] The joint solution module is configured to: use the vertically dominant element to solve the local vertical infiltration and interpolate to form a background field; combine the corrected spatial distribution data of actual evapotranspiration to establish a joint equation system for the laterally dominant element with the background field as the initial value and the correlation coefficient and phase lag time as constraints, and solve for the external lateral inflow.

[0022] The aggregation module is configured to: statistically aggregate the external lateral inflows located within the inter-basin groundwater exchange zone and crossing the surface watershed boundary to obtain the inter-basin groundwater circulation flux.

[0023] The above one or more technical solutions have the following beneficial effects:

[0024] This invention provides a remote sensing estimation method and system for groundwater flux that separates infiltration and lateral inflow. By introducing a gridded unit division and a groundwater dynamic response characteristic parameter identification mechanism, the study area is finely divided into vertically dominant units and laterally dominant units, and the spatial location of the cross-basin groundwater exchange zone is located accordingly. By constructing a joint equation system with the local vertical infiltration background field as the initial value and the correlation coefficient and phase lag time as proportional constraints, the quantitative decoupling of vertical infiltration and lateral inflow is achieved, thereby realizing the quantitative estimation of cross-basin groundwater circulation flux.

[0025] Advantages of additional aspects of the invention will be set forth in part in the description which follows, and in part will be obvious from the description, or may be learned by practice of the invention. Attached Figure Description

[0026] The accompanying drawings, which form part of this invention, are used to provide a further understanding of the invention. The illustrative embodiments of the invention and their descriptions are used to explain the invention and do not constitute an improper limitation of the invention.

[0027] Figure 1 This is a flowchart of the method in Example 1. Detailed Implementation

[0028] It should be noted that the following detailed descriptions are exemplary and intended to provide further illustration of the invention. Unless otherwise specified, all technical and scientific terms used herein have the same meaning as commonly understood by one of ordinary skill in the art to which this invention pertains.

[0029] It should be noted that the terminology used herein is for the purpose of describing particular implementations only and is not intended to limit the exemplary implementations of the present invention.

[0030] Where there is no conflict, the embodiments and features in the embodiments of the present invention can be combined with each other.

[0031] Terminology Explanation:

[0032] GRACE gravity satellite products: This is a dual-satellite system jointly launched by NASA and the German Aerospace Center. It inverts the spatiotemporal changes of the Earth's gravity field by measuring the changes in the distance between the two satellites, thereby calculating information on changes in land water storage.

[0033] MOD16 remote sensing evapotranspiration product: This is a dataset of global land surface evapotranspiration calculated based on medium resolution imaging spectrometer data combined with meteorological reanalysis data.

[0034] The overall concept proposed in this invention is as follows:

[0035] Existing technologies suffer from the difficulty of effectively separating local vertical infiltration from external lateral inflow, thus making it impossible to directly obtain cross-basin groundwater circulation flux. Therefore, a method is needed that can effectively separate local vertical infiltration from external lateral inflow, thereby quantitatively estimating cross-basin groundwater circulation flux.

[0036] Based on this, the present invention provides a remote sensing estimation method and system for groundwater flux that separates infiltration and lateral inflow. The method generally includes: acquiring data on changes in terrestrial water storage, remote sensing precipitation data, and land use data of the study area; dividing the study area into vegetated and bare soil zones, calculating evapotranspiration correction coefficients for each zone, and obtaining the corrected actual evapotranspiration; dividing the study area into multiple units of equal area, extracting water storage changes and precipitation time series from each unit, calculating correlation coefficients and phase lag times, and classifying the units into vertically dominant and laterally dominant types; determining the inter-basin groundwater exchange zone based on the laterally dominant units, establishing a joint equation system with the background field as the initial value and response characteristics as constraints, and solving for the external lateral inflow; summarizing the external lateral inflows located within the exchange zone and crossing the basin boundary to obtain the inter-basin groundwater circulation flux.

[0037] This invention divides the study area into units and calculates the response characteristics of land water storage changes and precipitation within each unit. Based on this, it identifies vertically dominant units and laterally dominant units. Then, based on the vertically dominant units, it solves the local vertical infiltration to form a spatial background field. For the laterally dominant units, it constructs a joint equation system with the background field as the initial value and the response characteristics as constraints to invert the external lateral inflow. This achieves effective separation of local vertical infiltration and external lateral inflow, solving the technical problem in the prior art that it is difficult to effectively separate local vertical infiltration and external lateral inflow, thus making it impossible to directly obtain the cross-basin groundwater circulation flux.

[0038] Example 1

[0039] This embodiment discloses a remote sensing estimation method for groundwater flux that separates infiltration and lateral inflow.

[0040] like Figure 1 As shown, the remote sensing estimation method for groundwater flux, which separates infiltration and lateral inflow, includes the following steps:

[0041] S1. Obtain a multi-source remote sensing dataset for the study area. The multi-source remote sensing dataset includes land water storage change data, remote sensing precipitation data, remote sensing evapotranspiration data, digital elevation model data, and land use type data. Resample all data to a uniform spatial resolution to obtain a spatially matched multi-source remote sensing dataset.

[0042] S2. Based on land use type data, the study area is divided into vegetated area and bare soil area. Based on the principle of water balance, the evapotranspiration correction coefficients of vegetated area and bare soil area are calculated respectively. The evapotranspiration correction coefficients are applied to the remote sensing evapotranspiration data of the corresponding area to obtain the corrected spatial distribution data of actual evapotranspiration.

[0043] S3. Divide the study area into units of equal area, extract the time series of land water storage change and precipitation for each unit, calculate the correlation coefficient and phase lag time of the two, and divide the units into vertically dominant units and laterally dominant units based on the correlation coefficient and phase lag time of each unit.

[0044] S4. Use digital elevation models to extract the boundaries of surface sub-basins and identify the strip-shaped areas composed of laterally dominant units belonging to different sub-basins as cross-basin groundwater exchange zones.

[0045] The local vertical infiltration rate is solved using vertically dominant elements and interpolated to form a background field. A joint equation system is established for laterally dominant elements with the background field as the initial value and the correlation coefficient and phase lag time as constraints, and the external lateral inflow rate is obtained by solving it.

[0046] The external lateral inflows located within the inter-basin groundwater exchange zone and crossing the surface watershed boundary are statistically summarized to obtain the inter-basin groundwater circulation flux.

[0047] This invention introduces a unit zoning mechanism driven by groundwater dynamic characteristics to identify vertically dominant and laterally dominant units, and constructs a decoupling mechanism for joint vertical and lateral inversion to establish water balance equations and additional constraint equations for the two types of units respectively. This achieves effective separation of local vertical infiltration and external lateral inflow, thereby quantitatively estimating the cross-basin groundwater circulation flux.

[0048] In the above technical solution:

[0049] (1) The study area was divided into units of equal area instead of the traditional method of taking the entire watershed as the whole unit because the watershed-scale overall correction cannot identify the spatial location of the cross-watershed exchange zone, while grid division can capture the spatial differences in groundwater storage changes and precipitation response in different regions.

[0050] (2) The correlation coefficient and phase lag time are used as the basis for classifying unit types because these two parameters characterize the degree of control of local precipitation on groundwater storage from the two dimensions of correlation and timeliness, respectively. Using either parameter alone may lead to misjudgment. For example, a unit with high correlation but long phase lag may actually be controlled by lateral recharge but is mistakenly classified as vertically dominant due to its high correlation.

[0051] (3) The local vertical infiltration obtained from solving the vertical dominant unit is interpolated to form the background field and then used for the constraint solution of the lateral dominant unit. This is because if the local vertical infiltration of the lateral dominant unit is directly assumed to be zero or a uniform empirical value is used, the continuous and gradual characteristics of vertical infiltration in space will be ignored, resulting in the inversion result of the external lateral inflow being systematically overestimated.

[0052] (4) The reason for establishing a joint set of equations constrained by correlation coefficient and phase lag time instead of a single water balance equation is that the water balance equation alone cannot separate the two unknowns of local vertical infiltration and external lateral inflow. Introducing response characteristics as proportional constraints fills this information gap and decouples the two.

[0053] (5) Finally, the external lateral inflows within the cross-basin groundwater exchange zone are statistically summarized instead of the inflows of all laterally dominant units. This is because only lateral flows that cross the surface watershed boundary constitute a true cross-basin circulation. If boundary screening is not performed, the lateral circulation within the sub-basin will be mistakenly included in the cross-basin flux, thus overestimating the cross-basin exchange volume.

[0054] Based on this, the evapotranspiration correction coefficients for vegetated and bare soil areas are calculated using the principle of local water balance. The ratio of the multi-year average precipitation to the original MOD16 evapotranspiration in each area is used as the correction coefficient, and the actual evapotranspiration in water areas is replaced by the water surface evaporation observed at meteorological stations.

[0055] Furthermore, the corrected spatial distribution data of actual evapotranspiration is obtained by multiplying the original MOD16 evapotranspiration of each pixel by the correction coefficient of the corresponding land cover type, with vegetated area pixels multiplied by the vegetated area correction coefficient and bare soil area pixels multiplied by the bare soil area correction coefficient.

[0056] The technical solution of this embodiment divides the study area into three categories: vegetated area, bare soil area, and water body area, and performs evapotranspiration correction separately for each category, instead of using a globally uniform correction coefficient. This is because there are fundamental differences in the energy balance and water exchange mechanisms of different land cover types. Water bodies are mainly subject to free evaporation, vegetated areas include both canopy transpiration and soil evaporation, and bare soil areas are dominated by soil evaporation. If a single coefficient is used for global correction, evapotranspiration will be underestimated in densely vegetated areas and overestimated in water body areas, resulting in a systematic shift in the spatial distribution of net recharge in subsequent water balance calculations.

[0057] The ratio of the multi-year average precipitation to the original MOD16 evapotranspiration in each region is used as the correction coefficient for vegetated and bare soil areas because in areas lacking measured evapotranspiration data, the multi-year average precipitation can be approximated as the sum of the actual evapotranspiration and runoff. For vegetated and bare soil areas, the surface runoff generation capacity is weak and negligible. Therefore, the precipitation itself is a reasonable estimate of the actual evapotranspiration. Using the ratio of the precipitation to the original remotely sensed evapotranspiration as the correction coefficient is consistent with the water balance principle and avoids introducing additional parameterization assumptions.

[0058] The reason for directly using the water surface evaporation observed by meteorological stations to replace the remote sensing evapotranspiration product in the water body area is that the MOD16 evapotranspiration product did not consider the water evaporation process in its algorithm design. The evapotranspiration value of the water body pixel is actually the interpolation result of the adjacent land surface pixel rather than the actual water surface evaporation. If used directly, the evapotranspiration of the water body area will be seriously underestimated, which will lead to a large error in the calculation of the net recharge of the area near the water body.

[0059] By using land cover type zoning correction, the causes of different types of evapotranspiration errors can be corrected separately, so that the corrected actual evapotranspiration more accurately reflects the water consumption characteristics under different land cover conditions.

[0060] Step S1 of this invention involves multi-source remote sensing data acquisition and preprocessing, aiming to acquire multi-source remote sensing and geographic information data of the study area within the target time period, and to perform spatiotemporal matching and standardization processing on all data. The multi-source remote sensing dataset specifically includes the following five data types, which are then standardized to ensure the accuracy and comparability of subsequent calculations:

[0061] 1. Land water storage change data: Data products are from the GRACE and GRACE-FO gravity satellites; specific parameters are equivalent water height in centimeters, representing the monthly variation of land water storage anomalies; time parameters are monthly temporal resolution, and spatial parameters are the original spatial resolution of 0.25º×0.25º; data sources are NASA's Jet Propulsion Laboratory (JPL) or the University of Texas Center for Space Research (CSR).

[0062] 2. Remote sensing precipitation data: The data product is a high spatiotemporal resolution surface meteorological element driven dataset for China. The specific parameters are monthly total precipitation in millimeters; monthly temporal resolution and 1 km spatial resolution; the data source is the National Tibetan Plateau Data Science Center (TPDC).

[0063] 3. Remote sensing evapotranspiration data: The data product is the MOD16 series evapotranspiration product, with specific parameters representing actual evapotranspiration in millimeters per year; the time parameter is annual, and the spatial parameter is 500 meters; the data source is NASA LP DAAC.

[0064] 4. Land use type data: The specific parameters are land cover classification type codes; the time resolution is yearly, and the spatial resolution is 500 meters; the data source is USGS LP DAAC.

[0065] 5. Digital Elevation Model Data: The Digital Elevation Model (DEM) uses SRTM data, with specific parameters including surface elevation in meters and a spatial parameter of 90 meters. The data source is NASA CGIAR-CSI.

[0066] It can be understood that the aforementioned NASA LP DAAC specifically refers to the NASA Land Processes Distributed Active Archive Center. It is part of NASA's Earth Observation System Data and Information System, specifically responsible for processing and distributing land remote sensing products such as MODIS (Medium Resolution Imaging Spectroradiometer).

[0067] USGS LP DAAC specifically refers to the US Geological Survey Land Processes Distributed Active Archive Center.

[0068] NASA CGIAR-CSI specifically refers to the NASA Consultative Group on International Agricultural Research - Consortium for Spatial Information. This is a collaborative project between NASA and the Spatial Information Consortium, a branch of CGIAR (Consortium for International Agricultural Research), specifically responsible for processing and distributing SRTM (Resolution-Based Digital Elevation Model) data.

[0069] Among them, the GRACE gravity satellite land water storage change data is monthly grid data, which is post-processed and corrected using decorrelation filtering and scale factor recovery methods; the remote sensing precipitation data is high-resolution precipitation spatial distribution grid data with a spatial resolution of 1 km and a temporal resolution of monthly, which is clipped and projected to keep it consistent with the spatial coordinate system of the land water storage change data.

[0070] Monthly terrestrial water storage variation grid data for the study area within the target time series were obtained from GRACE gravity satellite products. This data reflects the monthly variations in the total storage of surface water, soil water, and groundwater in the study area and is a key input for subsequent water balance calculations. To address issues such as coarse spatial resolution and signal leakage in the terrestrial water storage variation data, decorrelation filtering and scale factor recovery methods were employed for post-processing correction to improve data reliability.

[0071] From the high-resolution precipitation dataset for the Chinese region released by the Qinghai-Tibet Plateau Data Center, monthly spatial distribution grid data of precipitation in the study area within the target time series were obtained. This dataset has a spatial resolution of 1 km and a monthly temporal resolution, which can well reflect the spatial distribution characteristics of precipitation at the regional scale. After acquisition, the data was cropped and projected to be consistent with the spatial coordinate system of the land water storage change data.

[0072] The annual spatial distribution grid data of evapotranspiration in the study area within the target time series was obtained from the MOD16 series evapotranspiration products. This product has a spatial resolution of 500 meters, and its formula is as follows:

[0073]

[0074] In the formula, Evaporation amount;

[0075] The slope of the saturated water vapor pressure curve;

[0076] Net surface radiation;

[0077] Soil heat flux;

[0078] air density;

[0079] The specific heat of air at constant pressure;

[0080] It is the saturated vapor pressure;

[0081] This is the actual water vapor pressure;

[0082] Aerodynamic impedance;

[0083] The surface impedance of the canopy;

[0084] This is the constant of the wet / dry meter.

[0085] After acquisition, the evapotranspiration bands are extracted and mosaicked to obtain a complete evapotranspiration data layer covering the study area.

[0086] Furthermore, digital elevation model data for the study area, with a spatial resolution of no less than 90 meters, will be collected for subsequent topographic analysis and sub-basin delineation. Land use type data for the study area will be collected to distinguish between vegetated areas, bare soil areas, and water bodies, providing a basis for zonal correction of subsequent evapotranspiration data. Soil texture spatial distribution data for the study area will be collected for parameter setting in the vadose zone numerical simulation.

[0087] All the remote sensing data and auxiliary geographic information data acquired above underwent a unified projection transformation to the same geodetic coordinate system. All raster data were resampled to the same spatial resolution of 500 meters, consistent with the original resolution of the MOD16 evapotranspiration data, ensuring that each data layer could achieve pixel-level point-by-point calculations in space. Missing values ​​after resampling were filled using nearest-neighbor interpolation to ensure complete data coverage for all spatial locations within the study area.

[0088] Step S2 of the present invention is zonal correction of evapotranspiration data based on water balance. It aims to correct the spatial differences of MOD16 remote sensing evapotranspiration products by dividing land use types and applying the principle of local water balance, so as to obtain spatial distribution data of actual evapotranspiration that is more consistent with the actual situation of the study area.

[0089] In step S1, land use type data for the study area was acquired, which divides land cover types into three basic units: vegetated areas, bare soil areas, and water bodies. Since the surface energy balance and water exchange processes differ significantly among different land cover types, the systematic error characteristics of their remote sensing evapotranspiration products also vary. Therefore, corrections need to be performed separately for each type of unit.

[0090] For water body units, the actual evapotranspiration is the water surface evaporation. Water surface evaporation is obtained by collecting evaporation pan observation data from meteorological stations near the study area and converting it to standard water surface evaporation using the evaporation pan coefficient. Since the area of ​​the water body usually accounts for a small proportion of the total area of ​​the study area, and its evaporation process is relatively stable, the water surface evaporation observed from meteorological stations can be directly used as the actual evapotranspiration of the area without the need for conversion using a correction factor.

[0091] For vegetated and bare soil units, evapotranspiration correction coefficients are calculated separately using a method based on the principle of local water balance. Taking the vegetated area as an example, all pixels covered by vegetation are first extracted from the land use type data. These pixels are then spatially arranged into a complete set, which can be considered a locally balanced area in terms of water balance. For this locally balanced area, assuming that the water storage change is approximately zero on a multi-year average scale, its water balance relationship can be expressed as:

[0092]

[0093] In the formula, The multi-year average precipitation of this local equilibrium zone is obtained by statistically averaging the remote sensing precipitation data acquired in step S1 within this region.

[0094] This represents the multi-year average actual evapotranspiration of the region.

[0095] This represents the surface runoff output of the region.

[0096] For vegetated areas, due to their typically rough surface and well-developed root systems, surface runoff capacity is weak. In the absence of measured runoff data, the surface runoff term can be approximated as zero. Therefore, it can be deduced that the actual evapotranspiration in vegetated areas is approximately equal to the precipitation in that region.

[0097] On the other hand, the original MOD16 evapotranspiration of the vegetation zone was obtained by statistically averaging the evapotranspiration data acquired in step S1 within this region, denoted as . Since there is a systematic proportional relationship between MOD16 evapotranspiration and actual evapotranspiration, an evapotranspiration correction coefficient for vegetation zones can be defined. satisfy:

[0098]

[0099] Substituting the relationship between actual evapotranspiration and precipitation into the above formula, the evapotranspiration correction coefficient for the vegetation zone can be obtained:

[0100]

[0101] Using the same method, the evapotranspiration correction factor was calculated for the bare soil area unit. The actual evapotranspiration in bare soil areas is also approximately equal to the precipitation in that area, and the formula for calculating the correction factor is:

[0102]

[0103] In the formula, This represents the multi-year average precipitation in the locally balanced zone of the bare soil area.

[0104] The statistical average value of the original MOD16 evapotranspiration in this region.

[0105] After obtaining the evapotranspiration correction coefficients for each type of unit, they are applied to every pixel of the corresponding type throughout the entire study area. For any pixel located in a vegetated area... ,in The row index (vertical numbering) of the cell. The column index (horizontal number) of the cell, its corrected actual evapotranspiration. Equals the original MOD16 evapotranspiration of that pixel multiplied by the vegetation zone correction factor. For pixels in bare soil areas, multiply by the bare soil area correction factor. For pixels in the water body area, the water surface evaporation observed at the meteorological station is directly assigned. Through the above zonal correction, the spatial distribution data of the actual evapotranspiration in the entire study area after correction is obtained. This data fully considers the differences in evapotranspiration error characteristics between different land cover types, providing a more reliable input for subsequent steps.

[0106] Step S3 of this invention is the unit partitioning and type identification driven by groundwater dynamic characteristics. It aims to divide the study area into several interconnected units and identify the dominant water source type of each unit based on the response characteristics between groundwater storage changes and precipitation, laying the foundation for the subsequent decoupling of vertical and lateral components.

[0107] The study area was divided into equal-area units according to the unified spatial resolution in step S1, with each unit measuring 500 meters by 500 meters, consistent with the resolution of each remote sensing data layer. For each unit, monthly land water storage variation data and monthly precipitation data were extracted over the entire time series, forming two time series curves respectively. The land water storage variation data was provided by GRACE gravity satellite products, reflecting the monthly fluctuations in the total storage of groundwater, soil water, and surface water within the unit. The precipitation data was provided by remote sensing precipitation products acquired in step S1, reflecting the monthly variation characteristics of atmospheric precipitation input within the unit.

[0108] For each unit, the correlation coefficient between its land water storage change time series and precipitation time series was calculated. This correlation coefficient measures the degree to which the unit's groundwater storage changes respond to precipitation input. A higher correlation coefficient indicates that the water volume changes in the unit are mainly controlled by local precipitation, and that precipitation infiltration can quickly translate into changes in groundwater storage. A lower correlation coefficient indicates that the water volume changes in the unit may be significantly affected by other factors, such as lateral groundwater inflow or outflow from adjacent units, which weaken the direct correlation between local precipitation and groundwater storage.

[0109] For each grid cell, extract its values ​​over the entire study period (total). Monthly time series of changes in terrestrial water storage (monthly) and monthly precipitation time series The correlation coefficient is calculated using the Pearson correlation coefficient, and the specific formula is as follows:

[0110]

[0111] In the formula: The correlation coefficient between changes in terrestrial water storage and precipitation in the grid cell is given, with a value range of [value range missing]. ;

[0112] For the first The change in land water storage of this unit over a month, i.e., equivalent water height, in mm (millimeters). This represents the change in terrestrial water storage in this unit for the first month. This represents the change in terrestrial water storage in this unit for the second month. For the first The change in terrestrial water storage in this unit over the past month;

[0113] For the first The monthly precipitation in this unit, in mm; This represents the rainfall in this unit for the first month. This represents the rainfall in this unit for the second month. For the first The monthly precipitation in this unit;

[0114] This represents the multi-year monthly average of the changes in terrestrial water storage in this unit during the study period;

[0115] This represents the multi-year monthly average of precipitation in this unit during the study period;

[0116] The total number of months in the time series. Indicates the first Months.

[0117] Simultaneously, the phase lag time of the land water storage change time series relative to the precipitation time series for each unit was calculated. The phase lag time reflects the time delay between the occurrence of a precipitation event and the generation of a groundwater storage response. In regions dominated by vertical infiltration, precipitation recharges groundwater through the vadose zone; this process typically has a relatively stable time lag, with short lag times and small interannual fluctuations. In regions dominated by lateral recharge, changes in groundwater storage depend more on groundwater flow in adjacent units, resulting in a more ambiguous temporal relationship with local precipitation, and longer and unstable phase lag times.

[0118] Specifically, the calculation steps for phase lag time are as follows: Construct a time series of changes in terrestrial water storage. With precipitation time series cross-correlation function Its formula is:

[0119]

[0120] In the formula, The number of months is the lag. Lagging Changes in terrestrial water storage over a month.

[0121] Calculate different lag months Below Value, find Number of months lag before reaching the maximum value The phase lag time is defined as follows: The unit is month. When the value is less than the second preset threshold, it indicates that the process of precipitation infiltration replenishing groundwater is relatively rapid. When the value is large, it indicates that the water volume of this unit is significantly affected by the delay effect of the long-range lateral flow.

[0122] Based on the two characteristic parameters mentioned above, all units within the study area are divided into two basic types. The first type is the vertically dominant unit, which meets the following two conditions: a correlation coefficient greater than a first preset threshold and a phase lag time less than a second preset threshold. The water volume variation in the vertically dominant unit mainly originates from the vertical infiltration recharge of local precipitation. There is a close synchronous relationship between its groundwater storage and precipitation. The proportion of lateral groundwater flow from adjacent units in the total water volume of this type of unit is negligible.

[0123] The second type is the laterally dominant unit. This type of unit does not meet the above conditions, namely, the correlation coefficient is less than the first preset threshold or the phase lag time is greater than the second preset threshold. The water volume change of the laterally dominant unit is significantly affected by the lateral flow of groundwater in adjacent units. Local vertical infiltration and external lateral inflow both play a role in its water volume budget. The proportion of contributions from the two sources needs to be further solved in subsequent steps.

[0124] The first preset threshold is used to distinguish whether the correlation between changes in groundwater storage and precipitation is significant. It is determined as follows: the correlation coefficients of all grid cells within the study area are calculated, and a frequency distribution histogram of the correlation coefficients is plotted. Since vertically dominant and laterally dominant cells statistically exhibit a bimodal distribution, the valley value between the two peaks is taken as the first preset threshold.

[0125] The second preset threshold is used to determine whether the response rate of precipitation infiltration to replenish groundwater is within the normal physical time range of vertical infiltration. The determination method is as follows: the grid cells with the highest correlation coefficient (such as the top 20%) are initially screened as quasi-vertical cells, the phase lag time distribution of these cells is statistically analyzed, and their upper quartiles are calculated as the second preset threshold.

[0126] Building upon the identification of unit types, the spatial location of the inter-basin exchange zone is further identified. First, digital elevation model (DEM) data is used to extract the surface sub-basin boundaries of the study area, determining the spatial extent of each sub-basin. Then, the unit types of each unit are overlaid with the surface sub-basin boundaries to identify geographically adjacent laterally dominant units belonging to different surface sub-basins. These units, located on either side of the surface watershed and exhibiting laterally dominant water source characteristics, represent potential inter-basin groundwater exchange zones. The strip-shaped area formed by these continuously distributed laterally dominant units is defined as the precise spatial location of the inter-basin groundwater exchange zone, providing spatial constraints for the quantitative calculation of inter-basin flux in subsequent steps.

[0127] Step S4 of this invention is a decoupled calculation of cross-basin flux jointly inverted by vertical and lateral directions. It aims to accurately invert flux for laterally dominant units within the cross-basin groundwater exchange zone based on unit partitioning and the determined spatial location of the exchange zone. The overall process is as follows: Based on the input unit classification results, the paths of vertically dominant units and laterally dominant units are obtained; the local vertical infiltration is calculated based on the vertically dominant units, and a background field of local vertical infiltration covering the entire area is formed using spatial interpolation; then, for the laterally dominant units within the cross-basin groundwater exchange zone identified in step S3, a joint equation set including water balance equations and groundwater dynamic characteristic constraints is established. Using the background field interpolation results as the initial iterative values ​​of the local vertical infiltration, the external lateral inflow of each laterally dominant unit is solved. Finally, determine whether the obtained external lateral inflow vector crosses the surface watershed boundary: if so, accumulate and statistically analyze this part of the flux, and finally summarize and output the cross-watershed groundwater circulation flux; if not, treat it only as the internal circulation of the sub-watershed and do not include it in the cross-watershed flux.

[0128] For each unit, the water balance equation on a multi-year average scale is first established. Under multi-year average conditions, the change in total water storage within the unit equals the unit's precipitation minus actual evapotranspiration minus surface runoff minus lateral groundwater outflow plus the lateral groundwater inflow from adjacent units. This equation can be expressed as:

[0129]

[0130] In the formula, For the first The total water storage changes in each unit are provided by GRACE satellite data;

[0131] This refers to precipitation.

[0132] This represents the actual evaporation rate.

[0133] Surface runoff;

[0134] The lateral output of groundwater in unit i;

[0135] For from adjacent units The lateral input of groundwater.

[0136] Precipitation and actual evapotranspiration are obtained through steps S1 and S2, respectively. Surface runoff can be approximated in areas lacking measured data by using the runoff coefficient method based on digital elevation models and land use data. Lateral groundwater output and lateral groundwater input are unknowns and need to be obtained by solving a system of equations.

[0137] For the vertically dominant units identified in step S3, since their water volume changes mainly originate from the vertical infiltration of local precipitation, the impact of lateral groundwater flow from adjacent units on the water volume budget of these units is negligible. Therefore, for vertically dominant units, their lateral groundwater input can be approximated as zero relative to local vertical infiltration, and their lateral groundwater output can also be ignored accordingly.

[0138] Under these conditions, the water balance equation for this type of unit simplifies to: precipitation minus actual evapotranspiration minus surface runoff equals the change in total water storage within the unit. All variables in this simplified equation are known quantities or have been obtained through previous steps; therefore, the local vertical infiltration rate of this type of unit can be directly solved. This local vertical infiltration rate represents the net amount of groundwater recharged by precipitation in a vertically dominant unit and is an important parameter known as input in subsequent steps.

[0139] The local vertical infiltration amounts obtained from the above solutions for all vertically dominant units are extended to the entire study area using spatial interpolation methods to form a spatial background field for local vertical infiltration. This background field reflects the potential vertical recharge at various locations within the study area, considering only local precipitation infiltration. For laterally dominant units, the local vertical infiltration amount in their water budget can be initially estimated using spatial interpolation based on this background field, while the external lateral inflow amount is the unknown quantity that needs to be solved.

[0140] For laterally dominant units, a joint set of equations is established, including two unknowns: local vertical infiltration and external lateral inflow. The first equation is the complete water balance equation for this unit. Approximating the lateral groundwater output as zero, the equation can be rewritten as:

[0141]

[0142] In the formula, For flow into the unit from the adjacent unit External lateral inflow;

[0143] Let i be the local vertical infiltration rate of unit i;

[0144] The equation expresses the external lateral inflow as a functional relationship between precipitation, actual evapotranspiration, surface runoff, changes in total water storage, and local vertical infiltration.

[0145] The second equation, constructed based on the groundwater dynamic characteristic parameters calculated in step S3, constrains the relative proportions of local vertical infiltration and external lateral inflow in the total water volume of the unit. For laterally dominant units, the correlation coefficient between groundwater storage changes and precipitation is low or the phase lag is long. This characteristic reflects the degree of contribution of external lateral inflow to the water volume changes of the unit. The lower the correlation coefficient or the longer the phase lag, the greater the contribution of external lateral inflow. The constraint equations are constructed accordingly as follows:

[0146]

[0147] In the formula, For unit Correlation coefficient between changes in groundwater storage and precipitation;

[0148] This refers to the phase lag time.

[0149] It is a monotonically increasing function, and its specific form is calibrated using prior information from a small number of known points within the region.

[0150] The equation expresses the ratio of external lateral inflow to local vertical infiltration as a function of correlation coefficient and phase lag time; the lower the correlation coefficient or the longer the phase lag, the larger the ratio.

[0151] By simultaneously solving the two equations above, a system of equations with two unknowns is formed for each laterally dominant unit. The local vertical infiltration is given an initial value based on the spatial background field, while the external lateral inflow is the variable to be solved. Solving this system of equations yields the external lateral inflow for each laterally dominant unit. This external lateral inflow represents the lateral recharge of groundwater from adjacent units into that unit and is a core component of the inter-basin groundwater circulation flux.

[0152] After obtaining the external lateral inflows for all laterally dominant units, they are correlated according to their spatial location. For laterally dominant units located on both sides of the surface sub-basin boundary, the cross-basin groundwater flow path is identified based on the direction of the external lateral inflows.

[0153] The external lateral inflows crossing the surface watershed boundaries along these pathways are statistically summarized to obtain the total lateral groundwater exchange between each sub-watershed. This total exchange is used as the final inter-watershed groundwater circulation flux output, which is then used for regional water resource assessment and environmental impact analysis of inter-watershed water transfer projects.

[0154] The present invention also includes step S5, which is a result verification based on the vadose zone simulation. The aim is to obtain the simulated value of local vertical infiltration at typical points through numerical simulation of moisture transport in the vadose zone, and compare it with the estimated value of local vertical infiltration at the corresponding points obtained by inversion in step S4, thereby verifying the reliability of the estimation result of the method of the present invention.

[0155] In this embodiment, the local vertical infiltration is simulated using the Hydrus-1D vadose zone moisture transport numerical model. The model is calibrated with measured soil moisture data and driven by meteorological data of the target time series. It outputs the cumulative water flux at the bottom boundary of the vadose zone as the simulated value, which is used to compare and verify with the estimated local vertical infiltration obtained by the solution.

[0156] Hydrus-1D Vadose Zone Moisture Transport Numerical Model refers to a numerical simulation software developed by the USDA's US Salinity Laboratory, specifically designed to simulate the transport of moisture, heat, and solutes in the one-dimensional vertical vadose zone (unsaturated zone). It is a finite element numerical model based on physical processes, specifically designed to simulate the entire process of water entering the soil from the surface (infiltration), moving downwards in the soil (redistribution), being absorbed by plant roots, and ultimately replenishing groundwater (deep seepage).

[0157] Soil moisture monitoring equipment should be deployed at at least two typical locations within the study area, with at least one location located within a vertically dominant unit identified in step S3 and at least one location located within a laterally dominant unit. When selecting locations, factors such as accessibility, geomorphological representativeness, and the typicality of land use types should be comprehensively considered to ensure that the selected locations represent the hydrogeological characteristics of their respective unit types.

[0158] Soil moisture sensors are installed at each selected location. The burial depth of the sensors is determined based on the thickness of the vadose zone and the groundwater depth. Typically, no fewer than five monitoring depths are set, located at 10 cm, 40 cm, 70 cm, 100 cm, and 120 cm below the surface, to obtain complete vadose zone soil moisture profile data. Simultaneously, precipitation measurement equipment and evaporation pan measurement equipment are deployed at these locations to obtain atmospheric boundary input conditions. The reading interval for all equipment is set to 2 hours to ensure the rapid response of soil moisture after precipitation events can be captured.

[0159] The soil moisture monitoring equipment will be operated continuously for at least one full hydrological year to acquire dynamic data on soil moisture changes covering different seasons and precipitation conditions. Simultaneously, meteorological data such as air temperature, wind speed, and relative humidity will be collected during the monitoring period as input for numerical simulation. The collected measured data will be divided into two parts: the observation data from the first two-thirds of the time period will be used for model calibration, and the observation data from the latter one-third of the time period will be used for model validation.

[0160] Numerical models of moisture transport in the vadose zone were established using Hydrus-1D software at various typical locations. The upper boundary of the model was set as a water-accumulating atmospheric boundary, and the actual infiltration and evaporation fluxes were calculated based on measured precipitation and pan evaporation data. The lower boundary of the model was determined based on the groundwater depth at each location. When the groundwater depth was greater than 5 meters, the lower boundary was set as a free drainage boundary, assuming that moisture at the bottom of the vadose zone could freely drain into the groundwater system under gravity.

[0161] The model depth is determined based on the actual vadose zone thickness, typically set to 150-200 cm to ensure complete coverage of the soil moisture sensor's monitoring profile. The model is vertically divided into several layers with a spatial step size of 1 cm to guarantee the stability and accuracy of numerical calculations. Soil hydraulic parameters are adopted... The model description includes five parameters: residual moisture content, saturated moisture content, reciprocal of the inlet air value, pore connectivity parameter, and saturated permeability coefficient.

[0162] The model is a numerical simulation software developed by the US Salinity Laboratory of the U.S. Department of Agriculture (USDA) specifically for simulating the transport of moisture, heat, and solutes in a one-dimensional vertical vadose zone (unsaturated zone).

[0163] The model parameters were calibrated based on measured soil moisture data. The objective function was to minimize the root mean square error between the measured moisture content and the simulated moisture content at each depth during the monitoring period. An automatic inversion algorithm was then used to calibrate the model parameters. The five parameters in the model were optimized and adjusted. During the calibration process, the measured moisture content data for the first two-thirds of the monitoring period were used as the fitting target. The model parameters were continuously adjusted to achieve the best fit between the simulated and measured moisture content curves. After calibration, the calibrated model was verified using the measured moisture content data for the remaining one-third of the monitoring period, comparing the degree of agreement between the simulated and measured values ​​at each depth. If the correlation coefficient between the simulated and measured moisture content at each depth is greater than the preset threshold, it indicates that the calibrated model can reliably reflect the moisture transport pattern in the vadose zone at that location.

[0164] The calibrated and validated model is applied to the simulation calculation of the target time series. The year consistent with the target time series in steps S1 to S4 is selected as the simulation period. Daily precipitation, daily pan evaporation, and daily meteorological data for that year are input into the model to drive the simulation of moisture transport in the vadose zone. After the model runs, the cumulative moisture flux at the bottom boundary of the vadose zone within that year is output. This cumulative moisture flux is the simulated value of local vertical infiltration at that location on a multi-year average scale.

[0165] The simulated values ​​of local vertical infiltration at each typical location obtained from the above simulation are compared with the estimated values ​​of local vertical infiltration at the corresponding location obtained in step S4. If the relative error between the two is within a preset threshold range, for example, the relative error is less than 15%, it indicates that the local vertical infiltration estimated by the method of the present invention has high reliability, thereby verifying the effectiveness of the external lateral inflow inversion results in step S4.

[0166] If the relative error between the two exceeds the preset threshold range, it is necessary to return to step S4 to adjust the function parameters in the constraint equations and resolve the equations for the laterally dominant element until the verification result meets the accuracy requirements. Through the above closed-loop verification mechanism, the cross-basin groundwater circulation flux estimated by the method of this invention is reliably accurate.

[0167] Compared with the prior art, the embodiments of the present invention have the following technical advantages:

[0168] 1. This invention, by introducing a gridded unit division and a groundwater dynamic response characteristic parameter identification mechanism, finely divides the study area into vertically dominant units and laterally dominant units, thereby locating the spatial position of the inter-basin groundwater exchange zone. Compared with existing technologies that treat the entire basin as a whole unit and simply regard net recharge as the residual amount of vertical precipitation infiltration, this invention can effectively distinguish between the vertical and lateral contributions of water sources, avoiding the systematic bias of erroneously attributing external lateral inflows to local vertical infiltration, and providing a clear diagnostic basis for the spatial location of inter-basin groundwater circulation flux.

[0169] 2. This invention achieves quantitative decoupling of vertical infiltration and lateral inflow by constructing a joint equation system with the local vertical infiltration background field as the initial value and correlation coefficient and phase lag time as proportional constraints. Compared with the limitations of existing technologies that rely solely on a single water balance equation and cannot separate two unknowns, this invention introduces groundwater dynamic characteristics as an additional constraint to fill the information gap, making the inversion of external lateral inflow possible. Simultaneously, by setting a boundary screening mechanism within the inter-basin groundwater exchange zone, only external lateral inflow crossing the surface watershed boundary is statistically analyzed, avoiding interference from lateral circulation within sub-basins on inter-basin flux estimation, thus obtaining a more physically meaningful inter-basin groundwater circulation flux.

[0170] Example 2

[0171] This embodiment discloses a remote sensing estimation system for groundwater flux that separates infiltration and lateral inflow.

[0172] A remote sensing estimation system for groundwater flux that separates infiltration and lateral inflow includes:

[0173] The data acquisition module is configured to acquire a multi-source remote sensing dataset of the study area, which includes land water storage change data, remote sensing precipitation data, remote sensing evapotranspiration data, digital elevation model data, and land use type data.

[0174] The zoning correction module is configured to: divide the study area into vegetated areas and bare soil areas according to land use type data; calculate the evapotranspiration correction coefficients for vegetated areas and bare soil areas respectively based on the water balance principle; apply the evapotranspiration correction coefficients to the remote sensing evapotranspiration data of the corresponding areas; and obtain the corrected spatial distribution data of actual evapotranspiration.

[0175] The unit division module is configured to: divide the study area into units of equal area, extract the time series of changes in terrestrial water storage and precipitation for each unit, and calculate the correlation coefficient and phase lag time between the two.

[0176] The type identification module is configured to: classify the corresponding unit into vertically dominant units or laterally dominant units based on the correlation coefficient and phase lag time of each unit; extract the surface sub-basin boundary using the digital elevation model; and identify the strip-shaped area composed of laterally dominant units belonging to different surface sub-basins as the cross-basin groundwater exchange zone.

[0177] The joint solution module is configured to: use the vertically dominant element to solve the local vertical infiltration and interpolate to form a background field; combine the corrected spatial distribution data of actual evapotranspiration to establish a joint equation system for the laterally dominant element with the background field as the initial value and the correlation coefficient and phase lag time as constraints, and solve for the external lateral inflow.

[0178] The aggregation module is configured to: statistically aggregate the external lateral inflows located within the inter-basin groundwater exchange zone and crossing the surface watershed boundary to obtain the inter-basin groundwater circulation flux.

[0179] Those skilled in the art will understand that the modules or steps of the present invention described above can be implemented using general-purpose computer devices. Optionally, they can be implemented using computer-executable program code, thereby allowing them to be stored in a storage device for execution by a computer device, or they can be fabricated as separate integrated circuit modules, or multiple modules or steps can be fabricated as a single integrated circuit module. The present invention is not limited to any particular combination of hardware and software.

[0180] While the specific embodiments of the present invention have been described above in conjunction with the accompanying drawings, this is not intended to limit the scope of protection of the present invention. Those skilled in the art should understand that various modifications or variations that can be made by those skilled in the art without creative effort based on the technical solutions of the present invention are still within the scope of protection of the present invention.

Claims

1. A remote sensing estimation method for groundwater flux that separates infiltration and lateral inflow, characterized in that, Includes the following steps: Acquire a multi-source remote sensing dataset for the study area, which includes land water storage change data, remote sensing precipitation data, remote sensing evapotranspiration data, digital elevation model data, and land use type data; Based on land use type data, the study area was divided into vegetated area and bare soil area. Based on the principle of water balance, the evapotranspiration correction coefficients of vegetated area and bare soil area were calculated respectively. The evapotranspiration correction coefficients were applied to the remote sensing evapotranspiration data of the corresponding area to obtain the corrected spatial distribution data of actual evapotranspiration. The study area was divided into units of equal area. For each unit, time series of land water storage changes and precipitation were extracted, and the correlation coefficient and phase lag time were calculated. The correlation coefficient is the Pearson correlation coefficient between the land water storage change time series and the precipitation time series within the unit. The phase lag time is the time delay of the land water storage change time series relative to the precipitation time series, which is calculated through cross-correlation analysis and is used to reflect the time interval between precipitation events and groundwater storage response. Based on the correlation coefficient and phase lag time of each unit, the corresponding units are divided into vertically dominant units or laterally dominant units. The boundaries of surface sub-basins are extracted using digital elevation models, and the strip-shaped area composed of laterally dominant units belonging to different surface sub-basins is identified as the cross-basin groundwater exchange zone. The local vertical infiltration rate is solved using vertically dominant units and interpolated to form a background field. Combined with the corrected spatial distribution data of actual evapotranspiration, a joint equation system is established for laterally dominant units with the background field as the initial value and the correlation coefficient and phase lag time as constraints. The external lateral inflow rate is then obtained by solving the equation system. The external lateral inflows located within the inter-basin groundwater exchange zone and crossing the surface watershed boundary are statistically summarized to obtain the inter-basin groundwater circulation flux.

2. The remote sensing estimation method for groundwater flux separation of infiltration and lateral inflow as described in claim 1, characterized in that, In the multi-source remote sensing dataset: The terrestrial water storage change data are from the GRACE gravity satellite. The specific parameters are equivalent water height, monthly temporal resolution, and 0.25º×0.25º spatial resolution. They are used to reflect the monthly changes in the total storage of surface water, soil water, and groundwater in the study area. They are a key input in the subsequent water balance calculation and provide a data foundation for the time series of terrestrial water storage changes. The specific parameters of the remote sensing precipitation data are the total monthly precipitation, the time resolution is monthly, and the spatial resolution is 1 kilometer, which are used to provide a data basis for precipitation time series. The remote sensing evapotranspiration data are from the MOD16 series evapotranspiration products. The specific parameters are the actual evapotranspiration, the temporal resolution is yearly, and the spatial resolution is 500 meters. The specific parameters of the digital elevation model data are the ground elevation and the spatial resolution is 90 meters. The specific parameters of the land use type data are the land cover classification type code and the spatial resolution is 500 meters. It is used to distinguish between vegetated areas, bare soil areas and water areas, and provides a basis for the zoning of the study area.

3. The remote sensing estimation method for groundwater flux separation of infiltration and lateral inflow as described in claim 2, characterized in that, For vegetated and bare soil areas, the ratio of multi-year average precipitation to original MOD16 evapotranspiration in each area is used as the evapotranspiration correction factor, where: Evapotranspiration correction factor for vegetated areas Calculation formula: in, This refers to the multi-year average precipitation in the locally balanced vegetation zone. The original MOD16 evapotranspiration is the statistical average value of the vegetation area. The original MOD16 evapotranspiration is the remote sensing evapotranspiration data before correction.

4. The remote sensing estimation method for groundwater flux separation of infiltration and lateral inflow as described in claim 3, characterized in that, Evapotranspiration correction factor for bare soil area The calculation formula is: In the formula, This represents the multi-year average precipitation in the locally balanced zone of the bare soil area. The average value of the original MOD16 evapotranspiration in the bare soil area.

5. The remote sensing estimation method for groundwater flux separation of infiltration and lateral inflow as described in claim 4, characterized in that, Based on the correlation coefficient and phase lag time of each element, the corresponding elements are divided into vertically dominant elements or laterally dominant elements, specifically including: If the correlation coefficient of a certain unit is greater than the first preset threshold and the phase lag time is less than the second preset threshold, then it is a vertically dominant unit. If the correlation coefficient of a certain unit is less than the first preset threshold or the phase lag time is greater than the second preset threshold, then it is a laterally dominant unit. The first preset threshold is determined as follows: Calculate the correlation coefficients of all units within the study area and plot a frequency distribution histogram of the correlation coefficients. Since vertically dominant units and laterally dominant units exhibit a bimodal distribution characteristic in statistics, the valley value between the two peaks is taken as the first preset threshold. The second preset threshold is determined as follows: The top 20% of the units with the highest correlation coefficients were selected as quasi-vertical units; The phase lag time distribution of the quasi-vertical unit is statistically analyzed, and the upper quartile is calculated as the second preset threshold.

6. The remote sensing estimation method for groundwater flux separation of infiltration and lateral inflow as described in claim 1, characterized in that, The correlation coefficient is used to measure the degree of response of a unit's water quantity changes to local precipitation input. The specific formula for calculating the correlation coefficient is as follows: ； In the formula: The correlation coefficient between changes in terrestrial water storage and precipitation in the grid cell; For the first The change in terrestrial water storage in the current unit over the past month; For the first The current month's precipitation in the current unit; This represents the monthly average change in land water storage in the current unit during the research period; This represents the monthly average precipitation for the current unit during the research period. The total number of months in the time series. Indicates the first Months.

7. The remote sensing estimation method for groundwater flux separation of infiltration and lateral inflow as described in claim 6, characterized in that, The specific formula for calculating phase lag time is as follows: Constructing a time series of changes in terrestrial water storage With precipitation time series cross-correlation function The formula is: ； In the formula, The number of months is the lag. For lag Monthly change in terrestrial water storage; Calculate different lag months Below Value, find Number of months lag before reaching the maximum value Then the phase lag time Defined as .

8. The remote sensing estimation method for groundwater flux separation of infiltration and lateral inflow as described in claim 1, characterized in that, The local vertical infiltration rate is solved using vertically dominant elements and interpolated to form the background field, specifically including: The lateral groundwater input of the vertically dominant unit is set to zero; Ignore the lateral groundwater output of vertically dominant units; The water balance equation for a vertically dominant unit is: precipitation minus actual evapotranspiration minus surface runoff equals the change in total water storage within the unit, specifically: ； in, This refers to precipitation. This represents the actual amount of evaporation. Surface runoff; For the first Changes in total water storage in each unit; Finally, the local vertical infiltration rate of the vertically dominant unit is directly calculated.

9. The remote sensing estimation method for groundwater flux separation of infiltration and lateral inflow as described in claim 8, characterized in that, The combined equations include the water balance equation and the proportionality constraint equation, where: The water balance equation is as follows: In the formula, For flow into the unit from the adjacent unit External lateral inflow; Let i be the local vertical infiltration rate of unit i; The proportional constraint equation is as follows: In the formula, For unit Correlation coefficient between changes in groundwater storage and precipitation; This refers to the phase lag time. It is a monotonically increasing function.

10. A remote sensing estimation system for groundwater flux that separates infiltration and lateral inflow, characterized in that, include: The data acquisition module is configured to acquire a multi-source remote sensing dataset of the study area, which includes land water storage change data, remote sensing precipitation data, remote sensing evapotranspiration data, digital elevation model data, and land use type data. The zoning correction module is configured to: divide the study area into vegetated areas and bare soil areas according to land use type data; calculate the evapotranspiration correction coefficients for vegetated areas and bare soil areas respectively based on the water balance principle; apply the evapotranspiration correction coefficients to the remote sensing evapotranspiration data of the corresponding areas; and obtain the corrected spatial distribution data of actual evapotranspiration. The unit division module is configured to: divide the study area into units of equal area; extract the land water storage change time series and precipitation time series for each unit; and calculate the correlation coefficient and phase lag time between the two. The correlation coefficient is the Pearson correlation coefficient between the land water storage change time series and the precipitation time series within the unit; the phase lag time is the time delay of the land water storage change time series relative to the precipitation time series, which is calculated through cross-correlation analysis and is used to reflect the time interval between precipitation events and groundwater storage responses. The type identification module is configured to: classify the corresponding unit into vertically dominant units or laterally dominant units based on the correlation coefficient and phase lag time of each unit; extract the surface sub-basin boundary using the digital elevation model; and identify the strip-shaped area composed of laterally dominant units belonging to different surface sub-basins as the cross-basin groundwater exchange zone. The joint solution module is configured to: use the vertically dominant element to solve the local vertical infiltration and interpolate to form a background field; combine the corrected spatial distribution data of actual evapotranspiration to establish a joint equation system for the laterally dominant element with the background field as the initial value and the correlation coefficient and phase lag time as constraints, and solve for the external lateral inflow. The aggregation module is configured to: statistically aggregate the external lateral inflows located within the inter-basin groundwater exchange zone and crossing the surface watershed boundary to obtain the inter-basin groundwater circulation flux.

Images