A method for constructing water-land seamless DEM based on multi-source satellite data fusion

Abstract

This invention discloses a method for constructing a seamless land-water DEM based on multi-source satellite data fusion. The study area is divided into bare land, drawdown zone, and perennial waterlogged area based on the annual average water frequency of Sentinel-2 imagery. Using the ICESat-2 ATL08 ground point cloud during the dry season and contemporaneous Sentinel-2 imagery, the underwater topography of the drawdown zone is constructed using an XGBoost model. Using the ATL13 water surface elevation and depth during the wet season and contemporaneous Sentinel-2 imagery, the underwater topography of the perennial waterlogged area is constructed. Based on the water frequency, the results from the drawdown zone and the perennial waterlogged area are fused to generate a complete and continuous underwater topography. The elevation of the bare land ALOS DEM is then corrected using the ATL08 point cloud and fused with the topography of the drawdown zone and the perennial waterlogged area, achieving a seamless connection between land and water topography and constructing a high-precision integrated land-water DEM, thus improving the accuracy of the water area and shoreline topography.

Description

Technical Field

[0001] This invention relates to the field of remote sensing information technology, and in particular to a method for constructing a seamless land and water DEM based on multi-source satellite data fusion. Background Technology

[0002] With the increasing demand for high-precision topographic data in fields such as water resource management, flood simulation, and ecological environment monitoring, the construction of digital elevation models (DEMs) for water bodies and riparian zones has become a key technical challenge. However, traditional DEM data, such as ALOS, often suffer from insufficient accuracy in water-covered areas and significant elevation jumps in the riparian zone, making it difficult to meet the needs of refined hydrological analysis and underwater topographic research. Furthermore, topographic measurements of water-covered areas are significantly affected by seasonal variations, data sparsity, and optical conditions, making it difficult for single-source satellite data to provide continuous, high-precision integrated land and water topographic information.

[0003] ICESat-2 ATL08 data can provide elevation point clouds during the dry season, suitable for topographic construction of drawdown zones, but its water coverage is limited. ATL13 data provides water surface elevation and depth information, which can be used for underwater topographic construction in perennially flooded areas, but the data distribution is sparse and greatly affected by observation conditions. Using either data source alone makes it difficult to simultaneously guarantee land area accuracy and underwater topographic continuity. Optical remote sensing data, such as Sentinel-2 imagery, can provide rich temporal spectral features, offering auxiliary information for water body and shoreline topographic modeling, but it still has limitations in underwater and complex shoreline topographic construction. Summary of the Invention

[0004] To address the above technical problems, this invention provides a method for constructing a seamless land and water DEM based on multi-source satellite data fusion, comprising the following steps: Step S1: Acquire ICESat-2 ATL08 point cloud during the dry season, ICESat-2 ATL13 water surface elevation and depth data during the wet season, and corresponding temporal Sentinel-2 optical images, and perform preprocessing operations to eliminate differences between different data sources, including unified coordinate system transformation, data filtering, and cloud removal; calculate the annual average water frequency of Sentinel-2 images to divide the pixel inundation probability, and divide the study area into drawdown zone, perennial water accumulation zone, and bare land based on water frequency; Step S2: Taking the drawdown zone as the research object, the ICESat-2 ATL08 point cloud during the dry season is used as the true elevation label. Combined with the band reflectance, ratio and neighborhood of the Sentinel-2 image at the same time, the XGBoost machine learning model is used to construct the nonlinear relationship between spectral features and surface elevation to obtain the underwater topography of the drawdown zone. Step S3: Taking the perennial waterlogged area as the research object, the difference between the water surface elevation and water depth of ICESat-2 ATL13 during the high water season is used as the underwater elevation training sample. The Sentinel-2 multispectral image of the same period is matched and the band reflectance, ratio and neighborhood are extracted. The XGBoost regression model is used to establish a nonlinear mapping between the water body spectrum and the underwater elevation to obtain the underwater topography of the perennial waterlogged area. Step S4: Use the annual average water frequency obtained in step S1 to perform weighted fusion of underwater topography; in low-frequency water areas where the water frequency is closer to 0, give higher weight to the underwater topography of the drawdown zone and use the drawdown zone topography obtained in step S2; in high-frequency water areas where the water frequency is closer to 1, give higher weight to the underwater topography of the perennial water accumulation area and use the perennial water accumulation area topography obtained in step S3, and generate a continuous and consistent fused underwater topography after smoothing. Step S5: Using the ICESat-2 ATL08 altimeter as control data, the deviation of the bare area ALOS 12.5-meter DEM is corrected; the elevation error is calculated by spatial matching, and combined with spectral and neighborhood features, the elevation of the bare area is corrected pixel by pixel using XGBoost nonlinear fitting, so that the bare area DEM benchmark is consistent with the satellite altimeter data. Step S6: Perform terrain fusion of drawdown zone, perennial water accumulation area and bare area. The bare area adopts the ALOS DEM obtained in step S5 after correction by ATL08, and the water area adopts the fused underwater terrain obtained in step S4. The elevation jump between land and water is eliminated by spatial stitching and boundary smoothing, and finally the integrated land and water DEM product is constructed.

[0005] The beneficial effects of this invention are: This invention integrates ICESat-2 ATL08 ground point cloud data during the dry season, ATL13 water surface elevation and depth data during the wet season, and Sentinel-2 multispectral imagery from the same period to reconstruct bare land and underwater topography. First, the multi-source data undergoes preprocessing and outlier removal, and water body frequencies are calculated to divide the study area into drawdown zones, perennial waterlogged areas, and bare land. An XGBoost regression model is used to establish a nonlinear mapping between spectral features and the elevations of the drawdown zones and perennial waterlogged areas, thus constructing the underwater topography. Through weighted fusion of annual average water body frequencies, the topography of perennial waterlogged areas is preserved in high-frequency water areas, while the topography of the drawdown zones is preserved in low-frequency water areas, resulting in a smooth and continuous land-water boundary. Simultaneously, ATL08 control points are used to correct deviations and optimize the elevation of the ALOS DEM, eliminating anomalies in bare land areas. After boundary smoothing and local filtering, a high-precision DEM with spatial continuity, uniform accuracy, and natural land-water transition is generated. This fully utilizes the complementarity of multi-source data and is suitable for high-resolution monitoring and research in various types of water bodies. Attached Figure Description

[0006] Figure 1This is a schematic diagram of the overall process of the present invention; Figure 2 The images shown are Sentinel-2 B11 band images of the sample area during the dry and wet seasons in the example. Figure 3 This is a schematic diagram of the annual average water frequency of the Sentinel-2 sample area in the example; Figure 4 This is a schematic diagram of the ATL08 drawdown zone topography construction results in the sample area of ​​the example; Figure 5 This is a schematic diagram of the topographic reconstruction results of the ATL13 perennial waterlogged area in the sample area of ​​the example; Figure 6 This is a schematic diagram of the underwater topography results after the sample area is fused in the example. Figure 7 This is a schematic diagram of the ALOS 12.5-meter DEM of the sample area in the example; Figure 8 This is a schematic diagram of the corrected ALOS DEM of the bare land in the sample area in the example; Figure 9 This is a schematic diagram of the seamless land-water DEM result after fusion of the sample area in the example. Detailed Implementation

[0007] This embodiment provides a method for constructing a seamless land and water DEM based on multi-source satellite data fusion, such as... Figure 1 As shown, it includes the following steps: Step S1: Acquire ICESat-2 ATL08 point cloud during the dry season, ICESat-2 ATL13 water surface elevation and depth data during the wet season, and corresponding temporal Sentinel-2 optical images. Perform preprocessing work such as unified coordinate system transformation, data filtering, and cloud removal to eliminate differences from different data sources. Calculate the annual average water frequency of the Sentinel-2 images to divide the pixel inundation probability and divide the study area into drawdown zones, perennial waterlogged areas, and bare land.

[0008] The specific process of multi-source data preprocessing and data filtering in step S1 is as follows: acquire ICESat-2 ATL08 point cloud data during the dry season, where the variable h_te_best_fit_20m is used to reflect the elevation of the drawdown zone as the high-precision true value of the drawdown zone; at the same time, acquire the water surface elevation (ht_water_surf) and water depth (water_depth) data from the ICESat-2 ATL13 product during the wet season, perform difference processing to obtain the underwater elevation, which is used to construct the underwater topography of the perennial waterlogged area; at the same time, acquire the Sentinel-2 multispectral image of the corresponding time phase; in order to ensure the spatial consistency and analysis accuracy of the multi-source data, all data are uniformly converted to the WGS84 coordinate system and elevation datum to make the point cloud, DEM and imagery accurately aligned in space.

[0009] In the quality control phase, outlier removal was performed on the ATL 08 and ATL 13 point cloud and underwater data. The principle of noise reduction is as follows: ; In the formula, This refers to the elevation or water depth value at a single point. and These are the sample mean and standard deviation, respectively.

[0010] To delineate the study area, water bodies and non-water bodies were identified using Sentinel-2 imagery, and the Modified Normalized Difference Water Index (MNDWI) was used to distinguish between them. ; In the formula, It is the Sentinel-2 green band (B3). It is the shortwave infrared band 1 (B11).

[0011] A water body frequency dataset was constructed based on multi-year Sentinel-2 time-series imagery. The long-term hydrological inundation characteristics were characterized by statistically analyzing the proportion of times each pixel was identified as a water body out of the total number of valid observations throughout the observation period. A water body frequency value closer to 1 indicates that the area has been in a state of long-term water accumulation; conversely, a frequency value closer to 0 indicates that the area has remained exposed for a long time.

[0012] The formula for calculating water body frequency is as follows: ; In the formula, This represents the annual average water body frequency of pixel (x,y); The water mask value at time t is assigned a value of 1 when a pixel is identified as water, and a value of 0 otherwise. To effectively observe the mask value, a value of 1 is assigned when the pixel is not affected by factors such as clouds or cloud shadows and can participate in the classification; otherwise, a value of 0 is assigned. T is the total number of Sentinel-2 images participating in the statistics.

[0013] Based on the fact that water body frequency can effectively characterize the lake's water-land transition process and long-term inundation characteristics, the study area was further divided into bare areas ( ), drawdown zone ( ) and perennial waterlogged areas ( The three regions provide spatial constraints for subsequent zonal terrain reconstruction and seamless integration of land and water.

[0014] Step S2: Taking the drawdown zone as the research object, the ICESat-2 ATL08 high-precision point cloud during the dry season is used as the true elevation label. Combined with the band reflectance, ratio and neighborhood features of the Sentinel-2 image at the same time, the XGBoost machine learning model is used to construct the nonlinear relationship between spectral features and surface elevation to obtain the underwater topography of the drawdown zone.

[0015] The variable h_te_best_fit_20m from the ICESat-2 ATL08 point cloud was used as the ground truth label for elevation, providing accurate measurement data for the elevation of the drawdown zone. Spectral features were constructed by combining visible light (R, G, B) and near-infrared (NIR) bands from concurrently acquired Sentinel-2 multispectral imagery data. Original band reflectance, band ratios, water-related indices (NDWI, Stumpf index, Lyzenga index), and local neighborhood mean and standard deviation features were extracted to reflect surface type, humidity status, and other characteristics. Using the processed point cloud elevation and spectral feature data, an XGBoost machine learning model was constructed to fit the nonlinear relationship between spectral features and surface elevation. Based on the trained XGBoost model, the study area's topography was predicted, enabling the construction of underwater topography within the drawdown zone.

[0016] Step S3: Taking the perennial waterlogged area as the research object, the difference between the water surface elevation and water depth of ICESat-2 ATL13 during the high-water season is used as the underwater elevation training sample. The same Sentinel-2 multispectral image is matched, and features such as band reflectance, ratio and neighborhood are extracted. The XGBoost regression model is used to establish a nonlinear mapping between the water body spectrum and the underwater elevation to obtain the underwater topography of the perennial waterlogged area.

[0017] In underwater terrain construction, the underwater elevation difference is calculated as a training label based on the water surface elevation (ht_water_surf) and water depth (water_depth) data provided by the ICESat-2 ATL13 product: ; In the formula, ht_water_surf is the water surface elevation in the ICESat-2 ATL13 product, and water_depth is the water depth in the ICESat-2 ATL13 product.

[0018] Spectral features were extracted using Sentinel-2 multispectral and SWIR band images, including reflectance, band ratio, normalized water index (NDMI), Lyzenge index, and Stumpf index in the visible, near-infrared, and shortwave infrared bands. Simultaneously, local multi-scale statistical features (mean and standard deviation) were calculated to reflect spectral variations in the spatial neighborhood, improving the stability and modelability of the features. An XGBoost regression model was constructed to establish a nonlinear mapping relationship between spectral features and underwater elevation differences. The trained model was then used to predict underwater topography in the entire perennial waterlogged area, thus constructing the underwater topography of the perennial waterlogged area.

[0019] Step S4: Use the annual average water frequency obtained in step S1 to perform weighted fusion of underwater topography; in low-frequency water areas, the underwater topography of the drawdown zone is the main feature, and the drawdown zone topography obtained in step S2 is used; in high-frequency water areas, the underwater topography of the perennial water accumulation area is the main feature, and the perennial water accumulation area topography obtained in step S3 is used, and after smoothing, a continuous and consistent fused underwater topography is generated.

[0020] When performing underwater topography fusion, the annual average water frequency is calculated to incorporate the drawdown zone topography. Topography of perennial waterlogged areas The water body frequency-weighted fusion is calculated using the following formula: ; In the formula, For water body frequency, The terrain is characterized by perennial water accumulation. The underwater terrain is then processed by mean smoothing to eliminate elevation jumps and discontinuities, resulting in a continuous and clearly defined underwater base terrain product covering the entire area.

[0021] Step S5: Using the annual average ICESat-2 ATL08 high-precision altimeter points as control data, the deviation of the bare area ALOS12.5-meter DEM is corrected; the elevation error is calculated by spatial matching, and combined with spectral and neighborhood features, XGBoost nonlinear fitting is used to correct the elevation of the bare area pixel by pixel, so that the bare area DEM benchmark is consistent with the satellite altimeter data.

[0022] To ensure the consistency of elevation benchmarks in bare areas, annual average ICESat-2 ATL08 measured elevation points were used as high-precision control data to perform systematic bias correction and optimization on the original ALOS 12.5-meter resolution DEM. Spatial matching was performed between the ATL08 ground point cloud and the DEM raster, and the elevation difference of each matched pixel was calculated as a sample of the fitting error. ; In the formula, For ATL08 point cloud elevation data, This is ALOS DEM elevation data.

[0023] Subsequently, using the multispectral features of Sentinel-2 optical images, the original reflectance, logarithmic spectrum, band ratio, and vegetation-related indices (NDVI, RI) were constructed. The spatial features of the neighborhood mean and standard deviation were extracted as feature inputs. A nonlinear mapping relationship between elevation error and spectral features was established through an XGBoost regression model to achieve spatial prediction and correction of DEM system deviations.

[0024] During model training, by dividing the model into training and testing sets, the predicted error raster is applied to the original DEM to perform pixel-by-pixel optimization of the bare area elevation: ; In the formula, The corrected DEM, To correct the predicted elevation errors, the elevation anomalies in bare areas are effectively corrected, ensuring that the bare land DEM benchmark is highly consistent with satellite altimetry data, thus laying a high-precision foundation for the final integrated land and water DEM construction.

[0025] Step S6: Perform terrain fusion of drawdown zone, perennial water accumulation area and bare area. The bare area adopts the ALOS DEM obtained in step S5 after correction by ATL08. The water area adopts the adaptive fusion of underwater terrain obtained in step S4. Through spatial stitching and boundary smoothing, the elevation jump between land and water is eliminated. Finally, a high-precision integrated land and water DEM product with spatial continuity, uniform accuracy and natural and uninterrupted connection between land and water is constructed.

[0026] For bare land areas, DEMs corrected and optimized using ATL08 high-precision control points are prioritized to ensure elevation accuracy. For water areas, underwater topographic results obtained through adaptive fusion are used to combine the topographic information of the drawdown zone and the perennial water accumulation area. During the stitching process, boundary smoothing and local filtering are used to eliminate possible elevation jumps at the water-land interface, making the water-land transition zone natural and continuous. At the same time, outliers and extreme values ​​are removed to ensure that the entire DEM product is spatially continuous and consistent, with uniform accuracy and a natural, uninterrupted connection between water and land.

[0027] In this embodiment, a typical area in the center of a lake is selected as the research object, and the following operation steps are performed: Step 1: Acquiring and preprocessing multi-source remote sensing and altimetry data based on hydrological conditions during the dry and wet seasons. First, high-precision ICESat-2 ATL08 point cloud data was acquired. This product contains various variables, among which h_te_best_fit_20m reflects the actual elevation of the drawdown zone and is suitable as a high-precision ground truth for topographic reconstruction of the drawdown zone. Based on the hydrological time series characteristics of the lake area, ATL08 data from December 2024 to March 2025, during the period of exposed mudflats, were selected to ensure the temporal matching and representativeness of the point cloud data. After acquisition, [further processing / processing] was performed. Outlier removal ensures the accuracy and stability of the point cloud. The processed point cloud is exported as a TXT file in longitude, latitude, and elevation formats, providing standard input for subsequent modeling and training.

[0028] In the topographic reconstruction of the perennially flooded area, water surface elevation (ht_water_surf) and water depth (water_depth) information from the ICESat-2 ATL13 product were acquired, covering the high-water inundation period of the lake area. Considering the sparseness of ATL13 data points, observation data from the high-water season (June-September) of each year from 2019 to 2025 were selected to ensure that underwater topographic features cover typical years. The ATL13 point cloud and water depth data were also processed... Outlier removal was performed, abnormal observation points were cleaned up, and the data was exported as a txt file containing longitude, latitude, water surface elevation, and water depth. Sentinel-2 L2A imagery time-matched to ATL08 and ATL13 was acquired simultaneously. Low cloud cover (<30%) images from January 10, 2025 (dry season) and August 28, 2025 (high water season) were selected. The short near-infrared wavelengths of both images are as follows: Figure 2 As shown, the selected images cover the visible light, near-infrared and short-wave infrared bands (B, G, R, NIR, SWIR) and are exported as txt files (longitude, latitude, blue band, green band, red band, near-infrared band, short-infrared band).

[0029] In Sentinel-2 imagery, water body identification and analysis were performed on multi-temporal Sentinel-2 images covering the entire year of 2025 using the GEE platform. The Improved Normalized Difference Water Index (MNDWI) was calculated using selected image bands, distinguishing water bodies from non-water bodies by specifying the green (G) and shortwave infrared (SWIR) bands. A binary water body mask was generated based on an empirical threshold (MNDWI>0). Then, pixel-level statistics were performed on the entire time-series imagery, calculating the frequency of water body occurrence and the number of valid observations for each pixel to obtain the water body frequency (0–100%). This frequency reflects the perennial spatial distribution of water bodies; the high-frequency region (75–100%) represents perennial water accumulation areas, the mid-frequency region (10–75%) corresponds to drawdown zones, and the low-frequency region (0–10%) represents bare land. The water body frequency results were exported in GeoTIFF format, such as... Figure 3 As shown, the data was converted into a txt file (longitude, latitude, and water body frequency) using Python for subsequent weighted fusion analysis. To ensure consistency of multi-source data in spatial and elevation dimensions, all point clouds, DEMs, and remote sensing images were uniformly projected onto the WGS84 geodetic coordinate system and ellipsoidal elevation datum, achieving unification of planar position and elevation.

[0030] Step 2: Within the drawdown zone, a spectral feature set is constructed using the original blue (B), green (G), red (R), and near-infrared (NIR) bands of Sentinel-2. To address the differences between water bodies and tidal flat substrates, various indices, including NDWI, Stumpf index, and Lyzenga index, are calculated to enhance the spectral response to both water bodies and bare land. Simultaneously, the local mean and standard deviation of the spectra and indices are calculated using a 5×5 neighborhood to reflect spatial neighborhood characteristics. The spectral features are spatially matched one-to-one with the ATL08 high-precision point cloud. Spatial index matching is performed using Python to align the latitude and longitude of the point cloud with the image pixel positions, forming a training dataset. The training dataset contains feature vectors and a target value (ATL08 elevation h_te_best_fit_20m) for subsequent machine learning modeling.

[0031] The XGBoost nonlinear regression model was selected for terrain construction. The model parameters were set as follows: 2000 iteration trees, learning rate 0.05, maximum tree depth 12, subsample ratio 0.8, feature subsample ratio 0.8, and L1 / L2 regularization coefficients (0.5 / 2.0). The training set was divided into an 80% training set and a 20% validation set to evaluate model performance and overfitting, as shown in Table 1.

[0032] Table 1. Accuracy Evaluation of XGBoost Model in the Drawdown Zone

[0033] After training, predictions are made for each pixel within the study area to generate elevation data of the drawdown zone, such as... Figure 4 As shown in the figure, the data was exported as a txt file (longitude, latitude, and low-water season elevation). In verifying accuracy, the surface elevation exposed in the drawdown zone was extracted using ATL03 photon point clouds with a signal confidence level (signal_conf_ph>3) as an independent validation set. The root mean square error (RMSE) and mean absolute error (MAE) were calculated, as shown in Table 2. The results demonstrate that the method in this embodiment can accurately construct the drawdown zone terrain, providing a reliable elevation basis for subsequent terrain fusion.

[0034] Table 2. Accuracy Evaluation of the Dip-off Zone Topography and ATL03 Data

[0035] Step 3: Within the perennially flooded water body area, underwater topographic spectral features are constructed using blue (B), green (G), red (R), near-infrared (NIR), and short-wave infrared (SWIR) bands from Sentinel-2 multispectral imagery. Feature extraction includes raw band reflectance, band ratios (such as G / R, G / B, SWIR / R, SWIR / G), NDMI, and Lyzenga and Stumpf indices to characterize the optical properties of the water body and differences in bottom sediment. Simultaneously, the local mean and standard deviation of the spectra and indices are calculated using a 5×5 neighborhood to reflect the spectral heterogeneity of the spatial neighborhood and improve feature stability.

[0036] The training labels use the water surface elevation (ht_water_surf) and water depth (water_depth) provided by the ATL13 product, and the underwater elevation difference is calculated as the prediction target. The training set and feature set are also matched one-to-one using Python spatial matching. An XGBoost regression model is used for modeling, with parameters consistent with the drawdown zone modeling: 2000 iteration trees, a learning rate of 0.05, a maximum tree depth of 12, a subsample ratio of 0.8, a feature subsample ratio of 0.8, and L1 / L2 regularization coefficients of 0.5 / 2.0, achieving nonlinear prediction of underwater elevation, such as... Figure 5 As shown in Table 3, after the prediction is completed, the underwater elevation is exported as a txt file (longitude, latitude, and elevation during the high-water season).

[0037] Table 3. Accuracy Evaluation of XGBoost Model in Perennial Flooded Areas

[0038] Simultaneously, the accuracy was verified using ATL03 high-signal-confidence photon point clouds, as shown in Table 4. The verification results show that the method in this embodiment can accurately construct the terrain of perennial waterlogged areas, providing a solid foundation for terrain fusion.

[0039] Table 4. Accuracy Evaluation of Topography in Perennially Waterlogged Areas Compared with ATL03 Data

[0040] Step four involves underwater topographic fusion processing. An adaptive weighted fusion strategy is constructed based on water body frequency. In high-frequency water areas, underwater elevation topographic features of perennial water accumulation areas are prioritized to ensure the accuracy of elevation data within the water body. In low-frequency water areas, the elevation of the drawdown zone is prioritized to ensure a natural transition between the shoreline and seasonal water areas. Specifically, the drawdown zone elevation txt file and the perennial water accumulation area elevation txt file are matched pixel-by-pixel according to latitude and longitude, and weighted calculations are performed based on water body frequency to obtain the fused elevation data. The fused elevation data is also exported as a txt file (longitude, latitude, fused elevation), and a 3×3 window neighborhood mean smoothing process is applied to the fused elevation data. Figure 6 As shown, this provides input for subsequent DEM generation.

[0041] Step 5: To ensure consistent elevation datum in the bare area, ICESat-2 ATL08 high-precision altimeter points throughout 2025 were selected as control points. Systematic deviation correction was then performed on the original 12.5-meter resolution ALOS DEM. Figure 7 As shown in the figure. First, spatial matching is performed between the ATL08 point cloud and the DEM, and the elevation difference is calculated for each pixel, with the elevation difference used as the prediction target. Python is used to extract multispectral image features corresponding to the DEM raster, including original band reflectance, band ratio, and vegetation indices (such as NDVI and RI), and spatial features such as neighborhood mean and standard deviation are calculated using a multi-scale sliding window to form the elevation error prediction input. An XGBoost regression model is established to predict the DEM elevation error using spectral and spatial features. The model parameters are consistent with those of the aforementioned underwater topography construction (2000 iteration trees, learning rate 0.05, maximum tree depth 12, subsample ratio 0.8, feature subsample ratio 0.8, L1 / L2 regularization 0.5 / 2.0) to achieve nonlinear elevation error fitting. After the model prediction is completed, the elevation deviation is exported as a txt file (longitude, latitude, elevation deviation), and applied pixel-by-pixel to the original DEM for correction. The accuracy is verified with ATL08 data, as shown in Table 5. The results are as follows. Figure 8 As shown.

[0042] Table 5. Accuracy Evaluation of DEM Before and After Correction

[0043] Step Six: In the final DEM generation stage, for bare land areas, a DEM optimized with high-precision ATL08 point clouds is prioritized to ensure the elevation accuracy of the bare land. For water areas, the fused underwater topographic elevation is used to fully preserve the underwater features of perennial water accumulation areas while also considering drawdown zone information. Elevation data from water and non-water bodies are stitched pixel-by-pixel using a mask to form a preliminary complete DEM. Then, local filtering is performed, using a 3×3 window mean filtering method to eliminate outliers. The final generated DEM has spatial continuity, a unified elevation benchmark, and retains detailed features within water bodies, along shorelines, and on land, without elevation faults or anomalies. The final generated DEM product can be exported as GeoTIFF, such as... Figure 9 As shown.

[0044] Constructing a high-precision, seamless land-water DEM requires fully integrating the advantages of multi-source satellite data to achieve a seamless connection between underwater topography and land. By combining ATL08 drawdown zone data during the dry season, ATL13 perennial water accumulation zone elevation and depth data during the wet season, and corresponding temporal Sentinel-2 optical features, underwater topography can be constructed under different water cover conditions, and fused using the annual average water body frequency. Simultaneously, by combining bare land DEM correction processing, elevation jumps at the land-water junction can be eliminated, achieving a continuous, high-precision integrated land-water DEM, providing reliable basic data for hydrological analysis, environmental monitoring, and water resource management. Based on the above, this embodiment provides a method for constructing a seamless land-water DEM based on multi-source satellite data fusion, which can achieve a seamless connection between land and water topography and construct a high-precision integrated land-water DEM. It integrates the advantages of multi-source satellite temporal and optical features to improve the accuracy of water area and shoreline topography.

[0045] In addition to the embodiments described above, the present invention may have other implementations. All technical solutions formed by equivalent substitution or equivalent transformation fall within the protection scope claimed by the present invention.

Claims

1. A method for constructing a seamless land and water DEM based on multi-source satellite data fusion, characterized in that: Includes the following steps: Step S1: Acquire ICESat-2 ATL08 point cloud during the dry season, ICESat-2 ATL13 water surface elevation and depth data during the wet season, and corresponding temporal Sentinel-2 optical images, and perform preprocessing operations to eliminate differences between different data sources, including unified coordinate system transformation, data filtering, and cloud removal; calculate the annual average water frequency of Sentinel-2 images to divide the pixel inundation probability, and divide the study area into drawdown zone, perennial water accumulation zone, and bare land based on water frequency; Step S2: Taking the drawdown zone as the research object, the ICESat-2 ATL08 point cloud during the dry season is used as the true elevation label. Combined with the band reflectance, ratio and neighborhood of the Sentinel-2 image at the same time, the XGBoost machine learning model is used to construct the nonlinear relationship between spectral features and surface elevation to obtain the underwater topography of the drawdown zone. Step S3: Taking the perennial waterlogged area as the research object, the difference between the water surface elevation and water depth of ICESat-2 ATL13 during the high water season is used as the underwater elevation training sample. The Sentinel-2 multispectral image of the same period is matched and the band reflectance, ratio and neighborhood are extracted. The XGBoost regression model is used to establish a nonlinear mapping between the water body spectrum and the underwater elevation to obtain the underwater topography of the perennial waterlogged area. Step S4: Use the annual average water frequency obtained in step S1 to perform weighted fusion of underwater topography; in low-frequency water areas where the water frequency is closer to 0, give higher weight to the underwater topography of the drawdown zone and use the drawdown zone topography obtained in step S2; in high-frequency water areas where the water frequency is closer to 1, give higher weight to the underwater topography of the perennial water accumulation area and use the perennial water accumulation area topography obtained in step S3, and generate a continuous and consistent fused underwater topography after smoothing. Step S5: Using the ICESat-2 ATL08 altimeter as control data, the deviation of the bare area ALOS 12.5-meter DEM is corrected; the elevation error is calculated by spatial matching, and combined with spectral and neighborhood features, the elevation of the bare area is corrected pixel by pixel using XGBoost nonlinear fitting, so that the bare area DEM benchmark is consistent with the satellite altimeter data. Step S6: Perform terrain fusion of drawdown zone, perennial water accumulation area and bare area. The bare area adopts the ALOS DEM obtained in step S5 after correction by ATL08, and the water area adopts the fused underwater terrain obtained in step S4. The elevation jump between land and water is eliminated by spatial stitching and boundary smoothing, and finally the integrated land and water DEM product is constructed.

2. The method for constructing a seamless land and water DEM based on multi-source satellite data fusion according to claim 1, characterized in that: In step S1, ICESat-2 ATL08 point cloud data during the dry season is acquired, where the variable h_te_best_fit_20m is used to reflect the elevation of the drawdown zone and serves as the true value of the drawdown zone. Simultaneously, water surface elevation and water depth data from ICESat-2 ATL13 products during the wet season are acquired, and the difference is processed to obtain the underwater elevation, which is used to construct the underwater topography of the perennial waterlogged area. At the same time, Sentinel-2 multispectral images of the corresponding time phase are acquired. All data were uniformly converted to the WGS84 coordinate system and elevation datum to align point clouds, DEMs, and images in space; outliers were removed from ATL 08 and ATL 13 point clouds and underwater data. Water bodies and non-water bodies are identified in Sentinel-2 imagery. A water body frequency dataset is constructed based on multi-year Sentinel-2 time-series imagery. The water body frequency is obtained by statistically analyzing the proportion of times each pixel is identified as a water body out of the total number of valid observations throughout the observation period. At that time, the study area was divided into bare areas; water body frequency At that time, the study area was divided into drawdown zones; water body frequency At that time, the study area was divided into areas with perennial water accumulation.

3. The method for constructing a seamless land and water DEM based on multi-source satellite data fusion according to claim 2, characterized in that: In step S1, outlier removal is performed on the ATL 08 and ATL 13 point cloud and underwater data. The principle of noise reduction is as follows: ； In the formula, This is the original data. and These are the sample mean and standard deviation, respectively. The value is the standardized value; Water bodies and non-water bodies were identified in Sentinel-2 imagery, and the Modified Normalized Differential Water Index (MNDWI) was used to distinguish between them. ； In the formula, This is the Sentinel-2 green band. It is the shortwave infrared band 1; A water body frequency dataset was constructed based on multi-year Sentinel-2 time-series imagery. The long-term hydrological inundation characteristics were characterized by statistically analyzing the proportion of times each pixel was identified as a water body out of the total number of valid observations throughout the observation period. The formula for calculating the water body frequency is as follows: ； In the formula, This represents the annual average water body frequency of pixel (x,y); The water mask value at time t is assigned a value of 1 when a pixel is identified as water, and a value of 0 otherwise. To effectively observe the mask value, a value of 1 is assigned when the pixel is not affected by cloud shadows and can participate in the classification; otherwise, a value of 0 is assigned. T is the total number of Sentinel-2 images participating in the statistics.

4. The method for constructing a seamless land and water DEM based on multi-source satellite data fusion according to claim 1, characterized in that: In step S2, the terrain of the drawdown zone is constructed. The h_te_best_fit_20m variable in the ICESat-2 ATL08 point cloud is used as the true elevation label to provide a measurement basis for the elevation of the drawdown zone. Spectral features are constructed by combining the visible and near-infrared bands of the Sentinel-2 multispectral image data acquired at the same time. The original band reflectance, band ratio, water-related index, and local neighborhood mean and standard deviation are extracted. The water-related index includes NDWI, Stumpf index, and Lyzenga index. Using the processed point cloud elevation and spectral feature data, an XGBoost machine learning model is constructed to fit the nonlinear relationship between spectral features and surface elevation. Based on the trained XGBoost model, the range terrain of the study area is predicted to complete the construction of the underwater terrain of the drawdown zone.

5. The method for constructing a seamless land and water DEM based on multi-source satellite data fusion according to claim 1, characterized in that: In step S3, the terrain of the perennial waterlogged area is constructed. In the underwater terrain construction, the underwater elevation difference is calculated as a training label based on the water surface elevation and water depth data provided by the ICESat-2 ATL13 product. ； In the formula, ht_water_surf is the water surface elevation in the ICESat-2 ATL13 product, and water_depth is the water depth in the ICESat-2 ATL13 product; Spectral features were extracted using Sentinel-2 multispectral and SWIR band images, including reflectance, band ratio, normalized water index, Lyzenge index, and Stumpf index in the visible, near-infrared, and shortwave infrared bands. Simultaneously, local multi-scale statistical features reflecting spectral variations in the spatial neighborhood were calculated, including mean and standard deviation. An XGBoost regression model was constructed to establish a nonlinear mapping relationship between spectral features and underwater elevation differences. The trained model was then used to predict the underwater topography of the entire perennial waterlogged area, thus completing the underwater topography construction of the perennial waterlogged area.

6. The method for constructing a seamless land and water DEM based on multi-source satellite data fusion according to claim 1, characterized in that: In step S4, underwater topography fusion is performed. By calculating the annual average water frequency, the topography of the drawdown zone and the topography of the perennial water accumulation area are weighted and fused according to the water frequency, as shown in the following formula: ； In the formula, For water body frequency, The terrain is characterized by perennial water accumulation. The underwater topography is represented by the drawdown zone, and the merged underwater topography is then smoothed using mean value processing.

7. The method for constructing a seamless land and water DEM based on multi-source satellite data fusion according to claim 1, characterized in that: In step S5, DEM deviation correction is performed. Using the annual average ICESat-2 ATL08 measured elevation points as control data, the original ALOS 12.5-meter resolution DEM is systematically corrected and optimized. The ATL08 ground point cloud and the DEM raster are spatially matched, and the elevation difference of each matched pixel is calculated as a sample to be fitted for error. ； In the formula, For ATL08 point cloud elevation data, This is ALOS DEM elevation data.

8. The method for constructing a seamless land and water DEM based on multi-source satellite data fusion according to claim 6, characterized in that: In step S5, the original reflectance, logarithmic spectrum, band ratio and vegetation correlation index are constructed using the multispectral features of Sentinel-2 optical image. The vegetation correlation index includes NDVI and RI. The spatial features of neighborhood mean and standard deviation are extracted as feature input. The nonlinear mapping relationship between elevation error and spectral features is established through XGBoost regression model to perform spatial prediction and correction of DEM system deviation.

9. The method for constructing a seamless land and water DEM based on multi-source satellite data fusion according to claim 8, characterized in that: In step S5, during model training, by dividing the training set and the test set, the predicted error raster is applied to the original DEM to perform pixel-by-pixel optimization of the bare area elevation: ； In the formula, The corrected DEM, For the predicted elevation error, This is ALOSDEM elevation data; After the elevation anomaly in the bare area was corrected, the bare area DEM benchmark was consistent with the satellite altimetry data.

10. The method for constructing a seamless land and water DEM based on multi-source satellite data fusion according to claim 1, characterized in that: In step S6, a seamless land and water DEM is constructed. The bare land area uses a DEM that has been corrected and optimized by ATL08 control points; the water area uses the underwater topographic results obtained by adaptive fusion, in combination with the topographic information of the drawdown zone and the perennial water accumulation area. During the stitching process, boundary smoothing and local filtering are used to eliminate possible elevation jumps at the land-water interface, while outliers and extreme values ​​are removed.

Images