Water body and water level cooperative extraction method and system

Abstract

The present application relates to the field of hydrological remote sensing and geographic information technology, and more particularly to a water body and water level cooperative extraction method and system. The method comprises the following steps: collecting multi-source remote sensing data of a target area within a preset time period; extracting water body information of optical image data in the multi-source remote sensing data to generate a water surface mask, and quantitatively extracting water body morphological characteristic parameters based on the water surface mask; extracting observation data within the spatio-temporal range of the water surface mask in the spaceborne radar altimeter planar water level data that meets the water body morphological characteristic parameters to construct an initial planar water level field; performing spatio-temporal matching and spatial correction on the spaceborne laser radar point elevation data in the multi-source remote sensing data and the initial planar water level field; and performing spatio-temporal fusion according to the spatially corrected planar water level field and the spatio-temporal change information of the water surface mask. The present application can break through the limitations of a single data source in terms of spatio-temporal resolution and continuous monitoring capability, while improving the water level inversion accuracy and reliability.

Description

Technical Field

[0001] This invention relates to the field of hydrological remote sensing and geographic information technology, and in particular to a method and system for the collaborative extraction of water bodies and water levels. Background Technology

[0002] Accurate and continuous water level information is a crucial foundation for hydrological research, water resource management, flood and drought early warning, and water conservancy project scheduling. While traditional hydrological stations offer high monitoring accuracy, their spatial representativeness is limited, and their distribution globally (especially in remote areas) is sparse. Satellite remote sensing technology provides a transformative tool for conducting large-scale, dynamic water level monitoring.

[0003] Currently, mainstream satellite water level extraction technologies mainly rely on a single data source. However, no single satellite data source can simultaneously meet the requirements for high-precision, high spatiotemporal continuity, and high spatial coverage water level monitoring. In existing technologies, although some studies have attempted to combine different data, most are limited to simple data comparison or empirical interpolation, lacking a systematic, adaptive collaborative extraction framework that can fully leverage the advantages of each data source and overcome its disadvantages. Summary of the Invention

[0004] Therefore, it is necessary for the present invention to provide a method and system for the coordinated extraction of water body and water level to solve at least one of the above-mentioned technical problems.

[0005] To achieve the above objectives, a method for the coordinated extraction of water body and water level includes the following steps: Step S1: Collect multi-source remote sensing data of the target area within a preset time period, including optical image data, surface water level data from spaceborne radar altimeter, and point elevation data from spaceborne lidar. Step S2: Extract water body information from optical image data in multi-source remote sensing data to generate a water surface mask, and extract water body morphology feature parameters based on the water surface mask; Step S3: Extract observation data that satisfy the water body morphology characteristic parameters within the corresponding spatiotemporal range of the water surface mask from the spaceborne radar altimeter water level data, in order to construct the initial water level field; Step S4: Perform spatiotemporal matching between the point elevation data of the spaceborne lidar in the multi-source remote sensing data and the initial areal water level field, and perform spatial correction on the initial areal water level field; Step S5: Perform spatiotemporal fusion based on the spatially corrected planar water level field and the spatiotemporal variation information of the water surface mask to generate a high spatiotemporal resolution water surface elevation model.

[0006] This application achieves the synergistic utilization of multi-source information by fusing three types of remote sensing data: optical imagery, spaceborne radar altimeters, and spaceborne lidar. This overcomes the performance limitations of single data sources in terms of spatial coverage, observation accuracy, and temporal continuity. By leveraging the high-precision characteristics of lidar point elevation data, spatial local correction is performed on the areal water level data acquired by the radar altimeter. This effectively eliminates systematic biases generated by the radar altimeter in near-shore areas and complex terrain conditions, significantly improving the accuracy and reliability of water level retrieval results. Simultaneously, by combining the spatiotemporal variation information of water bodies acquired from optical imagery, spatiotemporal fusion and expansion of low-frequency water level observation results enable continuous estimation of discrete water level observation data during periods without observation, thereby achieving high spatiotemporal resolution dynamic water level monitoring. Furthermore, by constructing a hierarchical processing flow and introducing a spatially adaptable regression correction mechanism, the method can adapt to water environments of rivers, lakes, and reservoirs with different widths and shapes. It maintains stable processing performance and high applicability even in complex water scenarios, significantly improving the automation and robustness of the overall water level extraction method.

[0007] Optionally, this application also provides a water body and water level co-extraction system for performing the water body and water level co-extraction method described above, the water body and water level co-extraction system comprising: The data acquisition module is used to collect multi-source remote sensing data of the target area within a preset time period. The multi-source remote sensing data includes optical image data, surface water level data from the spaceborne radar altimeter, and point elevation data from the spaceborne lidar. The water body information extraction module is used to extract water body information from optical image data in multi-source remote sensing data, generate a water surface mask, and extract water body morphology feature parameters based on the water surface mask. The water level field construction module is used to extract observation data that satisfy the water body morphology characteristic parameters within the corresponding spatiotemporal range of the water surface mask in the planar water level data of the spaceborne radar altimeter, so as to construct the initial planar water level field; The water level field correction module is used to perform spatiotemporal matching between the point elevation data of the spaceborne lidar in the multi-source remote sensing data and the initial surface water level field, and to perform spatial correction on the initial surface water level field. The spatiotemporal fusion module is used to perform spatiotemporal fusion based on the spatially corrected isometric water level field and the spatiotemporal variation information of the water surface mask, thereby generating a high spatiotemporal resolution water surface elevation model.

[0008] The water body and water level collaborative extraction system of this application can realize any of the water body and water level collaborative extraction methods of this invention. It is used as a medium for combining the operation and signal transmission between various modules to complete the water body and water level collaborative extraction method. The internal modules of the system cooperate with each other, thereby improving the accuracy and reliability of water level inversion, while breaking through the limitations of single data sources in terms of spatiotemporal resolution and continuous monitoring capabilities, so that the obtained water level products have higher spatiotemporal continuity and adaptability. Attached Figure Description

[0009] Other features, objects, and advantages of the invention will become more apparent from the following detailed description of non-limiting embodiments with reference to the accompanying drawings: Figure 1 This is a schematic diagram of the steps of the water body and water level co-extraction method of the present invention; Figure 2 This is a comparison chart showing the before and after correction of the initial areal water level field using a geographically weighted regression model in an embodiment of the present invention. Figure 3 This is a block diagram of the water body and water level collaborative extraction system in this invention; The realization of the objective, functional features and advantages of the present invention will be further explained in conjunction with the embodiments and with reference to the accompanying drawings. Detailed Implementation

[0010] The technical method of the present invention will now be clearly and completely described with reference to the accompanying drawings. Obviously, the described embodiments are only some, not all, of the embodiments of the present invention. All other embodiments obtained by those skilled in the art based on the embodiments of the present invention without inventive effort are within the scope of protection of the present invention.

[0011] Furthermore, the accompanying drawings are merely illustrative of the invention and are not necessarily drawn to scale. The same reference numerals in the drawings denote the same or similar parts, and therefore repeated descriptions of them will be omitted. Some block diagrams shown in the drawings are functional entities and do not necessarily correspond to physically or logically independent entities. These functional entities can be implemented in software, in one or more hardware modules or integrated circuits, or in different network and / or processor methods and / or microcontroller methods.

[0012] To achieve the above objectives, please refer to Figures 1 to 3 This invention provides a method for the coordinated extraction of water body and water level, the method comprising the following steps: Step S1: Collect multi-source remote sensing data of the target area within a preset time period, including optical image data, surface water level data from spaceborne radar altimeter, and point elevation data from spaceborne lidar. In one embodiment, optical remote sensing images, spaceborne radar altimeter water level products, and spaceborne lidar elevation data covering the target area for 1 to 3 months are retrieved from a database or public database. The optical images can be multispectral images with a spatial resolution of 10m, and radiometric calibration and atmospheric correction are performed on the images to uniformly convert pixel reflectance to normalized reflectance values ​​between 0 and 1. The spaceborne radar altimeter areal water level data uses preprocessed two-dimensional raster products, with each grid containing the water level value and quality label field at the corresponding time. The spaceborne lidar point elevation data consists of discrete observation points along the orbit, with each observation point containing latitude and longitude coordinates, elevation value, and data quality label, and is uniformly converted to the same geodetic datum (e.g., EGM2008). Simultaneously, all data undergoes time normalization processing, converting the observation time to UTC time format.

[0013] Step S2: Extract water body information from optical image data in multi-source remote sensing data to generate a water surface mask, and extract water body morphology feature parameters based on the water surface mask; In a further embodiment, the Normalized Difference Water Index (NDWI) is calculated on the preprocessed optical images. This is achieved by constructing a water index using the difference in reflectance between the green and near-infrared bands, and setting a threshold of 0.2 to classify image pixels. Pixels with an NDWI value greater than 0.2 are identified as water pixels, thus forming a binary water map. Subsequently, time-series statistics are performed on the binary water maps of all phases within a preset time period. Probability statistics are performed using a sliding window of 5 images on the time axis, and median filtering is used to eliminate short-term misclassifications. Pixels are used to obtain a water surface probability map; when the probability of a certain pixel in the time series is greater than 0.6, it is determined to be a stable water body, thus generating a water surface mask; based on the water surface mask, the Zhang-Suen thinning method is used to extract the center line of the water body with the width of a single pixel, and a sampling point is set every 30m along the center line. At each sampling point, a normal line perpendicular to the center line is established to calculate the distance between the water body boundaries on both sides as the water body width. At the same time, the shoreline tortuosity is obtained by calculating the ratio of the shoreline length to the equivalent circumference, so as to form a set of water body morphology characteristic parameters.

[0014] Step S3: Extract observation data that satisfy the water body morphology characteristic parameters within the corresponding spatiotemporal range of the water surface mask from the spaceborne radar altimeter water level data, in order to construct the initial water level field; In a further embodiment, the surface water level data of the water surface mask and the satellite radar altimeter are spatially overlaid. It is determined whether the center coordinates of each water level observation pixel fall within the water body pixel corresponding to the surface water mask. Only water level observation pixels that meet the spatial location conditions are retained. Then, in the time dimension, the number of times each observation pixel is marked as a water body by the surface water mask within a preset time tolerance of ±1 day is counted, and the frequency of water body occurrence is calculated. When the frequency is greater than or equal to 80%, the observation pixel is determined as a valid SWOT observation point. All SWOT observation points are processed by regular gridding according to their spatial location to establish a regular grid with a resolution of 100m×100m. The water level observation values ​​falling within the same grid cell are statistically averaged. If the number of observation points in a certain grid cell is not less than 3, the average water level value is taken as the representative water level value of the grid cell. The grid cells are arranged according to their spatial location to form continuous grid data, thereby constructing an initial surface water level field covering the target area.

[0015] Most importantly, the method for determining time tolerance is as follows: ,in This refers to the optical image revisit cycle.

[0016] Step S4: Perform spatiotemporal matching between the point elevation data of the spaceborne lidar in the multi-source remote sensing data and the initial areal water level field, and perform spatial correction on the initial areal water level field; In a further embodiment, the quality identifier field in the point elevation data of the spaceborne lidar is used to filter the lidar points, retaining only valid observation points with a quality level of 0 or 1. Then, the filtered lidar points are spatially matched with the corresponding water level values ​​in the initial areal water level field, and the elevation difference between the two is calculated as a residual sequence. When the absolute value of the residual exceeds the median plus three times the absolute deviation of the median, it is identified as an outlier and removed. For the remaining sample points, a geographic weighted regression relationship is established based on their spatial coordinates, water width, and shoreline tortuosity. The regression coefficients are solved using the weighted least squares method, and a spatial weight matrix is ​​constructed using a double-squared kernel function. Local regression coefficient surfaces are calculated within a neighborhood with a bandwidth parameter of 5 km. The initial areal water level field is then corrected pixel-by-pixel based on the spatial variation trend of the regression coefficient surfaces. The corrected pixel water level values ​​are then subjected to spatial interpolation smoothing with a neighborhood radius of 3 grid units to generate a spatially corrected areal water level field.

[0017] Step S5: Perform spatiotemporal fusion based on the spatially corrected planar water level field and the spatiotemporal variation information of the water surface mask to generate a high spatiotemporal resolution water surface elevation model.

[0018] In a further embodiment, a Gaussian process regression relationship is established using the spatially corrected planar water level field and the temporal variation information of the water surface mask. The Gaussian process regression model is constructed using a zero-mean function and a squared exponential kernel function, with the kernel function taking the form of… ,in Indicates the spatiotemporal distance between samples. For length scale parameters, The signal variance parameter is used; the spatial coordinates (longitude and latitude) of the water level observation and the observation time are used to form a three-dimensional input feature vector, and the spatially corrected water level value is used as the training sample output. The kernel function parameter is determined by maximum likelihood estimation so that the model can describe the continuous change of water level in time and space. After the model training is completed, any target time and target position are input into the model for prediction, and the prediction results are constrained and corrected by combining the expansion or contraction information of the water surface mask in the time series, thereby generating a high spatiotemporal resolution water surface elevation model.

[0019] In another embodiment, the Gaussian process regression model specifically includes: Define water level For a Gaussian process: Mean function Composed of trend terms: ,in This is a time series of water body area extracted from a water surface mask. The parameter to be estimated. Covariance function. Employ the space-time separable Matern kernel function: in, For smoothing parameters, the preferred option is... Corresponding The function is in the form of ; For spatial distance, For time intervals; and These are parameters for spatial and temporal scales, respectively. The standard deviation of the signal; The standard deviation of noise; This refers to the Kroneckerdelta function. Hyperparameters Estimate by maximizing the log-marginal likelihood function: ;in, To observe the water level vector, For the corresponding spatiotemporal coordinate matrix, , Let represent the identity matrix. The L-BFGS-B optimization algorithm is used to solve the problem, with reasonable constraints on the hyperparameters. The trained Gaussian process model can be used to predict arbitrary spatiotemporal locations. The water level value and its uncertainty.

[0020] Optionally, generating the water surface mask in step S2 includes: The optical image data is preprocessed, the water index of the preprocessed optical image data is calculated, and water pixels are extracted to generate a binary water map. In one embodiment, radiometric calibration and atmospheric correction are performed on optical image data from each time phase to eliminate sensor response differences and atmospheric scattering effects, and the data is uniformly converted into surface reflectance data. Subsequently, resampling and geometric registration are performed at a uniform spatial resolution (preferably 10 m) to ensure strict alignment of pixels in each band within the same coordinate grid. Based on this, the Normalized Difference Water Index (NDWI) is calculated using the near-infrared and short-wave infrared bands in the image, and its expression is: A water body index raster matrix is ​​generated across the entire image area. A discrimination threshold (preferably between 0 and 0.1, 0.05 in this embodiment) is set for the water body index raster matrix to distinguish pixels. When the pixel index value is greater than the discrimination threshold, it is marked as a water body pixel, and when it is less than the discrimination threshold, it is marked as a non-water body pixel, thereby forming a binary water body map containing only 0 and 1, where "1" represents water body and "0" represents non-water body.

[0021] Time series analysis is performed on the binary water body maps of all time phases within a preset time period. Through time series filtering and synthesis, a spatiotemporally continuous water surface probability map is generated. Pixels in the water surface probability map with a probability greater than a preset stable water body determination threshold are identified as water bodies, forming a water surface mask.

[0022] In a further embodiment, the binary images of each time phase are stacked in a time series according to a unified grid coordinate system to construct a three-dimensional binary raster data volume. ,in Represents spatial coordinates, The observation time is indicated; subsequently, the frequency of water body occurrence is calculated for each spatial pixel along the time dimension. ,in This represents the number of times the pixel is labeled as a water body in the time series. This represents the total number of observations for that pixel in the statistics. To reduce random noise caused by cloud cover or misjudgment, a median filter with a sliding time window of length 3 is applied to the frequency sequence, and probabilistic smoothing is performed in conjunction with the spatial neighborhood (3×3 pixel window) to generate a water surface probability map. When the value exceeds the preset stable water body determination threshold (preferably between 0.6 and 0.8, and 0.7 in this embodiment), the corresponding pixel is determined to be a stable water body region, thereby forming a water surface mask for the target region.

[0023] Optionally, the water body morphology parameters in step S2 include water body width and shoreline tortuosity; wherein the extracted water body morphology parameters include: The Zhang-Suen thinning algorithm is used, and the center line of the water body with a single pixel width is extracted based on the water surface mask; In this embodiment, connected component identification processing is performed on the water surface mask to identify each independent water body unit, and a water body pixel matrix is ​​constructed in a grid with a uniform spatial resolution (e.g., 10 m). Subsequently, Zhang-Suen thinning processing is performed on this matrix. This method is a skeleton extraction method based on iterative neighborhood rules. It performs a two-stage iterative deletion judgment on the topological relationships of each pixel's eight neighbors, gradually stripping boundary pixels while maintaining connectivity, thereby shrinking the water body region into a skeleton structure with a single pixel width. The maximum number of iterations is set to 50, and pixel topological stability is used as the termination condition. After thinning, path tracing is performed on the obtained skeleton pixel set, and the pixels are connected in order of their center coordinates to form a water body centerline vector sequence.

[0024] Sampling points are set at fixed intervals along the centerline of the water body; at each sampling point, the normal line of each sampling point is drawn, and the distance between the intersection of the normal line of each sampling point and the boundary of the water body is calculated as the width of the water body at that sampling point; the widths of all sampling points are statistically analyzed to obtain the average width of the water body. In a further embodiment, sampling points are arranged at fixed intervals along the centerline of the water body, for example, one sampling point is set every 30 m, and the spatial coordinates of the sampling points are recorded. For each sampling point, its tangent direction is calculated based on the local direction of the centerline, and a normal vector perpendicular to the tangent direction is constructed. Then, the water body boundary is searched pixel by pixel along the normal direction. When the pixel value changes from water body (1) to non-water body (0), the intersection position is recorded, and the Euclidean distance between the two intersection points is calculated as the local water body width at that sampling point. To reduce the influence of local anomalies, the width values ​​of 5 consecutive sampling points are processed by moving average, and the width of all sampling points is statistically calculated to obtain the average width parameter, that is, the average width of the water body.

[0025] The shoreline length of each independent water body unit is calculated based on the water surface mask, and the circumference of the equivalent circle with the same area as each independent water body unit is calculated. The shoreline tortuosity is obtained by dividing the shoreline length of each independent water body unit by the circumference of the equivalent circle.

[0026] In a further embodiment, boundary extraction processing is performed on the water surface mask. The shoreline pixel set is identified by detecting the adjacency relationships between water pixels and non-water pixels, and a shoreline polygonal structure is constructed based on the pixel center coordinates. The total length of the shoreline polygonal structure is then calculated. During calculation, a pixel scale correction factor is used to eliminate raster orientation errors. For example, the horizontal and vertical adjacency lengths are taken as 1 pixel length, and the diagonal adjacency lengths are taken as... The pixel length is then multiplied. Subsequently, the area A of the corresponding water body unit is calculated based on the water surface mask, and an equivalent circle with the same area is constructed, with a radius... according to Calculate the corresponding perimeter as follows: The actual shoreline length The shoreline tortuosity index is obtained by calculating the ratio of the shoreline to the equivalent circumference. , among which when A value close to 1 indicates a relatively smooth shoreline. A value significantly greater than 1 indicates a more complex shoreline morphology.

[0027] Optionally, constructing the initial areal water level field in step S3 includes: The spatial position of each water level observation pixel in the surface water level data of the water surface mask and the spaceborne radar altimeter is superimposed and judged, and the water level observation pixels whose center coordinates fall within the water body pixels of the water surface mask are retained. In one embodiment, the water surface mask resolution is unified to a 30m grid using bilinear resampling. Subsequently, the coordinates of water level observation pixels and their corresponding water level values ​​are read from the radar altimeter product to construct a set of water level observation points. For the center coordinates of each observed pixel Spatial overlay judgment is performed, which involves querying the pixel value of the coordinate in the water surface mask grid. If the corresponding pixel value is a water body identifier (value 1), it is retained; otherwise, it is discarded. To ensure spatial matching accuracy, a maximum positioning error threshold of 15m can be set. The final result is a set of effective water level observation pixels falling within the water body area.

[0028] The water level observation pixels that are marked as water bodies at a frequency greater than or equal to 80% within the preset time tolerance are used as SWOT observation points. In a further embodiment, a time-series record is established for each retained water level observation cell, the structure of which includes the observation time. Spatial coordinates And the status of water body markings. With a preset time tolerance. As a time window, among which Statistical analysis was performed on the water body determination results at the same spatial location within a time window. The statistical analysis of the pixel's position within the time window was conducted. The number of times the area is identified as a water body by the water surface mask, and the frequency of water body occurrence is calculated. ,in Number of times for water body determination, The number of observations within the window. When the frequency... When the value is ≥0.8, the observed pixel is identified as a stable water body observation point and recorded in the SWOT observation point set. .

[0029] The SWOT observation points are rasterized according to their spatial location, and a regular grid is established with a preset grid resolution. The water level observations in each grid cell of the regular grid are statistically calculated, and the average value is taken as the representative water level value of the grid cell. An initial planar water level field is generated according to the spatial location of the grid cells.

[0030] In a further embodiment, a regular grid structure is established based on the spatial distribution of SWOT observation points. First, a two-dimensional regular grid is generated according to the target area, for example, setting the grid resolution to 500m × 500m, with each grid cell representing its spatial location using its center coordinates. Then, the SWOT observation points are assigned to their corresponding grid cells according to their spatial locations, and a grid observation matrix is ​​constructed, where each cell records all water level observation values ​​within that grid. Statistical calculations are performed on each grid cell; for example, if the number of observation points within a cell is not less than 3, its average water level is calculated as the representative water level; if there are insufficient observation points, it is marked as a null value. Finally, a two-dimensional water level matrix is ​​constructed according to the center coordinates of the grid cells and the corresponding representative water level values, thus forming an initial planar water level field.

[0031] Optionally, step S4 includes: Step S41: Filter lidar points based on the quality identifier field in the point elevation data of the spaceborne lidar; In one embodiment, the quality identification field in the point elevation data of the spaceborne lidar is parsed. The quality identification field includes parameters such as echo confidence marker, cloud interference marker, and signal saturation marker, where the echo confidence value ranges from 0 to 1. Lidar points with an echo confidence value of not less than 0.7 and not marked as cloud interference or signal anomalies are determined as valid observation points, and a point sequence index table is established according to the orbit scanning order. Simultaneously, based on the average trajectory spacing of approximately 20-30m between points, the spatial uniformity of the point data is checked. If the point density in a local area is more than twice that of the surrounding area, sparsification is performed using a neighborhood distance threshold (e.g., 15m) to form a lidar elevation point set.

[0032] Step S42: Match the selected lidar points with the values ​​of the initial planar water level field at the corresponding locations, and calculate the residual sequence based on the matching results; remove abnormal lidar points and their corresponding data points based on the robustness of the residual sequence. In a further embodiment, the selected lidar points are spatially matched with the initial areal water level field according to their geographic coordinates. Specifically, the initial areal water level field can be uniformly resampled to a 30m×30m regular water level grid, and the grid position of each lidar point is determined according to its latitude and longitude coordinates, thereby extracting the water level value in the corresponding grid. Then, the difference between the elevation value of each lidar point and the water level value of the corresponding water level grid is calculated to form a residual sequence (unit: m). To ensure the stability of the residual sequence, a robustness test is performed: the median and absolute deviation (MAD) of the residual sequence are calculated, and an outlier threshold of 3×MAD is set. When a residual value exceeds this threshold, the lidar point is determined to be an outlier, and the corresponding grid data is simultaneously removed, ultimately forming a stable residual sequence after outlier removal.

[0033] Step S43: Construct a geographic weighted regression model based on the residual sequence remaining after removing abnormal lidar points and water morphology parameters; In a further embodiment, the residual sequence is jointly modeled with water body morphology characteristic parameters. These parameters are derived from the spatial analysis results of water body extent data, including indicators such as water surface width, shoreline curvature, and local water area change rate. Subsequently, a geographically weighted regression model is constructed using the residual sequence as the dependent variable and water surface width, shoreline curvature, and spatial coordinates as independent variables, with an adaptive kernel function used to determine the spatial weights.

[0034] Step S44: Perform model fitting on the geographic weighted regression model to obtain the regression coefficient surface at each spatial location; In a further embodiment, the center of each water level grid is used as the local regression calculation location, and weighted regression calculations are performed on neighboring sample points based on kernel function weights. During the regression process, the spatial weight function is set to a Gaussian kernel function, whose weights decay exponentially with increasing spatial distance. After calculating each local regression coefficient, the regression coefficients are rasterized according to their corresponding spatial coordinates, thereby generating multiple regression coefficient distribution surfaces, such as width coefficient distribution surfaces and curvature coefficient distribution surfaces. Subsequently, smoothing processing is performed on each coefficient distribution surface, and a unified output is a regression coefficient surface with the same resolution as the initial planar water level field.

[0035] Step S45: Select the correction algorithm corresponding to the initial isometric water level field based on the data density, spatial heterogeneity index and width variation index of the initial isometric water level field, and perform spatial correction.

[0036] In a further embodiment, the number of water level grids per unit area is counted to obtain a data density index; then, the variance of the water level difference between adjacent grids is calculated to form a spatial heterogeneity index; simultaneously, the width variation index is calculated based on the change in the water body boundary spacing. When the data density is less than 100 grids per square kilometer and the spatial heterogeneity index is high, a local interpolation correction method based on regression coefficient weighting is adopted; when the data density is high and the width variation index is less than 0.2, a regional smoothing correction method is adopted. Subsequently, the water level values ​​of each grid are corrected according to the selected method, thereby obtaining a spatially continuous corrected water level field consistent with actual observations.

[0037] Optionally, constructing a geographically weighted regression model includes: Sample points are selected based on the Euclidean distance between lidar points; In one embodiment, each lidar point is used as a center point, and its Euclidean distance to surrounding points is calculated and sorted in ascending order of distance. A threshold of 80-120 neighborhood samples is set, or a maximum spatial distance threshold of 5 km is set. When the number of samples within this distance range for a given center point meets the threshold condition, these neighborhood points are identified as the local regression sample set for that center point, and their corresponding horizontal and vertical coordinates are recorded, forming a local sample data table containing spatial location, residual values, and corresponding water morphology parameters. To avoid overly concentrated sample distribution, a distance uniformity test is performed after neighborhood sample selection. When the distance between any two sample points is less than 50 m, the sample point with higher residual stability is retained to ensure the spatial representativeness of the local samples.

[0038] Based on the corresponding values ​​of each sample point in the residual sequence and the corresponding spatial coordinates of the water body morphology parameters, a geographically weighted regression model expression is established to construct the geographically weighted regression model; the specific expression of the geographically weighted regression model is as follows: ; in For the first The x-axis coordinates of each sample point For the first The ordinate of each sample point on the y-axis For correction amount, For the first Initial SWOT water level at each sample point; For the first The width of the water body at each sample point For the first shoreline tortuosity at each sample point For random error parameters, For the first The local regression coefficients for each sample point, where Including local intercept coefficient Local initial water level regression coefficient Water body width regression coefficient and the regression coefficient of shoreline tortuosity .

[0039] In a further embodiment, a local regression expression relationship is established based on the values ​​of each sample point in the residual sequence and the corresponding spatial coordinates of the water morphology characteristic parameters. Specifically, the correction value corresponding to each sample point is extracted from the sample data table. And obtain the initial SWOT water level value at the same location, and the water width calculated through the water body boundary line. And the shoreline tortuosity calculated by the continuous angular variation of the shoreline polygonal line. Then the coordinates of the sample points were... As spatial location identifiers, local regression relationships are constructed for the aforementioned variables to establish a linear correlation between the correction values ​​and the initial water level and water morphology parameters. When constructing the expressions, a random error term is set as a disturbance term with a mean of 0 and a standard deviation of approximately 0.05m to reflect the influence of observation errors and unmodeled factors. The data for each variable are stored in a sample matrix format, where each row corresponds to a sample point and each column corresponds to a variable parameter, thus forming the basic data structure for the geographic weighted regression model.

[0040] It is worth noting that, It can be used as a correction factor, or it can be obtained from the difference between the elevation of the lidar point in the residual sequence and the initial water level of SWOT.

[0041] Optionally, methods for obtaining local regression coefficients include: With the first Using the spatial coordinates of the nth sample point as the local regression center, the calculation of the nth... Euclidean distance between each sample point and other sample points; The kernel function with the first square is used to obtain the kernel with the first square as the kernel. A spatial weight set centered at each sample point is used to construct a spatial weight matrix, where the double-square kernel function is specifically: ; In the formula, For the first The sample point and the first Spatial weights of each sample point For the first The sample point and the first Euclidean distance between sample points For bandwidth parameters; The local regression coefficient vector is solved using the weighted least squares method and based on the weight matrix. The formula for solving the local regression coefficient vector is: ; The expression for solving the local regression coefficient vector is: In the formula For the first The local regression coefficient vector of each sample point To observe the independent variable matrix, This is a sequence of water level differences. This is the spatial weight matrix.

[0042] It is worth noting that, It is the product of weighted design matrices, representing the weighted covariance of the data; It is the inverse matrix of the product, ensuring that the regression coefficients can be estimated using the least squares method; It is the product of the weighted observations and the design matrix, representing the error correction term of the regression model.

[0043] In this embodiment, the spatial weight matrix It is a diagonal matrix, and its diagonal elements are the spatial weights of the corresponding sample points. This matrix is ​​used to reflect the position of each sample point pair. The degree of influence of the regression estimate. Due to the weight matrix In each spatial location Since they are all different, the regression coefficients obtained by the above weighted least squares calculation also vary with spatial location, thus forming a spatially adaptive local regression coefficient distribution, which is used to describe the local relationship between SWOT water level error and water body morphology characteristics in different regions, and to provide parameter basis for subsequent spatial correction of water level field.

[0044] The water level difference sequence can take the form of: ; The form of the observation independent variable matrix can be: ; The first column is a constant term used to estimate the intercept coefficient, and the other three columns correspond to the initial SWOT water level, water width, and shoreline tortuosity, respectively.

[0045] It is worth noting that the correction amount for each sample point The observed dependent variable vector is formed by arranging the samples in order. Simultaneously construct the observation independent variable matrix Its structure is an n×4 matrix, where the first column is all 1s to represent constant terms, and the other three columns are filled with the SWOT initial water levels of the corresponding sample points. The width of the water body calculated using the normal to the centerline of the water body. And the shoreline tortuosity calculated based on the shoreline length and the equivalent circumference. Subsequently, a diagonal weight matrix W is constructed based on the local spatial weights. The main diagonal elements of this matrix represent the spatial correlation between each sample point and the location to be estimated. For example, the weights can be calculated based on the Euclidean distance between the sample point and the target location, and the bandwidth parameter is set to approximately 3 km. The regression coefficients are then solved using weighted least squares to obtain the coefficient vector. The first term of the β vector is the estimated intercept coefficient. This intercept coefficient is used to represent the baseline correction when the initial SWOT water level, water width, and shoreline tortuosity are all baseline values; however, in the geographic weighted regression framework, because the weight matrix varies with spatial location, each spatial location will receive a different value. ,Should This is the local intercept coefficient, which describes the spatial difference in basic offset caused by differences in topography, hydrodynamic conditions and observation system errors in different river sections or lake areas, so that water level correction can reflect spatial non-uniformity.

[0046] Of particular importance are the methods for obtaining bandwidth parameters, including: A set of candidate bandwidths is set based on the spatial distribution density of the sample points; In one embodiment, a candidate bandwidth range is set based on the spatial coordinate distribution density of the sample points. Specifically, the Euclidean distances between sample points are first statistically analyzed, and the average distance between the 8th and 12th nearest neighbors is calculated to reflect the typical spatial spacing of the sample points. For example, when the average nearest neighbor distance of the sample points is approximately 1.2 km, a candidate bandwidth sequence can be set within the range of 2 km to 10 km, and a set of candidate bandwidths {2 km, 3 km, 4 km, ..., 10 km} is generated at 1 km intervals. The aforementioned candidate bandwidths are used to describe the spatial scale of the local neighborhood in geographically weighted regression. Each bandwidth corresponds to a local spatial weight range, and its function is to control the number of sample points participating in the local regression calculation.

[0047] Using each candidate bandwidth as a neighborhood scale, local weighted regression fitting is performed on the sample points to calculate the predicted value of each sample point under the corresponding candidate bandwidth, thus obtaining the bandwidth residual sequence. In a further embodiment, local weighted regression fitting is performed on the sample points sequentially using each candidate bandwidth as a neighborhood scale. Specifically, a local neighborhood is established centered on the spatial coordinates of the sample point, and spatial weights are calculated based on the Euclidean distance between the sample point and the center point, for example, using a Gaussian weighting function. ,in This represents the distance from the sample point to the center position. This represents the current candidate bandwidth. Then, weighted least squares is performed by combining the observed independent variable matrix and the dependent variable vector to obtain the local regression coefficients under this candidate bandwidth, and the predicted water level correction for each sample point is calculated accordingly. By subtracting the predicted values ​​from the actual corrections in the residual sequence, a new residual sequence, i.e., the bandwidth residual sequence, is obtained.

[0048] The error standard deviation estimate is calculated based on the bandwidth residual sequence, and the corresponding smoothing matrix is ​​constructed based on the local regression fitting relationship obtained by local weighted regression fitting, and then the trace of the smoothing matrix is ​​calculated. In a further embodiment, the error standard deviation estimate is calculated and a smoothing matrix is ​​constructed. Specifically, the error standard deviation estimate is first calculated based on the residual sequence. ,in The residual value, The number of sample points. To account for the number of regression parameters, in this embodiment... The value is set to 4. Then, a smoothing matrix S is constructed based on the fitting relationship of the locally weighted regression. This matrix describes the linear mapping relationship between the observed and predicted values, and its structure is an n×n matrix, where the nth... The rows represent the weighted contribution of the predicted value of the i-th sample point to each observation. The trace of the matrix can be obtained by summing the diagonal elements of the smoothing matrix. This parameter reflects the effective degrees of freedom of the model. For example, when the sample size is 300, It is usually between 15 and 40.

[0049] The modified Akaike Information Criterion value is calculated for each candidate bandwidth. The specific formula for calculating the modified Akaike Information Criterion is as follows: ; In the formula, To correct the Akaike Information Criterion value, The number of sample points. This is the estimated standard deviation of the error. The trace of the smoothing matrix; Compare the modified Akaike Information Criterion values ​​of each candidate bandwidth, and select the candidate bandwidth corresponding to the smallest modified Akaike Information Criterion value as the bandwidth parameter.

[0050] Optionally, the correction algorithm selected in step S45 corresponding to the initial areal water level field includes: If the data density is less than 0.1 Then, the preset random forest algorithm is selected to perform spatial correction on the entire initial areal water level field; In one embodiment, if the statistically obtained data density of effective observation points in the initial areal water level field is less than 0.1... The entire initial isometric water level field is then subjected to spatial correction processing based on random forest. Specifically, the water level value of the pixels in the initial isometric water level field, the water width at the corresponding location, and the shoreline tortuosity are used as input features, and the residual values ​​between the selected spaceborne lidar points and the initial isometric water level field are used as training labels to construct training samples. The random forest consists of multiple decision trees. In this embodiment, the number of decision trees is set to 200, and the maximum depth of each tree is 12. Three features are randomly selected for node splitting in a single split. Multiple sub-training sets are constructed by bootstrapping the training samples. After training each decision tree, the average value of its output is taken to obtain the predicted correction amount for each pixel position. This correction amount is then superimposed on the corresponding initial isometric water level value to complete the spatial correction of the entire initial isometric water level field.

[0051] If the data density is greater than or equal to 0.1 If the spatial heterogeneity index is greater than 15% or the width variation index is greater than 40%, then the geographic weighted regression model is selected to perform pixel-by-pixel correction on the entire initial isometric water level field.

[0052] In another embodiment, if the statistically obtained data density of effective observation points in the initial areal water level field is greater than or equal to 0.1... If the calculated spatial heterogeneity index is greater than 15% or the width variation index is greater than 40%, then a pixel-by-pixel spatial correction based on geographic weighted regression is performed on the entire initial isometric water level field. Specifically, the spatial coordinates of each pixel in the initial isometric water level field are used as the regression center point, and sample points within the neighborhood bandwidth participate in the local regression calculation. In this embodiment, the bandwidth is set to 6km. Different weights are assigned to the neighborhood sample points through a spatial weighting function, and the local regression coefficients are solved by combining the three variables of SWOT initial water level, water width, and shoreline tortuosity, forming a regression coefficient surface corresponding to the spatial location. Then, the corresponding local regression coefficients are read according to the spatial location of each pixel, and substituted into the regression expression to calculate the water level correction amount for that pixel. The calculated correction amount is superimposed on the initial isometric water level value, thereby realizing the pixel-by-pixel spatial correction of the entire initial isometric water level field.

[0053] Of particular importance is the use of a geographically weighted regression model to perform pixel-by-pixel correction of the entire initial areal water level field, specifically: The regression coefficient surface obtained by the geographic weighted regression model is mapped to the spatial location corresponding to the initial planar water level field. In this embodiment, an inverse distance weighting method can be used to spatially expand the local intercept coefficient, local initial water level regression coefficient, water width regression coefficient, and shoreline tortuosity regression coefficient at the sample point, resulting in four corresponding continuous coefficient surfaces. Subsequently, based on the center coordinates of each grid cell, the four regression coefficient values ​​are read at the corresponding positions and mapped to the same cell position of the initial planar water level field, thereby forming a regression coefficient distribution structure corresponding to the initial planar water level field pixel by pixel.

[0054] Based on the mapping results and the spatial variation trend of the regression coefficient surface, the entire initial planar water level field is corrected by pixels. In a further embodiment, water level correction calculation is performed for each pixel. Specifically, the original water level value of the pixel is used as the input variable, and the water body width and shoreline tortuosity parameters corresponding to that pixel are extracted and substituted into the geographic weighted regression expression to calculate the water level correction amount for that pixel. In this embodiment, a reasonable constraint range of ±2m is set for the correction amount calculated for each pixel to avoid overcorrection caused by local anomaly coefficients. Subsequently, the correction amount is superimposed on the initial water level value of the corresponding pixel to obtain the water level result updated pixel by pixel. The spatial correction processing of the entire initial areal water level field is completed by traversing the entire regular grid, so that the water level distribution can reflect the spatial change trend caused by the water body morphology characteristics.

[0055] Spatial interpolation smoothing is performed based on the pixel correction results to construct a spatially corrected planar water level field.

[0056] In a further embodiment, spatial interpolation smoothing is performed on the corrected water level results to eliminate discontinuous variations between local pixels. Specifically, the corrected water level field is used as input data, a spatial neighborhood structure is constructed within the water surface mask area, and a Gaussian weighted kernel function is used to smooth the pixels within the neighborhood. In this embodiment, the spatial smoothing window is set to 5×5 pixels, and the neighborhood weights are calculated from the Euclidean distance between pixel centers. The water level value of each pixel is weighted and averaged with the water level values ​​of its neighboring pixels to generate a spatially continuous water level distribution result. After smoothing, a spatially corrected planar water level field is formed, enabling the overall water level field to maintain good spatial continuity while preserving local variation characteristics.

[0057] Optionally, step S5 includes: Step S51: Establish the regression relationship between the spatially corrected planar water level field and the spatiotemporal variation information of the water surface mask in order to construct a Gaussian process regression model; In this embodiment, the spatially corrected planar water level field and the water surface mask are temporally aligned within the same spatial coordinate system, and the water area change information corresponding to each observation time is extracted. This embodiment uses a 30m×30m regular grid as the unified spatial unit, and constructs a time-series data table with a one-day time interval. For each grid cell, its water level value and the corresponding water area change ratio at that time are recorded, thus forming a data sample set containing spatial coordinates, time index, water level value, and water area change. Based on this, a Gaussian process regression model structure is constructed, where the model consists of two parts: a mean function and a covariance function. The mean function describes the overall water level change trend, and the covariance function characterizes the influence of spatial location and temporal changes on the correlation of water level, thereby establishing a regression relationship between the water level field and the spatiotemporal changes of the water surface mask.

[0058] Step S52: Use the spatial location and temporal information of the planar water level field as input features for the Gaussian process regression model, and input the model into the Gaussian process regression model; In a further embodiment, the grid center coordinates and corresponding time information in the spatially corrected isometric water level field are used to construct an input feature matrix. In this embodiment, the input feature matrix consists of three columns: the horizontal axis coordinate, the vertical axis coordinate, and the time series number, where the time series number increments daily based on the observation start time. This feature matrix is ​​used as the model input, and the water level values ​​of the corresponding grids are used as the observation output vector. Subsequently, the spatial scale parameters, time scale parameters, and noise parameters in the model are estimated by maximizing the log-marginal likelihood function, and the L-BFGS-B numerical optimization method is used to solve for the parameters. The spatial scale parameter range is set to 1 km to 50 km, and the time scale parameter range is set to 1 day to 30 days to obtain the trained Gaussian process regression model.

[0059] Step S53: The output of the Gaussian process regression model is fused with the spatiotemporal variation information of the water surface mask to generate the predicted water level field values ​​at each time and location, thereby constructing a high spatiotemporal resolution water surface elevation model.

[0060] In a further embodiment, the coordinate information of any spatial location and time node is input into a Gaussian process regression model for water level prediction. In this embodiment, the spatial extent of the water body at each time point is first extracted based on a water surface mask, and the predicted location is restricted within the water surface mask range to ensure consistency between the prediction result and the actual water body distribution. Subsequently, corresponding spatial coordinates and time indices are generated for each predicted location and input into the trained Gaussian process regression model to obtain the predicted water level value. The predicted water level result output by the model is constrained and fused with the water body area change information at the corresponding time point. Pixels located in the water body expansion area have their predicted water level values ​​adjusted using neighborhood weights, thereby generating a continuous water level distribution corresponding to each time point and spatial location. Based on this, a high spatiotemporal resolution water surface elevation model is constructed by reconstructing the water level using regular gridding and performing time series smoothing on the water levels at adjacent time nodes.

[0061] Optionally, this application also provides a water body and water level co-extraction system for performing the water body and water level co-extraction method described above, the water body and water level co-extraction system comprising: Data acquisition module 101 is used to acquire multi-source remote sensing data of the target area within a preset time period, including optical image data, surface water level data of spaceborne radar altimeter and point elevation data of spaceborne lidar. The water body information extraction module 102 is used to extract water body information from optical image data in multi-source remote sensing data, so as to generate a water surface mask and extract water body morphology feature parameters based on the water surface mask. The water level field construction module 103 is used to extract observation data that satisfy the water body morphology characteristic parameters within the corresponding spatiotemporal range of the water surface mask in the planar water level data of the spaceborne radar altimeter, so as to construct the initial planar water level field. The water level field correction module 104 is used to perform spatiotemporal matching between the point elevation data of the spaceborne lidar in the multi-source remote sensing data and the initial surface water level field, and to perform spatial correction on the initial surface water level field. The spatiotemporal fusion module 105 is used to perform spatiotemporal fusion based on the spatially corrected planar water level field and the spatiotemporal variation information of the water surface mask, thereby generating a high spatiotemporal resolution water surface elevation model.

[0062] Therefore, the embodiments should be considered as exemplary and non-limiting in all respects, and the scope of the invention is defined by the appended claims rather than the foregoing description. Thus, all variations falling within the meaning and scope of the equivalents of the application are intended to be included within the invention.

[0063] The above description is merely a specific embodiment of the present invention, enabling those skilled in the art to understand or implement the invention. Various modifications to these embodiments will be readily apparent to those skilled in the art, and the general principles defined herein may be implemented in other embodiments without departing from the spirit or scope of the invention. Therefore, the present invention is not to be limited to the embodiments shown herein, but is to be accorded the widest scope consistent with the principles and novel features of the invention herein.

Claims

1. A method for the coordinated extraction of water body and water level, characterized in that, Includes the following steps: Step S1: Collect multi-source remote sensing data of the target area within a preset time period, including optical image data, surface water level data from spaceborne radar altimeter, and point elevation data from spaceborne lidar. Step S2: Extract water body information from optical image data in multi-source remote sensing data to generate a water surface mask, and extract water body morphology feature parameters based on the water surface mask; Step S3: Extract observation data that satisfy the water body morphology characteristic parameters within the corresponding spatiotemporal range of the water surface mask from the spaceborne radar altimeter water level data, in order to construct the initial water level field; Step S4: Perform spatiotemporal matching between the point elevation data of the spaceborne lidar in the multi-source remote sensing data and the initial areal water level field, and perform spatial correction on the initial areal water level field; Step S5: Perform spatiotemporal fusion based on the spatially corrected planar water level field and the spatiotemporal variation information of the water surface mask to generate a high spatiotemporal resolution water surface elevation model.

2. The method for synergistic extraction of water body and water level according to claim 1, characterized in that, Step S2, generating the water surface mask, includes: The optical image data is preprocessed, the water index of the preprocessed optical image data is calculated, and water pixels are extracted to generate a binary water map. Time series analysis is performed on the binary water body maps of all time phases within a preset time period. Through time series filtering and synthesis, a spatiotemporally continuous water surface probability map is generated. Pixels in the water surface probability map with a probability greater than a preset stable water body determination threshold are identified as water bodies, forming a water surface mask.

3. The method for synergistic extraction of water body and water level according to claim 1, characterized in that, The water morphology parameters in step S2 include water body width and shoreline tortuosity; the extracted water morphology parameters include: The Zhang-Suen thinning algorithm is used, and the center line of the water body with a single pixel width is extracted based on the water surface mask; Sampling points are set at fixed intervals along the centerline of the water body; at each sampling point, the normal line of each sampling point is drawn, and the distance between the intersection of the normal line of each sampling point and the boundary of the water body is calculated as the width of the water body at that sampling point; the widths of all sampling points are statistically analyzed to obtain the average width of the water body. The shoreline length of each independent water body unit is calculated based on the water surface mask, and the circumference of the equivalent circle with the same area as each independent water body unit is calculated. The shoreline tortuosity is obtained by dividing the shoreline length of each independent water body unit by the circumference of the equivalent circle.

4. The method for synergistic extraction of water body and water level according to claim 1, characterized in that, Step S3, which involves constructing the initial planar water level field, includes: The spatial position of each water level observation pixel in the surface water level data of the water surface mask and the spaceborne radar altimeter is superimposed and judged, and the water level observation pixels whose center coordinates fall within the water body pixels of the water surface mask are retained. The water level observation pixels that are marked as water bodies at a frequency greater than or equal to 80% within the preset time tolerance are used as SWOT observation points. The SWOT observation points are rasterized according to their spatial location, and a regular grid is established with a preset grid resolution. The water level observation values ​​in each grid cell of the regular grid are statistically calculated, and the average value is taken as the representative water level value of the grid cell. An initial planar water level field is generated according to the spatial location of the grid cells.

5. The method for synergistic extraction of water body and water level according to claim 1, characterized in that, Step S4 includes: Step S41: Filter lidar points based on the quality identifier field in the point elevation data of the spaceborne lidar; Step S42: Match the selected lidar points with the values ​​of the initial planar water level field at the corresponding locations, and calculate the residual sequence based on the matching results; remove abnormal lidar points and their corresponding data points based on the robustness of the residual sequence. Step S43: Construct a geographic weighted regression model based on the residual sequence remaining after removing abnormal lidar points and water morphology parameters; Step S44: Perform model fitting on the geographic weighted regression model to obtain the regression coefficient surface at each spatial location; Step S45: Select the correction algorithm corresponding to the initial isometric water level field based on the data density, spatial heterogeneity index and width variation index of the initial isometric water level field, and perform spatial correction.

6. The method for synergistic extraction of water body and water level according to claim 5, characterized in that, Constructing a geographically weighted regression model includes: Sample points are selected based on the Euclidean distance between lidar points; Based on the corresponding values ​​of each sample point in the residual sequence and the corresponding spatial coordinates of the water body morphology parameters, a geographically weighted regression model expression is established to construct the geographically weighted regression model; the specific expression of the geographically weighted regression model is as follows: ； in For the first The x-axis coordinates of each sample point For the first The ordinate of each sample point For correction amount, For the first Initial SWOT water level at each sample point; For the first The width of the water body at each sample point For the first shoreline tortuosity at each sample point For random error parameters, For the first Local regression coefficients for each sample point, where Including local intercept coefficient Local initial water level regression coefficient Water body width regression coefficient and the regression coefficient of shoreline tortuosity .

7. The method for synergistic extraction of water body and water level according to claim 6, characterized in that, Methods for obtaining local regression coefficients include: With the first Using the spatial coordinates of the nth sample point as the local regression center, the calculation of the nth... Euclidean distance between each sample point and other sample points; The kernel function with the first square is used to obtain the kernel with the first square as the kernel. A spatial weight set centered at each sample point is used to construct a spatial weight matrix, where the double-square kernel function is specifically: ； In the formula, For the first The sample point and the first Spatial weights of each sample point For the first The sample point and the first Euclidean distance between sample points For bandwidth parameters; The local regression coefficient vector is solved using the weighted least squares method and based on the weight matrix. The formula for solving the local regression coefficient vector is: ； The expression for solving the local regression coefficient vector is: In the formula For the first The local regression coefficient vector of each sample point To observe the independent variable matrix, This is a sequence of water level differences. This is the spatial weight matrix.

8. The method for synergistic extraction of water body and water level according to claim 5, characterized in that, Step S45 selects the correction algorithm corresponding to the initial planar water level field, including: If the data density is less than 0.1 Then, the preset random forest algorithm is selected to perform spatial correction on the entire initial areal water level field; If the data density is greater than or equal to 0.1 If the spatial heterogeneity index is greater than 15% or the width variation index is greater than 40%, then the geographic weighted regression model is selected to perform pixel-by-pixel correction on the entire initial isometric water level field.

9. The method for synergistic extraction of water body and water level according to claim 1, characterized in that, Step S5 includes: Step S51: Establish the regression relationship between the spatially corrected planar water level field and the spatiotemporal variation information of the water surface mask in order to construct a Gaussian process regression model; Step S52: Use the spatial location and temporal information of the planar water level field as input features for the Gaussian process regression model, and input the model into the Gaussian process regression model; Step S53: The output of the Gaussian process regression model is fused with the spatiotemporal variation information of the water surface mask to generate the predicted water level field values ​​at each time and location, thereby constructing a high spatiotemporal resolution water surface elevation model.

10. A water body and water level co-extraction system, characterized in that, For performing the water body and water level co-extraction method as described in claim 1, the water body and water level co-extraction system comprises: The data acquisition module is used to collect multi-source remote sensing data of the target area within a preset time period. The multi-source remote sensing data includes optical image data, surface water level data from spaceborne radar altimeter, and point elevation data from spaceborne lidar. The water body information extraction module is used to extract water body information from optical image data in multi-source remote sensing data, generate a water surface mask, and extract water body morphology feature parameters based on the water surface mask. The water level field construction module is used to extract observation data that satisfy the water body morphology characteristic parameters within the corresponding spatiotemporal range of the water surface mask in the planar water level data of the spaceborne radar altimeter, so as to construct the initial planar water level field; The water level field correction module is used to perform spatiotemporal matching between the point elevation data of the spaceborne lidar in the multi-source remote sensing data and the initial surface water level field, and to perform spatial correction on the initial surface water level field. The spatiotemporal fusion module is used to perform spatiotemporal fusion based on the spatially corrected isometric water level field and the spatiotemporal variation information of the water surface mask, thereby generating a high spatiotemporal resolution water surface elevation model.

Images