A method for constructing a lake hydroperiod zone terrain based on SWOT satellite data

Abstract

The application discloses a kind of lake drawdown zone terrain construction methods based on SWOT satellite data, including extracting monthly scale water body range using lake optical remote sensing image data, monthly drawdown zone area is constructed based on the difference set of monthly water body range and the maximum water body range;With monthly drawdown zone area as spatial constraint, corresponding area is cut out from height measurement data;The spatial filtering of the multi-temporal SWOT satellite height measurement data of cutting area is carried out, with satellite laser height measurement data as true value, the accuracy of corresponding phase SWOT satellite height measurement data is evaluated and fusion weight is determined, then the weighted average of time-phase SWOT satellite height measurement data is obtained Fusion elevation information is obtained, and the drawdown zone terrain is constructed based on the fusion elevation information.The method of the application is particularly suitable for natural lakes and artificial reservoirs with large interannual shrinkage or seasonal variation.The lake drawdown zone terrain constructed based on the method can support the monitoring of large-area long-term lake water level and water volume change.

Description

TECHNICAL FIELD

[0001] The present application relates to the field of hydrology and water resources remote sensing, and particularly relates to a method for constructing a lake hydroperiod terrain based on SWOT satellite data. BACKGROUND

[0002] Although lakes (including artificial reservoirs) only account for 1.8% of the earth's land surface, they contain a variety of ecological services such as water resource supply, climate regulation, biodiversity maintenance, and clean energy production. In recent decades, under the influence of climate change and human activities, the water storage of global lakes is undergoing significant changes, which has important impacts on the surrounding ecosystems and socio-economic development. The changes in lake water storage show significant regional differences on the one hand, and complex interannual and seasonal variation characteristics on the other hand. In the face of the complex changes of global lakes, accurate dynamic monitoring of large-scale lake water storage is crucial for understanding the response characteristics of lakes to climate change and developing sustainable water resource management strategies.

[0003] The rapid development of remote sensing technology has made it an important way to carry out large-scale lake water resource monitoring, and the key is to capture the dynamic changes of lake water level. Traditional satellite altimetry techniques, such as ICESat, ICESat-2, and Sentinel-3, use laser or radar to measure water surface elevation, which has high accuracy. However, the effectiveness of this method is largely limited by the temporal and spatial coverage of altimetry satellites, and it is usually only applicable to large lakes and specific time periods. Another method is based on the lake hydroperiod terrain of low water level measurement, combined with the lake water area range, to calculate the water level elevation corresponding to different water edge lines. This method can complement the coverage of traditional satellite altimetry water level monitoring, but how to obtain the lake hydroperiod terrain is the key. Although there are global public DEM data such as SRTM DEM and TanDEM-X DEM that can provide part of the lake hydroperiod terrain, due to the fixed data collection time, they do not have a multi-temporal update mechanism, and a single data has limited expression ability for the hydroperiod terrain. In addition, through field measurement, complete terrain information of the underwater and hydroperiod area can be obtained, but considering its efficiency and cost, it is not suitable for large-scale application. In summary, there are still bottlenecks in the observation data and key methods for carrying out large-area, long-time series lake water storage change monitoring at the present stage.

[0004] The launch of the SWOT satellite makes up for the shortcomings of traditional satellite altimetry in monitoring lake water levels, and many studies have conducted high-time-frequency dynamic monitoring of lake water levels based on the SWOT satellite. In addition to providing high-time-frequency monitoring data in water areas, the SWOT satellite can also cover wet or rough surfaces with strong backscatter signals, especially in the drawdown zone around the lake, thus having the potential to construct the drawdown zone terrain. However, due to the uncertainty of the quality of the single-time-phase data of the SWOT data and the difference in coverage capacity, the accuracy and continuity of the drawdown zone terrain constructed directly using the original SWOT data cannot be guaranteed. SUMMARY

[0005] In view of the problem that the lake drawdown zone terrain data is scarce and the actual measurement cost is high, the present application provides a lake drawdown zone terrain construction method based on a new generation of surface water resource satellite SWOT, which takes SWOT data as the core, fully utilizes the advantages of multi-time-phase data, eliminates the shortcomings of single-time-phase data in monitoring range and data error, and combines high-precision laser altimetry satellite data to apply the SWOT satellite to the construction of the drawdown zone terrain and verify the effectiveness and feasibility of the drawdown zone terrain construction. The drawdown zone terrain constructed based on this method can provide important basic data for lake (especially continuously shrinking and seasonally rapidly changing lakes) water level change monitoring, especially in combination with high-time-frequency optical image to carry out lake water level and water volume change monitoring, which can effectively make up for the lack of time and space coverage of altimetry satellites and provide technical support for large-scale fine lake water resource change monitoring. With the advantages of high-time-frequency push-broom imaging of the SWOT satellite and the strategy of spatiotemporal fusion of multi-time-phase data, high-precision construction of the lake drawdown zone terrain is realized. The construction result combined with multi-source optical satellite image data with high-time-frequency observation advantages in the target period can realize long-term, high-time-frequency and large-scale water level change monitoring of the lake in the past, present and future through the water edge line method.

[0006] To achieve the above technical purposes, the present application adopts the following technical solutions:

[0007] A lake drawdown zone terrain construction method based on SWOT satellite data, the method comprising:

[0008] extracting the monthly scale water body range using lake optical remote sensing image data, and constructing the monthly drawdown zone area based on the difference set of the monthly water body range and the maximum water body range;

[0009] cutting out the corresponding area from the lake SWOT satellite altimetry data and satellite laser altimetry data with the monthly drawdown zone area as the spatial constraint;

[0010] The multi-temporal SWOT satellite altimetry data of the cropped area is spatially filtered, satellite laser altimetry data is used as true value to evaluate the accuracy of the corresponding temporal SWOT satellite altimetry data, the fusion weight of the temporal SWOT satellite altimetry data is determined based on the accuracy, then the weighted average of the temporal SWOT satellite altimetry data is performed to obtain the fusion elevation information, and the drawdown zone terrain is constructed based on the fusion elevation information.

[0011] In some embodiments of the present application, the lake optical remote sensing image data, the lake SWOT satellite altimetry data and the satellite laser altimetry data cover the maximum water body range of the lake for multiple years.

[0012] In some embodiments of the present application, the method for extracting the monthly scale water body range of the lake by using the optical remote sensing image data comprises:

[0013] The water body is extracted by using different indexes on the optical remote sensing image respectively, and the synthetic water body frequency of the water body extracted by using different indexes is calculated.

[0014] The monthly water body range is extracted based on the synthetic water body frequency.

[0015] In some embodiments of the present application, the SWOT satellite altimetry data is subjected to quality control and multiple rounds of data preprocessing, which comprises:

[0016] The data of the PGC processing version of the SWOT satellite altimetry data is preferentially used, and the data of the PIC processing version of the SWOT satellite altimetry data is used to supplement the missing temporal phase.

[0017] The wse and wse_qual variables of the SWOT satellite altimetry data are extracted and converted into a raster format, and the raster format wse data pixels are screened and matched with wse_qual <a.

[0018] The obtained raster results are converted to a unified reference system.

[0019] Based on the SWOT scientific orbit data, the pixels in the b km wide gap area along the track direction are removed.

[0020] Wherein a and b are preset threshold values.

[0021] In some embodiments of the present application, the method further comprises determining the maximum constructible area of the lake drawdown zone based on the monthly drawdown zone area, and cropping the SWOT satellite altimetry data based on the maximum constructible area, and the pixels exceeding the monthly drawdown zone area are not involved in the fusion calculation.

[0022] In some embodiments of the present application, the spatial filtering manner is:

[0023] The SWOT satellite altimetry data is sequentially stacked in time sequence to construct a spatial cubic data structure.

[0024] construct a spatial filtering block, the z-axis of the spatial filtering block is the number of SWOT satellite altimetry data images, and the initial size of the x-y plane is smaller than the minimum outer envelope rectangle of the lake;

[0025] The spatial filtering block is used to traverse the spatial cube data structure, and invalidation marking is performed on the pixels in the sliding window that exceed the double standard deviation confidence interval.

[0026] After the spatial cube data structure is traversed, the x-y plane size of the spatial filtering block is reduced for the next round of traversal processing, and the process is repeated until the x-y plane size of the spatial filtering block is smaller than a preset value.

[0027] In some embodiments of the application, the initial edge length of the x-y plane of the spatial filtering block is half of the length of the short side of the minimum outer envelope rectangle, and after each traversal is completed, the x-y plane edge length of the spatial filtering block is halved before the next traversal.

[0028] In some embodiments of the application, the calculation of the fusion weight comprises:

[0029] The precision index of the per-time-phase SWOT satellite altimetry data is calculated based on the satellite laser altimetry data as the true value, the initial weight is obtained, and the initial weight of all pixels in each time phase is the same.

[0030] When fusing the per-time-phase SWOT satellite altimetry data corresponding to a pixel, the normalized weight is calculated based on the initial weight of each to-be-fused pixel, and the to-be-fused pixels are weighted and fused based on the normalized weight.

[0031] In some embodiments of the application, the method further comprises using a multi-scale sliding window to remove discrete pixels and fill data holes for the constructed drawdown zone terrain, comprising:

[0032] Pixels with a number of valid pixels in the detection range of the sliding window lower than a minimum threshold are determined as discrete pixels and are removed.

[0033] When the number of valid pixels in the detection range of the sliding window is higher than a preset threshold and there are missing pixels in the window, the mean value of all valid pixels in the detection range of the sliding window is used to fill the pixel value of the missing pixels.

[0034] In some embodiments of the application, the satellite laser altimetry data is ICESat-2 ATL08 product data, and the optical remote sensing image data is Landsat-8 / 9 optical image data.

[0035] The application has the following two advantages:

[0036] (1) The present application is based on SWOT satellite observation data to construct the terrain of lake drawdown zone, which has the characteristics of low cost and high efficiency compared with traditional in-situ measurement under the premise of ensuring data quality. Traditional in-situ measurement is mainly based on lake water level, and the lake drawdown zone terrain is constructed based on ship-borne sonar and unmanned vehicle in high and low water period. The present application makes full use of the high spatio-temporal resolution observation mode of SWOT satellite, combines SIF-WAF spatio-temporal filtering fusion algorithm, realizes efficient construction of lake drawdown zone terrain under large-scale and complex hydrological conditions, effectively breaks through the spatial coverage and cost limit of traditional observation means, significantly improves the global lake sounding database, especially fills the key data gap of the drawdown zone terrain of the lake with continuous shrinkage and seasonal change.

[0037] (2) Another advantage of the present application is that the algorithm is robust, and the constructed drawdown zone terrain has high precision. By fully combining multi-temporal SWOT data, using time and space dimension filtering strategy, completing data self-validation and abnormal value removal, and combining high-precision reference elevation to realize weighted fusion, the elevation error and error of single temporal data are eliminated, and the precision and robustness of the drawdown zone construction are improved. In addition, the combination of multi-temporal SWOT data can continuously monitor the rapid change process of the lake, and realize the capture of the drawdown zone terrain with the largest range at the minimum water level and the minimum water area. The present application breaks through the data and technical bottlenecks faced by existing remote sensing methods in the construction of lake drawdown zone terrain, and is especially suitable for natural lakes and artificial reservoirs with large interannual shrinkage or seasonal change. The lake drawdown zone terrain constructed based on the method can support the monitoring of lake water level and water volume change in a large area and a long time sequence, and provide a scientific basis for the management of lake water resources. At the same time, the high-precision lake drawdown zone terrain constructed based on the present application can also provide basic data for lake water environment governance and improvement of lake ecological service function. BRIEF DESCRIPTION OF DRAWINGS

[0038] Figure 1 is the water body frequency map of reservoir A in embodiment 1 of the present application.

[0039] Figure 2 is the main process flowchart of the inventive method in embodiment 1 of the present application.

[0040] Figure 3 is the complete processing flow of spatial iterative filtering in embodiment 1 of the present application.

[0041] Figure 4 is the processing schematic diagram of weighted average fusion in embodiment 1 of the present application.

[0042] Figure 5 is the drawdown zone terrain construction result of reservoir A based on SWOT satellite data in embodiment 1 of the present application and the error distribution map of ICESat-2 ATL08 reference point data.

[0043] Figure 6 is a reservoir A drawdown zone terrain elevation distribution map constructed based on SWOT satellite data in embodiment 1 of the present application.

[0044] Figure 7 is a reservoir A drawdown zone terrain multi-section elevation result constructed based on SWOT satellite data in embodiment 1 of the present application.

[0045] Figure 8 is a reservoir A drawdown zone terrain correlation analysis with ICESat-2 ATL08 reference points constructed based on SWOT satellite data in embodiment 1 of the present application.

[0046] Figure 9 is a reservoir water level sequence from 2021 to 2023 generated by combining Sentinel-2 optical images using ICESat-2 ATL08 and the reservoir A drawdown zone terrain constructed based on SWOT satellite data in embodiment 1 of the present application. DETAILED DESCRIPTION

[0047] The specific embodiments of the present application will be further described in detail below in conjunction with the accompanying drawings and examples. The following implementation cases are used to illustrate the present application, but not to limit the scope of the present application.

[0048] Embodiment 1

[0049] This embodiment combines the multi-year water body frequency synthesis results of the global surface water dataset (JRC Global Surface Water) and the global surface water dynamics dataset (GLAD Global Surface Water Dynamics), selects a reservoir A as the research object, and takes the construction of the drawdown zone terrain in reservoir A as an example to further illustrate the technical solutions of the present application.

[0050] The geographic spatial distribution of reservoir A is shown in Figure 1 The reservoir is a large reservoir with an area of about 350 km², and its water body range shows obvious seasonal changes, with a large area of drawdown zone range. In this area, due to the rapid change of water and land, it still has high soil moisture when it is exposed, which can produce strong enough backscatter to be observed by SWOT satellite.

[0051] The main method flow of this embodiment is shown in Figure 2 , which specifically includes the following steps:

[0052] Step 1, 1984-2020 maximum water body range as spatial constraint. Download ICESat-2 ATL08 and SWOT Water Mask Raster Image Level-2 data products with 100 m spatial resolution from the data website (https: / / podaac.jpl.nasa.gov / dataset) according to the selected reservoir A maximum water body range as spatial boundary. Filter and organize Landsat-8 / 9 optical image data in the GEE platform that matches the time period and covers the reservoir A area.

[0053] Overlay the obtained multi-source remote sensing data (ICESat-2 ATL08, SWOT data, Landsat-8 / 9 optical image data) within the maximum water body range of reservoir A.

[0054] Step 2, multi-round preprocessing of SWOT data for reservoir A, format conversion and quality control:

[0055] 1) Preferentially download the PGC processing version of SWOT data with better quality, and use the PIC processing version of data for supplement of missing time phases to ensure fusion accuracy and sufficient data volume;

[0056] 2) Extract and convert the wse (water surface elevation) and wse_qual (uncertainty classification) variables of the downloaded SWOT data stored in NetCDF format to raster format;

[0057] 3) Perform quality control using the uncertainty level data variable of the data itself, and use wse_qual < 2 threshold to filter the matched high-quality pixel of raster format wse data;

[0058] 4) Convert the obtained raster results to a unified EGM2008 / WGS84 reference system;

[0059] 5) According to the SWOT scientific orbit data, remove the pixels located in the 20 km wide gap area along the track direction to ensure consistency with the actual use and theory.

[0060] Step 3, Firstly, the cloud-affected images and pixels were removed by cloud cover ratio and quality control band on the GEE platform, and the image cloud cover screening threshold was set to 10%. Based on the Landsat-8 / 9 image data matching the SWOT satellite operation period, the enhanced vegetation index (EVI), normalized difference vegetation index (NDVI), and modified normalized difference water index (MDNWI) were calculated, and the three-index discrimination rule (Formula (1)) was used to extract the water body range of each image. The water body frequency was calculated at the monthly scale based on the three indexes (see Wang, X., Xiao, X., Zou, Z. et al. Gainers and losers of surface and terrestrial water resources in China during 1989-2016. Nat Commun 11, 3471 (2020)). The monthly water body frequency data was extracted by setting a threshold of 25%, and the extracted results were quality controlled and corrected by manual visual inspection to ensure the accuracy of the water-land range division. The monthly drawdown zone range was obtained by removing the monthly water body range from the maximum water body range obtained in step 1.

[0061] (1)

[0062] In the formula, Water refers to the identified water body.

[0063] Step 4, Based on the monthly drawdown zone range obtained in step 3, the SWOT data obtained in step 2 and the ICESat-2 ATL08 data obtained in step 1 were effectively cropped to ensure that the data used for fusion and evaluation were located in the non-water body coverage area, avoiding the influence of water body elevation on the construction of drawdown zone terrain. Since the drawdown zone range varies from month to month, the cropped SWOT data range varies from month to month, but subsequent fusion needs to convert monthly data into a completely overlapping spatial cube data structure, and different range of SWOT data cannot be overlaid. To eliminate the influence of inconsistent drawdown zone range in the same study area on the fusion process, all monthly drawdown zone ranges were synthesized (i.e., taking the union) into the maximum constructable range of the drawdown zone, and all cropped SWOT data were expanded and unified based on this range. The expanded pixels were marked as invalid values not included in subsequent fusion.

[0064] Step 5, Based on the spatial autocorrelation theory and the principle of terrain continuity, a spatial iterative filtering model was constructed to remove outliers from monthly observation data, achieve dynamic convergence of data differences, and improve the quality of fused data.

[0065] AsFigure 3 As shown in FIG. 4, the SWOT data set of the unified geographical boundary in step 4 is serialized and stacked in the time dimension to construct a spatial cube data structure. The initial size of the spatial filtering block is set according to the minimum bounding rectangle range of the reservoir A, and the x-y plane length s is half of the shorter side of the minimum bounding rectangle (s = 150 m) ), and the z-axis is the number of all SWOT data of the reservoir. The minimum bounding range is continuously shifted in the longitude and latitude directions, and the non-overlapping coverage characteristics are maintained during the process. In each geographical plane traversal process, the pixels outside the double standard deviation confidence interval in the sliding window are marked as invalid. The invalidation marking refers to setting the pixel value to an invalid value (such as 9999) that does not participate in calculation. During the iterative filtering process, when the filtering block contains an invalidated pixel, the pixel value is not included in the calculation.

[0066] After each traversal, the geographical plane length of the spatial filtering block is halved (the z-axis length remains unchanged), and the traversal algorithm is re-implemented, and the iteration is continued until the set window length threshold condition (≤ 300 meters / 3 pixel level, in this embodiment, the image spatial resolution is 100 m, and 3 pixels are ≤ 300 m) is met, and the internal self-validation and quality control of the multi-temporal SWOT data are completed.

[0067] Based on the ICESat-2 ATL08 reference elevation data set obtained in step 1, the data set is divided into a training group and a verification group in a ratio of 7:3. The ICESat-2 ATL08 elevation data is taken as the true value, and the MAE of each geographical plane of the filtered SWOT data cube is calculated, so as to determine the weight during the fusion of the SWOT data.

[0068] As shown in FIG. 5, first, the initial weight of the weighted average fusion is calculated according to formula (2) Figure 4 . Then, according to formula (3), the time series data at each pixel position is taken as the calculation unit, and the dynamic normalized weight of each valid pixel in the time series data is calculated . By matching the corresponding valid elevation value, the multi-time sequence fusion result of this pixel position is calculated by using the weighted average method (formula (4)). By reducing the dimension of the SWOT data spatial cube along the z-axis direction, the terrain of the drawdown zone of the reservoir A is preliminarily constructed.

[0069] (2)

[0070] (3)

[0071] (4)

[0072] In the formula,​ MAE for each geographic plane of the filtered SWOT data cube. It is the common weight of each geographic plane, that is, the weight of all valid pixels on the same geographic plane. They are the same. It's important to note that invalid pixels... It is set to zero and excluded from the fusion process. Among them... The weights are dynamically normalized and range from 0 to 1. The variable n represents the number of all valid Z-axis temporal pixels used for fusion (for example, when fusing the Z-axis temporal elevation value corresponding to a certain pixel, n is the number of all valid pixels participating in the fusion on that Z-axis temporal axis). To integrate elevation values, This represents the elevation value of a specific pixel.

[0073] Step 6: Based on the preliminary topography of the drawdown zone of Reservoir A constructed in Step 5, a multi-scale sliding window optimization framework is built to implement a two-step optimization strategy for discrete pixel removal and data hole filling:

[0074] First, discrete pixels are removed. 3×3 and 5×5 sliding windows are applied sequentially to remove discrete pixels that meet the condition that the number of valid values ​​in the window is less than the conservative judgment threshold of 2. This process avoids the generation of new data holes.

[0075] Further data gap filling was performed on the results after discrete pixel removal. A 5×5 sliding window was applied to fill in missing pixels that met the condition of having more than 4 valid pixels in the corresponding window, using the average elevation values ​​of all valid pixels within the window to ensure the continuity of the terrain data. Through this two-step optimization strategy, the initial construction results of Reservoir A were optimized, and the final drawdown zone terrain was obtained.

[0076] like Figure 5 , 6 As shown, the accuracy of the constructed A drawdown zone of the reservoir was evaluated using the remaining 30% of the ICESat-2 ATL08 reference points. The error distribution of the constructed results was uniform, with elevations ranging from 480.12 to 503.39 m, and the elevation distribution pattern basically conformed to the width trend of the upstream and downstream of the reservoir. Figure 7 As shown, the constructed terrain can effectively reflect the changing trends of complete or partial cross-sections in different areas of the reservoir. For cross-section 4 ( Figure 7 In the last image, only fixed elevations can be obtained for the sections passing through permanent water bodies. For example... Figure 8 As shown, MAE=0.92 m, RMSE=1.48 m, R²=0.77, and the constructed drawdown zone topography fits well with ICESat-2ATL08 at the reference point. In summary, this demonstrates the feasibility and effectiveness of this invention in constructing lake drawdown zone topography.

[0077] The corresponding lake drawdown zone terrain can be constructed by the above method, which can be further popularized to large-scale research and can supplement the global lake depth information gap. For example, Figure 9 As shown in the figure, combined with high-frequency remote sensing satellite optical images (Sentinel-2) of a specific period, stable, high-frequency, long-term lake dynamic monitoring of the past, present and future can be realized. In addition, the lake drawdown zone terrain also provides an important data basis for carbon cycle, water pollution prevention and control, and biogeochemical research.

Claims

1. A method for constructing a lake hydroperiod zone terrain based on SWOT satellite data, characterized in that, The method comprises: extracting the water body range on a monthly scale using lake optical remote sensing image data, constructing the monthly drawdown zone area based on the difference set of the monthly water body range and the maximum water body range; cutting out the corresponding area from the lake SWOT satellite altimetry data and satellite laser altimetry data as spatial constraints using the monthly drawdown zone area as the spatial constraint; spatial filtering of the multi-temporal SWOT satellite altimetry data of the cutting area, using satellite laser altimetry data as the true value to evaluate the accuracy of the corresponding temporal SWOT satellite altimetry data, determining the fusion weight of the temporal SWOT satellite altimetry data based on the accuracy, and then performing weighted averaging on the temporal SWOT satellite altimetry data to obtain fused elevation information, and constructing the drawdown zone terrain based on the fused elevation information.

2. The method of claim 1, wherein, The lake optical remote sensing image data, lake SWOT satellite altimetry data and satellite laser altimetry data cover the maximum water body range of the lake for many years.

3. The method of claim 1, wherein, The method for extracting the water body range on a monthly scale using lake optical remote sensing image data comprises: extracting water bodies using different indexes on optical remote sensing images respectively, and calculating the synthetic water body frequency of water bodies extracted by different indexes; extracting the monthly water body range based on the synthetic water body frequency.

4. The method of claim 1, wherein, The SWOT satellite altimetry data has undergone quality control and multiple rounds of data preprocessing, including: preferentially using the PGC processed version of the SWOT satellite altimetry data, and using the PIC processed version of the data to supplement the missing time phases; extracting the wse and wse_qual variables of the SWOT satellite altimetry data and converting them into raster format, and screening the matched raster format wse data pixels with wse_qual < a; converting the obtained raster results to a unified reference system; based on the SWOT scientific orbit data, removing the pixels in the b km wide gap area along the track direction; wherein a and b are preset threshold values.

5. The method of claim 1, wherein, The method further comprises determining the maximum constructible area of the lake drawdown zone based on the monthly drawdown zone area, and performing cutting of the SWOT satellite altimetry data based on the maximum constructible area, and the pixels outside the monthly drawdown zone area are not involved in the fusion calculation.

6. The method according to claim 1 or 5, characterized in that, The spatial filtering method is: serializing and stacking the SWOT satellite altimetry data in chronological order to construct a spatial cubic data structure; constructing a spatial filtering block, the z-axis of the spatial filtering block being the number of SWOT satellite altimetry image data, and the initial size of the x-y plane being smaller than the minimum envelope rectangle of the lake; using the spatial filtering block to traverse the spatial cubic data structure, and performing invalidation marking on the pixels outside the double standard deviation confidence interval in the sliding window; after the spatial cubic data structure is traversed, reducing the x-y plane size of the spatial filtering block for the next round of traversal processing, and repeating until the x-y plane size of the spatial filtering block is smaller than a preset value.

7. The method of claim 6, wherein, The initial length of the x-y plane of the spatial filtering block is half the length of the short side of the minimum envelope rectangle, and after each traversal is completed, the x-y plane length of the spatial filtering block is halved before the next traversal.

8. The method of claim 1, wherein, The calculation of the fusion weight comprises: The precision index of the per-time-phase SWOT satellite altimetry data is calculated by taking satellite laser altimetry data as true value to obtain initial weight, and the initial weight of each pixel in each time phase is the same; When the per-time-phase SWOT satellite altimetry data corresponding to a pixel is fused, the normalized weight is calculated based on the initial weight of each to-be-fused pixel, and the to-be-fused pixel is weighted and fused based on the normalized weight.

9. The method of claim 1, wherein, The method further comprises removing discrete pixels and filling data holes for the constructed water-level-fluctuation belt terrain by using a multi-scale sliding window, comprising: Pixels with an effective pixel number in the detection range of the sliding window lower than a minimum threshold value are determined as discrete pixels and removed; When the effective pixel number in the detection range of the sliding window is higher than a preset threshold value and there are missing pixels in the window, the pixel value of the missing pixel is filled by using the mean value of all effective pixels in the detection range of the sliding window.

10. The method of claim 1, wherein, The satellite laser altimetry data is selected from ICESat-2 ATL08 product data, and the optical remote sensing image data is selected from Landsat-8 / 9 optical image data.

Images