High-precision DEM fusion construction method based on LiDAR photon statistics

Abstract

The invention discloses a high-precision DEM fusion construction method based on LiDAR photon statistics, and the method comprises the steps: collecting ICESat-2 L3A and GEDI L2A satellite-borne LiDAR data in a region, and carrying out the extraction and strip filtering processing of photon data; integrating a plurality of open source DEM products, and carrying out geodetic reference calibration, spatial offset correction and pixel-level alignment; constructing a photon spatial index by using a CK-Tree algorithm; retrieval of at most 50 nearest photon points is carried out on each DEM pixel position, and errors of each DEM relative to the photon points are calculated; the DEM pixel precision is sorted, DEM data with the minimum error is selected to allocate an elevation value for each pixel, and therefore the optimal DEM is generated; and a conservative smooth filtering technology is adopted for the final optimal DEM, and spatial discontinuity is eliminated. According to the method, the optimal DEM data can be automatically integrated, and a rapid and efficient solution is provided.

Description

Technical Field

[0001] The present invention relates to the field of geographic information, and in particular to a high-precision DEM fusion construction method based on LiDAR photon statistics. Background Art

[0002] DEMs are essential data in the field of geographic information. They store surface elevation information in a grid format and serve fields such as hydrology, geomorphology, urban planning, and environmental monitoring. Although a variety of open-source DEM products exist globally, such as ASTERGDEM, SRTM V3, and NASA DEM, these products exhibit significant variations in accuracy and quality across regions, primarily due to the diverse sensor technologies and data processing methods they rely on. To overcome the limitations of a single DEM product, researchers often employ fusion techniques to integrate the strengths of different DEMs and improve the accuracy and reliability of terrain data.

[0003] However, existing DEM fusion methods, such as simple averaging, weighted averaging, and variational models, often lack sufficient flexibility and interpretability and fail to fully utilize the latest remote sensing technologies. Furthermore, these methods fail to consider the dynamic updating of DEM data, limiting their application in monitoring rapidly changing terrain.

[0004] The emergence of spaceborne LiDAR technology, particularly ICESat-2 L3A and GEDI L2A data, has provided a new approach to acquiring high-precision terrain data. These sensors can penetrate vegetation layers and obtain precise three-dimensional information about the subsurface structure. Although spaceborne LiDAR data is spatially discrete, its high-precision measurements offer new possibilities for DEM fusion. Summary of the Invention

[0005] The purpose of this invention is to provide a high-precision DEM fusion construction method based on LiDAR photon statistics. The technical problem it aims to solve is how to achieve high-precision fusion generation of DEM through photon-based statistical technology driven by space-borne LiDAR data.

[0006] The technical solution to achieve the purpose of the present invention is: a high-precision DEM fusion construction method based on LiDAR photon statistics, comprising the following steps:

[0007] Collect ICESat-2 L3A and GEDI L2A space-borne LiDAR data in the area, and perform photon data extraction and strip filtering;

[0008] Integrate multiple open source DEM products and perform geodetic calibration, spatial offset correction, and pixel-level alignment;

[0009] The CK-Tree algorithm is used to construct a photon spatial index based on ICESat-2 L3A and GEDI L2A data;

[0010] Retrieve up to 50 nearest neighbor photon points for each DEM pixel location and calculate the error of each DEM relative to these photon points;

[0011] The accuracy of DEM pixels is sorted according to the error estimation value, and the DEM data with the smallest error is selected to assign elevation values ​​to each pixel, thereby generating the optimal DEM;

[0012] The final optimal DEM is smoothed and filtered to eliminate spatial discontinuities.

[0013] A computer device comprises a memory, a processor and a computer program stored in the memory and executable on the processor, wherein the steps of the above method are implemented when the processor executes the program.

[0014] A computer-readable storage medium stores a computer program, which implements the steps of the above method when executed by a processor.

[0015] A computer program product comprises a computer program, which implements the steps of the above method when executed by a processor.

[0016] Compared with the prior art, the present invention has the following beneficial effects:

[0017] The present invention is mainly applicable to the automatic fusion scenario of open source DEM. It can automatically, quickly and efficiently generate DEM products with the highest vertical accuracy in the region based on open source DEM, providing a new method for spatial data fusion.

[0018] This method not only improves the accuracy and interpretability of DEM fusion but also avoids arbitrary assumptions about DEM degradation factors by utilizing actual measured data as a priori conditions for fusion. Another significant advantage of this method is its scalability. It can adapt to continuously updated LiDAR and open-source DEM data, providing an effective solution for monitoring dynamic terrain changes. The application potential of this method will further expand as the amount of open-source DEM and space-borne LiDAR data increases, making it particularly suitable for geoscience research and applications that require high-precision terrain information. BRIEF DESCRIPTION OF THE DRAWINGS

[0019] The present invention will be further described below with reference to the accompanying drawings.

[0020] Figure 1 This is a flow chart of the high-precision DEM fusion construction method based on LiDAR photon statistics.

[0021] Figure 2It is the nearest photon search, photon weight allocation, error sorting and pixel-by-pixel optimal terrain judgment.

[0022] Figure 3 The ICESat-2 L3A metric is used to fused DEM and other open source products for accuracy.

[0023] Figure 4 The GEDI L2A metric is used to measure the accuracy of fused DEM and other open source products. DETAILED DESCRIPTION

[0024] This paper proposes a method for fusion construction of a high-precision digital elevation model (DEM) based on LiDAR photon statistics. The specific steps are as follows:

[0025] Step 1. LiDAR data collection: Collect ICESat-2 L3A and GEDI L2A satellite LiDAR data for the area. For areas with stable terrain, use the entire five-year data; for areas with unstable terrain, use only the current year's data.

[0026] Step 2. ICESat-2 L3A photon information processing: Extract key parameters from the ICESat-2 L3A data, including latitude and longitude, best-fit altitude, altitude uncertainty, number of photons, interpolated altitude, median altitude, atmospheric cloud flag, and night flag. Perform high-confidence photon screening, with the following criteria: a night flag of 1, a photon count of at least 50, a difference between the interpolated and median altitudes of less than 2 meters, an atmospheric cloud flag of 0 or 1, and an altitude uncertainty of less than 10 meters.

[0027] Step 3. ICESat-2 L3A strip filtering: compare the MAE of all photons in a strip with the ETOPO2022 error, and exclude data from entire strips where the error exceeds the 95th percentile;

[0028] Step 4. GEDI L2A photon information processing: Extract key parameters from the GEDI L2A data, including timestamp, lowest mode latitude and longitude, degradation flag, altitude deviation flag, quality flag, lowest mode altitude, and reception evaluation flag; perform high-confidence photon screening, with conditions including night time, quality flag of 1, degradation flag not equal to 0, altitude deviation flag of 0, and reception evaluation flag of 0;

[0029] Step 5. GEDI L2A strip filtering: compare the MAE of all photons in a strip with the ETOPO2022 error, and exclude data from entire strips where the error exceeds the 90th percentile;

[0030] Step 6. DEM data collection: Obtain all open source DEM products in the region, including GDEM V3, AW3D30, COP30, FABDEM, NASADEM, SRTM V3, and TAN30;

[0031] Step 7. Geodetic datum coordination: Convert the geoid heights of GDEM V3, AW3D30, NASADEM, and SRTM V3 from EGM96 to EGM08 to be consistent with COP30, FABDEM, and TAN30. The conversion formula is:

[0032]

[0033] Where, and Represent the geoid heights based on EGM96 and EGM08, respectively. and represent the geoid intervals between the EGM96 and EGM08 quasi-geoids and the WGS84 ellipsoid, respectively;

[0034] Step 8. Spatial offset correction and pixel-level alignment: Combine ICESat-2 L3A and GEDI L2A to create a fused LiDAR product , using bilinear interpolation to The elevation of the intersection of the mid-photon and the DEM is sampled and the MAE is calculated. The position with the smallest MAE is used as the optimal offset position. The DEM is moved in a 30m×30m grid with an increment of 0.5m. The sampling and error calculation are repeated to determine the offset. The DEM after the offset is aligned to all pixels using bilinear resampling.

[0035] Step 9. Build photon index: Use CK-Tree algorithm to Create a spatial index;

[0036] Step 10. Nearest neighbor photon search: for a specific pixel position , searching for no more than 50 nearest neighbor photons within a maximum radius of 5 km, that is:

[0037]

[0038] Where, represents the set of searched photons, represents photons, represents the search radius, represents the number of photons;

[0039] Step 11. Photon weight distribution: Figure 2As shown in (a), the size of the photon shows its distance from the central pixel, and the weight is assigned according to the type of photon and its proximity to the pixel. :

[0040]

[0041] Where, Photons The weight of is the distance between the photon and the pixel location, is the photon source function, defined as:

[0042]

[0043] Step 12. Calculate the error of each DEM at the pixel location :

[0044]

[0045] Where, Photons Compared to The mean absolute error.

[0046] Step 13. Error sorting and pixel-by-pixel optimal terrain determination: For each pixel's DEM, select the DEM with the lowest MAE as the elevation value of the current pixel;

[0047] Step 14. Smoothing Filtering: Apply a conservative smoothing filter to the optimal DEM to eliminate potential spatial discontinuities.

[0048] The present invention is described in detail below based on the accompanying drawings and embodiments to make the technical route and operation steps of the present invention clearer.

[0049] Example

[0050] The Loess Plateau covers a wide range of landform types and can comprehensively demonstrate the results and accuracy of DEM fusion. Therefore, the Loess Plateau was selected as the experimental area in this embodiment. The process of the high-precision DEM fusion construction method based on LiDAR photon statistics proposed in this embodiment is as follows: Figure 1 As shown, the following steps are included:

[0051] The first step, LiDAR data collection, involves collecting ICESat-2 L3A and GEDI L2A LiDAR data for the region. This step includes 1,897 strips of ICESat-2 L3A data from 2018 to 2022 and 1,883 strips of GEDI L2A data from 2019 to 2023.

[0052] The second step involves processing ICESat-2 L3A photon information: Key parameters are extracted from the ICESat-2 L3A data, including latitude and longitude, best-fit altitude, altitude uncertainty, number of photons, interpolated altitude, median altitude, atmospheric cloud flag, and night flag. High-confidence photon screening is performed, with criteria including a night flag of 1, a photon count of at least 50, a difference between the interpolated and median altitudes of less than 2 meters, an atmospheric cloud flag of 0 or 1, and an altitude uncertainty of less than 10 meters.

[0053] The third step is ICESat-2 L3A band filtering: comparing the MAE of all photons in a band with the ETOPO2022 error, and excluding data from entire bands where the error exceeds the 95th percentile. After this step, a total of 2,454,380 reliable ICESat-2 L3A photons were acquired.

[0054] Step 4: GEDI L2A photon information processing: Extract key parameters from the GEDI L2A data, including timestamp, lowest mode latitude and longitude, degradation flag, altitude deviation flag, quality flag, lowest mode altitude, and reception evaluation flag. High-confidence photon screening is performed, based on conditions including nighttime, quality flag being 1, degradation flag not being 0, altitude deviation flag being 0, and reception evaluation flag being 0.

[0055] The fifth step is GEDI L2A band filtering: the MAE of all photons in a band is compared to the ETOPO2022 error, and the data of the entire band with an error exceeding the 90th percentile is excluded. After this step, a total of 10,187,997 reliable GEDI L2A photons were obtained.

[0056] Step 6: DEM data collection: Obtain all open source DEM products in the area, including GDEM V3, AW3D30, COP30, FABDEM, NASADEM, SRTM V3 and TAN30.

[0057] Step 7: Geodetic datum coordination: Convert the geoid heights of GDEM V3, AW3D30, NASADEM, and SRTM V3 from EGM96 to EGM08 to be consistent with COP30, FABDEM, and TAN30. The conversion formula is:

[0058]

[0059] Where, and Represent the geoid heights based on EGM96 and EGM08, respectively. and Represents the geoid intervals between the EGM96 and EGM08 quasi-geoids and the WGS84 ellipsoid, respectively.

[0060] Step 8: Spatial offset correction and pixel-level alignment: Combine ICESat-2 L3A and GEDI L2A to create a fused LiDAR product , using bilinear interpolation to The elevation of the intersection of the mid-photon and the DEM is sampled and the MAE is calculated. The DEM is shifted in 0.5m increments within a 30m x 30m grid, and the sampling and error calculation are repeated to determine the offset. The shifted DEM is then aligned using bilinear resampling.

[0061] The ninth step is to build a photon index: using the CK-Tree algorithm, Create a spatial index;

[0062] Step 10: Nearest neighbor photon search: for a specific pixel position , searching for no more than 50 nearest neighbor photons within a maximum radius of 5 km, that is:

[0063]

[0064] Where, represents the set of searched photons, represents photons, represents the search radius, Represents the number of photons. Figure 2 As shown in (a), during the pixel search process, the farthest photon distance is 452.15m, which is much smaller than the set maximum photon search radius.

[0065] Step 11: Photon weight distribution: Figure 2 As shown in (a), the size of the photon shows its distance from the central pixel, and then a weight is assigned based on the type of photon and its proximity to the pixel. :

[0066]

[0067] Where, Photons The weight of is the distance between the photon and the pixel location, is the number of photons, is the photon source function, defined as:

[0068]

[0069] Step 12: Calculate the error of each DEM in the pixel location :

[0070]

[0071] Where, Photons Compared to The mean absolute error.

[0072] Step 13: Error sorting and pixel-by-pixel optimal terrain determination: Figure 2 As shown in (b), the accuracy of the current pixel in the open source DEM is ranked from high to low as follows: COPDEM, TAN30, FABDEM...; Figure 2 As shown in (c), for each pixel's DEM, the DEM with the lowest MAE is selected as the elevation value of the current pixel;

[0073] Step 14: Smoothing filter: Apply a conservative smoothing filter to the optimal DEM to eliminate potential spatial discontinuities. Figure 3 As shown in Figure 2, the accuracy of ICESat-2 L3A fusion DEM and other open source products is measured. Figure 4 Figure 2 shows the accuracy of the fused DEM and other open source products using GEDI L2A. It can be seen that regardless of the validation metric, the accuracy of the fused product is much higher than that of existing open source DEMs, demonstrating the effectiveness of this method.

[0074] The present invention is mainly applicable to the automatic fusion scenario of open source DEM. It can automatically, quickly and efficiently generate DEM products with the highest vertical accuracy in the region based on open source DEM, providing a new method for spatial data fusion.

[0075] The above descriptions are merely embodiments of the present invention and are not intended to limit the patent scope of the present invention. Any equivalent structure or equivalent process transformation made using the contents of the present invention's description and drawings, or directly or indirectly applied in other related technical fields, are also included in the patent protection scope of the present invention.

Claims

1. A high-precision DEM fusion construction method based on LiDAR photon statistics, characterized by: The following steps are involved: Collect ICESat-2 L3A and GEDI L2A space-borne LiDAR data in the area, and perform photon data extraction and strip filtering; Integrate multiple open source DEM products and perform geodetic calibration, spatial offset correction, and pixel-level alignment; The CK-Tree algorithm is used to construct a photon spatial index based on ICESat-2 L3A and GEDI L2A data; Retrieve up to 50 nearest neighbor photon points for each DEM pixel location and calculate the error of each DEM relative to these photon points; The accuracy of DEM pixels is sorted according to the error estimation value, and the DEM data with the smallest error is selected to assign elevation values ​​to each pixel, thereby generating the optimal DEM; The final optimal DEM is smoothed and filtered to eliminate spatial discontinuities.

2. The high-precision DEM fusion construction method based on LiDAR photon statistics according to claim 1 is characterized in that: Perform photon information extraction and stripe filtering on ICESat-2 L3A data, specifically: Extract key parameters from ICESat-2 L3A data, including latitude and longitude, best-fit altitude, altitude uncertainty, number of photons, interpolated altitude, median altitude, atmospheric cloud flag, and night flag; and perform high-confidence photon screening; Compare the MAE of all photons in a stripe to the ETOPO2022 error, excluding data from entire stripes where the error exceeds the 95th percentile.

3. The high-precision DEM fusion construction method based on LiDAR photon statistics according to claim 2 is characterized in that: The high-confidence photon screening is performed under the following conditions: the night flag is 1, the number of photons is not less than 50, the difference between the interpolation and the median height is less than 2 meters, the atmospheric cloud flag is 0 or 1, and the height uncertainty is less than 10 meters.

4. The high-precision DEM fusion construction method based on LiDAR photon statistics according to claim 1 is characterized in that: Photon information extraction and stripe filtering are performed on the GEDI L2A data, specifically: Extract key parameters from GEDI L2A data, including timestamp, lowest mode latitude and longitude, degradation flag, altitude deviation flag, quality flag, lowest mode altitude, and reception assessment flag; perform high-confidence photon screening; Compare the MAE of all photons in a stripe to the ETOPO2022 error, excluding data from entire stripes where the error exceeds the 90th percentile.

5. The high-precision DEM fusion construction method based on LiDAR photon statistics according to claim 4 is characterized in that: The conditions for executing high-confidence photon screening include night time, quality flag being 1, degradation flag being not 0, height deviation flag being 0, and reception evaluation flag being 0.

6. The high-precision DEM fusion construction method based on LiDAR photon statistics according to claim 1 is characterized in that: Integrate multiple open source DEM products and perform geodetic datum calibration, spatial offset correction, and pixel-level alignment, specifically: Access all open source DEM products in the region, including GDEM V3, AW3D30, COP30, FABDEM, NASADEM, SRTM V3, and TAN30; Converted geoid heights for GDEM V3, AW3D30, NASADEM, and SRTM V3 from EGM96 to EGM08 to align with COP30, FABDEM, and TAN30; ICESat-2 L3A and GEDI L2A are combined to create a fused LiDAR product S. The elevation of the intersection of photons and DEM in S is sampled using bilinear interpolation and the MAE is calculated.

7. The high-precision DEM fusion construction method based on LiDAR photon statistics according to claim 1 is characterized in that: For each DEM pixel location, a maximum of 50 nearest neighbor photon points are retrieved, and the error of each DEM relative to these photon points is calculated, specifically: For the DEM pixel location, no more than 50 nearest neighbor photons are searched within a maximum radius of 5 km; The size of the photon indicates its distance from the central pixel, with weights assigned based on the type of photon and its proximity to the pixel. : ； Where, Photons The weight of is the distance between the photon and the pixel location, is the number of photons, is the photon source function, defined as: ； Calculate the error of each DEM in the pixel location : ； Where, Photons Compared to The mean absolute error.

8. A computer device comprising a memory, a processor, and a computer program stored in the memory and executable on the processor, wherein: When the processor executes the program, the steps of the method according to any one of claims 1 to 7 are implemented.

9. A computer-readable storage medium having a computer program stored thereon, characterized in that: When the program is executed by a processor, the steps of the method according to any one of claims 1 to 7 are implemented.

10. A computer program product comprising a computer program, characterized in that When the computer program is executed by a processor, the steps of the method according to any one of claims 1 to 7 are implemented.

Images