A method for extracting elevation control points from spaceborne photon counting lidar suitable for islands

Abstract

The present invention provides a method for extracting elevation control points of satellite-borne photon counting laser radars suitable for islands, and belongs to the field of satellite-borne photon counting laser radar data processing technology. The method of the present invention can extract a certain number of elevation control points in island (reef) areas with relatively small land areas and relatively complex terrains, and existing experimental results show that the overall mean error of the elevation control points extracted by this method meets the required accuracy threshold for elevation control of the land part of large-scale nautical charts of 1:10000 and above, and can provide strong data support for subsequent global nautical chart revision and quality control.

Description

Technical Field

[0001] The invention belongs to the technical field of satellite-borne photon counting laser radar data processing and relates to a satellite-borne photon counting laser radar elevation control point extraction method suitable for islands. Background Art

[0002] Island (reef) topographic survey is an important part of marine surveying and mapping, and elevation control points are important basic data for topographic surveying and precision verification. Due to realistic factors such as being far away from the mainland and scattered in distribution, obtaining elevation control points through traditional field methods is inefficient, costly, and risky, and it is difficult to ensure the timeliness of the data. In addition, for some islands located in sensitive areas, it is almost impossible to obtain measured control information. With its characteristics of no geographical restrictions, all-day, and long-distance earth observation, satellite-borne laser radar (Light Detection and Ranging, LiDAR) can quickly and directly obtain three-dimensional spatial information of land objects around the world, and can assist in satellite image regional network adjustment and improve image positioning accuracy without ground control points. It has become one of the reliable sources of elevation control data and can provide a solid data foundation for marine battlefield environment assurance applications such as chart revision and quality control.

[0003] At present, many earth observation satellites are equipped with LiDAR altimetry systems. According to their measurement systems, they can be divided into two categories: full-waveform LiDAR and photon counting LiDAR. Among them, satellite-borne full-waveform LiDAR technology is relatively mature and has many platforms, and has been widely used in the extraction of elevation control points. However, the first-generation full-waveform LiDAR satellite ICESat (Ice, Cloud and land Elevation Satellite) failed in 2009, and the data is less current; and the footprint spacing of the laser altimetry data of the Gaofen-7 satellite independently launched by my country is large (about 2.4km), and the distribution is relatively sparse after quality control, which limits its application in the acquisition of island elevation control points. As the only earth observation satellite equipped with photon counting LiDAR, ICESat-2 (Ice, Cloud and land Elevation Satellite-2) was launched in September 2018, providing a new means for obtaining more detailed three-dimensional surface information. Due to its measurement system, ICESat-2 transmits and receives weak signals, which are inevitably interfered by atmospheric scattering, solar background, instrument noise and other factors during the data acquisition process. Therefore, the data needs to be screened before it can be used as island elevation control information. Summary of the invention

[0004] In order to achieve rapid acquisition of island elevation control information and provide data support for global nautical chart revision and quality control, the present invention proposes a method for extracting elevation control points of space-borne photon counting lidar suitable for islands.

[0005] The technical solution adopted by the present invention to achieve the above-mentioned purpose is:

[0006] A method for extracting elevation control points of a spaceborne photon counting laser radar applicable to an island mainly comprises the following steps:

[0007] a. Read various attribute information of elevation points in ICESat-2 ATL08 data, and extract sea area mask based on the land cover type parameter (segment_landcover) in the read data and the land area layer in the publicly released electronic nautical chart, quickly remove the sea area points in the data, and construct a land area elevation point dataset;

[0008] b. Based on the median surface elevation (h_te_median), interpolated surface elevation (h_te_interp) and reference surface elevation (dem_h) in the read data, calculate the elevation difference and remove gross error points whose elevation difference does not meet the conditions;

[0009] c. Comprehensively analyze the influence of atmospheric environment, signal-to-noise ratio, fitting accuracy and other factors on elevation values, and construct multi-parameter constraint preliminary screening data based on the multiple scattering warning flag (msw_flag), cloud cover flag (cloud_flag_atm), signal-to-noise ratio parameter (snr), number of surface photons in statistical unit (n_te_photons) and best fitting elevation of surface photons in 20m sub-statistical unit (h_te_best_fit_20m) in the read data;

[0010] d. Make full use of the slope parameter (terrain_slope) in the data, combine the elevation error estimate (sigma_atlas_land), elevation standard deviation (h_te_std), surface photon elevation skewness (h_te_skew) and statistical unit surface photon number and other parameter information, and design an adaptive terrain threshold based on the accuracy index of elevation control points for large-scale nautical chart revision to extract elevation control points in a refined manner.

[0011] In step a, the sea area points are removed in the following manner:

[0012] (1) Based on the land cover type parameter (segment_landcover) in the read data, as shown in formula (1), remove the elevation points with parameter values ​​of 80 or 200;

[0013]

[0014] (2) Extract the land area (LNDARE) layer from the publicly available international standard electronic nautical chart, perform an overlay analysis on all elevation points, and remove elevation points that do not intersect with the layer.

[0015] In step b, the judgment condition of the gross error point is:

[0016]

[0017] In formula (2), σ ref The absolute elevation accuracy of the reference elevation is 3σ in the island (reef) area without loss of generality. ref The value can be taken as 25m. Based on the above judgment, only the elevation points whose elevation difference meets the above conditions are retained.

[0018] In step c, the method for implementing the preliminary screening of elevation control points with multi-parameter constraints is:

[0019] (1) Based on the multiple scattering warning flag (msw_flag) in the read data, determine whether it is affected by scattering, and remove the elevation points affected by complex scattering; for the elevation points only affected by clouds, further use the cloud cover flag (cloud_flag_atm) to retain the elevation points with cloud cover less than the threshold. The specific judgment conditions for multiple scattering and cloud cover flags are:

[0020]

[0021] Based on the above judgment, only the elevation points that meet the above multiple scattering and cloud cover conditions are retained.

[0022] (2) Based on the signal-to-noise ratio parameter (SNR) in the read data, retain the elevation points whose SNR parameter is greater than the threshold. The specific judgment conditions of the SNR parameter are:

[0023]

[0024] Based on the above judgment, only the elevation points that meet the above signal-to-noise ratio conditions are retained.

[0025] (3) Based on the number of surface photons in the statistical unit (n_te_photons) in the read data, remove the elevation points where the number of surface photons is less than the threshold. The specific judgment conditions for the number of surface photons in the statistical unit are:

[0026] n_te_photons ≥ 50 (5)

[0027] Based on the above judgment, only the elevation points whose surface photon numbers within the statistical unit meet the above conditions are retained.

[0028] (4) Based on the best fitting elevation of the surface photons in the 20m sub-statistical unit in the read data (h_te_best_fit_20m), remove the elevation points with invalid values. The specific judgment conditions are:

[0029] h_te_best_fit_20m≠inf (6)

[0030] Based on the above judgment, only the elevation points whose surface photon numbers within the statistical unit meet the above conditions are retained.

[0031] In step d, the method for extracting elevation control points based on the adaptive terrain threshold is:

[0032] The discrete degree of the distribution of surface photons along the track and in the elevation direction is analyzed. The slope θ of the obtained elevation points after executing step c is preliminarily screened based on the slope parameter (terrain_slope) in the data. The elevation points under three terrain types, namely flat land (θ<2°), hills (2°≤θ<6°), and mountains (6°≤θ<25°), are obtained according to their slope conditions. Combined with the elevation error estimate (sigma_atlas_land), elevation standard deviation (h_te_std), elevation skewness (h_te_skew) of surface photons in the statistical unit in the ICESat-2ATL08 data fitting process, and the number of surface photons in the statistical unit, the refined extraction criteria shown in the following formula are adaptively set according to the required accuracy of the elevation control points:

[0033]

[0034]

[0035] In the formula, σ est represents the estimated standard deviation based on the slope parameter, and the calculation form is shown in formula (8), N te , ΔH are the number of surface photons in the statistical unit and the terrain height difference estimated based on the slope parameter; T is the required elevation control point accuracy threshold. In the compilation of nautical charts, the land terrain elements are mainly transferred from the topographic map, and the revision and quality control of nautical charts with a scale of 1:10000 and above mainly use 1:50000 mapping data as a reference. Based on this, without loss of generality, the present invention takes the values ​​of T in flat land, hilly and mountainous terrain types as 0.8m, 1.0m and 1.2m respectively.

[0036] After the above extraction process, the island elevation control points are finally obtained.

[0037] Beneficial effects of the invention: The satellite-borne photon counting lidar elevation control point extraction method proposed in the invention can realize the extraction of a certain number of elevation control points in island (reef) areas with relatively small land area and relatively complex terrain, and existing experimental results show that the overall mean error of the elevation control points extracted by this method meets the accuracy threshold required for land elevation control of large-scale nautical charts of 1:10000 and above, and can provide strong data support for subsequent global nautical chart revision and quality control. BRIEF DESCRIPTION OF THE DRAWINGS

[0038] Figure 1 This is a main flow chart of the method for extracting elevation control points of a space-borne photon counting laser radar applicable to an island, which is proposed in the present invention. DETAILED DESCRIPTION

[0039] The present invention is further described below in conjunction with the accompanying drawings and specific examples. This example is implemented based on the technical solution of the present invention, and a specific implementation plan and specific operation process are given, but the protection scope of the present invention is not limited to the following examples. Based on ICESat-2 ATL08 data, the present invention is used to extract island elevation control points, including the following steps: Figure 1 As shown:

[0040] Step a, read various attribute information of elevation points in ICESat-2 ATL08 data, and extract the sea area mask based on the surface cover type parameter (segment_landcover) in the read data, combined with the land area layer in the publicly released electronic nautical chart, quickly remove the sea area points in the data, and construct a land area elevation point dataset. The sea area point removal method used in the present invention is:

[0041] (1) Based on the land cover type parameter (segment_landcover) in the read data, as shown in formula (1), remove the elevation points with parameter values ​​of 80 or 200;

[0042]

[0043] (2) Extract the land area (LNDARE) layer from the publicly available international standard electronic nautical chart, perform an overlay analysis on all elevation points, and remove elevation points that do not intersect with the layer.

[0044] Step b, based on the median surface elevation (h_te_median), the interpolated surface elevation value (h_te_interp) and the reference surface elevation value (dem_h) in the read data, calculate the elevation difference and remove the gross error points whose elevation difference is greater than the threshold. The gross error point judgment condition used in the present invention is:

[0045]

[0046] In formula (10), σ ref The absolute elevation accuracy of the reference elevation is 3σ in the island (reef) area without loss of generality. ref The value can be taken as 25m. Based on the above judgment, only the elevation points whose elevation difference meets the above conditions are retained.

[0047] Step c, comprehensively analyzing the influence of atmospheric environment, signal-to-noise ratio, fitting accuracy and other factors on the elevation value, based on the multiple scattering warning flag (msw_flag), cloud cover flag (cloud_flag_atm), signal-to-noise ratio parameter (snr), number of surface photons in statistical unit (n_te_photons) and best fitting elevation of surface photons in 20m sub-statistical unit (h_te_best_fit_20m) in the read data, construct multi-parameter constraint preliminary screening data. The multi-parameter constraint screening criteria constructed by the present invention are specifically expressed as:

[0048] (1) Based on the multiple scattering warning flag (msw_flag) in the read data, determine whether it is affected by scattering, and remove the elevation points affected by complex scattering; for the elevation points only affected by clouds, further use the cloud cover flag (cloud_flag_atm) to retain the elevation points with cloud cover less than the threshold. The specific judgment conditions for multiple scattering and cloud cover flags are:

[0049]

[0050] Based on the above judgment, only the elevation points that meet the above multiple scattering and cloud cover conditions are retained.

[0051] (2) Based on the signal-to-noise ratio parameter (SNR) in the read data, retain the elevation points whose SNR parameter is greater than the threshold. The specific judgment conditions of the SNR parameter are:

[0052]

[0053] Based on the above judgment, only the elevation points that meet the above signal-to-noise ratio conditions are retained.

[0054] (3) Based on the number of surface photons in the statistical unit (n_te_photons) in the read data, remove the elevation points where the number of surface photons is less than the threshold. The specific judgment conditions for the number of surface photons in the statistical unit are:

[0055] n_te_photons ≥ 50 (13)

[0056] Based on the above judgment, only the elevation points whose surface photon numbers within the statistical unit meet the above conditions are retained.

[0057] (4) Based on the best fitting elevation of the surface photons in the 20m sub-statistical unit in the read data (h_te_best_fit_20m), remove the elevation points with invalid values ​​in the value. The specific judgment conditions are:

[0058] h_te_best_fit_20m≠inf (14)

[0059] Based on the above judgment, only the elevation points whose surface photon numbers within the statistical unit meet the above conditions are retained.

[0060] Step d, make full use of the slope parameter (terrain_slope) in the data, combine the elevation error estimate (sigma_atlas_land), elevation standard deviation (h_te_std), surface photon elevation skewness (h_te_skew) and statistical unit surface photon number and other parameter information, and design an adaptive terrain threshold based on the accuracy index of large-scale nautical chart revision elevation control points to extract elevation control points in a refined manner:

[0061] The discrete degree of the distribution of surface photons in the along-track and elevation directions is comprehensively analyzed. The slope θ of the elevation points obtained after executing step c is preliminarily screened based on the slope parameter (terrain_slope) in the data. The elevation points under three terrain types, namely flat land (θ<2°), hills (2°≤θ<6°), and mountains (6°≤θ<25°), are obtained according to their slope conditions. Combined with the elevation error estimate (sigma_atlas_land), elevation standard deviation (h_te_std), surface photon elevation skewness (h_te_skew), and number of surface photons in the statistical unit in the ICESat-2ATL08 data fitting process, the refined extraction criteria shown in the following formula are adaptively set according to the required accuracy of the elevation control points:

[0062]

[0063]

[0064] In the formula, σ est represents the estimated standard deviation based on the slope parameter, and the calculation form is shown in formula (8), N te , ΔH are the number of surface photons in the statistical unit and the terrain height difference estimated based on the slope parameter; T is the required elevation control point accuracy threshold. In the compilation of nautical charts, the land terrain elements are mainly transferred from the topographic map, and the revision and quality control of nautical charts with a scale of 1:10000 and above mainly use 1:50000 mapping data as a reference. Based on this, without loss of generality, the present invention takes the values ​​of T in flat land, hilly and mountainous terrain types as 0.8m, 1.0m and 1.2m respectively.

[0065] After the above extraction process, the island elevation control points are finally obtained.

[0066] The following is an example of the present invention using some of the outlying islands of Puerto Rico and the U.S. Virgin Islands on the northern edge of the Caribbean Sea. In the example area, the outlying islands of Puerto Rico include Vieques Island, Culebra Island and some of its surrounding islands (reefs), with a total land area of ​​about 163.92 km 2 The terrain is mainly hilly with steep coasts. The Virgin Islands include St. Croix, St. Thomas, St. John and some of their surrounding islands (reefs), with a total land area of ​​about 350.58 km 2 , St. Croix Island is high in the north and low in the south, and its area is relatively large, while the other two islands have large undulating terrain, many hills and mountains, and relatively small areas. The embodiment area is a typical archipelago area with complex terrain features and rich surface cover types, which is highly representative. This embodiment first collects all ATL08 data from October 2018 to August 2022 in the embodiment area for elevation control point extraction. Furthermore, in order to verify the accuracy of the elevation control points extracted by the method of the present invention, this embodiment also collects the 1m resolution DEM data products publicly released by the United States Geological Survey as high-precision verification data. The valid land points after the sea area points are removed are regarded as original data, and the number of elevation points, mean absolute error (MAE), root mean square error (RMSE) and data retention rate obtained under different screening conditions are shown in Table 1.

[0067] Table 1. Verification results of elevation control point accuracy

[0068]

[0069]

[0070] After removing the sea area points, there are a total of 15,388 land elevation points in the study area. The MAE and RMSE of the original elevation points are 3.86m and 7.89m respectively, which cannot be used directly for elevation control. After extraction by the method of the present invention, a total of 1,990 elevation control points are retained, the data retention rate is 12.93%, the MAE and RMSE are 0.36m and 0.69m respectively, and the accuracy has been greatly improved. It is not difficult to see that the elevation control points extracted by the method of the present invention have a high elevation accuracy as a whole. In terms of the spatial distribution of elevation control points, except for St. John Island, which has a sparse distribution of elevation control points due to complex terrain, dense vegetation and poor data integrity, the other four islands have certain elevation control points after extraction by the method of the present invention, which can meet the needs of island elevation control. This shows that even in island (reef) areas with relatively small land areas and relatively complex terrain, the method of the present invention can extract a certain number of elevation control points from ICESat-2 ATL08 data, which is suitable for the extraction of island elevation control points.

[0071] The extraction accuracy of elevation control points under different terrain types was further analyzed, and the number of elevation points, MAE, RMSE and data retention rate obtained under different terrain types and different screening conditions were statistically analyzed. The results are shown in Table 2.

[0072] Table 2 Verification results of elevation control point accuracy under different terrain types

[0073]

[0074] Analyzing the above statistical results, since most of the islands in the embodiment area are volcanic islands, with many mountains and hills, steep coasts and complex terrain, after slope classification, the land elevation points under the types of flat land, hills and mountains are 5618, 4301 and 5138 respectively, which are roughly the same in number. The elevation residuals of the original elevation points of the three terrain types are relatively scattered, and a considerable number of elevation point residuals are large. The MAE and RMSE are 2.65m / 7.23m, 3.92m / 7.65m, and 4.93m / 8.29m respectively, showing a trend of increasing step by step. The accuracy compliance rates are 61.65%, 43.87% and 32.21% respectively, showing a trend of decreasing step by step, indicating that the terrain undulation does have a significant impact on the elevation accuracy of the ATL08 data.

[0075] After removing gross error points based on the reference elevation values, the elevation points of the three terrain types decreased by 18.58%, 30.64% and 35.36%, respectively, and the accuracy compliance rates increased by 12.27%, 16.07% and 13.05%, respectively. This shows that there is indeed a certain amount of gross error data in the ICESat-2 ATL08 data. By simply comparing with the reference elevation values, the number of elevation points involved in multi-parameter constraint screening is reduced, which can improve the computational efficiency of subsequent elevation point extraction.

[0076] After multi-parameter constraints, the MAE of elevation points under the three terrain types increased by 57.14%, 56.29% and 52.79%, respectively, and the RMSE increased by 49.76%, 47.25% and 39.72%, respectively. This shows that the multi-parameter constraints proposed in the present invention can comprehensively consider the influence of atmospheric environment, signal-to-noise ratio and data integrity on the accuracy of elevation points, and preliminarily obtain data with relatively reliable accuracy. Moreover, in the embodiment area, the constraints are not limited by the terrain environment, and play a greater role in the three terrain types.

[0077] The amount of data removed by the adaptive terrain threshold is small, indicating that the elevation points obtained by the preliminary screening based on multi-parameter constraints have good elevation accuracy, but the MAE and RMSE of the elevation control points after fine extraction are still improved by 28.20% / 42.45%, 30.30% / 45.14% and 42.72 / 58.99% respectively. This shows that there are still some elevation gross error points that do not meet the accuracy threshold in the preliminary extracted elevation points. The present invention refines the extraction of elevation points by setting the adaptive terrain threshold, so that the overall RMSE of the elevation control points extracted under different types can meet the accuracy threshold requirements, thereby improving the robustness of the method.

[0078] In summary, the method of the present invention can extract a certain number of elevation control points from ICESat-2 ATL08 data even in island (reef) areas with relatively small land areas and relatively complex terrain. The overall mean square error of the elevation control points extracted under different terrain types can meet the elevation control accuracy threshold requirements of land elements of nautical maps above 1:10000. It is suitable for the extraction of island elevation control points and can provide strong data support for subsequent global nautical chart revision and quality control.

[0079] The present invention is described by way of embodiments, and those skilled in the art will appreciate that various changes or equivalent substitutions may be made to these features and embodiments without departing from the spirit and scope of the present invention. In addition, under the teachings of the present invention, these features and embodiments may be modified to adapt to specific circumstances and materials without departing from the spirit and scope of the present invention. Therefore, the present invention is not limited to the specific embodiments disclosed herein, and all embodiments falling within the scope of the claims of this application belong to the protection scope of the present invention.

Claims

1. A method for extracting elevation control points of a spaceborne photon counting laser radar suitable for an island, characterized in that: The method comprises the following steps: a. Read various attribute information of elevation points in ICESat-2 ATL08 data, and extract sea area mask based on the land cover type parameter segment_landcover in the read data and the land area layer in the publicly released electronic nautical chart, quickly remove the sea area points in the data, and construct a land area elevation point dataset; b. Based on the median surface elevation h_te_median, the interpolated surface elevation h_te_interp and the reference surface elevation dem_h in the read data, the elevation difference is calculated, and the gross error points whose elevation difference does not meet the conditions are removed; c. Comprehensively analyze the influence of atmospheric environment, signal-to-noise ratio and fitting accuracy on elevation value, and construct preliminary screening data with multi-parameter constraints based on the multiple scattering warning flag msw_flag, cloud cover flag cloud_flag_atm, signal-to-noise ratio parameter snr, number of surface photons in statistical unit n_te_photons and best fitting elevation of surface photons in 20m sub-statistical unit h_te_best_fit_20m in the read data; d. Make full use of the slope parameter terrain_slope in the data, combine the elevation error estimate sigma_atlas_land, the elevation standard deviation h_te_std, the surface photon elevation skewness h_te_skew and the surface photon number information of the statistical unit, and design an adaptive terrain threshold based on the accuracy index of the elevation control points measured by large-scale nautical chart revision to extract the elevation control points in a refined manner.

2. The method for extracting elevation control points of a space-borne photon counting laser radar suitable for an island according to claim 1, characterized in that: In step a, the sea area points are removed in the following manner: (1) Based on the land cover type parameter segment_landcover in the read data, as shown in formula (1), remove the elevation points with parameter values ​​of 80 or 200; (2) Extract the land area LNDARE layer from the publicly available international standard electronic nautical chart, perform overlay analysis on all elevation points, and remove elevation points that do not intersect with the layer.

3. The method for extracting elevation control points of a space-borne photon counting laser radar applicable to an island according to claim 1 or 2, characterized in that: In step b, the judgment condition of the gross error point is: In formula (2), σ ref Indicates the absolute elevation accuracy of the reference elevation, 3σ in the island / reef area ref The value is 25m; Based on the above judgment, only the elevation points whose elevation difference meets the above conditions are retained.

4. The method for extracting elevation control points of a space-borne photon counting laser radar applicable to an island according to claim 1 or 2, characterized in that: In step c, the method for implementing the preliminary screening of elevation control points with multi-parameter constraints is: (1) Based on the multiple scattering warning flag msw_flag in the read data, determine whether it is affected by scattering, and remove the elevation points affected by complex scattering; for the elevation points only affected by clouds, further use the cloud cover mark cloud_flag_atm to retain the elevation points with cloud cover less than the threshold; the specific judgment conditions for multiple scattering and cloud cover marks are: Only the elevation points that satisfy the above multiple scattering and cloud cover conditions are retained; (2) Based on the signal-to-noise ratio parameter snr in the read data, retain the elevation points where the signal-to-noise ratio parameter is greater than the threshold; the specific judgment condition of the signal-to-noise ratio parameter is: Only the elevation points that meet the above signal-to-noise ratio conditions are retained; (3) Based on the surface photon number n_te_photons of the statistical unit in the read data, remove the elevation points where the surface photon number is less than the threshold; the specific judgment condition of the surface photon number in the statistical unit is: n_te_photons ≥ 50 (5) Only the elevation points whose surface photon counts within the statistical unit meet the above conditions are retained; (4) Based on the best fitting elevation of the surface photons in the 20m sub-statistical unit in the read data, remove the elevation points with invalid values. The specific judgment conditions are: h_te_best_fit_20m≠inf (6) Only the elevation points whose surface photon counts within the statistical unit meet the above conditions are retained.

5. The method for extracting elevation control points of a space-borne photon counting laser radar applicable to an island according to claim 3, characterized in that: In step c, the method for implementing the preliminary screening of elevation control points with multi-parameter constraints is: (1) Based on the multiple scattering warning flag msw_flag in the read data, determine whether it is affected by scattering, and remove the elevation points affected by complex scattering; for the elevation points only affected by clouds, further use the cloud cover mark cloud_flag_atm to retain the elevation points with cloud cover less than the threshold; the specific judgment conditions for multiple scattering and cloud cover marks are: Only the elevation points that satisfy the above multiple scattering and cloud cover conditions are retained; (2) Based on the signal-to-noise ratio parameter snr in the read data, retain the elevation points where the signal-to-noise ratio parameter is greater than the threshold; the specific judgment condition of the signal-to-noise ratio parameter is: Only the elevation points that meet the above signal-to-noise ratio conditions are retained; (3) Based on the surface photon number n_te_photons of the statistical unit in the read data, remove the elevation points where the surface photon number is less than the threshold; the specific judgment condition of the surface photon number in the statistical unit is: n_te_photons ≥ 50 (5) Only the elevation points whose surface photon counts within the statistical unit meet the above conditions are retained; (4) Based on the best fitting elevation of the surface photons in the 20m sub-statistical unit in the read data, remove the elevation points with invalid values. The specific judgment conditions are: h_te_best_fit_20m≠inf (6) Only the elevation points whose surface photon counts within the statistical unit meet the above conditions are retained.

6. A method for extracting elevation control points of a space-borne photon counting laser radar applicable to an island according to claim 1, 2 or 5, characterized in that: In step d, the method for extracting elevation control points based on the adaptive terrain threshold is: The discrete degree of the distribution of surface photons along the track and in the elevation direction is analyzed. The slope θ of the elevation points obtained after executing step c is preliminarily screened based on the slope parameter terrain_slope in the data. The elevation points of three terrain types, namely flat land, hills and mountains, are obtained according to their slope conditions. Combined with the elevation error estimate sigma_atlas_land, elevation standard deviation h_te_std, surface photon elevation skew h_te_skew and the number of surface photons in the statistical unit in the ICESat-2 ATL08 data fitting process, the refined extraction criteria shown in the following formula are adaptively set according to the required accuracy of the elevation control point: In the formula, σ est represents the estimated standard deviation based on the slope parameter, and the calculation form is shown in formula (8), N te , ΔH are the number of surface photons in the statistical unit and the terrain height difference estimated based on the slope parameters; T is the required accuracy threshold of the elevation control point.

7. The method for extracting elevation control points of a space-borne photon counting laser radar applicable to an island according to claim 6, characterized in that: The slope of flat land is θ<2°, the slope of hills is 2°≤θ<6°, and the slope of mountains is 6°≤θ<25°.

8. The method for extracting elevation control points of a space-borne photon counting laser radar applicable to an island according to claim 6, characterized in that: The values ​​of T in flat, hilly and mountainous terrain types are taken as 0.8m, 1.0m and 1.2m respectively.

Images