IGG iii anti-difference estimation-based polar surface elevation change estimation method and system

Abstract

The application provides a polar surface elevation change estimation method and system based on IGG III anti-difference estimation, satellite altimetry data extraction; satellite altimetry data processing, including extracting qualified data points, performing geophysical correction, performing along-track median filtering based on a sliding window; separating the data into ascending and descending orbits and storing according to the orbit number; performing elevation change calculation based on IGG III weight function, including extracting the point with the largest number of points in the orbit as a reference orbit, interpolating the points on other repeat orbits to the points with consistent latitude of the reference orbit through inverse distance weighting; assuming that the elevation value is linearly related to the distance under the condition that a group of repeat orbit points are at the same latitude, performing least squares based on an elevation change function model, introducing the standard deviation of the least squares as a priori accuracy into the IGG III weight function for calculation, and iteratively obtaining elevation change scattered points; performing Kriging interpolation method to generate an elevation change result and output.

Description

Technical Field

[0001] This invention belongs to the field of geodesy and surveying engineering, and specifically relates to a technical solution for estimating surface elevation changes in polar regions based on IGGⅢ robust estimation. Background Technology

[0002] Antarctica and Greenland are two major cold sources in the global climate, playing a crucial role in regulating global climate and environmental change. Changes in the atmosphere, oceans, and glaciers in the polar regions directly or indirectly affect global atmospheric and oceanic circulation. The Antarctic ice sheet regulates the Earth's energy balance by reflecting sunlight and absorbing heat, and participates in atmospheric circulation through evaporation and precipitation. Therefore, monitoring the mass balance of the Antarctic ice sheet can help understand and predict global climate change trends. However, with the intensification of climate change caused by human activities, the melting rate of the Antarctic ice sheet is accelerating. This is leading to a faster rate of sea-level rise. The melting of the Antarctic ice sheet has already caused sea levels to rise by approximately 15 millimeters, threatening low-lying areas such as coastal cities and islands. Therefore, special attention needs to be paid to changes in the mass balance of the Antarctic ice sheet to predict potential sea-level rise and address its risks.

[0003] Antarctic elevation changes are a key parameter for retrieving polar mass balance. By calculating the elevation changes of the Antarctic ice sheet and combining this with surface area and snow density, the mass balance of the Antarctic region can be determined. The accuracy of this calculation significantly impacts the assessment of polar mass balance. To better understand polar mass balance and explore the coupling effects of atmosphere-sea ice-ocean, accurate assessment of Antarctic elevation changes is crucial.

[0004] Existing satellite altimetry methods for obtaining elevation changes on the polar surface include the intersection method, the repeat orbit method, and the DEM difference method. The intersection method obtains the intersection points between the satellite's ascending and descending orbits, and calculates the elevation change by the ratio of the time and elevation difference between two satellite passages through the intersection point. This method has low data utilization, especially in the low-latitude Antarctic region where orbits are not densely packed. The DEM difference method generates a DEM from the satellite's annual altimetry data and compares the DEMs of adjacent years to obtain the elevation change. This method requires generating a DEM before comparison, which can easily introduce other errors. The repeat orbit method uses the satellite's revisit period in the polar region to obtain multiple periods of data along the same orbit, and combines these multiple periods of data for a comprehensive solution. While this method makes full use of the satellite's nadir data, combining multiple periods of data significantly increases the uncertainty of the calculation process. Traditional methods use the Laida criterion to iterate the data, calculating the residuals and standard deviations for each iteration, and discarding values ​​with residuals greater than three times the standard deviation, but do not process residuals within the three-standard-deviation range. Traditional methods are greatly affected by individual geoid anomalies, significantly impacting the uncertainty of elevation change results. Therefore, this invention aims to propose a novel method for estimating polar ice surface elevation changes based on IGIII robust estimation. While maintaining high spatial resolution through repeated orbit calculations, this method reduces the uncertainty in estimating ice shelf elevation changes by refining data within a three-standard-deviation range. Summary of the Invention

[0005] To address the shortcomings of the existing technologies, this invention proposes a method and system for estimating the elevation changes of polar ice surfaces based on IGGⅢ robust estimation.

[0006] To achieve the above objectives, this invention proposes a method for estimating the surface elevation variation in polar regions based on IGGⅢ robust estimation, comprising the following processing steps:

[0007] Satellite altimetry data extraction includes extracting longitude, latitude, elevation values, time, mass labels, tidal corrections, and counter-pressure correction parameters;

[0008] Satellite altimetry data processing includes extracting qualified data points based on data quality labels, performing geophysical corrections on the data, performing median filtering along the orbit based on a sliding window to remove elevation outliers along the orbit, and separating the data into ascending and descending orbits and storing it according to the orbit number.

[0009] The elevation change calculation is performed based on the IGGⅢ weight function, and the implementation method is as follows.

[0010] For each set of orbit numbers, the point with the most points in the orbit is extracted as the reference orbit, and the points on other repeating orbits are interpolated to the point with the same latitude as the reference orbit using an inverse distance weighting method.

[0011] Establish an elevation variation function model for points at the same latitude;

[0012] Assuming that for a set of repeating orbit points at the same latitude, the elevation value is linearly related to the distance, the least squares model based on the elevation change function is used, and the standard deviation of the least squares is introduced as the prior accuracy into the IGGⅢ weight function for solution, and the elevation change scatter results are obtained iteratively.

[0013] Perform kriging interpolation on the scatter plot of elevation changes to generate and output the elevation change results.

[0014] Furthermore, after extracting qualified data points based on data quality labels, data points with an elevation greater than 5000 meters are removed.

[0015] Furthermore, the method for separating the ascending and descending rails is to calculate the slope based on the latitude and longitude of the ascending and descending rails. If the slope is positive, it is determined to be a descending rail; if the slope is negative, it is determined to be an ascending rail.

[0016] Furthermore, for points at the same latitude, the following elevation variation function model can be established:

[0017]

[0018] In the formula, h0 represents the elevation value at the initial moment. α represents the annual average elevation change, A and B represent the amplitude of seasonal variation, t represents the time of the nadir measurement relative to the initial moment, D represents the distance between the point on the repeating orbit and the reference orbit, and tanα represents the slope along the latitudinal direction.

[0019] Furthermore, the IGGⅢ equivalent weight function adopts the following model.

[0020]

[0021] In the formula, k0 and k1 are preset parameters, σ represents the mean square error of each iteration, and v i Represents the residual of each iteration, |v i / | represents the ratio of the residual to the mean error after each iteration, p i This represents the weight after each iteration.

[0022] On the other hand, the present invention provides a polar surface elevation change estimation system based on IGGⅢ robust estimation, which is used to implement the polar surface elevation change estimation method based on IGGⅢ robust estimation as described above.

[0023] Moreover, it includes the following modules,

[0024] The first module is used for satellite altimetry data extraction, including the extraction of longitude, latitude, elevation values, time, mass labels, tidal correction, and counter-pressure correction parameters;

[0025] The second module is used for satellite altimetry data processing, including extracting qualified data points based on data quality labels, performing geophysical corrections on the data, performing median filtering along the orbit based on a sliding window to remove elevation anomalies along the orbit, and separating the data into ascending and descending orbits and storing it according to the orbit number.

[0026] The third module is used to calculate elevation changes based on the IGGⅢ weight function. The implementation method is as follows: for each group of track numbers, the point with the most points in the track is extracted as the reference track, and the points on other repeated tracks are interpolated to the point with the same latitude as the reference track by inverse distance weighting.

[0027] Establish an elevation variation function model for points at the same latitude;

[0028] Assuming that for a set of repeating orbit points at the same latitude, the elevation value is linearly related to the distance, the least squares model based on the elevation change function is used, and the standard deviation of the least squares is introduced as the prior accuracy into the IGGⅢ weight function for solution, and the elevation change scatter results are obtained iteratively.

[0029] The fourth module is used to perform kriging interpolation on the scatter plot of elevation changes, generate elevation change results, and output them.

[0030] Alternatively, it may include a processor and a memory, the memory being used to store program instructions, and the processor being used to call the stored instructions in the memory to execute, as described above, a polar surface elevation change estimation method based on IGGⅢ robust estimation.

[0031] Alternatively, it may include a readable storage medium storing a computer program that, when executed, implements a polar surface elevation variation estimation method based on IGGⅢ robust estimation as described above.

[0032] The ice shelf surface elevation change estimation method based on repeated orbits proposed in this invention can effectively reduce the uncertainty caused by outliers and improve the reliability of the ice shelf surface elevation change results.

[0033] The present invention is simple and convenient to implement, highly practical, and solves the problems of low practicality and inconvenience in actual application of related technologies. It can improve the accuracy of ice sheet elevation change calculation and has important scientific significance. Attached Figure Description

[0034] Figure 1 This is a flowchart of an embodiment of the present invention;

[0035] Figure 2 This is a schematic diagram illustrating the inverse distance weighting method used in this embodiment of the invention to obtain interpolation values ​​along a repeating trajectory.

[0036] Figure 3 The results of the ice shelf surface elevation change distribution of the Ross Ice Shelf in Antarctica from December 2019 to March 2021 are from an embodiment of the present invention. Detailed Implementation

[0037] The technical solution of the present invention will be described in detail below with reference to the accompanying drawings and embodiments.

[0038] This invention provides a method for estimating the elevation changes of polar ice shelves based on robust IGGIII estimation. Current methods for estimating elevation changes using satellite altimetry require obtaining parameters such as the latitude, longitude, elevation, and geophysical corrections of the satellite's nadir. Appropriate data preprocessing is then performed based on regional characteristics, primarily including data quality monitoring, geophysical correction, and along-orbit filtering. Multiple elevation observations of the same area collected at different times are obtained using the satellite's designed revisit period, and an elevation change function model is established to calculate the surface elevation changes of the ice shelf. In the least-squares calculation of elevation changes, this invention does not use the traditional Laida criterion as a threshold for data deletion; instead, it employs the IGGIII (Institute of Geodesy and Geophysics III) weighting function for evaluation. Based on the IGGⅢ weight function, the technical solution of this invention assumes that the main body of elevation changes still follows a normal distribution, with doubtful and rejected values ​​accounting for only a small portion. In the robust estimation process, reliable information from the main observation area is fully utilized, retaining its original weights. For observations falling within the doubtful segment, weighting is reduced according to their reliability. Significant anomalies are disregarded, and their weights are zero. Compared to the traditional Laida criterion, the IGGⅢ robustness provides more refined processing for doubtful segments, effectively improving the reliability of elevation change calculations.

[0039] Current research papers have proposed the application of IGG robust estimation in elevation fitting, but these are still far from addressing the technical problems that this invention aims to solve. This invention proposes several innovative designs: the data used is derived from laser satellite altimetry data, rather than GPS static measurement data. The accuracy of the two types of data differs significantly, with GPS static measurement data being significantly more accurate than laser altimetry satellite data. Therefore, applying IGG robust estimation to laser satellite altimetry data will optimize the results better than applying it to GPS data processing. This invention focuses on elevation variation, not elevation fitting; therefore, different fitting models are used. Unlike the paper that uses surface-based quadratic fitting over a range of regions, this invention uses a set of repeating orbital points at the same latitude for seasonal variation-based fitting. Furthermore, due to the use of different models, different ranges need to be selected for IGG III robust estimation for different models. After experiments, this invention optimized the research ranges for k0 = 1.5 and k1 = 2.5.

[0040] See Figure 1 The embodiment proposes a method for estimating the surface elevation change in polar regions based on IGGⅢ robust estimation, which includes the following steps:

[0041] Step 1, Extraction of satellite altimetry data

[0042] This embodiment is based on the ATL06 laser altimeter product released by ICESat-2. Taking the Ross Ice Shelf as the study area, it extracts parameters such as longitude, latitude, elevation (h_li), time (delta_time), quality summary (atl06_quality_summary), tide correction (tide_ocean), and reverse pressure correction (dac) for all subsatellite points within the ice shelf. Longitude, latitude, elevation, and time are used to construct an elevation change function to calculate elevation change values. Quality summary, tide correction, reverse pressure correction, and other quality summaries, along with geophysical parameters, are used for data preprocessing and gross error removal. It can also record duplicate orbit numbers to facilitate subsequent operations.

[0043] Step 2, Satellite altimetry data processing

[0044] The embodiment further proposes a preferred implementation of step 2, which includes the following sub-steps:

[0045] Step 2.1: Extract data points that meet the required quality standards based on the data quality labels.

[0046] The implementation example, based on the evaluation in the NSIDC ATL06 user manual, identifies points with a quality label of 1 as those with poor acquisition conditions or significant errors in the fitting process. Furthermore, since the elevation across Antarctica is less than 5000m, points with elevations greater than 5000m are considered to be photon points reflecting incorrect cloud layers. To ensure calculation accuracy, points with a data quality label of 1 and points with elevations greater than 5000m are removed.

[0047] Step 2.2, perform geophysical corrections on the data, such as tidal corrections and inverse pressure corrections:

[0048] Because ice shelves float on the sea surface, geophysical parameter corrections are necessary. ICESat-2 already implemented some corrections when releasing the ATL06 product; therefore, tidal and reverse pressure corrections are also required to obtain the corrected elevation.

[0049]

[0050] Among them, h i The elevation is before correction; tide_ocean is for tidal correction, and dac is for reverse pressure correction.

[0051] Step 2.3, median filtering along the track based on a sliding window: According to the preset sliding window size, median filtering is performed on each track to remove elevation outliers along the track.

[0052] In practice, the size of the sliding window can be set by the user according to the situation. In this embodiment of the invention, the empirical value of 15 km obtained from polar region experiments is preferably used.

[0053] Step 2.4: Separate the data into ascending and descending tracks and store them according to track number: Compare the start and end points of each track to distinguish between ascending and descending tracks, and store them separately by period for subsequent calculations. The method for separating ascending and descending tracks uses the latitude and longitude of the ascending and descending tracks to calculate the slope. The specific formula is as follows:

[0054]

[0055] In the formula: lon1,lat1 represent the latitude and longitude of the starting point; lon2,lat2 represent the latitude and longitude of the ending point; k1 represents the calculated slope. If the slope is positive, it is judged as a descending orbit; if the slope is negative, it is judged as an ascending orbit.

[0056] In practice, to improve efficiency, it can be determined whether the number of repeated tracks meets the starting requirements. If yes, proceed to step 3; otherwise, remove the track.

[0057] Step 3: Estimate ice shelf elevation changes based on IGGⅡ weighting function

[0058] This step, based on the first law of geosciences, assumes that for a set of repeating orbit points at the same latitude, elevation values ​​are linearly correlated with distance, and that the least squares iterative process utilizes the IGGⅢ robust function. This improves the spatial resolution of the repeating orbits, avoids the unprocessed data segments encountered when using the Laida criterion, and enhances the reliability of elevation change calculations.

[0059] The embodiment further proposes a preferred implementation of step 2, which includes the following sub-steps:

[0060] Step 3.1: For each group of orbital numbers, extract the point with the most points within that orbit as the reference orbit, and interpolate the points on other repeating orbits using an inverse distance weighting method to the point with the same latitude as the reference orbit.

[0061] In this embodiment, the orbits of each repeating cycle are taken, with the point with the most points used as the reference orbit and the remaining orbits as repeating orbits. Inverse distance weighting is used to correct the longitude, elevation, and time of the points on the repeating orbits to the same latitude as the reference orbit. Based on the points on the reference orbit, the two points on a single repeating orbit closest to the reference orbit point are selected for inverse distance weighting. The formula for inverse distance weighting is as follows:

[0062]

[0063] In the above formula: d1, d2 represent the distances from the two nearest points on the repeating orbit to the point on the reference orbit; p1, p2 represent the weights of the two nearest points; lat1, lat2, h1, h2, and t1, t2 represent the longitude, elevation, and time of the two nearest points on the repeating orbit, respectively; lat inter ,h inter ,t inter These represent the interpolated longitude, elevation, and time of the point, respectively. The specific interpolation format is as follows: Figure 2 In the diagram, squares represent measured elevation points on the reference track, dots represent measured elevation points on the repeating track, and dots represent interpolated elevation points on the repeating track obtained by interpolating the elevation points on the reference track.

[0064] Step 3.2: Establish an elevation variation function model for points at the same latitude.

[0065] The example involves treating points at the same latitude as a group and substituting them into the elevation variation function model for least squares calculation. The preferred function model is:

[0066]

[0067] Where: h i (t) represents the elevation value of a point on the reference orbit at time t, and h0 represents the elevation value of the reference orbit point at the initial time. denoted by α, A represents the annual average elevation change, B represents the amplitude of seasonal variation, t represents the time measured at the nadir point relative to the initial moment, D represents the distance between the point on the repeating orbit and the reference orbit, tanα represents the slope along the latitudinal direction, and T represents the time of one year in seconds. Considering leap years, it is taken as 31,557,600 seconds here.

[0068] Step 3.3: Based on the elevation change function model, perform least squares calculations. Introduce the standard deviation of the least squares as prior accuracy into the IGGⅢ weight function to reweight the data. Then, use the changed weights to perform least squares calculations. Finally, introduce the recalculated residuals and standard deviations into the IGGⅢ weight function to reweight the data again, iterating to obtain the elevation change scatter plot results.

[0069] Traditional methods based on repeated orbits obtain results through least squares, then substitute these results back into the elevation variation function model to calculate residuals, discarding points whose residuals do not conform to the Laida criterion. However, suspicious values ​​within three standard deviations are not addressed. Therefore, this invention employs the IGGⅢ robust estimation model, which requires a least squares pass to provide prior accuracy. During subsequent iterations, the IGGⅢ equivalent weight function model is used to modify the weights. When the absolute value of the residual is greater than k1 times the standard error, it is considered to have a large error, and the IGGⅢ equivalent weight function discards values ​​with large errors by assigning them zero. When the absolute value of the residual falls between k0 and k1 times the standard error, it is considered to have a potential error, and suspicious values ​​are weighted using a weighting function. When the residual is less than k0 times the standard error, it is considered to contain only a small error, and small errors are weighted using least squares unit weights. This process is repeated iteratively until all residuals are less than k0 times the standard error, thus achieving robustness. The preferred embodiment provides the following IGGⅢ equivalent weight function model:

[0070]

[0071] In the formula, k0 and k1 are preset parameters, with k0 preferably taking the value 1.5 and k1 preferably taking the value 2.5. σ represents the mean square error of each iteration, and v i Represents the residual of each iteration, |v i / | represents the ratio of the residual to the mean error after each iteration, used to measure the magnitude of the residual, p i This represents the weight after each iteration.

[0072] That is, the following iterations are performed:

[0073] 1) Obtain the mean square error from least squares as the prior accuracy (i.e., σ in the IGGⅢ equivalent weight function model). The formula for calculating the mean square error is as follows. Proceed to step 2).

[0074]

[0075] In the formula: σ represents the mean square error of each iteration, x represents the parameter value obtained by least squares calculation, p represents the weight matrix, and n represents the number of points involved in the calculation.

[0076] 2) Substitute the mean error and residual calculated in step 1 into the IGGⅢ equivalent weight function model for reweighting, and then perform least squares calculation based on the new weight matrix;

[0077] Traditional methods use unit weights for calculation, directly discarding results that exceed three times the mean square error. This method breaks with convention by proposing a reweighting approach, assigning different weights to residuals of varying sizes. The weights in the least squares process of step 1 are unit weights, while subsequent steps 2 and 3 use weights calculated using the IGGⅢ weight function.

[0078] In other words, the mean square error is compared with the least squares residual. If the residual falls within the elimination region (i.e., greater than k1(2.5) times the mean square error), the weight is set to 0. If the residual falls within the doubt region (i.e., greater than k0(1.5) times and less than k1(2.5) times the mean square error), the weight is calculated using the IGGⅢ equivalent weight function. Perform weight reduction processing; if the residual falls within the normal region, that is, less than k0(1.5) times the mean error, then assign a weight of 1;

[0079] 3) If the residuals obtained from the least squares calculation are all less than k0(1.5) times the standard error, then end the iteration and proceed to step 4.

[0080] Otherwise, return to step 2) for iteration. Due to the increased precision during iteration, the residuals and mean square error will change, leading to changes in the weights of results falling within the region of doubt. This, in turn, results will differ. Furthermore, due to the increased precision, some results falling within the normal region will enter the region of doubt, causing their weights to change as well, similarly altering the results.

[0081] Step 4: Grid the estimated ice shelf elevation changes: Use the Kriging interpolation method on the scatter plot of elevation changes to generate an elevation change grid with a preset resolution.

[0082] For the scatter plot results of ice shelf elevation changes, points with an absolute elevation change greater than 10 m / yr were removed as outliers. The remaining points were then used to generate a 1 km resolution grid using Kriging interpolation. The resulting grid data is shown below. Figure 3 As shown.

[0083] In specific implementation, the method proposed in the technical solution of this invention can be automatically executed by those skilled in the art using computer software technology. System devices for implementing the method, such as computer-readable storage media storing the corresponding computer program of the technical solution of this invention and computer equipment including the computer program running the corresponding computer program, should also be within the protection scope of this invention.

[0084] In some possible embodiments, a polar surface elevation variation estimation system based on IGGⅢ robust estimation is provided, comprising the following modules:

[0085] The first module is used for satellite altimetry data extraction, including the extraction of longitude, latitude, elevation values, time, mass labels, tidal correction, and counter-pressure correction parameters;

[0086] The second module is used for satellite altimetry data processing, including extracting qualified data points based on data quality labels, performing geophysical corrections on the data, performing median filtering along the orbit based on a sliding window to remove elevation anomalies along the orbit, and separating the data into ascending and descending orbits and storing it according to the orbit number.

[0087] The third module is used to calculate elevation changes based on the IGGⅢ weight function. The implementation method is as follows: for each group of track numbers, the point with the most points in the track is extracted as the reference track, and the points on other repeated tracks are interpolated to the point with the same latitude as the reference track by inverse distance weighting.

[0088] Establish an elevation variation function model for points at the same latitude;

[0089] Assuming that for a set of repeating orbit points at the same latitude, the elevation value is linearly related to the distance, the least squares model based on the elevation change function is used, and the standard deviation of the least squares is introduced as the prior accuracy into the IGGⅢ weight function for solution, and the elevation change scatter results are obtained iteratively.

[0090] The fourth module is used to perform kriging interpolation on the scatter plot of elevation changes, generating and outputting the elevation change results.

[0091] In some possible embodiments, a polar surface elevation change estimation system based on IGGⅢ robust estimation is provided, including a processor and a memory. The memory is used to store program instructions, and the processor is used to call the stored instructions in the memory to execute the polar surface elevation change estimation method based on IGGⅢ robust estimation as described above.

[0092] In some possible embodiments, a polar surface elevation change estimation system based on IGGⅢ robust estimation is provided, including a readable storage medium storing a computer program, which, when executed, implements the polar surface elevation change estimation method based on IGGⅢ robust estimation as described above.

[0093] The specific embodiments described herein are merely illustrative of the spirit of the invention. Those skilled in the art to which this invention pertains may make various modifications or additions to the described specific embodiments or use similar methods to substitute them, without departing from the spirit of the invention or exceeding the scope defined by the appended claims.

Claims

1. A method for estimating surface elevation changes in polar regions based on IGGⅢ robust estimation, characterized in that, A set of repeating orbital points at the same latitude is selected for fitting based on seasonal variations, including the following processing: Satellite altimetry data extraction includes extracting longitude, latitude, elevation values, time, mass labels, tidal corrections, and counter-pressure correction parameters; Satellite altimetry data processing includes extracting qualified data points based on data quality labels and performing geophysical corrections on the data, including tidal corrections and inverted pressure corrections. Median filtering along the track is performed using a sliding window to remove elevation outliers; the slope is calculated based on the latitude and longitude of the ascending and descending tracks, and the data is separated into ascending and descending tracks and stored according to track number. The elevation change calculation is performed based on the IGGⅢ weight function, and the implementation method is as follows. For each set of orbit numbers, the point with the most points in the orbit is extracted as the reference orbit, and the points on other repeating orbits are interpolated to the point with the same latitude as the reference orbit using an inverse distance weighting method. For points at the same latitude, the following elevation variation function model is established based on the slope along the latitudinal direction. In the formula, This represents the elevation value at the initial moment. Indicates the annual average elevation change. and The amplitude representing seasonal variation. T represents the time measured at the nadir relative to the initial time, where T represents one year. This represents the distance between a point on the repeating track and the reference track. Indicates the slope along the latitudinal direction; Assuming that for a set of repeated orbit points at the same latitude, the elevation value is linearly related to the distance, the least squares are performed based on the elevation change function model. The standard deviation of the least squares is introduced into the IGGⅢ weight function as the prior accuracy to reweight, and then the changed weights are introduced to perform least squares solution. The recalculated residuals and standard deviations are then introduced into the IGGⅢ weight function to reweight again, and the elevation change scatter results are obtained iteratively. The iteration includes the following steps: 1) The mean square error obtained through least squares is used as the prior accuracy; 2) Substitute the mean error and residual calculated in step 1) into the IGGⅢ equivalent weight function model for reweighting, and then perform least squares calculation based on the new weight matrix; The reweighting method involves comparing the standard error with the least squares residual. If the residual falls within the elimination region, i.e., is greater than 1 / 3, then the error is considered weighted. If the mean square error is doubled, the weight is assigned to 0; if the residual falls within the region of doubt, i.e., greater than 1 / 2... times and less than The mean square error is calculated using the IGGⅢ equivalent weight function. Perform weight reduction processing; if the residual falls within the normal region, i.e., less than If the mean square error is doubled, then the weight is assigned to 1; where, , These are preset parameters. The value is 1.

5. The value is 2.

5. This represents the residual for each iteration. Represents the residual after each iteration. With mean error The ratio is used to measure the magnitude of the residual; 3) If the residuals obtained from least squares calculation are all less than If the error is not found, the iteration ends; otherwise, return to step 2) to continue the iteration. Perform kriging interpolation on the scatter plot of elevation changes to generate and output the elevation change results.

2. The method for estimating polar surface elevation changes based on IGGⅢ robust estimation according to claim 1, characterized in that: After extracting qualified data points based on data quality labels, data points with an elevation greater than 5000 meters are removed.

3. The method for estimating polar surface elevation changes based on IGGⅢ robust estimation according to claim 1, characterized in that: The method for separating the ascending and descending rails is to calculate the slope based on the latitude and longitude of the ascending and descending rails. If the slope is positive, it is determined to be a descending rail; if the slope is negative, it is determined to be an ascending rail.

4. A polar surface elevation variation estimation system based on IGGⅢ robust estimation, characterized in that: This method is used to implement the polar surface elevation variation estimation method based on IGGⅢ robust estimation as described in any one of claims 1-3.

5. The polar surface elevation change estimation system based on IGGⅢ robust estimation according to claim 4, characterized in that: Includes the following modules, The first module is used for satellite altimetry data extraction, including the extraction of longitude, latitude, elevation values, time, mass labels, tidal correction, and counter-pressure correction parameters; The second module is used for satellite altimetry data processing, including extracting qualified data points based on data quality labels, performing geophysical corrections on the data, performing median filtering along the orbit based on a sliding window to remove elevation anomalies along the orbit, and separating the data into ascending and descending orbits and storing it according to the orbit number. The third module is used to calculate elevation changes based on the IGGⅢ weight function, and is implemented as follows. For each set of orbit numbers, the point with the most points in the orbit is extracted as the reference orbit, and the points on other repeating orbits are interpolated to the point with the same latitude as the reference orbit using an inverse distance weighting method. Establish an elevation variation function model for points at the same latitude; Assuming that for a set of repeating orbit points at the same latitude, the elevation value is linearly related to the distance, the least squares model based on the elevation change function is used, and the standard deviation of the least squares is introduced as the prior accuracy into the IGGⅢ weight function for solution, and the elevation change scatter results are obtained iteratively. The fourth module is used to perform kriging interpolation on the scatter plot of elevation changes, generate elevation change results, and output them.

6. The polar surface elevation change estimation system based on IGGⅢ robust estimation according to claim 4, characterized in that: It includes a processor and a memory, the memory being used to store program instructions, and the processor being used to call the stored instructions in the memory to execute the polar surface elevation change estimation method based on IGGⅢ robust estimation as described in any one of claims 1-3.

7. The polar surface elevation change estimation system based on IGGⅢ robust estimation according to claim 4, characterized in that: The device includes a readable storage medium on which a computer program is stored, and when the computer program is executed, it implements the method for estimating polar surface elevation changes based on IGGⅢ robust estimation as described in any one of claims 1-3.

Images