A method for calculating water depth and storage of pit and pond by combining multi-source data

Abstract

The application discloses a kind of pit pond water depth estimation and storage capacity calculation method of combined multi-source data, comprising the following steps: collecting measured water depth data, data preprocessing is carried out, and is divided into training sample point set and test sample point set;Obtain the multispectral remote sensing image of the estimation area, and execute image preprocessing;Unified to the same coordinate system;The fractal dimension of each band and combined band is calculated, and the water area is segmented;Establish the mapping relationship between the water depth value of training sample point and standardization multispectral apparent reflectance;Establish water depth inversion model, and use test sample point set to verify water depth inversion model;The discrete water depth point of pit pond is calculated using optimal water depth inversion model;Establish minimum energy surface, obtain single pit pond water storage capacity, and accumulate to generate total water storage capacity.The application fuses measured and remote sensing data, constructs energy minimization water depth inversion model, realizes pit pond water storage capacity accurate calculation, has the advantages of high precision, strong applicability.

Description

Technical Field

[0001] This invention relates to the field of hydrological and water resources monitoring technology, and in particular to a method for estimating the depth and storage capacity of ponds using multi-source data. Background Technology

[0002] Currently, small water bodies such as ponds play an important role in agricultural irrigation, water resource allocation, and ecological environment management. Traditional methods for obtaining the depth and storage capacity of ponds mainly rely on manual measurement or estimation methods based on a single data source, such as using depth sounders to collect depth point data or indirectly estimating based on the water body extent determined by remote sensing imagery. These methods have certain limitations in practical applications: manual measurement is time-consuming and labor-intensive, making it difficult to carry out efficiently over large areas; while single remote sensing data processing methods are easily constrained by lighting conditions, atmospheric influences, and the complex spectral characteristics of water bodies, resulting in insufficient accuracy of the results.

[0003] Existing technologies for calculating the spatial distribution of water depth and storage capacity in ponds generally lack the fusion and coordination of multi-source data, making it difficult to simultaneously ensure measurement accuracy and regional coverage. Furthermore, most existing methods focus on simple regression modeling or empirical formula calculations, neglecting the precise identification of water body boundaries, the systematic construction of water depth factors, and spatial continuity constraints, which can easily lead to deviations in water depth estimation and errors in storage capacity estimation. Therefore, there is an urgent need for a comprehensive method that can integrate measured water depth data with remote sensing multi-source information to improve the accuracy and applicability of pond water depth inversion and storage capacity calculation. Summary of the Invention

[0004] One objective of this invention is to propose a method for estimating the water depth and storage capacity of ponds by combining multi-source data. This invention integrates measured and remote sensing data to construct an energy-minimizing water depth inversion model, thereby achieving accurate estimation of pond water storage capacity. It has the advantages of high accuracy and strong applicability.

[0005] A method for estimating pond water depth and calculating storage capacity based on combined multi-source data according to an embodiment of the present invention includes the following steps:

[0006] Collect measured water depth data of pits and ponds in the estimation area, perform data preprocessing, generate standardized water depth data, and divide it into training sample point set and test sample point set;

[0007] Acquire multispectral remote sensing images of the extrapolated area and perform image preprocessing to generate standardized multispectral apparent reflectance images;

[0008] Unify standardized water depth data and standardized multispectral apparent reflectance images into the same coordinate system;

[0009] Based on standardized multispectral apparent reflectance images, the fractal dimension of each band and combined bands is calculated. The fractal thresholds of water bodies and non-water bodies are determined using training sample point sets. Water body regions are segmented according to the fractal thresholds to generate pond surface boundary vectors.

[0010] Establish the mapping relationship between the water depth values ​​of the training sample points and the standardized multispectral apparent reflectance, and construct the water depth inversion factor vector;

[0011] Using the water depth inversion factor vector and the water depth values ​​of the training sample points as input, a water depth inversion model is established, and the water depth inversion model is validated using a set of test sample points to generate the optimal water depth inversion model.

[0012] The discrete water depth points of the pit are calculated using the optimal water depth inversion model, and a discrete water depth point dataset is generated.

[0013] Using the boundary vector of the pond water surface as the boundary condition and the discrete water depth point dataset as the interior point constraint, a minimum energy surface is established. Volume integration is performed on the minimum energy surface to obtain the water storage of a single pond. The water storage of all ponds is then summed to generate the total water storage of the estimated area.

[0014] Optionally, the partitioning of the training sample point set and the test sample point set specifically includes:

[0015] Collect raw measured water depth data in the pits and ponds in the estimation area, record the plane coordinates, collection time and water depth value of each depth measurement point, and form a set of raw measured water depth data.

[0016] Outlier removal was performed on the original measured water depth data set;

[0017] Perform duplicate value removal on the water depth data set after outlier removal;

[0018] The deduplicated water depth data set is standardized by using the average and standard deviation of the water depth values ​​in the set as a benchmark to standardize each record and generate a standardized water depth data set, which includes plane coordinates, acquisition time and standardized water depth values.

[0019] The standardized water depth data set is divided into training sample sets and test sample set according to the training ratio.

[0020] Optionally, the image preprocessing includes radiometric correction, atmospheric correction, geometric correction, orthorectification, image fusion, and noise removal.

[0021] Optionally, the generation of the boundary vector of the pond water surface specifically includes:

[0022] A band set and a combined band set are established on a standardized multispectral apparent reflectance image. The apparent reflectance of each pixel in a single band is recorded directly. For each combined band, the weighted average of the apparent reflectance of each single band is used as the apparent reflectance of the combined band. The band set includes blue light band, green light band, red light band and near-infrared band.

[0023] A multi-scale window sequence is constructed with the pixel as the center. At each scale, the number of covering boxes is calculated using the differential box counting method, and the relationship between the number of covering boxes and the scale is recorded. By fitting the slope of the straight line between the logarithm of the number of covering boxes and the logarithm of the scale, the fractal dimension of the pixel is generated.

[0024] The fractal dimensions of all single-band and combined-band data are summarized to construct the fractal feature vector of the pixel. The fractal feature vector is then matched with the training sample point set at the same coordinate position. Based on the water body label and non-water body label of the training sample point set, the fractal feature vector is divided into the water body set and the non-water body set.

[0025] For each pixel's fractal feature component, calculate the mean and variance of the water body set and the non-water body set respectively. Take the ratio of the mean difference to the variance sum as the discrimination criterion, and select the largest value of the discrimination criterion as the fractal threshold of the fractal feature component to obtain the fractal threshold set.

[0026] Based on the fractal threshold set, a pixel comprehensive judgment value is constructed using the fusion coefficients corresponding to each fractal feature component. A water body region mask is generated based on the comparison result between the pixel comprehensive judgment value and the judgment threshold. If the pixel comprehensive judgment value is greater than or equal to the judgment threshold, it is judged as a water body pixel. If the pixel comprehensive judgment value is less than the judgment threshold, it is judged as a non-water body pixel. In the water body region mask, water body pixels are marked as 1 and non-water body pixels are marked as 0.

[0027] Perform connected component screening and morphological processing on the water body mask, delete connected components with an area smaller than the minimum area threshold and perform closing operation to obtain water body regions with satisfactory connectivity and smooth boundaries.

[0028] Boundary tracing is performed on the outer boundary of the water body area to obtain a sequence of boundary polygons. The boundary polylines are then simplified with a simplification tolerance to generate vertex coordinates. An ordered vector boundary is constructed according to the vertex order to generate the boundary vector of the pond water surface.

[0029] Optionally, the construction of the water depth inversion factor vector specifically includes:

[0030] Under a unified coordinate system, the standardized water depth data and the standardized multispectral apparent reflectance image are spatially correlated to determine the corresponding pixel position of each training sample point on the standardized multispectral apparent reflectance image.

[0031] At the corresponding pixel location of the training sample point, extract the apparent reflectance of each single band, the apparent reflectance of the combined band, the difference factor between the two bands, and the ratio factor between the two bands, and correlate them with the water depth value of the training sample point to establish the mapping relationship between the water depth value and the apparent reflectance of the training sample point.

[0032] The mapping results of each training sample point are organized, and for each training sample point, a factor set is formed consisting of single-band apparent reflectance, combined-band apparent reflectance, difference factor between two bands, and ratio factor between two bands.

[0033] For each training sample point, its factor set is arranged in a fixed order to form the water depth inversion factor vector of that training sample point. The water depth inversion factor vectors of all training sample points are summarized to form the water depth inversion factor vector set.

[0034] Optionally, the generation of the optimal water depth inversion model specifically includes:

[0035] Based on the set of water depth inversion factor vectors and the water depth values ​​of training sample points, the training sample point set and the test sample point set are distinguished. For each training sample point in the training sample point set, the water depth inversion factor vector is recorded as the input variable, and the water depth value is recorded as the output variable.

[0036] Arrange the water depth inversion factor vectors of all training sample points in a fixed order to form the training sample point design matrix, and arrange the water depth values ​​of the training sample points in the same order to form the training sample point target vector.

[0037] The water depth inversion model is defined as an energy-minimizing water depth inversion model based on variational constraints. The energy functional of the energy-minimizing water depth inversion model consists of three parts: the first part is the weighted sum of the absolute values ​​of the differences between the predicted water depth values ​​of the training sample points and the actual water depth values ​​of the training sample points, which serves as a data consistency term; the second part is the sum of the squares of the differences between the predicted water depth values ​​of adjacent training sample points, which serves as a spatial smoothing constraint term; and the third part is the product of the sum of the squares of the parameter vectors corresponding to the water depth inversion model and the regularization coefficient, which serves as a regularization constraint term. The three parts are added together to form the energy functional.

[0038] An iterative optimization method based on gradient descent is adopted to determine the parameter vector of the water depth inversion model by minimizing the energy functional, and the parameter estimation results under different regularization coefficients are obtained.

[0039] The parameter estimation results are verified using a set of test sample points. For each test sample point, its water depth inversion factor vector is input, the predicted water depth value is output, and it is compared with the water depth value of the test sample point. The sum of the absolute values ​​of the differences between the predicted values ​​of all test sample points and the water depth values ​​of the test sample points is calculated as the verification error measure.

[0040] Within a pre-defined set of candidate regularization coefficients, the verification error metric corresponding to each regularization coefficient is calculated one by one. The regularization coefficient with the smallest verification error metric is selected, and the corresponding parameter estimate is used as the parameter of the optimal water depth inversion model to generate the optimal water depth inversion model.

[0041] Optionally, the generation of the estimated total water storage in the region specifically includes:

[0042] Under a unified coordinate system, the calculation area of ​​a single pond is defined by the pond water surface boundary vector, and the calculation area is triangulated in a plane to obtain a grid set composed of the pond water surface boundary vector and the discrete water depth point dataset.

[0043] Apply boundary conditions to the boundary vector of the pond water surface, fix the value of the continuous water depth function at the boundary to zero, and make the pond water surface boundary a constraint boundary.

[0044] Apply interior point constraints to each discrete water depth point in the discrete water depth point dataset that is located inside the computational region, so that the value of the continuous water depth function at the discrete water depth point is equal to the predicted water depth value of the discrete water depth point.

[0045] The objective function for constructing the minimum energy surface is composed of a membrane energy term and a thin plate bending energy term. The membrane energy term measures the change of the first-order partial derivative of the continuous water depth function in the computational domain, and the thin plate bending energy term measures the change of the second-order partial derivative of the continuous water depth function in the computational domain. The objective function is a weighted sum of the two parts.

[0046] Using the minimum objective function as the criterion, and simultaneously satisfying boundary conditions and interior point constraints, the minimum energy surface for a single pit is obtained;

[0047] The minimum energy surface is discretized on a mesh set using the finite element method. The continuous water depth function is represented by a linear combination of nodal basis functions and assembled into a discrete system of equations. The water depth values ​​of the mesh nodes are obtained by solving the system of equations.

[0048] After obtaining the water depth values ​​of the grid nodes, volume integration is performed on the computational domain, and discrete integration is performed using the vertex averaging method of the grid cells to obtain the water storage capacity of a single pit.

[0049] Calculate the water storage capacity of all ponds within the estimated area, and sum them up to generate the total water storage capacity of the estimated area.

[0050] The beneficial effects of this invention are:

[0051] This invention establishes a method for estimating pond water depth and calculating storage capacity by integrating measured water depth data with multispectral remote sensing image data. This method overcomes the limitations of traditional methods that rely on a single data source or manual measurement. In terms of technical approach, this method fully leverages the complementary advantages of high accuracy of measured data and wide coverage of remote sensing data. Through data preprocessing, standardization, and coordinate unification, data from different sources are incorporated into the same analytical framework, ensuring the scientific validity and consistency of the estimation results. At the same time, fractal dimension discrimination and boundary vector generation mechanisms are introduced during the water body extraction process, enabling precise identification of pond water body boundaries and accurate expression of geometric morphology, effectively improving the accuracy and reliability of water body area delineation.

[0052] Secondly, regarding water depth inversion, this invention constructs a systematic water depth inversion factor vector, which not only includes the apparent reflectance of single-band and combined-band data but also introduces difference and ratio factors, thus forming a richer and more comprehensive feature set. Combining a variational constraint energy minimization modeling method, and introducing data consistency, spatial smoothing, and regularization constraints, the model can ensure accuracy while maintaining spatial continuity and complexity control, avoiding the overfitting and spatial abrupt change problems of traditional regression models. Through gradient descent iterative optimization and test sample validation, the optimal regularization coefficient can be automatically selected, thereby obtaining the optimal water depth inversion model, significantly improving the accuracy and stability of water depth prediction.

[0053] Finally, in the process of calculating water storage, this invention introduces a minimum energy surface and the finite element method to organically combine water body boundary conditions with discrete water depth point constraints, establishing a continuous water depth distribution function within the pit / pond area. This surface not only satisfies the physical constraint of zero water depth at the boundary but also rigorously approximates the predicted point values ​​and achieves a smooth transition between the boundary and interior points, effectively reproducing the spatial variation law of water depth within the pit / pond. Based on this, the water storage of a single pit / pond is accurately calculated through volume integration, and the results of all pits / ponds within the region are accumulated to obtain the total storage of the estimated area. This method avoids errors caused by boundary simplification or inaccurate water depth interpolation in traditional volume estimation, improving the scientific rigor and reliability of the calculation. Attached Figure Description

[0054] The accompanying drawings are provided to further illustrate the invention and form part of the specification. They are used in conjunction with embodiments of the invention to explain the invention and do not constitute a limitation thereof. In the drawings:

[0055] Figure 1 This is a flowchart of a method for estimating pond water depth and calculating storage volume using combined multi-source data, as proposed in this invention.

[0056] Figure 2This diagram illustrates a method for calculating the depth and storage capacity of ponds using multi-source data proposed in this invention. It involves discretizing the minimum energy surface using the finite element method to calculate the water storage capacity of the pond. Detailed Implementation

[0057] The present invention will now be described in further detail with reference to the accompanying drawings. These drawings are simplified schematic diagrams, illustrating only the basic structure of the invention, and therefore only show the components relevant to the invention.

[0058] refer to Figure 1-2 A method for estimating pond water depth and calculating storage capacity using multi-source data, comprising the following steps:

[0059] Collect measured water depth data of pits and ponds in the estimation area, perform data preprocessing, generate standardized water depth data, and divide it into training sample point set and test sample point set;

[0060] Acquire multispectral remote sensing images of the extrapolated area and perform image preprocessing to generate standardized multispectral apparent reflectance images;

[0061] Unify standardized water depth data and standardized multispectral apparent reflectance images into the same coordinate system;

[0062] Based on standardized multispectral apparent reflectance images, the fractal dimension of each band and combined bands is calculated. The fractal thresholds of water bodies and non-water bodies are determined using training sample point sets. Water body regions are segmented according to the fractal thresholds to generate pond surface boundary vectors.

[0063] Establish the mapping relationship between the water depth values ​​of the training sample points and the standardized multispectral apparent reflectance, and construct the water depth inversion factor vector;

[0064] Using the water depth inversion factor vector and the water depth values ​​of the training sample points as input, a water depth inversion model is established, and the water depth inversion model is validated using a set of test sample points to generate the optimal water depth inversion model.

[0065] The discrete water depth points of the pit are calculated using the optimal water depth inversion model, and a discrete water depth point dataset is generated.

[0066] Using the boundary vector of the pond water surface as the boundary condition and the discrete water depth point dataset as the interior point constraint, a minimum energy surface is established. Volume integration is performed on the minimum energy surface within the boundary range of the pond water surface to obtain the water storage of a single pond. The water storage of all ponds is then summed to generate the total water storage of the estimated area.

[0067] In this embodiment, the division of the training sample point set and the test sample point set specifically includes:

[0068] Collect raw measured water depth data in the pits and ponds in the estimation area, record the plane coordinates, collection time and water depth value of each depth measurement point, and form a set of raw measured water depth data.

[0069] Outlier removal was performed on the original measured water depth data set;

[0070] Perform duplicate value removal on the water depth data set after outlier removal;

[0071] The deduplicated water depth data set is standardized by using the average and standard deviation of the water depth values ​​in the set as a benchmark to standardize each record and generate a standardized water depth data set, which includes plane coordinates, acquisition time and standardized water depth values.

[0072] The standardized water depth data set is divided into training sample sets and test sample set according to the training ratio.

[0073] In this embodiment, the image preprocessing includes radiometric correction, atmospheric correction, geometric correction, orthorectification, image fusion, and noise removal.

[0074] In this embodiment, the generation of the boundary vector of the pond water surface specifically includes:

[0075] A band set and a combined band set are established on a standardized multispectral apparent reflectance image. The apparent reflectance of each pixel in a single band is recorded directly. For each combined band, the weighted average of the apparent reflectance of each single band is used as the apparent reflectance of the combined band. The band set includes blue light band, green light band, red light band and near-infrared band.

[0076] A multi-scale window sequence is constructed with the pixel as the center. At each scale, the number of covering boxes is calculated using the differential box counting method, and the relationship between the number of covering boxes and the scale is recorded. By fitting the slope of the straight line between the logarithm of the number of covering boxes and the logarithm of the scale, the fractal dimension of the pixel is generated.

[0077] The generation of the fractal dimension of the pixel specifically includes:

[0078] In standardized multispectral apparent reflectance images, a multi-scale window sequence centered on each pixel is established, with the window side length taking odd numbers of pixels and increasing step by step.

[0079] Within each scale window, the differential box counting method is applied to the apparent reflectance of single bands and the apparent reflectance of combined bands respectively. The apparent reflectance within the scale window is linearly normalized to a fixed gray level range according to the minimum and maximum apparent reflectance, and divided into several intensity layers.

[0080] The scale window is divided into equal-sized lattices. The maximum and minimum layer numbers of the intensity layers are read in each lattice, and the difference between the two plus one is taken as the number of cover boxes in that lattice. The number of cover boxes in all lattices under the scale window is summed to obtain the total number of cover boxes.

[0081] For all scales, the correspondence between the number of cover boxes and the scale is recorded sequentially. The least squares line is fitted with the logarithm of the number of cover boxes as the ordinate and the logarithm of the scale as the abscissa. The slope of the fitted line is used as the fractal dimension of the target pixel under the apparent reflectance of the single band or the apparent reflectance of the combined band.

[0082] When a multi-scale window goes out of bounds, it is filled by neighborhood filling. When the window contains invalid pixels with cloud or shadow masks, the current scale calculation of that pixel is skipped. When the number of valid scales is insufficient, the fractal dimension of that pixel is not output.

[0083] Finally, the fractal dimension was obtained for all single-band apparent reflectance and combined-band apparent reflectance of each pixel.

[0084] The fractal dimensions of all single-band and combined-band data are summarized to construct the fractal feature vector of each pixel. This vector is then mapped to the training sample point set at the same coordinate position. Based on the water body and non-water body markings of the training sample point set, the fractal feature vector is divided into a water body set and a non-water body set. The water body and non-water body markings of the training sample point set represent a binary classification of the sample points based on the measured water depth value, indicating whether the sample point belongs to a water body area. If the measured water depth value of the sample point is greater than zero, the sample point is in a water body area and is marked as a water body. If the measured water depth value of the sample point is equal to zero or the sample point is in a non-water body location, it is marked as a non-water body.

[0085] For each pixel's fractal feature component, the mean and variance of the water body set and the non-water body set are calculated respectively. The ratio of the mean difference to the variance sum is taken as the discrimination criterion, and the largest value of the discrimination criterion is selected as the fractal threshold of the fractal feature component, thus obtaining the fractal threshold set. The fractal feature component is an element of the fractal feature vector.

[0086] Based on the fractal threshold set, a pixel comprehensive judgment value is constructed using the fusion coefficients corresponding to each fractal feature component. A water body region mask is generated based on the comparison result between the pixel comprehensive judgment value and the judgment threshold. If the pixel comprehensive judgment value is greater than or equal to the judgment threshold, it is judged as a water body pixel. If the pixel comprehensive judgment value is less than the judgment threshold, it is judged as a non-water body pixel. In the water body region mask, water body pixels are marked as 1 and non-water body pixels are marked as 0.

[0087] Perform connected component screening and morphological processing on the water body mask, delete connected components with an area smaller than the minimum area threshold and perform closing operation to obtain water body regions with satisfactory connectivity and smooth boundaries.

[0088] Boundary tracing is performed on the outer boundary of the water body region to obtain a sequence of boundary polygons. The boundary polylines are then simplified with a simplification tolerance to generate vertex coordinates. An ordered vector boundary is constructed according to the vertex order to generate the pond water surface boundary vector. Specifically, boundary tracing is performed on the outer boundary of the water body region, and the boundary positions of adjacent water body pixels are recorded sequentially to obtain a continuous sequence of boundary polygons. Based on the boundary polygon sequence, a simplification tolerance is set to simplify the boundary polylines. Redundant vertices are reduced while ensuring accurate representation of the boundary geometry. Vertex coordinates are generated for the simplified boundary polyline sequence. An ordered vector boundary is constructed according to the vertex order. Finally, the ordered vector boundaries of each water body region are summarized to generate the pond water surface boundary vector. The simplification tolerance represents the maximum allowable geometric deviation of the water body boundary vector polylines during the simplification process, used to balance accuracy and data complexity.

[0089] In this embodiment, the construction of the water depth inversion factor vector specifically includes:

[0090] Under a unified coordinate system, the standardized water depth data and the standardized multispectral apparent reflectance image are spatially correlated to determine the corresponding pixel position of each training sample point on the standardized multispectral apparent reflectance image.

[0091] At the corresponding pixel location of the training sample point, extract the apparent reflectance of each single band, the apparent reflectance of the combined band, the difference factor between the two bands, and the ratio factor between the two bands, and correlate them with the water depth value of the training sample point to establish the mapping relationship between the water depth value and the apparent reflectance of the training sample point.

[0092] The mapping results of each training sample point are organized, and for each training sample point, a factor set is formed consisting of single-band apparent reflectance, combined-band apparent reflectance, difference factor between two bands, and ratio factor between two bands.

[0093] For each training sample point, its factor set is arranged in a fixed order to form the water depth inversion factor vector of that training sample point. The water depth inversion factor vectors of all training sample points are summarized to form the water depth inversion factor vector set.

[0094] In this embodiment, the generation of the optimal water depth inversion model specifically includes:

[0095] Based on the set of water depth inversion factor vectors and the water depth values ​​of training sample points, the training sample point set and the test sample point set are distinguished. For each training sample point in the training sample point set, the water depth inversion factor vector is recorded as the input variable, and the water depth value is recorded as the output variable.

[0096] Arrange the water depth inversion factor vectors of all training sample points in a fixed order to form the training sample point design matrix, and arrange the water depth values ​​of the training sample points in the same order to form the training sample point target vector.

[0097] The water depth inversion model is defined as an energy-minimizing water depth inversion model based on variational constraints. The energy functional of the energy-minimizing water depth inversion model consists of three parts: the first part is the weighted sum of the absolute values ​​of the differences between the predicted water depth values ​​of the training sample points and the actual water depth values ​​of the training sample points, which serves as a data consistency term; the second part is the sum of the squares of the differences between the predicted water depth values ​​of adjacent training sample points, which serves as a spatial smoothing constraint term; and the third part is the product of the sum of the squares of the parameter vectors corresponding to the water depth inversion model and the regularization coefficient, which serves as a regularization constraint term. The three parts are added together to form the energy functional.

[0098] The variational constraint-based energy-minimizing water depth inversion model refers to constructing an energy functional model that includes data consistency constraints, spatial smoothness constraints, and regularization constraints during the water depth inversion process, using the water depth inversion factor vector of the training sample point set and the water depth values ​​of the training sample points as input and output relationships. The parameter vector of the water depth inversion model is determined by minimizing this energy functional. The data consistency constraint is used to measure the difference between the model's predicted water depth value and the water depth value of the training sample points. The spatial smoothness constraint is used to maintain the continuity of the predicted water depth values ​​of adjacent training sample points. The regularization constraint is used to prevent the model parameter vector from being too large or overfitting. By minimizing this energy functional, the model can maintain the smoothness of the spatial distribution and control the model complexity while ensuring the consistency between the predicted water depth value and the water depth value of the training sample points, thereby obtaining the optimal water depth inversion model.

[0099] The predicted water depth value of the training sample points is the output value obtained by inputting the water depth inversion factor vector of the training sample points into the current water depth inversion model and calculating the model parameters.

[0100] An iterative optimization method based on gradient descent is adopted to determine the parameter vector of the water depth inversion model by minimizing the energy functional, and the parameter estimation results under different regularization coefficients are obtained.

[0101] The parameter estimation results are verified using a set of test sample points. For each test sample point, its water depth inversion factor vector is input, the predicted water depth value is output, and it is compared with the water depth value of the test sample point. The sum of the absolute values ​​of the differences between the predicted values ​​of all test sample points and the water depth values ​​of the test sample points is calculated as the verification error measure.

[0102] Within a pre-defined set of candidate regularization coefficients, the verification error metric corresponding to each regularization coefficient is calculated one by one. The regularization coefficient with the smallest verification error metric is selected, and the corresponding parameter estimate is used as the parameter of the optimal water depth inversion model to generate the optimal water depth inversion model.

[0103] In this embodiment, the generation of the estimated total water storage in the region specifically includes:

[0104] Under a unified coordinate system, the calculation area of ​​a single pond is defined by the pond water surface boundary vector, and the calculation area is triangulated in a plane to obtain a grid set composed of the pond water surface boundary vector and the discrete water depth point dataset.

[0105] Apply boundary conditions to the boundary vector of the pond water surface, fix the value of the continuous water depth function at the boundary to zero, and make the pond water surface boundary a constraint boundary.

[0106] The continuous water depth function refers to a two-dimensional function that covers the boundary range of the pond water surface under a unified coordinate system. The input is planar coordinates, and the output is the water depth value corresponding to that position. At discrete water depth points, the value of the continuous water depth function is consistent with the predicted water depth value. On the boundary vector of the pond water surface, the value of the continuous water depth function is zero. Between the boundary and the discrete points, the continuous water depth function forms a smooth transition through the energy minimization method, which is used to describe the spatial distribution of water depth inside the pond.

[0107] Apply interior point constraints to each discrete water depth point in the discrete water depth point dataset that is located inside the computational region, so that the value of the continuous water depth function at the discrete water depth point is equal to the predicted water depth value of the discrete water depth point.

[0108] An objective function for constructing a minimum energy surface is defined, consisting of a membrane energy term and a thin-plate bending energy term. The membrane energy term measures the change of the first-order partial derivative of the continuous water depth function over the computational domain, and the thin-plate bending energy term measures the change of the second-order partial derivative of the continuous water depth function over the computational domain. The objective function is a weighted sum of these two parts.

[0109] ;

[0110] in, Let represent the objective function of the minimum energy surface. This represents the computational region defined by the boundary vector of the pond's water surface. Represents a continuous water depth function. , Indicates the continuous water depth function in and The first-order partial derivative in the direction corresponds to the membrane energy term. , , , representing the second-order partial derivative of the continuous water depth function, corresponding to the bending energy term of the thin plate. , This represents a weighting coefficient used to balance the importance of the membrane energy term and the sheet bending energy term.

[0111] Using the minimum objective function as the criterion, and simultaneously satisfying boundary conditions and interior point constraints, the minimum energy surface for a single pit is obtained;

[0112] In the process of constructing the minimum energy surface, the boundary condition refers to the requirement that the value of the continuous water depth function must be equal to zero at the position corresponding to the boundary vector of the pond water surface, so as to ensure that the water depth at the water surface boundary is zero. The interior point constraint refers to the requirement that the value of the continuous water depth function must be equal to the predicted water depth value of the discrete water depth point at the position of the discrete water depth point dataset, so as to ensure that the model is strictly consistent with the prediction result at these positions. The requirement to simultaneously satisfy the boundary condition and the interior point constraint means that when solving the minimum energy surface, the optimization process not only requires the continuous water depth function to minimize the energy of the entire region, but also requires it to satisfy the boundary condition of zero water depth at the boundary, and to satisfy the interior point constraint consistent with the predicted water depth value at all discrete water depth points, so that the final continuous water depth function not only conforms to the physical characteristics of the boundary, but also maintains consistency with the discrete point data.

[0113] The minimum energy surface is generated by selecting a specific solution from all continuous water depth functions that simultaneously satisfies the boundary conditions and interior point constraints and minimizes the objective function. This solution is strictly zero at the boundary, strictly equal to the predicted water depth value at discrete water depth points, and filled with a smooth transition at other locations. This generates a water depth distribution surface that satisfies both physical boundary characteristics and observation point constraints.

[0114] The minimum energy surface is discretized on a mesh set using the finite element method. The continuous water depth function is represented by a linear combination of nodal basis functions and assembled into a discrete system of equations. The water depth values ​​of the mesh nodes are obtained by solving the system of equations.

[0115] Under a unified coordinate system, a grid set is constructed within the computational region defined by the boundary vector of the pond water surface. The continuous water depth function is represented at the nodes of each grid cell using nodal basis functions. Specifically, nodes are selected in each grid cell, and nodal basis functions are defined corresponding to the nodes. The continuous water depth function is represented by a linear combination of nodal basis functions. Based on this, the objective function of the minimum energy surface is discretized on the grid set to form an algebraic equation for the water depth value of each node. The algebraic equations corresponding to all grid cells are assembled in sequence to obtain a discrete equation set covering the entire computational region. Boundary conditions corresponding to the boundary vector of the pond water surface and interior point constraints corresponding to the discrete water depth point dataset are applied to the discrete equation set to fix the water depth values ​​of the boundary nodes and interior point nodes. Finally, a numerical solution method is used to solve the discrete equation set to obtain the water depth value distribution of all nodes in the grid set.

[0116] The nodal basis function refers to a local function defined on the mesh element in the finite element method. It takes a value of 1 at its own node and a value of 0 at other nodes. It is used to represent the continuous water depth function as a linear combination of the water depth values ​​at each node.

[0117] Within the computational domain defined by the boundary vector of the pond water surface, the continuous water depth function is represented on the grid set as a linear combination of nodal basis functions. That is, within each grid cell, the nodal water depth value is used as the coefficient, and the continuous water depth function is represented by nodal basis functions. Based on this, the integral terms involving the continuous water depth function in the objective function of the minimum energy surface are replaced with the form of a linear combination of nodal basis functions, and the integration operation is transformed into a finite summation of each grid cell in the grid set. For each grid cell, the local stiffness matrix and local load vector are calculated separately, and the local stiffness matrix and local load vector of all grid cells are assembled into a discrete equation system covering the entire grid set. Finally, an algebraic equation system that depends only on the nodal water depth value is obtained, thereby realizing the discretization of the objective function of the minimum energy surface.

[0118] The local stiffness matrix is ​​a coefficient matrix obtained by discretizing the membrane energy term and the thin plate bending energy term in the objective function within a single grid cell, and then differentiating and integrating with respect to the nodal basis functions.

[0119] The local load vector is a vector obtained by integrating the nodal basis functions after discretizing the terms related to external constraints or boundary conditions in the objective function within a single grid cell.

[0120] After obtaining the water depth values ​​of the grid nodes, volume integration is performed on the computational domain, and discrete integration is performed using the vertex averaging method of the grid cells to obtain the water storage capacity of a single pit.

[0121] Calculate the water storage capacity of all ponds within the estimated area, and sum them up to generate the total water storage capacity of the estimated area.

[0122] Example 1:

[0123] To verify the feasibility of this invention in practice, it was applied to a contiguous network of ponds and waterways in the suburbs of a city. The area contains numerous natural and artificial water bodies of varying shapes and sizes. Local water resource management departments have long faced the challenge of lacking a precise understanding of the distribution of pond water depths and actual water storage. Previous methods relied heavily on manual water depth measurements combined with empirical methods for water storage estimation, which was time-consuming, labor-intensive, and difficult to cover a wide area. Furthermore, during the rainy season when water levels fluctuate frequently, the timeliness and accuracy of the data were even more difficult to guarantee. Simultaneously, simple area estimation methods based on remote sensing water body identification are affected by factors such as lighting conditions, changes in water reflectivity, and shoreline vegetation obstruction, resulting in significant errors in water body boundary and storage calculations, which cannot meet the current needs for high-precision water resource management and scheduling.

[0124] In this scenario, field measurements were conducted on multiple ponds within the area. Portable depth sounders and high-precision positioning instruments were used to acquire raw water depth data at several sampling points within the ponds. Simultaneously, the coordinates and time information of the measurement points were accurately recorded. The collected measured data underwent preliminary screening to remove abnormal records and duplicate measurement points. After processing using standardized methods, all depth sounding data were unified under the same benchmark and randomly divided into training and testing sets according to a predetermined ratio, providing support for subsequent model construction and validation.

[0125] Meanwhile, high-resolution multispectral remote sensing imagery of the region on a cloudless day was used to collect data across the entire area. The imagery underwent radiometric, atmospheric, geometric, and orthorectification corrections to ensure that the reflectance of each pixel accurately reflected the actual surface conditions. Further noise removal and image fusion were then performed to create a standardized multispectral apparent reflectance image spatially corresponding to the measured points. To facilitate subsequent analysis, all data were projected using a unified coordinate system, achieving seamless integration between remote sensing and measured data.

[0126] Based on the distribution characteristics of pits and ponds in remote sensing imagery, a multi-scale window sequence including blue, green, red, and near-infrared bands was constructed. The number of cover boxes at each scale was calculated using the differential box counting method, and the fractal dimension of each pixel was obtained by fitting. The fractal dimensions of all single-band and combined-band data were organized. Combined with the water and non-water body markers of the measured water depth point set, the differences in fractal characteristics between water and non-water bodies were systematically analyzed, and the set of fractal thresholds with the most discriminative power was extracted. Based on these fractal thresholds, a comprehensive water body judgment quantity was constructed, and a water body mask was generated from the remote sensing imagery. After connected component filtering and morphological processing, noise and small-area noise were effectively removed, and finally, a water surface boundary vector that conforms to the actual shoreline morphology of the pits and ponds was formed. Through boundary tracking and simplification tolerance processing, the complex boundary was transformed into ordered vector data with clear structure and refined vertices, laying a solid foundation for subsequent spatial modeling.

[0127] In the water depth inversion stage, standardized water depth data is precisely spatially mapped to multispectral imagery. Centered on each training sample point, the apparent reflectance of the corresponding pixel in both single-band and combined-band indices is extracted. Difference and ratio factors are calculated to construct a multidimensional water depth inversion factor vector set. These factors fully reflect the potential correlation between different spectral information and water depth, providing rich representations for the model input layer. Based on the training set, a water depth inversion model with energy minimization variational constraints is established. This model not only constrains the difference between predicted and measured water depths but also incorporates spatial smoothness of adjacent measurement points and parameter regularization requirements, effectively preventing overfitting and extreme value fluctuations. The model is validated multiple times using test set data. By adjusting the regularization coefficients, the model achieves an optimal balance between accuracy and generalization ability, ultimately determining the optimal water depth inversion model parameters for the region.

[0128] By applying the optimal water depth inversion model, water depth values ​​are batch-calculated for all pixels within the water body area of ​​remote sensing images, resulting in a discrete water depth point dataset covering the entire pond area. These calculated values ​​not only retain the accuracy of measured points but also achieve high-precision extrapolation based on remote sensing characteristics in non-measured areas. Using the water surface boundary vector as the boundary condition and the discrete water depth points as interior point constraints, a minimum energy surface is constructed within the pond area using the finite element method. The continuous water depth function is discretized through nodal basis functions, and a stiffness matrix and load vector are assembled for each grid cell, forming a discrete system of equations. After applying zero water depth conditions at boundary nodes and water depth calculation constraints at interior node nodes, a numerical method is used to efficiently solve for the water depth distribution of all nodes.

[0129] Finally, for each pond, within the boundary-defined computational area, volume integration is performed using the vertex averaging method based on the solved node water depth values ​​to obtain the estimated water storage capacity. The water storage capacities of all ponds in the entire region are then summed to obtain the total water storage capacity of the target area, providing a scientific basis for practical applications such as water resource allocation, irrigation management, and drought early warning.

[0130] To verify the performance of the present invention in practice, it was compared with traditional methods, and the results are shown in Table 1.

[0131] Table 1. Comparison of the effectiveness of the method of the present invention and the traditional method in estimating the water depth and storage capacity of ponds.

[0132]

[0133] As can be clearly seen from Table 1, the multi-source data fusion and energy minimization model method proposed in this invention is superior to traditional manual measurement estimation and single remote sensing image estimation in the task of estimating pond water depth and storage capacity.

[0134] Regarding water depth error, traditional manual methods suffer from an average error as high as 0.32 meters due to sparse measurement points and human error. Single remote sensing methods are limited by water spectral variations and resolution, and can only reduce the error to 0.28 meters. However, this invention, by fusing measured water depth points and remote sensing multispectral information, and combining fractal feature boundary identification, spatial smoothing, and energy minimization inversion, reduces the water depth estimation error to 0.10 meters, achieving higher estimation accuracy.

[0135] Regarding the relative error of water storage, traditional methods suffer from a calculation error of 14.7% due to the difficulty in accurately obtaining the boundary; single remote sensing methods also have an error of 11.3% due to the limited extrapolation of water depth. The method of this invention reduces the cumulative error by using minimum energy surfaces, finite element discretization, and high-density discrete point extrapolation under spatial constraints. The relative error of water storage is only 3.8%, which can provide more reliable data support for actual scheduling and ecological management.

[0136] The advantages of this invention are also evident in its operational efficiency and spatial resolution. Traditional manual methods are extremely time-consuming, requiring 22 hours to measure 10 ponds at once, and the spatial resolution is insufficient to cover the entire area. While remote sensing methods improve speed, their resolution is limited by the sensor, and there are still shortcomings in handling the boundaries of some small water bodies. In contrast, this invention utilizes automated remote sensing interpretation, mathematical modeling, and limited manual assistance to significantly shorten the operation time to 2.5 hours and improve the spatial resolution to 2 meters, achieving comprehensive and precise coverage of small ponds.

[0137] In terms of human involvement, this invention minimizes on-site requirements, requiring only one person to complete data collection and verification, reducing the workload of management departments, minimizing errors caused by human interference, and improving the level of process automation and intelligence.

[0138] The above description is only a preferred embodiment of the present invention, but the scope of protection of the present invention is not limited thereto. Any equivalent substitutions or modifications made by those skilled in the art within the scope of the technology disclosed in the present invention, based on the technical solution and inventive concept of the present invention, should be covered within the scope of protection of the present invention.

Claims

1. A method for estimating pond water depth and calculating storage capacity using multi-source data, characterized in that, Includes the following steps: Collect measured water depth data of pits and ponds in the estimation area, perform data preprocessing, generate standardized water depth data, and divide it into training sample point set and test sample point set; Acquire multispectral remote sensing images of the extrapolated area and perform image preprocessing to generate standardized multispectral apparent reflectance images; Unify standardized water depth data and standardized multispectral apparent reflectance images into the same coordinate system; Based on standardized multispectral apparent reflectance images, the fractal dimension of each band and combined bands is calculated. The fractal thresholds of water bodies and non-water bodies are determined using training sample point sets. Water body regions are segmented according to the fractal thresholds to generate pond surface boundary vectors. Establish the mapping relationship between the water depth values ​​of the training sample points and the standardized multispectral apparent reflectance, and construct the water depth inversion factor vector; Using the water depth inversion factor vector and the water depth values ​​of the training sample points as input, a water depth inversion model is established, and the water depth inversion model is validated using a set of test sample points to generate the optimal water depth inversion model. The discrete water depth points of the pit are calculated using the optimal water depth inversion model, and a discrete water depth point dataset is generated. Using the boundary vector of the pond water surface as the boundary condition and the discrete water depth point dataset as the interior point constraint, a minimum energy surface is established. Volume integration is performed on the minimum energy surface to obtain the water storage of a single pond. The water storage of all ponds is then summed to generate the total water storage of the estimated area. The generation of the optimal water depth inversion model specifically includes: Based on the set of water depth inversion factor vectors and the water depth values ​​of training sample points, the training sample point set and the test sample point set are distinguished. For each training sample point in the training sample point set, the water depth inversion factor vector is recorded as the input variable, and the water depth value is recorded as the output variable. Arrange the water depth inversion factor vectors of all training sample points in a fixed order to form the training sample point design matrix, and arrange the water depth values ​​of the training sample points in the same order to form the training sample point target vector. The water depth inversion model is defined as an energy-minimizing water depth inversion model based on variational constraints. The energy functional of the energy-minimizing water depth inversion model consists of three parts: the first part is the weighted sum of the absolute values ​​of the differences between the predicted water depth values ​​of the training sample points and the actual water depth values ​​of the training sample points, which serves as a data consistency term; the second part is the sum of the squares of the differences between the predicted water depth values ​​of adjacent training sample points, which serves as a spatial smoothing constraint term; and the third part is the product of the sum of the squares of the parameter vectors corresponding to the water depth inversion model and the regularization coefficient, which serves as a regularization constraint term. The three parts are added together to form the energy functional. An iterative optimization method based on gradient descent is adopted to determine the parameter vector of the water depth inversion model by minimizing the energy functional, and the parameter estimation results under different regularization coefficients are obtained. The parameter estimation results are verified using a set of test sample points. For each test sample point, its water depth inversion factor vector is input, the predicted water depth value is output, and it is compared with the water depth value of the test sample point. The sum of the absolute values ​​of the differences between the predicted values ​​of all test sample points and the water depth values ​​of the test sample points is calculated as the verification error measure. Within a pre-defined set of regularization coefficient candidates, the verification error metric corresponding to each regularization coefficient is calculated one by one. The regularization coefficient with the smallest verification error metric is selected, and the corresponding parameter estimate is used as the parameter of the optimal water depth inversion model to generate the optimal water depth inversion model. The generation of the estimated total water storage in the region specifically includes: Under a unified coordinate system, the calculation area of ​​a single pond is defined by the pond water surface boundary vector, and the calculation area is triangulated in a plane to obtain a grid set composed of the pond water surface boundary vector and the discrete water depth point dataset. Apply boundary conditions to the boundary vector of the pond water surface, fix the value of the continuous water depth function at the boundary to zero, and make the pond water surface boundary a constraint boundary. Apply interior point constraints to each discrete water depth point in the discrete water depth point dataset that is located inside the computational region, so that the value of the continuous water depth function at the discrete water depth point is equal to the predicted water depth value of the discrete water depth point. The objective function for constructing the minimum energy surface is composed of a membrane energy term and a thin plate bending energy term. The membrane energy term measures the change of the first-order partial derivative of the continuous water depth function in the computational domain, and the thin plate bending energy term measures the change of the second-order partial derivative of the continuous water depth function in the computational domain. The objective function is a weighted sum of the two parts. Using the minimum objective function as the criterion, and simultaneously satisfying boundary conditions and interior point constraints, the minimum energy surface for a single pit is obtained; The minimum energy surface is discretized on a mesh set using the finite element method. The continuous water depth function is represented by a linear combination of nodal basis functions and assembled into a discrete system of equations. The water depth values ​​of the mesh nodes are obtained by solving the system of equations. After obtaining the water depth values ​​of the grid nodes, volume integration is performed on the computational domain, and discrete integration is performed using the vertex averaging method of the grid cells to obtain the water storage capacity of a single pit. Calculate the water storage capacity of all ponds within the estimated area, and sum them up to generate the total water storage capacity of the estimated area.

2. The method for estimating pond water depth and calculating storage capacity using combined multi-source data as described in claim 1, characterized in that, The specific division of the training sample point set and the test sample point set includes: Collect raw measured water depth data in the pits and ponds in the estimation area, record the plane coordinates, collection time and water depth value of each depth measurement point, and form a set of raw measured water depth data. Outlier removal was performed on the original measured water depth data set; Perform duplicate value removal on the water depth data set after outlier removal; The deduplicated water depth data set is standardized by using the average and standard deviation of the water depth values ​​in the set as a benchmark to standardize each record and generate a standardized water depth data set, which includes plane coordinates, acquisition time and standardized water depth values. The standardized water depth data set is divided into training sample sets and test sample set according to the training ratio.

3. The method for estimating pond water depth and calculating storage capacity using combined multi-source data as described in claim 1, characterized in that, The image preprocessing includes radiometric correction, atmospheric correction, geometric correction, orthorectification, image fusion, and noise removal.

4. The method for estimating pond water depth and calculating storage capacity using combined multi-source data as described in claim 1, characterized in that, The generation of the boundary vector of the pond water surface specifically includes: A band set and a combined band set are established on a standardized multispectral apparent reflectance image. The apparent reflectance of each pixel in a single band is recorded directly. For each combined band, the weighted average of the apparent reflectance of each single band is used as the apparent reflectance of the combined band. The band set includes blue light band, green light band, red light band and near-infrared band. A multi-scale window sequence is constructed with the pixel as the center. At each scale, the number of covering boxes is calculated using the differential box counting method, and the relationship between the number of covering boxes and the scale is recorded. By fitting the slope of the straight line between the logarithm of the number of covering boxes and the logarithm of the scale, the fractal dimension of the pixel is generated. The fractal dimensions of all single-band and combined-band data are summarized to construct the fractal feature vector of the pixel. The fractal feature vector is then matched with the training sample point set at the same coordinate position. Based on the water body label and non-water body label of the training sample point set, the fractal feature vector is divided into the water body set and the non-water body set. For each pixel's fractal feature component, calculate the mean and variance of the water body set and the non-water body set respectively. Take the ratio of the mean difference to the variance sum as the discrimination criterion, and select the largest value of the discrimination criterion as the fractal threshold of the fractal feature component to obtain the fractal threshold set. Based on the fractal threshold set, a pixel comprehensive judgment value is constructed using the fusion coefficients corresponding to each fractal feature component. A water body region mask is generated based on the comparison result between the pixel comprehensive judgment value and the judgment threshold. If the pixel comprehensive judgment value is greater than or equal to the judgment threshold, it is judged as a water body pixel. If the pixel comprehensive judgment value is less than the judgment threshold, it is judged as a non-water body pixel. In the water body region mask, water body pixels are marked as 1 and non-water body pixels are marked as 0. Perform connected component screening and morphological processing on the water body mask, delete connected components with an area smaller than the minimum area threshold and perform closing operation to obtain water body regions with satisfactory connectivity and smooth boundaries. Boundary tracing is performed on the outer boundary of the water body area to obtain a sequence of boundary polygons. The boundary polylines are then simplified with a simplification tolerance to generate vertex coordinates. An ordered vector boundary is constructed according to the vertex order to generate the boundary vector of the pond water surface.

5. The method for estimating pond water depth and calculating storage capacity using combined multi-source data as described in claim 1, characterized in that, The construction of the water depth inversion factor vector specifically includes: Under a unified coordinate system, the standardized water depth data and the standardized multispectral apparent reflectance image are spatially correlated to determine the corresponding pixel position of each training sample point on the standardized multispectral apparent reflectance image. At the corresponding pixel location of the training sample point, extract the apparent reflectance of each single band, the apparent reflectance of the combined band, the difference factor between the two bands, and the ratio factor between the two bands, and correlate them with the water depth value of the training sample point to establish the mapping relationship between the water depth value and the apparent reflectance of the training sample point. The mapping results of each training sample point are organized, and for each training sample point, a factor set is formed consisting of single-band apparent reflectance, combined-band apparent reflectance, difference factor between two bands, and ratio factor between two bands. For each training sample point, its factor set is arranged in a fixed order to form the water depth inversion factor vector of that training sample point. The water depth inversion factor vectors of all training sample points are summarized to form the water depth inversion factor vector set.

Images