A method for correction and fusion of optical stereo satellite DEM elevation

Abstract

The present invention discloses an optical stereo satellite DEM elevation correction and fusion method. It includes the following steps: step 1: elevation correction seed point selection; step 2: seed point DEM elevation extraction; step 3: elevation correction model comparison selection; step 4: adopting a determined optimal fitting model to perform DEM correction; step 5: calculating the elevation accuracy reliability index of the corrected DEM; step 6: adjacent DEM fusion and transfer correction; step 7: obtaining a large range of DEM data. The present invention overcomes the problem of inconsistent accuracy between different data sources and different DEMs, and the problem that simple arithmetic average fusion is prone to obvious inconsistency; it has the advantage of being able to improve the elevation accuracy of DEM produced by uncontrolled satellite photogrammetry through scarce external control data.

Description

Technical Field

[0001] The invention relates to the technical field of satellite remote sensing mapping, and more particularly to an optical stereo satellite DEM elevation correction and fusion method. Background Art

[0002] Digital elevation models (DEMs), as essential foundational geographic information data, are crucial for scientific research, environmental assessment, soil and water conservation, and engineering survey and design. Satellites with stereo mapping or elevation measurement capabilities, such as high-resolution optical stereo mapping satellites, interferometric radar satellites, and gravity satellites, are called mapping satellites. High-resolution mapping satellites using optical sensors are currently the primary data source for mapping applications. Optical remote sensing, due to its excellent visibility and robust technical feasibility, has rapidly developed. Optical stereo mapping satellites are also one of the primary data sources for the production of large-scale small and medium-scale DEMs worldwide.

[0003] In the conventional satellite photogrammetry process, most of the measurement errors affected by the geometric parameters of the satellite sensors can be better controlled by using ground control points, and then multiple satellite stereo DEMs are fused to obtain a large-scale high-precision DEM (the average of the overlapping areas of multiple DEMs is directly taken). However, in satellite photogrammetry without ground control points, the optical stereo satellite DEM obtained thereby will have large systematic errors, and there will be obvious discrepancies between the multiple DEMs obtained at different phases and by different sensors, which will hinder the high-precision fusion of DEMs. Invention patent CN104931022B can improve the stereo mapping accuracy of satellites such as Ziyuan-3 by using a stereo adjustment method that obtains plane and elevation control from different data sources and separates plane and elevation control. The essential idea is still the traditional control point regional network adjustment method, but this method is more dependent on sufficient control data and requires more complex control extraction and adjustment calculations. Invention patent CN106960174B uses lidar data to extract elevation control points of high-resolution stereo images, and uses a high-resolution image positioning model with elevation control point constraints to improve the positioning accuracy of high-resolution images. However, to achieve good results, this method relies on LiDAR and high-resolution imagery data from the same or similar platforms, and is less applicable to small amounts of external control data other than laser altimetry data. Invention patent CN114612806A proposes using accurate coordinate points measured by GNSS RTK to compensate for elevation anomalies using a quadratic polynomial surface model to improve the accuracy of consumer drone DEMs. This approach is primarily targeted at small-scale, high-precision drone product applications, and requires manual RTK acquisition of a large number of high-precision measurement points.

[0004] my country's three satellites, Ziyuan-3 (Ziyuan-3) 01, 02, and 03, have formed the country's first stereo mapping satellite constellation. On November 3, 2019, my country successfully launched the Gaofen-7 stereo mapping satellite. The satellite's operational full-waveform laser altimeter, equipped with a high-performance laser altimeter, provides crucial support for rapidly acquiring elevation control points and improving the accuracy of stereo mapping without ground control points for stereo images from the same platform. However, this approach is also limited by the limited variety of satellite data and the limited temporal richness of the data.

[0005] China is a vast country, with numerous uninhabited areas, including plateaus, deserts, and glaciers. For obtaining high-precision three-dimensional terrain data in these vast areas, satellite photogrammetry without ground control points remains the preferred choice for both feasibility and affordability. Exploring a versatile optical stereo satellite DEM elevation correction method, while simultaneously overcoming the barriers to high-precision fusion between different data sources and DEMs (i.e., significant discrepancies when fusion occurs due to inconsistent precision), is a pressing challenge. Summary of the Invention

[0006] The purpose of the present invention is to provide a method for correcting and merging DEM elevation using optical stereo satellites. The method is highly versatile and makes comprehensive use of various archived satellite data and a small amount of other reference data that can provide elevation control (such as scientific survey navigation trajectory points, laser altimetry data from other satellite platforms, etc.), thereby improving the accuracy and reliability of large-scale DEM acquisition using optical stereo satellites. The method also solves the problems of low elevation accuracy of DEMs acquired using uncontrolled optical stereo satellites and obstacles to high-precision fusion of DEMs from different data sources and different widths in the prior art.

[0007] In order to achieve the above object, the technical solution of the present invention is: an optical stereo satellite DEM elevation correction and fusion method, characterized in that it includes the following steps:

[0008] Step 1: Select seed points for elevation correction;

[0009] Determine the seed points used for elevation correction based on the extent of each DEM (i.e., optical stereo satellite DEM) produced by uncontrolled satellite photogrammetry and the distribution of auxiliary elevation correction reference data;

[0010] Step 2: Extract DEM elevation of seed points;

[0011] Extract the grid elevation value of the DEM produced by uncontrolled satellite photogrammetry corresponding to the location of the seed point, and calculate the difference dz between the grid elevation value and the reference elevation value of the seed point i:

[0012] Step 3: Compare and select elevation correction models;

[0013] Evaluate the effect of height difference fitting correction and select the optimal fitting model by comparison;

[0014] Step 4: Use the determined optimal fitting model to correct DEM;

[0015] Using the elevation correction reference data, the DEM elevation is corrected based on the selected optimal fitting model;

[0016] Step 5: Calculate the relative reliability index of the elevation accuracy of the corrected DEM;

[0017] Step 6: Fusion and transfer correction of adjacent DEMs;

[0018] Adjacent DEMs are weighted and fused according to the relative reliability index. For uncontrolled satellite stereo DEMs that have not been corrected with external reference data, the elevation values ​​of the corrected DEM that overlaps with the DEM are used as the elevation data source for extracting seed points for correcting the DEM. The elevation control correction is then transferred and the corrections are made in the order of the ratio of the area of ​​the polygon surrounding the seed point to the area of ​​the DEM.

[0019] Step 7: Acquire large-scale DEM data. "Large-scale" in this invention refers to tens of thousands of square kilometers. This method, based on different satellite data sources and correction data, can be used to acquire large-scale DEM data, with high versatility.

[0020] In the above technical solution, in step 1, the method for selecting the seed point is as follows:

[0021] The number of seed points is greater than or equal to 3. Under the existing data, they are distributed as evenly as possible in the range covered by the DEM. The slope map of the range where the DEM is located is generated using the uncontrolled DEM, and the slope values ​​of each candidate point are extracted. Under the condition of reasonable distribution, points with slope values ​​less than 15° are selected as seed points. Reference points other than the seed points can be used as check points.

[0022] In the above technical solution, in step 2, the difference dz between the grid elevation value and the reference elevation value of the seed point is i The calculation method is:

[0023] dz i =Zdem i -Z i (1)

[0024] In formula (1): Zdem i is the uncontrolled DEM elevation corresponding to the seed point location, Zi The reference elevation value of the seed point.

[0025] In the above technical solution, in step three, the elevation correction model is compared and selected, the elevation difference fitting correction effect is evaluated, and the optimal fitting model is selected by comparison. The elevation correction model used for comparison and selection in the present invention includes three elevation correction models, namely: constant correction, first plane fitting correction and quadratic surface fitting correction;

[0026] Constant correction, the specific method is: directly use the seed point dz i The mean of is used as the correction parameter;

[0027] A plane fitting correction, the specific method is:

[0028] dz i =a+b×X i +c×Y i (2)

[0029] In formula (2): (X i ,Y i ) is the plane coordinate of the seed point, a, b, c are the parameters to be solved for the fitted height correction plane; an error equation can be established for each seed point, and the least squares method can be used to solve the fitted plane parameters;

[0030] Quadratic surface fitting correction, the specific method is:

[0031] dz i =a+b×X i +c×Y i +d×X i ×Y i +e×X i 2 +f×Y i 2 (3)

[0032] In formula (3): (X i ,Y i ) is the plane coordinate of the seed point, a, b, c, d, e, and f are the parameters to be solved for the fitted height correction surface; an error equation can be established for each seed point, and the least squares method can be used to solve the fitting surface parameters.

[0033] In the above technical solution, when a plane fitting correction is adopted, the number of seed points is greater than or equal to 3;

[0034] When using quadratic surface fitting correction, the number of seed points is greater than or equal to 6.

[0035] In the above technical solution, in step 3, the difference Δdz between the fitted elevation correction value and the actual elevation difference is calculated. iEvaluate the effect of the corrected model and select the optimal fitting model:

[0036] ① Statistical Δdz i If the systematic error of the satellite RFM model is well corrected, the residual elevation error is mainly accidental error. The residual fitting elevation error Δdz i Dividing it by the calibrated elevation accuracy of the satellite sensor gives the normalized difference in elevation error:

[0037]

[0038] In formula (4): Δdz i N The standardization of elevation error is poor. The calibrated elevation accuracy of the satellite sensor can generally be found on the official website of the competent department or in the research results of relevant researchers. Normalization of elevation errors can weaken the impact of different satellite sensor accuracies on the accuracy analysis results of the elevation correction model.

[0039] Ideally, the normalized difference in elevation error is Δdz i N Satisfying the standard normal distribution N(0,1), the test statistic can be constructed as:

[0040]

[0041] In formula (5): U is the test statistic; is the mean of the standardized difference of the DEM elevation error; μ0 is the expectation of the normal distribution, which is 0 here; σ0 is the standard deviation of the normal distribution, which is 1 here; n is the number of points.

[0042] Perform a two-sided hypothesis test based on the significance level α; if U>|U α / 2 |Then the null hypothesis is rejected, and it is believed that the standardization of the elevation error after fitting the DEM elevation correction model is poor and does not meet the standard normal distribution, that is, the elevation error fitting model fails to fit the elevation error well; otherwise, it is believed that the fitting model is accurate and can better remove the systematic error of the uncontrolled DEM;

[0043] ② In the process of actual DEM elevation error fitting, if the number of reference points is small and the texture differences in different regions are too large, the matching accuracy is inconsistent. The above-mentioned fitting model accuracy analysis method based on hypothesis testing is not always effective. The direct statistical residual fitting elevation error Δdz can be used. i The mean and the root mean square error M rms :

[0044]

[0045] In formula (6): Δdz i is the residual fitting height error of a reference point, n is the number of points. The absolute value of the mean of the residual fitting height error The closer it is to 0, the smaller the residual root mean square error M of the fitting height error rms The smaller the value, the better the model fitting effect.

[0046] According to the optimal elevation correction model determined by the above comparison, the uncontrolled DEM is subjected to elevation correction to obtain the corrected DEM.

[0047] When the amount of data is sufficient, method ① and method ② can be used simultaneously for joint evaluation to ensure more rigorous theory. When the amount of data is insufficient and method ① is ineffective, method ② can be directly used for evaluation.

[0048] In the above technical solution, in step 5, there are k DEMs corrected with external reference elevation data in the entire survey area. The area of ​​the outer polygon of the seed point used to correct the jth (j∈[1,k]) DEM is The total coverage area of ​​this DEM is The root mean square error of the residual fitting elevation error of this DEM is The relative reliability index of the elevation accuracy of the corrected DEM is defined as:

[0049]

[0050] In formula (7): is the mean of the root mean square error of the residual fitting elevation error of k DEMs, R j R is the relative reliability index of the elevation accuracy of the DEM, which mainly measures the relative reliability of the elevation accuracy of the corrected DEM among the multiple DEMs uniformly produced in the survey area. j The larger it is, the higher the relative reliability of the DEM elevation accuracy.

[0051] In the above technical solution, in step 6, the fusion of adjacent DEMs is carried out as follows:

[0052] When the DEMs in the overlapping areas are corrected using external reference data, the fusion of the DEMs in the overlapping areas is obtained by weighted averaging the relative reliability indices of the multiple DEMs:

[0053]

[0054] In formula (8): Z fused is the DEM elevation after fusion of overlapping areas, R l Z is the relative reliability index of the DEM to be fused.l is the elevation corresponding to the DEM; m is the number of DEM slices;

[0055] When a DEM lacks external reference data correction, the elevation values ​​of the corrected DEM that overlaps with the DEM are used as the elevation data source for extracting seed points for correcting the DEM, and the elevation control correction is transferred. The correction is performed in priority order according to the ratio of the area of ​​the polygon surrounding the seed point to the area of ​​the DEM. The seed points are selected through step one, and the DEM is corrected through steps two to four. The DEM corrected in this step is fused with the DEM of the surrounding overlapping area by directly taking the arithmetic average of the elevations of multiple DEMs.

[0056] The present invention addresses the problems of low DEM elevation accuracy obtained by uncontrolled optical stereo satellites and the difficulty in high-precision fusion of DEMs from different data sources and different widths. A highly versatile optical stereo satellite DEM elevation correction and fusion method is proposed. Compared with the existing technology, the present invention has the following advantages:

[0057] The present invention comprehensively utilizes various archived satellite data and a small amount of other reference data that can provide elevation control. By constructing an elevation error fitting model and selecting the optimal fitting model through statistical analysis, the systematic error in the original DEM is removed, and the elevation accuracy of the DEM is improved. By defining the relative reliability index of the elevation accuracy of the DEM after correction with external reference data, a scientific DEM fusion method is proposed. The elevation control transfer between the corrected and uncorrected DEMs is sequentially carried out for cyclic correction and fusion, overcoming the problem of inconsistent accuracy between different data sources and different DEMs, and the obvious discrepancies easily found in simple arithmetic average fusion. The invention realizes the improvement of the elevation positioning accuracy of uncontrolled satellite photogrammetry by utilizing abundant archived satellite data from different sources and a small amount of external reference data. Through the scientific DEM fusion method, the acquisition of large-scale DEMs covering tens of thousands of square kilometers in vast uninhabited areas is achieved. BRIEF DESCRIPTION OF THE DRAWINGS

[0058] Figure 1 The present invention provides a technical flow chart of an optical stereo satellite DEM elevation correction and fusion method.

[0059] Figure 2 This is a statistical histogram of the elevation accuracy of the uncontrolled optical stereo satellite DEM according to an embodiment of the present invention.

[0060] Figure 3 This is a statistical histogram of the elevation accuracy of the corrected uncontrolled optical stereo satellite DEM according to an embodiment of the present invention.

[0061] exist Figure 2 In the figure, the horizontal axis represents the elevation difference, unit: meter; the vertical axis represents the number of points, unit: point.

[0062] exist Figure 3 In the figure, the horizontal axis represents the elevation difference, unit: meter; the vertical axis represents the number of points, unit: point. DETAILED DESCRIPTION

[0063] The following detailed description of the embodiments of the present invention is given in conjunction with the accompanying drawings, which do not limit the present invention but are merely examples. The description makes the advantages of the present invention clearer and easier to understand.

[0064] Refer to the attached figure: An optical stereo satellite DEM elevation correction and fusion method, the technical process is shown in Figure 1 , the specific contents are as follows:

[0065] (1) Based on the range of each DEM produced by uncontrolled satellite photogrammetry and the distribution of auxiliary elevation correction reference data, the seed points used for elevation correction are determined. The selection of seed points mainly considers two aspects:

[0066] First, under the existing data, the seed points should be distributed as evenly as possible in the area covered by the DEM, and the number of seed points should be no less than 4. Second, the slope map of the range where the DEM is located is generated using the uncontrolled DEM, and the slope values ​​of each candidate point are extracted. Under the condition of reasonable distribution, points in flat areas (with smaller slope values) are selected as seed points, and reference points other than the seed points can be used as check points.

[0067] (2) Extract the grid elevation value of the DEM produced by uncontrolled satellite photogrammetry corresponding to the location of the seed point, and calculate the difference dz between the grid elevation value and the reference elevation value of the seed point i :

[0068] dz i =Zdem i -Z i (1)

[0069] In formula (1): Zdem i is the uncontrolled DEM elevation corresponding to the seed point location, Z i The reference elevation value of the seed point.

[0070] (3) Fitting elevation difference dz i With the plane coordinates (X i ,Y i )

[0071] The Rational Function Model (RFM) is a universal sensor model widely used in high-resolution remote sensing satellites. my country's Ziyuan-3 and Gaofen series stereo mapping satellites use this model instead of strict geometric models, achieving accuracy better than 1% pixel. RFM models may contain systematic errors and require correction using control points. A common RFM correction method involves inserting control points into the imagery, then using an affine transformation to establish a correction relationship between the model's calculated image point coordinates and the true image point coordinates. The affine transformation correction parameters are calculated using least squares, and these parameters are used to optimize the original RPC parameters. In uncontrolled satellite photogrammetry, the RFM is directly used to construct a stereo model. After matching congruent points and performing triangulation adjustment, a large number of congruent points are obtained using dense point matching techniques. Based on the forward intersection of the RFM, the three-dimensional coordinates of the ground points are obtained. Finally, an irregular triangulated network is constructed and interpolated to obtain a DEM. This method primarily uses elevation correction reference data to correct the DEM elevation for object space. Three elevation correction models are provided:

[0072] ① Constant correction

[0073] Directly use seed point dz i The mean of is used as the correction parameter.

[0074] ②One-time plane fitting correction

[0075] dz i =a+b×X i +c×Y i (2)

[0076] In formula (2): (X i ,Y i ) are the plane coordinates of the seed point, and a, b, and c are the parameters to be solved for the height correction plane. An error equation can be established for each seed point, and the least squares method can be used to solve for the plane parameters. It should be noted that when using a single plane fitting correction, at least three seed points are required.

[0077] ③ Quadratic surface fitting correction

[0078] dz i =a+b×X i +c×Y i +d×X i ×Y i +e×X i 2 +f×Y i 2 (3)

[0079] In formula (3): (X i ,Yi ) are the plane coordinates of the seed point, and a, b, c, d, e, and f are the parameters to be solved for the height correction surface. An error equation can be established for each seed point, and the least squares method can be used to solve for the fitted surface parameters. It should be noted that when using quadratic surface fitting correction, at least six seed points are required.

[0080] (4) Evaluate the effect of height difference fitting correction and select the optimal fitting model

[0081] The difference between the fitted elevation correction value and the actual elevation difference is Δdz i Evaluate the effectiveness of the corrected model:

[0082] ① Statistical Δdz i If the systematic error of the satellite RFM model is well corrected, the residual elevation error is mainly accidental error. The residual fitting elevation error Δdz i Dividing it by the calibrated elevation accuracy of the satellite sensor gives the normalized difference in elevation error:

[0083]

[0084] In formula (4): Δdz i N The standardization of elevation error is poor. The calibrated elevation accuracy of the satellite sensor can generally be found on the relevant authority's website or in the research results of relevant researchers. Normalizing the elevation error can reduce the impact of different satellite sensor accuracies on the accuracy analysis results of the elevation correction model.

[0085] Ideally, the normalized difference in elevation error is Δdz i N Satisfying the standard normal distribution N(0,1), the test statistic can be constructed as:

[0086]

[0087] In formula (5): U is the test statistic; is the mean of the standardized difference of the DEM elevation error; μ0 is the expectation of the normal distribution, which is 0 here; σ0 is the standard deviation of the normal distribution, which is 1 here; n is the number of points.

[0088] Perform a two-sided hypothesis test based on the significance level α. α / 2 |Then the null hypothesis is rejected, and it is considered that the normalization of the elevation error after fitting the DEM elevation correction model is poor and does not meet the standard normal distribution, that is, the elevation error fitting model fails to fit the elevation error well. Otherwise, it is considered that the fitting model is accurate and can effectively remove the systematic error of the uncontrolled DEM.

[0089] ② In the process of actual DEM elevation error fitting, if the number of reference points is small and the texture differences in different regions are too large, the matching accuracy is inconsistent. The above-mentioned fitting model accuracy analysis method based on hypothesis testing is not always effective. The direct statistical residual fitting elevation error Δdz can be used. i The mean and the root mean square error M rms :

[0090]

[0091] In formula (6): Δdz i is the residual fitting height error of a reference point, n is the number of points. The absolute value of the mean of the residual fitting height error The closer it is to 0, the smaller the residual root mean square error M of the fitting height error rms The smaller the value, the better the model fitting effect.

[0092] According to the best elevation correction model determined by the above comparison, the uncontrolled DEM is corrected to obtain the corrected DEM.

[0093] (5) There are k DEMs corrected using external reference elevation data in the entire survey area. The area of ​​the outer polygon of the seed point used to correct the jth (j∈[1,k]) DEM is The total coverage area of ​​this DEM is The root mean square error of the residual fitting elevation error of this DEM is The relative reliability index of the elevation accuracy of the corrected DEM is defined as:

[0094]

[0095] In formula (7): is the mean of the root mean square error of the residual fitting elevation error of k DEMs, R j R is the relative reliability index of the elevation accuracy of the DEM, which mainly measures the relative reliability of the elevation accuracy of the corrected DEM among the multiple DEMs uniformly produced in the survey area. j The larger it is, the higher the relative reliability of the DEM elevation accuracy.

[0096] (6) The fusion process of adjacent DEMs can be divided into the following two cases:

[0097] ① If the DEMs in the overlapping areas are all corrected using external reference data, the fusion of the DEMs in the overlapping areas is obtained by weighted averaging the relative reliability indices of multiple DEMs.

[0098]

[0099] In formula (8): Z fused is the DEM elevation after fusion of overlapping areas, R l Z is the relative reliability index of the DEM to be fused. l is the elevation corresponding to the DEM.

[0100] ② If a DEM lacks external reference data correction, the elevation value of the corrected DEM with overlapping areas is used as the elevation data source for extracting the seed points for correcting the DEM, and the elevation control correction is transferred. The correction is performed in the order of priority based on the ratio of the area of ​​the polygon surrounding the seed point to the area of ​​the DEM. In the process of seed point selection, the basic principles in (1) are still followed, and then the DEM is corrected according to the steps in (2) to (4). The DEM corrected in this step is directly fused with the DEM of the surrounding overlapping areas by the arithmetic mean of the elevations of multiple DEMs. If multiple DEMs lack external reference data correction, the correction and fusion are performed cyclically according to this step.

[0101] (7) Obtain large-scale DEM data through the above steps (1) to (6).

[0102] Example:

[0103] The present invention is now described in detail by taking the trial application of the present invention to a certain research area for optical stereo satellite DEM elevation correction and fusion as an example, which also has a guiding role in the application of the present invention to other areas for optical stereo satellite DEM elevation correction and fusion.

[0104] In this example, a study area is located in the Grove Mountains region of East Antarctica, which is the main area for geological and geophysical surveys from Zhongshan Station to Dome A. Grove Mountains is 400-500 km away from Zhongshan Station. The entire Grove Mountains region is geographically bounded by 73°40′ to 76°00′ east longitude and 72°15′ to 73°15′ south latitude, with a total area of ​​approximately 8,000 km. 2 It is a mountainous area with exposed peaks on the inland ice sheet, and the terrain gradually rises from northwest to southeast.

[0105] This example uses two-view archived Resource III satellite data covering the Grove Mountain area to perform uncontrolled satellite photogrammetry to produce a DEM. ICESat altimetry data is then used as external reference data. The DEM is corrected for object space according to the method of the present invention. Other reference data outside the correction are used as checkpoints, and the accuracy information of the checkpoints is calculated.

[0106] like Figure 2 As shown, Figure 2The left picture is the elevation accuracy statistical histogram of the DEM of the Jingziyuan-3 uncontrolled optical stereo satellite. Figure 2 The right picture is the elevation accuracy statistical histogram of another adjacent scene, Ziyuan-3 uncontrolled optical stereo satellite DEM; Figure 3 As shown, Figure 3 The left picture is Figure 2 The left picture corresponds to the elevation accuracy statistical histogram of the corrected stereo satellite DEM. Figure 3 The picture on the right is Figure 2 The right figure corresponds to the elevation accuracy statistical histogram of the corrected stereo satellite DEM. The elevation accuracy of the two DEMs obtained by uncontrolled satellite photogrammetry was -5.1±4.1m and -15.5±6.1m respectively. After correction using the method of the present invention, the elevation accuracy of the DEM was -0.04±1.6m and 0.1±3.0m respectively. The elevation accuracy of the DEM has been significantly improved, and the mean value of the DEM elevation error after correction is close to 0 and is normally distributed, indicating that the systematic error in the DEM measured by uncontrolled satellite photogrammetry has been well removed. In this embodiment, the accuracy difference between the two DEMs before correction is large, and direct fusion will result in large deviations. The external reference data involved in the correction in this embodiment is widely distributed and rich in data volume. The two adjacent DEMs have been well corrected. DEM fusion is performed according to the method of the present invention. There are no abnormal values ​​in the DEM elevation, the transition is natural, and the fusion effect is good, which verifies the effectiveness of the method of the present invention.

[0107] The above description is merely an exemplary embodiment of the present invention and should not be construed as limiting the scope of patent protection of the present invention. The scope of patent protection of the present invention is determined by the appended claims.

[0108] Other parts not described belong to the prior art.

Claims

1. An optical stereo satellite DEM elevation correction and fusion method, characterized by: The following steps are included: Step 1: Select the seed point for elevation correction; Determine the seed points for elevation correction based on the range of each DEM produced by uncontrolled satellite photogrammetry and the distribution of auxiliary elevation correction reference data; The method for selecting seed points is as follows: the number of seed points is greater than or equal to 3 points; and the slope map of the range of the DEM is generated using the uncontrolled DEM, and the slope values ​​of each candidate point are extracted. In the case of reasonable distribution, points with slope values ​​less than 15° are selected as seed points; Step 2: Extract DEM elevation of seed points and calculate difference; Extract the grid elevation value of the DEM produced by uncontrolled satellite photogrammetry corresponding to the location of the seed point, and calculate the difference dz between the grid elevation value and the reference elevation value of the seed point i : Step 3: Compare and select elevation correction models; Evaluate the correction effect of height difference fitting, and compare and select the optimal height correction model; The elevation correction model is compared and selected, and the elevation difference fitting correction effect is evaluated, and the optimal elevation correction model is compared and selected, including three elevation correction models, namely: constant correction, first plane fitting correction and second surface fitting correction; Constant correction, the specific method is: directly use the seed point dz i The mean of the correction parameter dz i 0 ; The specific method of primary plane fitting correction is as follows: d i 1 =a1+b1×X i +c1×Y i (2) In formula (2): dz i 1 is the difference between the plane grid elevation value and the reference elevation value of the seed point; (X i ,Y i ) is the plane coordinate of the seed point, a1, b1, c1 are the parameters to be solved for the fitted height correction plane; an error equation is established for each seed point, and the least squares method is used to solve the fitted plane parameters; Quadratic surface fitting correction, the specific method is: <h2 style=";text-align:left;direction:ltr">dz<h2 style=";text-align:left;direction:ltr"> i <h2 style=";text-align:left;direction:ltr"> 2 <h2 style=";text-align:left;direction:ltr"> =a2+b2×X<h2 style=";text-align:left;direction:ltr"> i <h2 style=";text-align:left;direction:ltr"> +c2×Y<h2 style=";text-align:left;direction:ltr"> i <h2 style=";text-align:left;direction:ltr"> +d2×X<h2 style=";text-align:left;direction:ltr"> i <h2 style=";text-align:left;direction:ltr"> ×Y<h2 style=";text-align:left;direction:ltr"> i <h2 style=";text-align:left;direction:ltr"> +e2×X<h2 style=";text-align:left;direction:ltr"> i <h2 style=";text-align:left;direction:ltr"> 2 <h2 style=";text-align:left;direction:ltr"> +f2×Y<h2 style=";text-align:left;direction:ltr"> i <h2 style=";text-align:left;direction:ltr"> 2 <h2 style=";text-align:left;direction:ltr"> (3) In formula (3): dz i 2 is the difference between the quadratic surface grid elevation value and the seed point reference elevation value; (X i ,Y i ) is the plane coordinate of the seed point, a2, b2, c2, d2, e2, f2 are the parameters to be solved for the fitted height difference correction surface; an error equation is established for each seed point, and the least squares method is used to solve the fitting surface parameters; Step 4: Use the determined optimal elevation correction model to correct DEM; Using the elevation correction reference data, the DEM elevation is corrected based on the selected optimal elevation correction model; Step 5: Calculate the relative reliability index of the elevation accuracy of the corrected DEM; There are k DEMs corrected using external reference elevation data in the entire survey area. The area of ​​the outer polygon of the seed point used to correct the jth DEM is Where j∈[1,k], the total coverage area of ​​the DEM is The root mean square error of the residual fitting elevation error of this DEM is The relative reliability index of the elevation accuracy of the corrected DEM is defined as: In formula (7): is the mean of the root mean square error of the residual fitting elevation error of k DEMs, R j is the relative reliability index of the elevation accuracy of the DEM; Step 6: Adjacent DEM fusion and transfer correction; Adjacent DEMs are weighted and fused according to the relative reliability index; For the uncontrolled satellite stereo DEM without external reference data correction, the elevation value of the corrected DEM with overlapping area with the DEM is used as the elevation data source for extracting the seed points for correcting the DEM, and the elevation control correction is transferred. The correction is carried out in the priority order of the ratio of the area of ​​the polygon outside the seed point to the area of ​​the DEM. Step 7: Obtain large-scale DEM data using the data obtained in step 6.

2. The optical stereo satellite DEM elevation correction and fusion method according to claim 1, characterized in that: In step 2, the difference dz between the grid elevation value and the reference elevation value of the seed point i The calculation method is: d i =Zdem i -WITH i (1) In formula (1): Zdem i is the uncontrolled DEM elevation corresponding to the seed point location, Z i The reference elevation value of the seed point.

3. The optical stereo satellite DEM elevation correction and fusion method according to claim 1, characterized in that: When using a single plane fitting correction, the number of seed points is greater than or equal to 3; When using quadratic surface fitting correction, the number of seed points is greater than or equal to 6.

4. The optical stereo satellite DEM elevation correction and fusion method according to claim 3 is characterized in that: In step 3, the difference between the fitted elevation correction value and the actual elevation difference is △dz i Evaluate the effect of the elevation correction model and select the optimal elevation correction model: ①Statistics △dz i If the systematic error of the satellite RFM model is well corrected, the residual elevation error is mainly accidental error. i Dividing by the calibrated elevation accuracy of the satellite sensor, we get the normalized difference in elevation error: In formula (4): △dz i N The standardization of elevation error is poor. The calibrated altitude accuracy of the satellite sensor; Ideally, the standardization of elevation error is poor △dz i N Satisfying the standard normal distribution N(0,1), the test statistic can be constructed as: In formula (5): U is the test statistic; is the standardized mean of the DEM elevation error; μ0 is the expectation of the normal distribution, which is 0 here; σ0 is the standard deviation of the normal distribution, which is 1 here; n is the number of points; Perform a bilateral hypothesis test based on the significance level α; if U>|U α / 2 |Then the null hypothesis is rejected, and it is believed that the standardization of the elevation error after fitting the DEM elevation correction model is poor and does not meet the standard normal distribution, that is, the elevation correction model fails to fit the elevation error well; otherwise, it is believed that the elevation correction model fits accurately and can better remove the systematic error of the uncontrolled DEM; ② In the process of fitting the actual DEM elevation error, if the number of reference points is small and the texture difference in different regions is too large, the matching accuracy is inconsistent. The accuracy analysis method of the elevation correction model based on the above hypothesis test is not always effective. The direct statistical residual fitting elevation error △dz is used. i The mean and the root mean square error M rms : In formula (6): △dz i is the residual fitting elevation error of a reference point, and n is the number of points.

5. The optical stereo satellite DEM elevation correction and fusion method according to claim 1, characterized in that: In step six, adjacent DEMs are fused as follows: When the DEMs in the overlapping areas are corrected using external reference data, the fusion of the DEMs in the overlapping areas is obtained by weighted averaging the relative reliability indexes of multiple DEMs: In formula (8): Z fused is the DEM elevation after fusion of the overlapping area, R l is the relative reliability index of the DEM to be fused, Z l is the elevation corresponding to the DEM; When a DEM lacks external reference data correction, the elevation value of the corrected DEM with overlapping areas with the DEM is used as the elevation data source for extracting seed points for correcting the DEM, and the elevation control correction is transferred. The correction is made in priority order according to the ratio of the area of ​​the polygon surrounding the seed point to the area of ​​the DEM. The seed points are selected through step one, and the DEM is corrected through steps two to four. The DEM corrected in this step is fused with the DEM of the surrounding overlapping areas by directly taking the arithmetic mean of the elevations of multiple DEMs.

Images