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

Abstract

The application discloses a high-precision DEM fusion construction method based on LiDAR photon statistics, comprising the following steps: collecting ICESat-2 L3A and GEDI L2A satellite-borne LiDAR data in a region, and performing extraction and strip filtering processing on the photon data; integrating multiple open-source DEM products, and performing geodetic datum calibration, spatial offset correction and pixel-level alignment; constructing a photon space index by using a CK-Tree algorithm; searching for up to 50 nearest neighbor photon points for each DEM pixel position, and calculating the error of each DEM relative to the photon points; sorting the DEM pixel accuracy, selecting the DEM data with the smallest error to assign an elevation value to each pixel, thereby generating an optimal DEM; and adopting a conservative smoothing filtering technique on the final optimal DEM to eliminate spatial discontinuity. The application can automatically integrate optimal DEM data and provide a fast and efficient solution.

Description

TECHNICAL FIELD

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

[0002] DEM is an important data in the field of geographic information, which stores the terrain elevation information in grid form and serves the fields of hydrology, geomorphology, urban planning and environmental monitoring. Although there are various open-source DEM products globally, such as ASTER GDEM, SRTM V3, Nasa DEM, etc., these products have significant differences in precision and quality in different regions, which is mainly due to the diversity of sensor technology and data processing methods they rely on. In order to overcome the limitations of single DEM product, researchers usually use fusion technology to integrate the advantages of different DEMs, improve the accuracy and reliability of terrain data.

[0003] However, existing DEM fusion methods, such as simple average, weighted average, variational model, etc., often lack sufficient flexibility and interpretability, and fail to fully utilize the latest remote sensing technology. In addition, these methods do not consider the dynamic updating of DEM data, limiting their application in rapid change terrain monitoring.

[0004] The emergence of spaceborne LiDAR technology, especially ICESat-2 L3A and GEDI L2A data, provides a new way for high-precision terrain data acquisition. These sensors can penetrate the vegetation layer and obtain accurate three-dimensional structure information below the ground surface. Although spaceborne LiDAR data is discrete in spatial distribution, its high-precision measurement provides new possibilities for DEM fusion. SUMMARY

[0005] The purpose of the present application is to provide a high-precision DEM fusion construction method based on LiDAR photon statistics, aiming to solve the technical problem of how to realize high-precision fusion generation of DEM through photon-based statistical technology driven by spaceborne LiDAR data.

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

[0007] Collecting ICESat-2 L3A and GEDI L2A spaceborne LiDAR data within the region, and extracting and strip filtering the photon data;

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

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

[0010] For each DEM pixel position, search for up to 50 nearest neighbor photon points, and calculate the error of each DEM relative to these photon points;

[0011] According to the error estimate value, the accuracy of the DEM pixel is sorted, and the DEM data with the smallest error is selected to assign an elevation value to each pixel, thereby generating an optimal DEM;

[0012] The final optimal DEM is subjected to smoothing filtering to eliminate spatial discontinuity.

[0013] A computer device comprises a memory, a processor and a computer program stored on the memory and executable on the processor, and the processor executes the program to realize the steps of the above method.

[0014] A computer readable storage medium has a computer program stored thereon, and the program is executed by a processor to realize the steps of the above method.

[0015] A computer program product comprises a computer program, and the computer program is executed by a processor to realize the steps of the above method.

[0016] Compared with the prior art, the beneficial effects of the present application are:

[0017] The present application is mainly applicable to the automatic fusion scene of open source DEM, and can automatically, quickly and efficiently generate DEM products with the highest vertical precision in the region according to the open source DEM, thereby providing a new method for spatial data fusion.

[0018] The present application not only improves the precision and interpretability of DEM fusion, but also avoids arbitrary assumptions about DEM degradation factors by using actual measurement data as a prior condition for fusion. Another significant advantage of the present application is its scalability, which can adapt to the constantly updated LiDAR and open source DEM data, and provide an effective solution for monitoring dynamic topographic changes; the application potential of the method will be further expanded with the increase of open source DEM and spaceborne LiDAR data, and the method is particularly suitable for geoscience research and application that requires high-precision topographic information. BRIEF DESCRIPTION OF DRAWINGS

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

[0020] Figure 1 It is a high-precision DEM fusion construction method flowchart based on LiDAR photon statistics.

[0021] Figure 2is the nearest neighbor photon search, photon weight assignment, error ranking, and pixel-by-pixel optimal terrain determination.

[0022] Figure 3 is the precision of using ICESat-2 L3A metric fused DEM and other open source products.

[0023] Figure 4 is the precision of using GEDI L2A metric fused DEM and other open source products. DETAILED DESCRIPTION

[0024] The present application proposes a high-precision digital elevation model (DEM) fusion construction method based on LiDAR photon statistics, and the specific steps are as follows:

[0025] Step 1. LiDAR data collection: collect ICESat-2 L3A and GEDI L2A satellite LiDAR data in the region, use all data in the past 5 years for stable terrain areas, and only use data for the current year for unstable terrain areas;

[0026] Step 2. ICESat-2 L3A photon information processing: extract key parameters from ICESat-2 L3A data, including latitude, longitude, best-fit height, height uncertainty, number of photons, interpolated height, median height, atmospheric cloud flag, and night flag; perform high-confidence photon filtering, conditions include night flag is 1, number of photons is not less than 50, interpolated and median height difference is less than 2 meters, atmospheric cloud flag is 0 or 1, height uncertainty is less than 10 meters;

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

[0028] Step 4. GEDI L2A photon information processing: extract key parameters from GEDI L2A data, including timestamp, lowest mode latitude and longitude, degradation flag, height bias flag, quality flag, lowest mode height, and reception evaluation flag; perform high-confidence photon filtering, conditions include nighttime, quality flag is 1, degradation flag is not 0, height bias flag is 0, and reception evaluation flag is 0;

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

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

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

[0032]

[0033] where, and represent the geoid height based on EGM96 and EGM08, respectively, and represent the geodetic separation between the EGM96 and EGM08 quasi-geoid and the WGS84 ellipsoid;

[0034] Step 8. Spatial offset correction and pixel-level alignment: Create a fused LiDAR product combining ICESat-2 L3A and GEDI L2A , sample the height of the photon-DEM intersection part in using bilinear interpolation and calculate MAE, use the position with the latest smallest MAE as the best offset position; move the DEM in a 30 m x 30 m grid, with an increment of 0.5 m, repeat sampling and error calculation to determine the offset; the DEM after offset uses bilinear resampling to align all pixels;

[0035] Step 9. Photon index construction: Use the CK-Tree algorithm to establish a spatial index for ;

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

[0037]

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

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

[0040]

[0041] where, is the weight of the photon , is the distance between the photon and the pixel position, is the photon source function, defined as:

[0042]

[0043] Step 12. Calculate the error of each DEM in the pixel position :

[0044]

[0045] where, is the mean absolute error of the photon compared to .

[0046] Step 13. Error ranking and per-pixel optimal terrain determination: for the DEM of each pixel, select the DEM with the lowest MAE as the elevation value of the current pixel;

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

[0048] The technical route and operation steps of the present application are more clearly described below in detail according to the accompanying drawings and examples.

[0049] Examples

[0050] The Loess Plateau covers a variety of landforms, which can comprehensively demonstrate the results and accuracy of DEM fusion. Therefore, the experimental area of this example is selected in the Loess Plateau. The process of the high-precision DEM fusion construction method based on LiDAR photon statistics in this example is shown in Figure 1 , which includes the following steps:

[0051] First, collect LiDAR data: collect ICESat-2 L3A and GEDI L2A spaceborne LiDAR data in the region. A total of 1897 strips of ICESat-2 L3A data from 2018 to 2022 and 1883 strips of GEDI L2A data from 2019 to 2023 are collected in this step.

[0052] Second step, ICESat-2 L3A photon information processing: Key parameters are extracted from ICESat-2 L3A data, including latitude, longitude, best-fit height, height uncertainty, photon count, interpolated height, median height, atmospheric cloud flag, and night flag. High-confidence photon filtering is performed with the conditions including night flag = 1, photon count >= 50, interpolated height - median height < 2 m, atmospheric cloud flag = 0 or 1, and height uncertainty < 10 m.

[0053] Third step, ICESat-2 L3A strip filtering: MAE of all photons in a strip compared to ETOPO2022 error is calculated, and data from the entire strip with error exceeding 95th percentile is excluded. After this step, 2454380 reliable ICESat-2 L3A photons are obtained.

[0054] Fourth step, GEDI L2A photon information processing: Key parameters are extracted from GEDI L2A data, including timestamp, lowest mode latitude and longitude, degradation flag, height bias flag, quality flag, lowest mode height, and reception evaluation flag. High-confidence photon filtering is performed with the conditions including nighttime, quality flag = 1, degradation flag!= 0, height bias flag = 0, and reception evaluation flag = 0.

[0055] Fifth step, GEDI L2A strip filtering: MAE of all photons in a strip compared to ETOPO2022 error is calculated, and data from the entire strip with error exceeding 90th percentile is excluded. After this step, 10187997 reliable GEDI L2A photons are obtained.

[0056] Sixth step, DEM data collection: All open-source DEM products within the region are obtained, including GDEM V3, AW3D30, COP30, FABDEM, NASADEM, SRTM V3, and TAN30.

[0057] Seventh step, geodetic datum coordination: Geodetic datum heights of GDEM V3, AW3D30, NASADEM, and SRTM V3 are converted from EGM96 to EGM08 to be consistent with COP30, FABDEM, and TAN30. The conversion formula is:

[0058]

[0059] where, and denote geodetic height based on EGM96 and EGM08, respectively, and denote geodetic separation between EGM96 and EGM08 quasi-geoid and WGS84 ellipsoid, respectively.

[0060] Step 8, Spatial offset correction and pixel-level alignment: Create a fused LiDAR product combining ICESat-2 L3A and GEDI L2A , the elevation of the portion of the photon intersecting the DEM is sampled and the MAE is calculated. The DEM is moved in 30 m x 30 m grids with 0.5 m increments, and the sampling and error calculation are repeated to determine the offset. The offset DEM is aligned to all pixels using bilinear resampling.

[0061] Step 9, Build photon index: Utilize the CK-Tree algorithm to create a spatial index for all photons in the DEM;

[0062] Step 10, Nearest neighbor photon search: For a specific pixel location , search for no more than 50 nearest neighbors within a maximum 5 km radius, i.e.:

[0063]

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

[0065] Step 11, Photon weight assignment: As shown in (a), the size of the photon shows its distance compared to the center pixel, and then the weight is assigned according to the photon type and its proximity to the pixel Figure 2 :

[0066]

[0067] where, is the weight of the photon , 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] In the formula, is a photon Compared with The average absolute error.

[0072] Thirteenth step, error sorting and per-pixel optimal terrain judgment: as shown in (b) of the figure, the accuracy of the current pixel in the open source DEM is ranked from high to low: COPDEM, TAN30, FABDEM…; as shown in (c) of the figure, for the DEM of each pixel, the DEM with the lowest MAE is selected as the elevation value of the current pixel; Figure 2 Figure 2 Thirteenth step, error sorting and per-pixel optimal terrain judgment: as shown in (b) of the figure, the accuracy of the current pixel in the open source DEM is ranked from high to low: COPDEM, TAN30, FABDEM…; as shown in (c) of the figure, for the DEM of each pixel, the DEM with the lowest MAE is selected as the elevation value of the current pixel;

[0073] Fourteenth step, smoothing filter: a conservative smoothing filter is applied to the optimal DEM to eliminate potential spatial discontinuity. As shown in the figure, Figure 3 The accuracy of the fused DEM and other open source products is shown in the figure. Figure 4 The accuracy of the fused DEM and other open source products is shown in the figure.It can be seen that under any verification index, the accuracy of the fused product is much higher than that of the existing open source DEM, proving the effectiveness of the method.

[0074] The present application is mainly applicable to the automatic fusion scene of open source DEM, and can automatically, quickly and efficiently generate the DEM product with the highest vertical accuracy in the region according to the open source DEM, providing a new method for spatial data fusion.

[0075] The above description is only an embodiment of the present application, and does not limit the patent scope of the present application, and any equivalent structure or equivalent process transformation using the content of the specification and drawings, or direct or indirect application in other related technical fields, are also included in the patent protection scope of the present application.

Claims

1. A high-precision DEM fusion construction method based on LiDAR photon statistics, characterized by, The method comprises the following steps: ICESat-2 L3A and GEDI L2A spaceborne LiDAR data in the collection area are collected, and photon data is extracted and strip filtering is performed; a plurality of open source DEM products are integrated, and geodetic datum calibration, spatial offset correction and pixel level alignment are performed; a fusion LiDAR product S is created in combination with ICESat-2 L3A and GEDI L2A, bilinear interpolation is used to sample the height of the intersection part of the photon and DEM in S and calculate MAE; the DEM is moved in the grid, and the sampling and error calculation are repeated to determine the offset; the offset DEM uses bilinear resampling to align all pixels; CK-Tree algorithm is used to construct a photon space index based on ICESat-2 L3A and GEDI L2A data; For each DEM pixel position, a search for a maximum of 50 nearest neighbor photon points is performed, and the error of each DEM relative to the photon points is calculated; The accuracy of the DEM pixels is sorted according to the error estimates, and the DEM data with the smallest error is selected to assign an elevation value to each pixel, thereby generating an optimal DEM; The final optimal DEM is subjected to smoothing filtering to eliminate spatial discontinuity.

2. The LiDAR photon statistics based high-precision DEM fusion construction method according to claim 1, characterized in that, ICESat-2 L3A data is subjected to photon information extraction and strip filtering, specifically: Key parameters are extracted from ICESat-2 L3A data, including latitude and longitude, best-fit height, height uncertainty, photon number, interpolated height, median height, atmospheric cloud flag and night flag; high-confidence photon screening is performed; The MAE of the error of all photons in a strip compared with ETOPO2022 is compared, and data of the entire strip with an error exceeding 95 percentiles is excluded.

3. The LiDAR photon statistics based high-precision DEM fusion construction method according to claim 2, characterized in that, The high-confidence photon screening includes conditions such as night flag being 1, photon number being not less than 50, difference between interpolated height and median height being less than 2 meters, atmospheric cloud flag being 0 or 1, and height uncertainty being less than 10 meters.

4. The LiDAR photon statistics based high-precision DEM fusion construction method according to claim 1, characterized in that, GEDI L2A data is subjected to photon information extraction and strip filtering, specifically: Key parameters are extracted from GEDI L2A data, including timestamp, lowest mode latitude and longitude, degradation flag, height bias flag, quality flag, lowest mode height and reception evaluation flag; high-confidence photon screening is performed; The MAE of the error of all photons in a strip compared with ETOPO2022 is compared, and data of the entire strip with an error exceeding 90 percentiles is excluded.

5. The LiDAR photon statistics based high-precision DEM fusion construction method according to claim 4, characterized in that, The high-confidence photon screening includes conditions such as night time, quality flag being 1, degradation flag not being 0, height bias flag being 0, and reception evaluation flag being 0.

6. The LiDAR photon statistics based high-precision DEM fusion construction method according to claim 1, characterized in that, Geodetic datum calibration is performed, specifically: All open source DEM products in the region are obtained, including GDEM V3, AW3D30, COP30, FABDEM, NASADEM, SRTM V3 and TAN30; The geoid heights of GDEM V3, AW3D30, NASADEM and SRTM V3 were converted from EGM96 to EGM08 to be consistent with COP30, FABDEM and TAN30.

7. The LiDAR photon statistics based high-precision DEM fusion construction method according to claim 1, characterized in that, A search for up to 50 nearest neighbour photon points is performed for each DEM pixel location and the error of each DEM is calculated with respect to these photon points, specifically: For a DEM pixel location, search for up to 50 nearest neighbour photons within a maximum 5 km radius; The size of the photon shows its distance compared to the center pixel, assigning a weight according to the type of photon and its proximity to the pixel : ； In the formula, For photons The weight, It is the distance between the photon and the pixel position. For the number of photons, It is the photon source function, defined as: ； Calculate error for each DEM in pixel position : ； In the formula, is a photon Compared to the average absolute error.

8. A computer device comprising a memory, a processor, and a computer program stored on the memory and executable on the processor, characterized in that, The processor, when executing the program, implements the steps of the method of any one of claims 1-7.

9. A computer readable storage medium having stored thereon a computer program, characterized in that, The program, when executed by the processor, implements the steps of the method of any one of claims 1-7.

10. A computer program product comprising a computer program, characterized in that, The computer program, when executed by the processor, implements the steps of the method of any one of claims 1-7.

Images