An InSAR atmospheric delay phase correction method weighted by ground distance

Abstract

The application discloses an InSAR atmospheric delay phase correction method along a slant range direction, comprising the following steps: S1, acquiring two SAR image subsets covering a research area and acquiring two clean water vapor data; S2, performing differential interference processing on the two SAR image subsets covering the research area to acquire differential interference phases; S3, accurately calculating atmospheric delay phases by using an atmospheric delay phase correction model; and S4, accurately obtaining real deformation of the ground along a radar line-of-sight direction. The application gives the latest weight value of the zenith path delay difference of each ground pixel on the slant range direction of the ground by correcting the sampling mode of the water vapor product along the ground distance direction, finally proposes an accurate and efficient slant range phase total delay calculation method along the ground distance, further optimizes and improves the repeat orbit synthetic aperture radar differential interferometry method, so that the measurement result is more scientific and reasonable, and therefore the application has a wide application prospect in the field of repeat orbit InSAR atmospheric phase correction.

Description

Technical Field

[0001] This invention relates to the fields of obtaining surface subsidence information in space-based Earth observation, geodesy, engineering surveying, geophysics, geological disaster monitoring and early identification, particularly in the fields of seismic coseismic deformation monitoring, investigation of the distribution of fault system slip and closure degree over time, urban land subsidence monitoring, surface subsidence monitoring in mining areas caused by underground mining, volcanic eruption movement monitoring, and iceberg and glacier movement monitoring. In particular, it relates to an InSAR atmospheric delay phase correction method with weighted sampling along the distance direction. Background Technology

[0002] Over the past few decades, Synthetic Aperture Radar (InSAR) interferometry has been successfully applied in fields such as topographic mapping, land subsidence monitoring, and seismic coseismic deformation monitoring. However, atmospheric phase delay errors severely restrict the high-precision application of InSAR technology when using InSAR interferometry to obtain surface deformation information. To mitigate the impact of atmospheric effects on the accuracy of InSAR measurements, scholars both domestically and internationally have proposed the following four methods for removing atmospheric phase delay in interferograms:

[0003] Method 1: Pair-wise logic. Massonnet first proposed pair-wise logic in 1994 when studying post-Landers earthquake deformation. By comparing interferometric results at different time intervals, pair-wise logic can qualitatively identify atmospheric signals and separate them from geophysical signals. Pair-wise logic is a fundamental method for analyzing InSAR atmospheric effects. Generally, it can only be used qualitatively to determine the existence of atmospheric effects and cannot quantitatively calculate atmospheric phase delay values.

[0004] Method 2: Interferogram stacking and relay-interferogram stacking. Interferogram stacking reduces atmospheric effects by averaging N independent interferometric phase maps generated from different SAR images of the same study area. As one of the early methods for atmospheric correction, interferogram stacking requires the assumption of a constant deformation rate, making it unsuitable for regions with nonlinear deformation rates. Furthermore, it is extremely sensitive to sudden weather events (such as large-scale rainfall or lightning) and can produce significant errors.

[0005] Method 3: Permanent Scatterers InSAR (PS InSAR). The PS InSAR method estimates geophysical signals using the phase values ​​of discrete, time-stable natural reflectors. Based on the different characteristics of components in the differential interferometric phase, such as residual topographic phase, orbital error phase, decorrelation phase, atmospheric delay phase, and noise phase, in their respective time and spatial domains (atmospheric delay phase is uncorrelated in interferograms with a time baseline greater than one day, but exhibits high spatial correlation between interferograms with a time baseline greater than one day), the total phase of each highly coherent point is analyzed. Various error phases are gradually removed from the total phase to improve the accuracy of unknown quantities such as surface deformation phase and atmospheric delay phase, thus eliminating atmospheric effects from points with high scattering characteristics in the interferogram. The PS InSAR method requires a large SAR image dataset (generally more than 30 scenes), which is not conducive to the practical application of SAR technology. Furthermore, PS InSAR technology has high requirements for the study area, which must have a sufficient number of reliable PS points.

[0006] Method 4: External Data (Water Vapor Product) Space Radiometer Monitoring. In recent years, with the improvement of water vapor product accuracy and spatial resolution, methods for removing InSAR atmospheric effects have received considerable attention due to the use of global-scale water vapor products, such as the Medium Resolution Imaging Spectrometer (MERS), the Moderate Resolution Imaging Spectroradiometer (MODIS), and the European Centre for Medium-Range Weather Forecasts (ECMWF). Related research results demonstrate that using space radiometer monitoring methods to remove atmospheric phase delay is effective. External data (water vapor product) space radiometer monitoring provides a high-resolution water vapor field for InSAR atmospheric delay phase correction and has become one of the mainstream research methods in InSAR atmospheric correction studies. However, its computational principle does not consider the combined effects of horizontal atmospheric disturbances and water vapor heterogeneity on the atmospheric delay phase, leading to a significant reduction in the accuracy of atmospheric correction results.

[0007] In summary, the space radiometer monitoring method using external data (water vapor products) has significant advantages over the other methods mentioned above. However, overcoming its shortcomings, systematically analyzing the impact of horizontal atmospheric disturbances and water vapor heterogeneity on atmospheric phase delay, comprehensively considering the ground-distance sampling direction of water vapor, and further improving the accuracy of the correction results of the space radiometer monitoring method to make the improved method feasible and universal in practice are currently hot topics and core contents in space-based Earth observation, geodesy, and geological hazard identification research. Summary of the Invention

[0008] To address the problems existing in the prior art, the present invention aims to provide an InSAR atmospheric delay phase correction method with weighted sampling along the distance direction. As an improved atmospheric delay phase correction method, it can be used in repeated orbit synthetic aperture radar interferometry to accurately calculate the true deformation of the Earth's surface.

[0009] To achieve the above objectives, the present invention adopts the following technical solution:

[0010] A weighted sampling method for InSAR atmospheric delay phase correction along the distance dimension includes the following steps:

[0011] S1. Perform multi-view and cropping processing on two synthetic aperture radar (SAR) images on the same orbit to obtain a subset of the two SAR images covering the study area. Based on the location of the study area, crop the near-infrared precipitable water vapor products, and then remove and fill in the cloud pollution pixels in the precipitable water vapor products to obtain two clean water vapor data.

[0012] S2. Perform differential interferometric processing on two subsets of SAR images covering the study area to obtain the differential interferometric phase. And the differential interference phase This includes atmospheric delay phase error;

[0013] S3. Based on the acquired water vapor data from the two images, the atmospheric delay phase correction model is used to accurately calculate the atmospheric delay phase caused by the difference in atmospheric conditions during the observation of the two SAR images. ;

[0014] S4. Remove atmospheric delay phase error from the acquired differential interferometric phase to accurately determine the true surface deformation along the radar line of sight. .

[0015] Preferably, the specific steps of step S1 are as follows:

[0016] S11. Perform multi-view processing on two synthetic aperture radar (SAR) images on the same orbit to generate a SAR image multi-view intensity map. Based on the row and column number and the view-to-number ratio of the study area in the SAR image multi-view intensity map, calculate the row and column number of the study area in the single-view complex image. Then, crop the SAR image to obtain SAR data pairs covering the study area.

[0017] S12. In the remote sensing image processing software, perform geometric correction on the near-infrared water vapor products and give the water vapor products precise geographic coordinates in the WGS84 coordinate system; then strictly follow the latitude and longitude range of the circumscribed rectangle of the study area to crop the near-infrared water vapor products.

[0018] S13. Removal of cloud-polluted pixels from water vapor products: In remote sensing image processing software, cloud mask products are first used to remove cloud-polluted pixels from water vapor products. For a small number of cloud-polluted pixels that are not completely removed, a threshold is set and band operations are used to further remove them, replacing all cloud-polluted pixels with null values.

[0019] S14. Cloud pollution pixel filling: Convert the water vapor product raster data into corresponding vector point data files, and then use the water vapor value attribute of the data file as the Z value field. Use software such as ArcGIS to perform spatial interpolation using the Kriging interpolation algorithm to obtain the subset of precipitable water vapor products after removing cloud pollution pixels.

[0020] Preferably, the specific steps of step S2 are as follows:

[0021] S21. Register the auxiliary SAR image onto the main SAR image, with a registration accuracy of at least 0.1 pixels;

[0022] S22. Resampling of auxiliary SAR image: The auxiliary SAR image is resampled according to the fine registration polynomial and phase resampling method to obtain the resampled auxiliary SAR image.

[0023] S23. Perform interferometric processing on the main SAR image and the resampled auxiliary SAR image to obtain the differential interferometric phase. .

[0024] Preferably, the specific steps of step S3 are as follows:

[0025] S31. Calculate the zenith wet delay: Based on the near-infrared water vapor product (PWV) provided by the space radiometer, convert precipitable water vapor into zenith wet delay (ZWD) using the surface temperature measurement method.

[0026]

[0027] in, It is a dimensionless conversion factor, with a value between 6.0 and 6.5; This represents the density of liquid water. and These represent the universal water vapor constant and the molar mass of water vapor, respectively. and It is the atmospheric reflectivity constant; This represents the average atmospheric temperature of the troposphere;

[0028] S32. Calculate the zenith wet delay difference plot ( ZPDDM ):

[0029]

[0030] Where t1 and t2 represent the imaging times of the two PWV tests, respectively. This represents the zenith wet delay at time t1. This indicates the zenith wet delay at time t2;

[0031] S33. Considering the influence of SAR satellite orbital azimuth angle, water vapor product sampling must be carried out along the ground slant range direction. Calculate the number of pixels in the ground slant range direction between ground pixel point B and A respectively. n ), number of east-west pixels ( j ) and the number of offset pixels in the north-south direction ( m ):

[0032]

[0033]

[0034]

[0035] in, Indicates the angle of incidence of the radar signal. Indicates the azimuth angle of the SAR satellite orbit. R This indicates the spatial resolution of precipitable water vapor products, measured in kilometers.

[0036] S34. Calculate the zenith path delay difference in the direction of ground slope distance. :

[0037]

[0038] Where gsd represents the direction of ground slant distance;

[0039] S35. Based on the vertical distribution of water vapor in the troposphere, calculate the number of pixels occupied by each vertical water vapor layer in the direction of its oblique distance from the ground:

[0040]

[0041] S36. Further, the latest weight value of each pixel on the ground in the zenith path delay difference in the ground slant range direction is given. ):

[0042]

[0043] in, h Indicates the average thickness of the troposphere;

[0044] S37. By combining the zenith path delay difference in the ground slant range direction with the weighted value of the zenith path delay difference for each pixel in the ground slant range direction, a further optimized slant range-to-atmospheric delay phase correction model is obtained:

[0045]

[0046] in, For the atmospheric delay phase, Indicates the radar signal wavelength. and These represent the path delay difference from the ground distance to the zenith at pixel i and its weight value, respectively, which are calculated by steps S34 and S36.

[0047] S38. Further calculate the true surface deformation phase of the study area after removing the atmospheric delay phase. :

[0048]

[0049] in, The differential interferometric phase in the interferogram, after removing the flat terrain effect and topographic phase, is calculated by differential interferometry software. This represents the final desired true deformation phase of the Earth's surface.

[0050] Compared with the prior art, the present invention has the following beneficial effects:

[0051] (1) This invention provides an InSAR atmospheric delay phase correction method with weighted sampling along the distance. This method only requires two SAR images and two water vapor products corresponding to the SAR images to accurately and efficiently realize the atmospheric delay phase correction process. Therefore, it has broad application prospects in the field of repeated orbit InSAR atmospheric phase correction.

[0052] (2) Compared with the traditional atmospheric correction method of space radiometer, the method proposed in this invention fully considers the influence of atmospheric disturbance and water vapor heterogeneity in the horizontal direction, gives the latest weight of each pixel on the path delay difference from the ground to the zenith, and clarifies the process of removing and filling water vapor product cloud pollution pixels. Therefore, it has higher accuracy in removing atmospheric delay phase.

[0053] (3) The method proposed in this invention corrects the sampling direction of the zenith path delay difference map in the prior art, that is, it should be along the distance from the ground rather than the east-west direction, so the result is more scientific and reasonable.

[0054] (4) The method proposed in this invention utilizes a slant range-to-atmospheric delay phase correction model, which allows for the accurate acquisition of surface deformation information for almost the entire land area of ​​the Earth using only SAR image pairs acquired by a single SAR satellite. This has significant scientific and practical value for space-based Earth observation and geological disaster monitoring and identification.

[0055] This invention corrects the sampling method of water vapor products along the ground distance direction, provides the latest weight value of the zenith path delay difference of each ground pixel in the slant distance direction, and finally proposes a method for calculating the total phase delay in the slant distance direction along the ground distance. It further optimizes and improves the differential interferometry method of repeating orbit synthetic aperture radar, making the measurement results more scientific and reasonable.

[0056] Attached image description.

[0057] To more clearly illustrate the technical solutions in the embodiments of the present invention, the accompanying drawings used in the embodiments will be briefly introduced below. Obviously, the drawings described below are only some embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0058] Figure 1 This is a flowchart of an InSAR atmospheric delay phase correction method using weighted sampling along the ground distance;

[0059] Figure 2 A schematic diagram illustrating the geometric relationship between the satellite's azimuth direction and its distance from the Earth;

[0060] Figure 3 The difference in zenith path delay between water vapor images ( ZPDDM Two-dimensional layout plan view;

[0061] Figure 4 A schematic diagram of the east-west stepped resampling of the zenith path delay differential for water vapor products;

[0062] Figure 5 This is a schematic diagram illustrating the geometric relationship between vertical tropospheric stratification and zenith wet delay and slant-range wet delay. Detailed Implementation

[0063] In the following description, specific details such as particular system architectures and techniques are set forth for illustrative purposes and not for limitation, in order to provide a thorough understanding of the embodiments of the invention. However, those skilled in the art will understand that the invention may be implemented in other embodiments without these specific details.

[0064] Reference Figure 1 A weighted sampling method for InSAR atmospheric delay phase correction along the distance direction includes the following steps:

[0065] S1. Perform multi-view and cropping processing on two synthetic aperture radar (SAR) images on the same orbit to obtain a subset of the two SAR images covering the study area. The key is to crop the near-infrared precipitable water vapor products according to the location of the study area, and then remove and fill in cloud-polluted pixels in the precipitable water vapor products to obtain two clean water vapor data. Specifically, this includes the following steps:

[0066] S11. Acquire two co-orbiting synthetic aperture radar (SAR) images, including a single-view complex SAR auxiliary image and a single-view complex SAR main image. Download digital elevation model (DEM) data and near-infrared precipitable water vapor products that fully cover the SAR images, namely precipitable water vapor 1 and precipitable water vapor 2. Use synthetic aperture radar interferometry processing software (including GAMMA, SARscape, ROI_PAC, Doris, etc.) to perform multi-view processing on the two co-orbiting SAR images to generate a SAR image multi-view intensity map. Based on the row and column number and the view-to-number ratio of the study area in the SAR image multi-view intensity map, calculate the row and column number of the study area in the single-view complex image. Then, crop the SAR image to obtain SAR data pairs covering the study area.

[0067] S12. In the remote sensing image processing software, perform geometric correction on the near-infrared water vapor products and give the water vapor products precise geographic coordinates in the WGS84 coordinate system; then strictly follow the latitude and longitude range of the circumscribed rectangle of the study area to crop the near-infrared water vapor products.

[0068] S13. Removal of cloud-polluted pixels from water vapor products: The presence of clouds causes the retrieved water vapor value to be much smaller than the true value, and very close to 0. In remote sensing image processing software, cloud mask products (such as the cloud aldebo cloud albedo band) are first used to remove as many cloud-polluted pixels as possible from the water vapor product. For the small number of cloud-polluted pixels that are not completely removed, a threshold is set and band operations are used to further remove them, replacing all cloud-polluted pixels with null values.

[0069] S14. Cloud pollution pixel filling: Convert the water vapor product raster data into the corresponding vector point data file, and then use the water vapor value attribute of the data file as the Z value field. Use software such as ArcGIS to perform spatial interpolation using the Kriging interpolation algorithm to obtain the subset of precipitable water vapor products after removing cloud pollution pixels.

[0070] S2. Perform differential interferometric processing on two subsets of SAR images covering the study area to obtain the differential interferometric phase. And the differential interference phase It includes atmospheric delay phase error; specifically, it includes the following steps:

[0071] S21. Register the auxiliary SAR image onto the main SAR image, with a registration accuracy of at least 0.1 pixels;

[0072] S22. Resampling of auxiliary SAR image: The auxiliary SAR image is resampled according to the fine registration polynomial and phase resampling method to obtain the resampled auxiliary SAR image.

[0073] S23. Perform interferometric processing on the main SAR image and the resampled auxiliary SAR image to obtain the differential interferometric phase. .

[0074] S3. Based on the acquired water vapor data from the two images, the atmospheric delay phase correction model is used to accurately calculate the atmospheric delay phase caused by the difference in atmospheric conditions during the observation of the two SAR images. Specifically, it includes the following steps:

[0075] S31. Calculate the zenith wet delay: Based on the near-infrared water vapor product (PWV) provided by the space radiometer, convert precipitable water vapor into zenith wet delay (ZWD) using the surface temperature measurement method.

[0076]

[0077] in, It is a dimensionless conversion factor, with a value between 6.0 and 6.5; This represents the density of liquid water. and These represent the universal water vapor constant and the molar mass of water vapor, respectively. and It is the atmospheric reflectivity constant; This represents the average atmospheric temperature of the troposphere;

[0078] S32. Calculate the zenith wet delay difference plot ( ZPDDM ):

[0079]

[0080] Where t1 and t2 represent the imaging times of the two PWV tests, respectively. This represents the zenith wet delay at time t1. This indicates the zenith wet delay at time t2;

[0081] S33. Considering the influence of SAR satellite orbital azimuth angle, water vapor product sampling must be carried out along the ground slant range direction. Calculate the number of pixels in the ground slant range direction between ground pixel point B and A respectively. n ), number of east-west pixels ( j ) and the number of offset pixels in the north-south direction ( m ):

[0082]

[0083]

[0084]

[0085] in, Indicates the angle of incidence of the radar signal. Indicates the azimuth angle of the SAR satellite orbit. R This indicates the spatial resolution of precipitable water vapor products, measured in kilometers.

[0086] This invention proposes a more scientific and reasonable method to correct the sampling direction of the zenith path delay differential map in the prior art, that is, it should be along the ground distance direction instead of the east-west direction. Specifically, it proposes a method for calculating the total phase delay along the slant range direction of the ground distance, which is explained in detail below:

[0087] Because SAR satellites generally do not orbit strictly along polar orbits (see...) Figure 2 The angle between its orbit and polar orbit is the satellite orbit azimuth angle α, which is generally between -167° and -170°. SAR images are imaged along the slant range-azimuth direction. The sampling direction of the zenith path delay difference should be along the slant range direction rather than the east-west direction. The SAR image data rows are not along the east-west-north-south direction, which is different from general remote sensing images (including water vapor products). Therefore, when using water vapor products for atmospheric phase delay correction, it is necessary to further consider the correspondence between the SAR data row arrangement direction and the corresponding water vapor product data row arrangement direction.

[0088] In view of this, the present invention considers using the two-dimensional structural distribution of water vapor data to present this correspondence, and displays a two-dimensional water vapor image data of m rows and j columns (see...). Figure 3 The water vapor data is uniformly distributed in a two-dimensional pattern along an east-west-north-south direction. Additionally, the zenith path delay difference calculated from water vapor images (such as MERIS)... ZPDDM Resampling is performed along the ground distance direction.

[0089] Specifically: The Earth distance is now divided into m equal segments, meaning that along the Earth distance direction, each segment contains n / m pixels. According to the nearest neighbor method, the path delay difference from the Earth distance to the zenith can be approximated by the path delay difference of its nearest east-west zenith. That is, the total path delay along the slant distance of the Earth distance is approximated by the step-like east-west path delays of m segments (see...). Figure 4 ).

[0090] S34. Create a piecewise function to calculate the zenith path delay difference in the ground slope direction. ):

[0091]

[0092] Where gsd represents the ground slant range direction; n, j, and m represent the number of ground slant range pixels, east-west pixels, and north-south offset pixels between ground point B and A, respectively;

[0093] According to the water vapor distribution in the troposphere published by China Weather Network, "48.7% of the total water vapor content is below 850 hPa, 77.5% is below 700 hPa, and 92.5% is below 550 hPa, meaning that water vapor in the lower troposphere accounts for more than 3 / 4 of the total atmospheric water vapor content." This invention further updates the vertical stratification pattern of water vapor in the troposphere as shown in Table 1 below, namely, 48.7% of the water vapor content is concentrated in areas below 1.5 km altitude; 28.8% of the water vapor is concentrated between 1.5 km and 3 km altitude; 15% of the water vapor is concentrated between 3 km and 4.5 km altitude; and the remaining 7.5% of the water vapor is concentrated between 4.5 km and 12 km altitude.

[0094] Table 1. Vertical stratification patterns of tropospheric water vapor

[0095]

[0096] The average thicknesses of the four vertical water vapor layers are 1.5 km, 1.5 km, 1.5 km and 7.5 km, respectively.

[0097] Therefore, S35, based on the vertical distribution of water vapor in the troposphere, calculate the number of pixels occupied by each vertical water vapor layer in the direction of the slant distance from the ground:

[0098]

[0099] S36, Reference Figure 5 Furthermore, the latest weight values ​​for the zenith path delay difference of each pixel on the ground in the ground slant range direction are given. ):

[0100]

[0101] in, h Indicates the average thickness of the troposphere;

[0102] S37, By combining the zenith path delay difference in the ground slant range direction ( ) and the weight value of each cell in the zenith path delay in the ground slant range direction ( ), thus obtaining a further optimized slant range-direction atmospheric delay phase correction model:

[0103]

[0104] in, For the atmospheric delay phase, Indicates the radar signal wavelength. and These represent the path delay difference from the ground distance to the zenith at pixel i and its weight value, respectively, which are calculated by steps S34 and S36.

[0105] S38. Further calculate the true surface deformation phase of the study area after removing the atmospheric delay phase. :

[0106]

[0107] in, The differential interferometric phase in the interferogram, after removing the flat terrain effect and topographic phase, is calculated by differential interferometry software. This represents the final desired true deformation phase of the Earth's surface.

[0108] S4. Remove atmospheric delay phase error from the acquired differential interferometric phase to accurately determine the true surface deformation along the radar line of sight. .

[0109] This invention is not limited to the specific embodiments described above. Any modifications made by those skilled in the art based on the above concept without creative effort are within the scope of protection of this invention.

Claims

1. A weighted sampling method for InSAR atmospheric delay phase correction along the ground distance, characterized in that, Includes the following steps: S1. Perform multi-view and cropping processing on two synthetic aperture radar (SAR) images on the same orbit to obtain a subset of the two SAR images covering the study area. Based on the location of the study area, crop the near-infrared precipitable water vapor products, and then remove and fill in the cloud pollution pixels in the precipitable water vapor products to obtain two clean water vapor data. S2. Perform differential interferometric processing on two subsets of SAR images covering the study area to obtain the differential interferometric phase. And the differential interference phase This includes atmospheric delay phase error; S3. Based on the acquired water vapor data from the two images, the atmospheric delay phase correction model is used to accurately calculate the atmospheric delay phase caused by the difference in atmospheric conditions during the observation of the two SAR images. ; Specifically, it includes the following steps: S31. Calculate the zenith wet delay: Based on the near-infrared water vapor product (PWV) provided by the space radiometer, convert precipitable water vapor into zenith wet delay (ZWD) using the surface temperature measurement method. in, It is a dimensionless conversion factor, with a value between 6.0 and 6.5; This represents the density of liquid water. and These represent the universal water vapor constant and the molar mass of water vapor, respectively. and It is the atmospheric reflectivity constant; This represents the average atmospheric temperature of the troposphere; S32. Calculate the zenith wet delay difference plot. ZPDDM : in, and These represent the imaging times of the two PWV tests, respectively. express The zenith wet delay time express Moisture delay at the zenith; S33. Considering the influence of SAR satellite orbital azimuth angle, water vapor product sampling must be carried out along the ground slant range direction. Calculate the number of pixels along the ground slant range direction between ground pixel point B and A respectively. n Number of east-west pixels j and the number of offset pixels in the north-south direction m : in, Indicates the angle of incidence of the radar signal. Indicates the azimuth angle of the SAR satellite orbit. R This indicates the spatial resolution of precipitable water vapor products, measured in kilometers. S34. Calculate the zenith path delay difference in the direction of ground slope distance. : Where gsd represents the direction of ground slant distance; S35. Based on the vertical distribution of water vapor in the troposphere, calculate the number of pixels occupied by each vertical water vapor layer in the direction of its oblique distance from the ground: ； S36. Further, the latest weight values ​​for the zenith path delay difference of each pixel on the ground in the ground slant range direction are given. : in, h Indicates the average thickness of the troposphere; S37. By combining the zenith path delay difference in the ground slant range direction with the weighted value of the zenith path delay difference for each pixel in the ground slant range direction, a further optimized slant range-to-atmospheric delay phase correction model is obtained: in, For the atmospheric delay phase, Indicates the radar signal wavelength. and These represent the path delay difference from the ground distance to the zenith at pixel i and its weight value, respectively, which are calculated by steps S34 and S36. S38. Further calculate the true surface deformation phase of the study area after removing the atmospheric delay phase. : in, The differential interferometric phase in the interferogram, after removing the flat terrain effect and topographic phase, is calculated by differential interferometry software. This represents the final desired true surface deformation phase; S4. Remove atmospheric delay phase error from the acquired differential interferometric phase to accurately determine the true surface deformation along the radar line of sight. .

2. The InSAR atmospheric delay phase correction method for weighted sampling along the distance direction according to claim 1, characterized in that, The specific steps of step S1 are as follows: S11. Perform multi-view processing on two synthetic aperture radar (SAR) images on the same orbit to generate a SAR image multi-view intensity map. Based on the row and column number and the view-to-number ratio of the study area in the SAR image multi-view intensity map, calculate the row and column number of the study area in the single-view complex image. Then, crop the SAR image to obtain SAR data pairs covering the study area. S12. In the remote sensing image processing software, perform geometric correction on the near-infrared water vapor products and give the water vapor products precise geographic coordinates in the WGS84 coordinate system; then strictly follow the latitude and longitude range of the circumscribed rectangle of the study area to crop the near-infrared water vapor products. S13. Removal of cloud-polluted pixels from water vapor products: In remote sensing image processing software, cloud mask products are first used to remove cloud-polluted pixels from water vapor products. For a small number of cloud-polluted pixels that are not completely removed, a threshold is set and band operations are used to further remove them, replacing all cloud-polluted pixels with null values. S14. Cloud pollution pixel filling: Convert the water vapor product raster data into corresponding vector point data files, then use the water vapor value attribute of the data file as the Z value field, and use the Kriging interpolation algorithm in ArcGIS software to perform spatial interpolation to obtain the subset of precipitable water vapor products after removing cloud pollution pixels.

3. The InSAR atmospheric delay phase correction method for weighted sampling along the distance direction according to claim 2, characterized in that, The specific steps of step S2 are as follows: S21. Register the auxiliary SAR image onto the main SAR image, with a registration accuracy of at least 0.1 pixels; S22. Resampling of auxiliary SAR image: The auxiliary SAR image is resampled according to the fine registration polynomial and phase resampling method to obtain the resampled auxiliary SAR image. S23. Perform interferometric processing on the main SAR image and the resampled auxiliary SAR image to obtain the differential interferometric phase. .

Images