CN118584488A - Method, device, storage medium and electronic device for correcting atmospheric delay error - Google Patents

Abstract

The application relates to the technical field of topographic observation, and particularly provides a method, a device, a storage medium and electronic equipment for correcting atmospheric delay errors, wherein the method can comprise the following steps: acquiring window parameters of each window in an InSAR interferogram to be corrected, wherein the window parameters comprise: terrain correlation coefficients, spatial autocorrelation coefficients, and phase gradients; comparing the window parameters with preset conditions, and determining the size of each window in the InSAR interferogram to be corrected to obtain a window corrected interferogram; and correcting the initial interference phase in the interference pattern after window correction through a phase fitting formula to obtain a target InSAR interference pattern, wherein the phase fitting formula is related to the horizontal position parameter and the elevation parameter of the observation point. Some embodiments of the application may improve the accuracy of the InSAR interferogram.

Description

Method, device, storage medium and electronic device for correcting atmospheric delay error

Technical Field

The present application relates to the field of terrain observation technology, and in particular to a method, device, storage medium and electronic device for correcting atmospheric delay errors.

Background Art

With the rapid development of interferometric synthetic aperture radar (InSAR) technology, it has become an important observation method for detecting surface deformation with high spatial resolution. At present, this technology has been successfully applied to the monitoring of surface deformation related to hydrology, volcanoes and tectonic movements. However, when using InSAR to measure small-amplitude and long-wavelength tectonic deformation signals (for example, interseismic deformation, slow slip events in subduction zones, and creep, etc.), it is seriously affected by atmospheric disturbances, which in turn affects the accuracy of subsequent surface deformation monitoring.

At present, in order to improve the accuracy of InSAR technology in the scene of surface deformation monitoring, an empirical method is usually used when facing the atmospheric delay in the InSAR interferogram. This empirical method is to determine the empirical fitting parameters according to the characteristics of atmospheric disturbances in different research areas to weaken the long-wavelength atmospheric errors related to terrain on the interferogram. However, this empirical method is highly subjective. Although it can correct errors to a certain extent, the accuracy of the interferogram cannot be guaranteed.

Therefore, how to provide a technical solution for a method of correcting atmospheric delay errors with high accuracy has become a technical problem that needs to be solved urgently.

Summary of the invention

The purpose of some embodiments of the present application is to provide a method, device, storage medium and electronic device for atmospheric delay error correction. The technical solutions of the embodiments of the present application can improve the accuracy of atmospheric delay error correction, thereby improving the accuracy of terrain observation.

In a first aspect, some embodiments of the present application provide a method for correcting an atmospheric delay error, comprising: obtaining window parameters of each window in a synthetic aperture radar InSAR interferogram to be corrected, wherein the window parameters include: a terrain correlation coefficient, a spatial autocorrelation coefficient, and a phase gradient; determining the size of each window in the InSAR interferogram to be corrected by comparing the window parameters with preset conditions, and obtaining a window-corrected interferogram; correcting an initial interference phase in the window-corrected interferogram by a phase fitting formula, and obtaining a target InSAR interferogram, wherein the phase fitting formula is related to a horizontal position parameter and an elevation parameter of an observation point.

Some embodiments of the present application compare the window parameters in the InSAR interferogram to be corrected with the preset conditions, correct the window size to obtain the window-corrected interferogram; then correct the initial interference phase in the window-corrected interferogram by the phase fitting formula to obtain the target InSAR interferogram. The embodiments of the present application can improve the accuracy and effectiveness of atmospheric delay error correction by correcting the interference phase after correcting the window size, thereby improving the accuracy of terrain observation.

In some embodiments, obtaining the window parameters of each window in the synthetic aperture radar InSAR interferogram to be corrected includes: performing fitting processing on the original InSAR interferogram to obtain the InSAR interferogram to be corrected; sliding a set initial window on the InSAR interferogram to be corrected according to a sliding step size to obtain a sliding window and the initial window in an adjacent state; and respectively calculating the window parameters corresponding to the sliding window and the initial window.

Some embodiments of the present application can facilitate obtaining corresponding window parameters by processing and sliding the original InSAR interferogram, thereby providing an effective data basis for subsequent corrections.

In some embodiments, determining the size of each window in the InSAR interferogram to be corrected by comparing the window parameters with preset conditions includes: if it is confirmed that the window parameters meet the preset conditions, merging the sliding window and the initial window to obtain a merged window size, wherein the merged window size is a window in the InSAR interferogram to be corrected; if it is confirmed that the window parameters do not meet the preset conditions, using the sliding window as a window in the InSAR interferogram to be corrected.

Some embodiments of the present application compare window parameters with preset conditions and process adjacent windows, which can improve the accuracy of window division and provide higher-precision data for subsequent corrections.

In some embodiments, confirming that the window parameters satisfy the preset conditions includes: confirming that the terrain correlation coefficients corresponding to the sliding window and the initial window are both greater than a first threshold, the phase values of the spatial autocorrelation coefficients corresponding to the sliding window and the initial window have the same sign, and the difference in phase gradients corresponding to the sliding window and the initial window is less than a second threshold.

Some embodiments of the present application compare the parameters in the window parameters with the corresponding conditions, and confirm that the preset conditions are met when both are met, so as to achieve accurate division of the window size.

In some embodiments, the method of correcting the initial interference phase in the window-corrected interference graph by a phase fitting formula to obtain a target InSAR interference graph includes: calculating the differential interference phase observation value of the observation point in each window in the window-corrected interference graph by the phase fitting formula; and correcting the initial interference phase according to the differential interference phase observation value to obtain the target InSAR interference graph.

Some embodiments of the present application determine the differential interferometric phase observation value of the observation point through a phase fitting formula, and then correct the initial interferometric phase to obtain the target InSAR interferogram, which can achieve accurate correction of the atmospheric delay error in the interferogram.

In some embodiments, the phase fitting formula is obtained by determining the estimated parameters in the phase fitting formula using the least squares method based on the observed data; wherein the observed data includes: observation point coordinates and phase difference values, and the observation point coordinates include: horizontal position parameters and elevation parameters.

In a second aspect, some embodiments of the present application provide a device for correcting an atmospheric delay error, comprising: a window parameter acquisition module, used to obtain the window parameters of each window in a synthetic aperture radar InSAR interference map to be corrected, wherein the window parameters include: terrain correlation coefficient, spatial autocorrelation coefficient and phase gradient; a window size correction module, used to determine the size of each window in the InSAR interference map to be corrected by comparing the window parameters with preset conditions, and obtain a window-corrected interference map; a phase correction module, used to correct the initial interference phase in the window-corrected interference map by a phase fitting formula to obtain a target InSAR interference map, wherein the phase fitting formula is related to the horizontal position parameters and elevation parameters of the observation point.

In a third aspect, some embodiments of the present application provide a computer-readable storage medium having a computer program stored thereon, which, when executed by a processor, can implement the method described in any embodiment of the first aspect.

In a fourth aspect, some embodiments of the present application provide an electronic device comprising a memory, a processor, and a computer program stored in the memory and executable on the processor, wherein the processor can implement a method as described in any embodiment of the first aspect when executing the program.

In a fifth aspect, some embodiments of the present application provide a computer program product, wherein the computer program product comprises a computer program, wherein the computer program, when executed by a processor, can implement the method described in any embodiment of the first aspect.

BRIEF DESCRIPTION OF THE DRAWINGS

In order to more clearly illustrate the technical solutions of some embodiments of the present application, the drawings required for use in some embodiments of the present application will be briefly introduced below. It should be understood that the following drawings only show some embodiments of the present application and therefore should not be regarded as limiting the scope. For ordinary technicians in this field, other related drawings can be obtained based on these drawings without paying creative work.

FIG1 is a flow chart of a method for correcting atmospheric delay errors provided in some embodiments of the present application;

FIG2 is a schematic diagram of window sliding provided in some embodiments of the present application;

FIG3 is a second flow chart of a method for correcting atmospheric delay errors provided in some embodiments of the present application;

FIG4 is a schematic diagram of a comparison of GPS deformation velocity fields provided by some embodiments of the present application;

FIG5 is a block diagram of an atmospheric delay error correction device provided by some embodiments of the present application;

FIG6 is a schematic diagram of an electronic device provided in some embodiments of the present application.

DETAILED DESCRIPTION

The technical solutions in some embodiments of the present application will be described below in conjunction with the drawings in some embodiments of the present application.

It should be noted that similar reference numerals and letters represent similar items in the following drawings, so once an item is defined in one drawing, it does not need to be further defined and explained in the subsequent drawings. At the same time, in the description of this application, the terms "first", "second", etc. are only used to distinguish the description and cannot be understood as indicating or implying relative importance.

In the related art, InSAR atmospheric delay is caused by the temporal and spatial variation of the refractive index of the air, which causes the interference phase to change. It changes with the pressure, temperature and water vapor content and has complex spatial characteristics. Phase delay is the combined effect of the mixture of turbulent disturbances and stratified disturbances in the atmosphere. InSAR atmospheric delay includes ionospheric delay and tropospheric delay. The ionosphere usually affects radar signals on a large scale, and it has little effect on radar signals in the X or C band. The atmospheric delay involved in this application refers to tropospheric delay.

Specifically, the atmospheric delay in the InSAR interferogram is caused by the change in atmospheric refractive index between two SAR acquisitions. It mainly includes three parts: (1) terrain-related atmospheric disturbances; (2) spatial autocorrelation disturbances; (3) spatially heterogeneous atmospheric disturbances. Atmospheric phase changes or delays can introduce significant errors of several centimeters in InSAR observations and often cover the deformation signals of interest. Atmospheric error correction is particularly important for deformation monitoring of tectonic movements on the order of millimeters. Therefore, it is of great significance to improve the monitoring ability of InSAR data for surface deformation and effectively weaken the InSAR atmospheric delay effect.

The existing InSAR atmospheric delay correction methods can be divided into two categories, namely prediction methods and empirical methods. Among them, the prediction method refers to the method of correction based on external atmospheric delay data. For example: using atmospheric delay disturbance data based on GNSS and remote sensing satellites and model-estimated tropospheric delay, as well as the GACOS (Generic Atmospheric Correction Online Service) atmospheric correction model obtained by integrating GPS, ECMWF and terrain data. The above prediction methods have been successfully applied to weaken the long-wavelength atmospheric errors related to terrain on the interferogram. The empirical method is a method of using the InSAR data itself to estimate and correct the atmospheric phase. The key to the empirical method of atmospheric delay correction is to determine the empirical fitting parameters (linear or exponential) according to the characteristics of atmospheric disturbances in different study areas. The empirical method usually uses linear estimation of the terrain-related components of the delay according to the terrain characteristics of the study area. However, based on the assumption of a single relationship between the entire interferogram, turbulence and coherent short-scale components cannot be explained. In order to consider turbulence and atmospheric disturbances at short spatial scales, the prior art proposes atmospheric delay correction methods at different spatial scales. For example, the prior art proposes to use Gaussian filters of different fixed widths to estimate terrain-related atmospheric disturbances. However, this method ignores the spatial variability of tropospheric signals. In addition, the prior art also proposes to simultaneously consider the spatial correlation power law model, divide the study area into multiple rectangular windows, and estimate the relationship between the local phase and terrain on these windows. Then, the window is weighted by using the distance from the data point and the estimated window uncertainty, and interpolated to all data points. The power law reference height h0 and the power law coefficient α are constants estimated based on balloon sounding data or weather model data. Therefore, external atmospheric model correction data is required as auxiliary data.

Moreover, the key to the empirical method lies in the determination of the filter window, and the filter window length should reflect the degree of correlation. The window length affects the results obtained using the spatiotemporal filtering method, and it can vary in space. The above empirical method for estimating atmospheric delay usually only divides the correction window according to the terrain or subjectively divides the filter window of the study area, and then uses linear or exponential functions to weaken the atmospheric disturbances related to the terrain. Specifically, it is necessary to divide the atmospheric delay correction windows of different spatial scales according to the specific spatial characteristics of the residual phase in the interferogram.

However, the correction method based on external atmospheric delay data is usually limited by the spatial and temporal resolution and the accuracy of auxiliary data. Due to the low spatial and temporal resolution and accuracy, the above external model may underestimate or overestimate the correction, especially for InSAR interferograms with large terrain changes and strong atmospheric disturbances. The above empirical method only divides the filter window according to the terrain, and is only effective for correcting terrain-related atmospheric delay errors, without taking into account the correlation between other data and terrain, as well as the spatial autocorrelation and spatial anisotropy of atmospheric delay errors. As a result, the accuracy of atmospheric delay error correction remains to be studied.

It can be seen from the above-mentioned related technologies that the accuracy of atmospheric delay error correction in the existing technology needs to be improved.

In view of this, some embodiments of the present application provide a method for correcting an atmospheric delay error, which can calculate various types of window parameters in the InSAR interferogram to be corrected, and then correct the window size through the window parameters and preset conditions to obtain the window-corrected interferogram; finally, the initial interference phase in the window-corrected interferogram is corrected through the phase fitting formula to obtain the corrected target InSAR interferogram. Some embodiments of the present application simultaneously consider the correlation between the phase and the terrain as well as the spatial autocorrelation and spatial anisotropy of the atmospheric delay error, so as to improve the accuracy of the correction, and then improve the accuracy of surface deformation monitoring in the future.

The following is an illustrative description of the implementation process of atmospheric delay error correction provided by some embodiments of the present application in conjunction with FIG1 .

Please refer to FIG. 1 , which is a flow chart of a method for correcting an atmospheric delay error provided by some embodiments of the present application. The method for correcting an atmospheric delay error may include:

S110, obtaining window parameters of each window in the synthetic aperture radar InSAR interferogram to be corrected, wherein the window parameters include: terrain correlation coefficient, spatial autocorrelation coefficient and phase gradient.

For example, in some embodiments of the present application, the terrain correlation coefficient, spatial autocorrelation coefficient and phase gradient of each window in the interference pattern are calculated, and the relationship between the characteristic parameters is comprehensively considered to improve the accuracy of atmospheric delay error correction. Among them, the spatial autocorrelation coefficient is mainly used to describe the average correlation degree of all spatial units with the surrounding areas in the entire area, which can be obtained by a commonly used calculation formula.

In some embodiments of the present application, S110 may include: performing fitting processing on the original InSAR interferogram to obtain the InSAR interferogram to be corrected; sliding a set initial window on the InSAR interferogram to be corrected according to a sliding step size to obtain a sliding window and the initial window in an adjacent state; and respectively calculating the window parameters corresponding to the sliding window and the initial window.

For example, in some embodiments of the present application, a differential interferometric phase map (as a specific example of the original InSAR interferometric map) is obtained by using D-InSAR technology based on satellite SAR observation data. Then, a second-order linear fitting (as a specific example of fitting processing) is used to eliminate the phase slope caused by the orbit error in the differential interferometric phase map. Finally, the sliding step size and the size of the initial window are set, and the initial window with the sliding step size is slid on the InSAR interferometric map to be corrected, and the window parameters of the generated adjacent windows are calculated and analyzed. For example, as shown in Figure 2 (Figure 2 is already in a state of sliding d0 to the right). The initial window is an area composed of black thick lines, whose size is a0×a0, and the sliding step size is d0. When the initial window slides to the right by d0, a sliding window of a0×d0 (i.e., the area surrounded by the dotted line in Figure 2) and an initial window of a0×a0 (i.e., the area surrounded by the black thick line in Figure 2) can be obtained. At this time, the two windows are adjacent windows. The two windows are calculated separately through the InSAR interferometric map to be corrected to obtain the window parameters corresponding to each window. It is understandable that d0 and the size of the initial window can be set according to the actual terrain conditions, and the embodiment of the present application does not make any specific limitations here.

S120, determining the size of each window in the InSAR interferogram to be corrected by comparing the window parameters with preset conditions, and obtaining the interferogram after window correction.

For example, in some embodiments of the present application, the terrain correlation coefficient, spatial autocorrelation coefficient and phase gradient are respectively compared with the corresponding preset conditions to determine the size of each corrected window, and after the area processing of the entire InSAR interferogram to be corrected is completed, the window corrected interferogram is obtained. Among them, the preset conditions may contain multiple conditions related to the window parameters, which are not specifically limited in the embodiments of the present application.

In order to facilitate the explanation of the window size correction process, the following is an exemplary explanation with reference to two adjacent windows in FIG. 2 .

In some embodiments of the present application, S120 may include: if it is confirmed that the window parameters meet the preset conditions, merging the sliding window and the initial window to obtain a merged window size, wherein the merged window size is a window in the InSAR interferogram to be corrected. Confirming that the window parameters meet the preset conditions includes: confirming that the terrain correlation coefficients corresponding to the sliding window and the initial window are both greater than a first threshold, the phase values of the spatial autocorrelation coefficients corresponding to the sliding window and the initial window have the same sign, and the difference in phase gradients corresponding to the sliding window and the initial window is less than a second threshold.

For example, in some embodiments of the present application, the preset conditions include a total of three conditions, namely: the terrain correlation coefficients of adjacent windows are all greater than the first threshold r3; the phase values of adjacent windows have the same sign; the difference in the phase values of adjacent windows is less than the second threshold r1. When the window parameters of the adjacent sliding window and the initial window in Figure 2 meet the above three conditions, they need to be corrected. The specific correction method is: merge the two to obtain a window with a merged window size of (a0+d0)×a0. It should be understood that the values of r1 and r3 can be set according to the actual application situation, and the embodiments of the present application are not specifically limited here. It should be noted that the content of the above preset conditions is an embodiment of the present application. In actual situations, the preset conditions can be adaptively adjusted based on the collected window parameters, and the embodiments of the present application are not limited to this.

After obtaining the merged window size, the initial window slides to the right again by d0. At this time, the size of the window adjacent to it on the left is 2d0×a0. The window and the initial window are compared again in the above manner. If they are merged, they continue to move to the right by d0 for analysis in the same manner, and so on.

In some other embodiments of the present application, S120 may include: if it is confirmed that the window parameter does not satisfy the preset condition, using the sliding window as a window in the InSAR interferogram to be corrected.

For example, in some embodiments of the present application, if the window parameters of the adjacent sliding window and the initial window in Figure 2 do not satisfy one of the above three conditions, no correction is required, that is, no window merging is required, and the sliding window is treated as an independent window.

After completing the first size adjustment, the initial window slides rightward again by d0, and the newly generated sliding window of a0×d0 is compared with the initial window in the above manner to determine whether to merge or form a new window separately, and so on.

After the entire area of the InSAR interferogram to be corrected is searched, a multi-scale window-corrected interferogram with the window size corrected is determined.

S130, correcting the initial interference phase in the window-corrected interference pattern by a phase fitting formula to obtain a target InSAR interference pattern, wherein the phase fitting formula is related to a horizontal position parameter and an elevation parameter of an observation point.

For example, in some embodiments of the present application, the initial interference phase of the observation point in the figure is corrected by the phase fitting formula newly proposed in the present application to obtain a corrected target InSAR interferogram. Since the phase fitting formula in the present application takes into account the overall situation of the horizontal position parameter change and the elevation parameter change at the same time during fitting, it can improve the accuracy of atmospheric delay error correction to a certain extent and has better universality.

In some embodiments of the present application, the phase fitting formula is obtained by determining the estimated parameters in the phase fitting formula based on the observation data using the least squares method; wherein the observation data includes: observation point coordinates and phase difference values, and the observation point coordinates include: horizontal position parameters and elevation parameters.

For example, in some embodiments of the present application, the phase fitting formula is:

Among them, a , b , c , c0 , C and α are the parameters to be estimated, _and x , y and h are the horizontal position parameters and elevation parameters of the observation point respectively. is the differential interferometric phase observation value (i.e., the error value caused by terrain and position-related atmospheric delay disturbances).

The parameters to be estimated are obtained by fitting the phase using the least squares method through known observation data.

In some embodiments of the present application, S130 may include: calculating the differential interference phase observation value of the observation point in each window in the window-corrected interference diagram through the phase fitting formula; correcting the initial interference phase according to the differential interference phase observation value to obtain the target InSAR interference diagram.

For example, in some embodiments of the present application, by inputting the coordinates ( x , y , h ) of the observation point of each window into the phase fitting formula, the corresponding observation point can be obtained. . The initial interference phase is The difference is taken as the corrected phase value. After each window observation point is corrected, the corrected target InSAR interferogram is obtained.

Finally, based on the target InSAR interferogram, time-series InSAR technology (e.g., SBAS algorithm) is used to estimate the deformation amount and deformation rate.

The specific process of atmospheric delay error correction provided by some embodiments of the present application is exemplarily described below with reference to FIG. 3 .

Please refer to Figure 3, which is a flow chart of a method for correcting atmospheric delay errors provided in some embodiments of the present application.

The above process is explained below as an example.

S310, obtaining original InSAR observation data.

S320, using D-InSAR technology to obtain the original InSAR interferogram.

S330, a second-order linear fitting is performed on the original InSAR interferogram to obtain the InSAR interferogram to be corrected.

S340, sliding the set initial window on the InSAR interferogram to be corrected according to the sliding step length to obtain adjacent windows in an adjacent state.

S350, respectively calculating window parameters corresponding to the sliding window and the initial window.

S360, compare the window parameters with the preset conditions to confirm whether to modify the window size, if so, execute S361, otherwise return to S370.

S361, merge adjacent windows to obtain a merged window, and execute S370.

S370, determine whether the entire area search is completed, if so, execute S380, otherwise return to S340.

Among them, after the entire area search is completed, the window-corrected interference map can be obtained.

S380, calculating the differential interference phase observation value of the observation point in each window of the interference pattern after window correction by using a phase fitting formula.

S390, correcting the initial interference phase according to the differential interference phase observation value to obtain the target InSAR interference pattern.

S391, based on the target InSAR interferogram, the time series InSAR technology is used to estimate the deformation amount and deformation rate.

It should be noted that the specific implementation process of S310~S391 can refer to the method embodiment provided above. To avoid repetition, the detailed description is appropriately omitted here.

Through some of the above embodiments of the present application, it can be seen that the present application does not require prior atmospheric or water vapor model data as a priori conditions, and does not need to divide the correction area according to experience, that is, the correction window can be automatically divided according to the spatial characteristics of terrain, interference phase and atmospheric delay. Moreover, atmospheric delay correction at multiple spatial scales can be achieved, while taking into account the spatial variation characteristics of atmospheric delay disturbances, with high accuracy and high universality.

In order to verify the feasibility of the method for correcting the atmospheric delay error provided by the present application, the present application compares the method provided by the present application with the existing atmospheric delay correction model.

Specifically, the ascending and descending orbit observation data of Sentinel-1A in a certain Haiyuan area are taken as an example, and the above three atmospheric correction methods including the present invention are compared and tested.

Since the tectonic deformation in the study area does not change significantly at a spatial scale of tens of kilometers, the sliding window size should be tens of kilometers. The optimal initial sliding window size is determined through multiple experiments. A sliding window of 24-32km is set in this study area. Considering that the spatial scale of atmospheric delay disturbance is several kilometers to more than ten kilometers, and the small spatial scale is not sensitive to tropospheric signals, the sliding step is set to about 8km. Since this experiment is aimed at the differential interferometric phase of a small space-time baseline, the time baseline is set to <180 days, and the phase threshold is set to 2rad, which is about 4.5mm. The tectonic deformation in this study area is less than 18mm/yr, so setting 2rad is reasonable to ensure that the tectonic signal will not be weakened.

Method 1 is the method proposed by the present invention, method 2 is an exponential model method of a single spatial scale, and method 3 is the GACOS atmospheric delay correction model.

The accuracy of the solution results is analyzed from three aspects: the residual phase of the interferogram, the accuracy analysis of the InSAR deformation velocity field, and the deformation time series analysis. First, the solution results of 105 pairs of interferometric images formed by 29 selected Sentinel-1A ascending orbit (T55) images show that this method can effectively weaken the atmospheric delay disturbances caused by terrain correlation and spatial autocorrelation. The above three methods are used to correct the atmospheric delay of the unwrapped interferometric phase. By correcting 105 pairs of interferometric phase images generated by 29 images of Sentinel-1A ascending orbit, Table 1 shows the statistical situation of the standard deviation of the residual phase of 105 interferometric phase images of T55 orbit and 113 interferometric phase images of T62 orbit. The average values of the interferometric phase before and after the atmospheric delay correction are obtained by three methods.

Table 1

Method 2 and method 3 are the phase averages corrected by the single scale exponential model and the GACOS atmospheric delay correction model, respectively. From the phase average results, it can be seen that the phase of the fitting result using this method is the smallest, which is better than the result without atmospheric delay correction. This method reduces the phase average calculated by method 2 by 0.2-0.3 rad. In addition, there are still vertical stratified phases of small and medium spatial scales related to terrain in the interference phases of methods 2 and 3. At the same time, there are also spatially autocorrelated atmospheric delay phases. After using this method to weaken the atmospheric delay, the average value of the ascending orbit interference phase residual is reduced from 2.1 rad to 1.66 rad, which is weakened by about 21%; the average value of the descending orbit interference phase residual is reduced from 1.79 rad to 1.04 rad, which is weakened by about 41%. However, the phase residual result after correction by method 3 is smaller than the result without correction, and is close to the phase residual of atmospheric delay related to terrain only considering large spatial scales, indicating that method 3 can correct most of the long-wave phases related to terrain, but the correction effect of atmospheric delay disturbances on small spatial scales is not good.

In order to verify the effectiveness of different atmospheric delay correction methods, the above three methods were used in the experiment to correct the atmospheric delay of the InSAR interferometric phase, and the InSAR deformation velocity field was obtained. By comparing the GPS deformation velocity field results, the effectiveness of the three types of atmospheric delay correction methods can be compared. The GPS three-dimensional deformation velocity field is projected onto the InSAR line of sight and compared with the InSAR deformation velocity field. Figure 4 shows the comparison of the GPS deformation velocity field (RMS of the residual of the InSAR deformation velocity field obtained by different atmospheric delay correction methods: (a) is the ascending orbit, and (b) is the descending orbit), and the RMS of the residual of the InSAR deformation velocity field corrected by the three methods. It can be seen from the calculation results that the RMS of the residual of the method proposed in the present invention (i.e., the 3D RMS in Figure 4) is the smallest, indicating that the correction effect of this method is the best.

Please refer to Figure 5, which shows a block diagram of the composition of the atmospheric delay error correction device provided by some embodiments of the present application. It should be understood that the atmospheric delay error correction device corresponds to the above method embodiment and can perform each step involved in the above method embodiment. The specific functions of the atmospheric delay error correction device can be found in the description above. To avoid repetition, the detailed description is appropriately omitted here.

The device for correcting the atmospheric delay error of FIG5 includes at least one software function module that can be stored in a memory in the form of software or firmware or solidified in the device for correcting the atmospheric delay error. The device for correcting the atmospheric delay error includes: a window parameter acquisition module 410, used to obtain the window parameters of each window in the synthetic aperture radar InSAR interference map to be corrected, wherein the window parameters include: terrain correlation coefficient, spatial autocorrelation coefficient and phase gradient; a window size correction module 420, used to determine the size of each window in the InSAR interference map to be corrected by comparing the window parameters with preset conditions, and obtain a window-corrected interference map; a phase correction module 430, used to correct the initial interference phase in the window-corrected interference map by a phase fitting formula to obtain a target InSAR interference map, wherein the phase fitting formula is related to the horizontal position parameter and elevation parameter of the observation point.

In some embodiments of the present application, the window parameter acquisition module 410 is used to perform fitting processing on the original InSAR interferogram to obtain the InSAR interferogram to be corrected; slide the set initial window on the InSAR interferogram to be corrected according to the sliding step size to obtain the sliding window and the initial window in an adjacent state; and calculate the window parameters corresponding to the sliding window and the initial window respectively.

In some embodiments of the present application, the window size correction module 420 is used to merge the sliding window and the initial window to obtain a merged window size if it is confirmed that the window parameters meet the preset conditions, wherein the merged window size is a window in the InSAR interferogram to be corrected; if it is confirmed that the window parameters do not meet the preset conditions, the sliding window is used as a window in the InSAR interferogram to be corrected.

In some embodiments of the present application, the window size correction module 420 is used to confirm that the terrain correlation coefficients corresponding to the sliding window and the initial window are both greater than a first threshold, the phase values of the spatial autocorrelation coefficients corresponding to the sliding window and the initial window have the same sign, and the difference in phase gradients corresponding to the sliding window and the initial window is less than a second threshold.

In some embodiments of the present application, the phase correction module 430 is used to calculate the differential interference phase observation value of the observation point in each window in the window-corrected interference diagram through the phase fitting formula; correct the initial interference phase according to the differential interference phase observation value to obtain the target InSAR interference diagram.

In some embodiments of the present application, the phase fitting formula is obtained by determining the estimated parameters in the phase fitting formula based on the observation data using the least squares method; wherein the observation data includes: observation point coordinates and phase difference values, and the observation point coordinates include: horizontal position parameters and elevation parameters.

Those skilled in the art can clearly understand that, for the convenience and brevity of description, the specific working process of the device described above can refer to the corresponding process in the aforementioned method, and will not be described in detail here.

Some embodiments of the present application further provide a computer-readable storage medium having a computer program stored thereon, which, when executed by a processor, can implement the operations of the method corresponding to any of the above methods provided in the above embodiments.

Some embodiments of the present application further provide a computer program product, which includes a computer program, wherein when the computer program is executed by a processor, it can implement the operations corresponding to any of the above methods provided in the above embodiments.

As shown in Figure 6, some embodiments of the present application provide an electronic device 500, which includes: a memory 510, a processor 520, and a computer program stored in the memory 510 and executable on the processor 520, wherein the processor 520 can implement a method as described in any of the above embodiments when reading the program from the memory 510 through a bus 530 and executing the program.

Processor 520 can process digital signals and can include various computing structures, such as complex instruction set computer structure, reduced instruction set computer structure, or a structure that implements a combination of multiple instruction sets. In some examples, processor 520 can be a microprocessor.

The memory 510 may be used to store instructions executed by the processor 520 or data related to the execution of instructions. These instructions and/or data may include codes for implementing some or all functions of one or more modules described in the embodiments of the present application. The processor 520 of the disclosed embodiment may be used to execute instructions in the memory 510 to implement the method shown above. The memory 510 includes a dynamic random access memory, a static random access memory, a flash memory, an optical memory, or other memory known to those skilled in the art.

The above description is only an embodiment of the present application and is not intended to limit the scope of protection of the present application. For those skilled in the art, the present application may have various changes and variations. Any modifications, equivalent substitutions, improvements, etc. made within the spirit and principles of the present application should be included in the scope of protection of the present application. It should be noted that similar reference numerals and letters represent similar items in the following drawings, so once an item is defined in one drawing, it does not need to be further defined and explained in the subsequent drawings.

The above is only a specific implementation of the present application, but the protection scope of the present application is not limited thereto. Any person skilled in the art who is familiar with the present technical field can easily think of changes or substitutions within the technical scope disclosed in the present application, which should be included in the protection scope of the present application. Therefore, the protection scope of the present application should be based on the protection scope of the claims.

It should be noted that, in this article, relational terms such as first and second, etc. are only used to distinguish one entity or operation from another entity or operation, and do not necessarily require or imply any such actual relationship or order between these entities or operations. Moreover, the terms "include", "comprise" or any other variants thereof are intended to cover non-exclusive inclusion, so that a process, method, article or device including a series of elements includes not only those elements, but also other elements not explicitly listed, or also includes elements inherent to such process, method, article or device. In the absence of further restrictions, the elements defined by the sentence "comprise a ..." do not exclude the presence of other identical elements in the process, method, article or device including the elements.

Claims (10)

1. A method of atmospheric delay error correction, comprising:

Acquiring window parameters of each window in an InSAR interferogram to be corrected, wherein the window parameters comprise: terrain correlation coefficients, spatial autocorrelation coefficients, and phase gradients;

comparing the window parameters with preset conditions, and determining the size of each window in the InSAR interferogram to be corrected to obtain a window corrected interferogram;

and correcting the initial interference phase in the interference pattern after window correction through a phase fitting formula to obtain a target InSAR interference pattern, wherein the phase fitting formula is related to the horizontal position parameter and the elevation parameter of the observation point.

2. The method of claim 1, wherein the obtaining window parameters for each window in the synthetic aperture radar InSAR interferogram to be corrected comprises:

Fitting the original InSAR interferogram to obtain the InSAR interferogram to be corrected;

sliding the set initial window on the InSAR interferogram to be corrected according to the sliding step length to obtain a sliding window in an adjacent state and the initial window;

and respectively calculating the window parameters corresponding to the sliding window and the initial window.

3. The method of claim 2, wherein determining the size of each window in the InSAR interferogram to be corrected by comparing the window parameters with preset conditions comprises:

If the window parameters meet the preset conditions, combining the sliding window with the initial window to obtain a combined window size, wherein the combined window size is one window in the InSAR interferogram to be corrected;

And if the window parameters are confirmed not to meet the preset conditions, taking the sliding window as one window in the InSAR interferogram to be corrected.

4. A method according to claim 3, wherein said confirming that the window parameter meets the preset condition comprises:

confirming that the topographic correlation coefficients corresponding to the sliding window and the initial window are larger than a first threshold value, the phase value signs of the spatial autocorrelation coefficients corresponding to the sliding window and the initial window are the same, and the difference value of the phase gradients corresponding to the sliding window and the initial window is smaller than a second threshold value.

5. The method of any one of claims 1-4, wherein correcting the initial interferometric phase in the window-corrected interferogram by a phase fitting formula to obtain a target InSAR interferogram comprises:

calculating a differential interference phase observation value of an observation point in each window in the interference diagram after window correction through the phase fitting formula;

and correcting the initial interference phase according to the differential interference phase observation value to obtain the target InSAR interferogram.

6. The method of any one of claims 1-4, wherein the phase fitting equation is derived from observed data after determining parameters to be estimated in the phase fitting equation using a least squares method; wherein the observation data includes: the observation point coordinates and the phase difference value, the observation point coordinates include: a horizontal position parameter and an elevation parameter.

7. An apparatus for correcting an atmospheric delay error, comprising:

The system comprises a window parameter acquisition module, a window parameter correction module and a window parameter correction module, wherein the window parameter acquisition module is used for acquiring the window parameter of each window in an InSAR interferogram to be corrected, and the window parameter comprises: terrain correlation coefficients, spatial autocorrelation coefficients, and phase gradients;

The window size correction module is used for determining the size of each window in the InSAR interferogram to be corrected by comparing the window parameters with preset conditions to obtain a window corrected interferogram;

the phase correction module is used for correcting the initial interference phase in the interference pattern after window correction through a phase fitting formula to obtain a target InSAR interference pattern, wherein the phase fitting formula is related to the horizontal position parameter and the elevation parameter of the observation point.

8. A computer readable storage medium, characterized in that the computer readable storage medium has stored thereon a computer program, wherein the computer program when run by a processor performs the method according to any of claims 1-6.

9. An electronic device comprising a memory, a processor, and a computer program stored on the memory and running on the processor, wherein the computer program when run by the processor performs the method of any one of claims 1-6.

10. A computer program product, characterized in that the computer program product comprises a computer program, wherein the computer program, when run by a processor, performs the method according to any of claims 1-6.