Long-distance high-precision deformation measurement method based on corner reflector SAR (Synthetic Aperture Radar) phase information

Abstract

The invention belongs to the field of synthetic aperture radars, and discloses a long-distance high-precision deformation measurement method based on corner reflector SAR phase information, which comprises the following steps: S1, filtering an SAR image by adopting an improved Lee filtering algorithm to obtain a filtered SAR image; s2, performing registration processing on the filtered SAR image; s3, obtaining the position of each corner reflector in the SAR image after registration; s4, estimating the elevation of the corner reflector; s5, compensating the flat ground phase and the elevation phase to obtain an initial deformation phase; s6, unwrapping the initial deformation phase to obtain a real deformation phase of each corner reflector; and S7, the deformation quantity is obtained based on the real deformation phase. According to the method, the existing Lee filtering algorithm is improved, and the SAR image is filtered based on the adaptive filtering value, so that effective filtering can be carried out while excessive smoothness is avoided, and a more accurate deformation measurement result can be obtained.

Description

Technical Field

[0001] The invention relates to the field of synthetic aperture radar, and in particular to a long-distance high-precision deformation measurement method based on corner reflector SAR phase information. Background Art

[0002] SAR, or synthetic aperture radar, is an active earth observation system that can be installed on aircraft, satellites, spacecraft and other flying platforms to observe the earth around the clock and in all weather conditions, and has a certain surface penetration capability. Therefore, the SAR system has unique advantages in disaster monitoring, environmental monitoring, ocean monitoring, resource exploration, crop yield estimation, mapping and military applications, and can play a role that other remote sensing methods cannot play.

[0003] In the process of measuring the deformation of the device based on the SAR image, it is usually necessary to filter the SAR image first. In the prior art, the Lee filter algorithm is generally used for filtering. However, the Lee filter algorithm needs to specify the size of the filter window in advance during the noise reduction process. When the filter window is too large, it is easy to cause over-smoothing of the image, resulting in loss or blurring of details; when the filter window is too small, it is highly affected by local noise, resulting in poor noise reduction effect. Summary of the invention

[0004] The purpose of the present invention is to disclose a long-distance high-precision deformation measurement method based on corner reflector SAR phase information to solve the technical problems pointed out in the background technology.

[0005] In order to achieve the above object, the present invention adopts the following technical solution:

[0006] The present invention provides a long-distance high-precision deformation measurement method based on corner reflector SAR phase information, comprising:

[0007] S1, filtering the SAR image using an improved Lee filtering algorithm to obtain a filtered SAR image;

[0008] S2, performing registration processing on the filtered SAR image;

[0009] S3, obtaining the position of each corner reflector in the registered SAR image;

[0010] S4, estimate the elevation of the corner reflector;

[0011] S5, compensating the flat ground phase and the elevation phase to obtain the initial deformation phase;

[0012] S6, unwrapping the initial deformation phase to obtain the real deformation phase of each corner reflector;

[0013] S7, obtaining the deformation amount based on the real deformation phase;

[0014] Among them, S1 includes:

[0015] Each pixel in the SAR image is filtered in the following way to obtain a filtered SAR image:

[0016] For pixel point p in the SAR image, calculate the adaptive filtering value of p;

[0017] Filter p based on the adaptive filter value to obtain the pixel value of p after filtering.

[0018] Preferably, the calculation formula of the adaptive filtering value is:

[0019]

[0020] adpw p represents the adaptive filter value, pres represents the preset logarithmic parameter, wp represents the set of pixels in a 5×5 window centered on p, nwp represents the total number of pixels in wp, pxl i Represents the pixel value of pixel i, pxl p Represents the pixel value of pixel point p.

[0021] Preferably, S2 includes:

[0022] In the registration process, one SAR image is selected from multiple SAR images as a reference image, and the remaining images are registered to the grid where the reference image is located.

[0023] Preferably, S3 includes:

[0024] Performing interpolation processing on the registered SAR image to obtain an interpolated image;

[0025] Get the position of each corner reflector in the interpolated image.

[0026] Preferably, S4 includes:

[0027] Generate an interferogram based on two SAR images taken adjacently in time;

[0028] The correlation formula between the phase difference ΔΦ and the elevation in the interference pattern is obtained:

[0029]

[0030] h represents the elevation, λ is the wavelength of the SAR signal, and θ is the angle between the SAR signal and the ground;

[0031] Calculate h:

[0032]

[0033] Preferably, S5 includes:

[0034] The flat ground phase difference is calculated using the following formula:

[0035] ΔΦflat=Φactual-Φflat

[0036] ΔΦflat is the flat ground phase difference; Φactual is the actual phase of the corner reflector; Φflat is the theoretical flat ground phase of the corner reflector when the ground is completely flat;

[0037] Calculate the elevation phase difference:

[0038]

[0039] ΔΦheight is the elevation phase difference; hcorner is the actual elevation of the corner reflector; hflat is the assumed elevation of an ideal flat ground; λ is the wavelength of the SAR signal;

[0040] Calculate the atmospheric delay phase:

[0041]

[0042] represents the atmospheric delay phase, Δd represents the atmospheric delay from the zenith;

[0043] The initial deformation phase is obtained by subtracting the ground phase difference, elevation phase difference and atmospheric delay phase from the interference phase of the corner reflector.

[0044] Preferably, S6 includes:

[0045]

[0046] n is an integer used to represent the integer ambiguity, represents the true deformation phase, represents the initial deformation phase.

[0047] Preferably, S7 includes:

[0048]

[0049] wd represents the deformation amount.

[0050] Beneficial effects:

[0051] Compared with the prior art, the present invention improves the existing Lee filtering algorithm and filters the SAR image based on an adaptive filtering value. It can perform effective filtering while avoiding over-smoothing. It can not only effectively reduce the probability of events in which the filtering window is too large or too small, but also effectively improve the quality of filtering, which is conducive to obtaining more accurate deformation measurement results. BRIEF DESCRIPTION OF THE DRAWINGS

[0052] In order to more clearly illustrate the technical solutions of the embodiments of the present invention, the drawings required for use in the embodiments are briefly introduced below. It should be understood that the following drawings only show certain embodiments of the present invention 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 creative work.

[0053] Figure 1 The present invention is a schematic diagram of a long-distance high-precision deformation measurement method based on corner reflector SAR phase information.

[0054] Figure 2 This is the optical image distribution diagram of the corner reflector.

[0055] Figure 3 Schematic diagram of corner reflector signal in SAR image.

[0056] Figure 4 This is the distribution map of corner reflector SAR images.

[0057] Figure 5 This is the distribution map of SAR images of power tower corner reflectors.

[0058] Figure 6 This is a schematic diagram of GACOS atmospheric data.

[0059] Figure 7 It is a Delaunay triangulation diagram.

[0060] Figure 8a Schematic diagram of the time-series deformation of the corner reflector of Tower 1.

[0061] Figure 8b Schematic diagram of the time-series deformation of the corner reflector of Tower 2.

[0062] Figure 8c Schematic diagram of the time-series deformation of the corner reflector of Tower 3. DETAILED DESCRIPTION

[0063] The technical solutions in the embodiments of the present invention will be clearly and completely described below in conjunction with the drawings in the embodiments of the present invention. Obviously, the described embodiments are only a part of the embodiments of the present invention, rather than all the embodiments. The components of the embodiments of the present invention generally described and shown in the drawings here can be arranged and designed in various different configurations. Therefore, the following detailed description of the embodiments of the present invention provided in the drawings is not intended to limit the scope of the claimed invention, but merely represents the selected embodiments of the present invention. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without making creative work belong to the scope of protection of the present invention.

[0064] like Figure 1 In one embodiment shown, the present invention provides a long-distance high-precision deformation measurement method based on corner reflector SAR phase information, comprising:

[0065] S1, filtering the SAR image using an improved Lee filtering algorithm to obtain a filtered SAR image;

[0066] S2, performing registration processing on the filtered SAR image;

[0067] S3, obtaining the position of each corner reflector in the registered SAR image;

[0068] S4, estimate the elevation of the corner reflector;

[0069] S5, compensating the flat ground phase and the elevation phase to obtain the initial deformation phase;

[0070] S6, unwrapping the initial deformation phase to obtain the real deformation phase of each corner reflector;

[0071] S7, obtaining the deformation amount based on the real deformation phase;

[0072] Among them, S1 includes:

[0073] Each pixel in the SAR image is filtered in the following way to obtain a filtered SAR image:

[0074] For pixel point p in the SAR image, calculate the adaptive filtering value of p;

[0075] Filter p based on the adaptive filter value to obtain the pixel value of p after filtering.

[0076] The present invention improves the existing Lee filtering algorithm and filters the SAR image based on an adaptive filtering value. It can perform effective filtering while avoiding over-smoothing. It can not only effectively reduce the probability of events in which the filtering window is too large or too small, but also effectively improve the quality of filtering, which is conducive to obtaining more accurate deformation measurement results.

[0077] Preferably, the calculation formula of the adaptive filtering value is:

[0078]

[0079] adpw p represents the adaptive filter value, pres represents the preset logarithmic parameter, wp represents the set of pixels in a 5×5 window centered on p, nwp represents the total number of pixels in wp, pxl i Represents the pixel value of pixel i, pxl p Represents the pixel value of pixel point p.

[0080] The adaptive filtering value of the present invention is calculated from the ratio between the pixel value of a pixel point and the variance of the pixel values ​​of pixels in the neighborhood. When the local signal-to-noise ratio of the neighborhood is higher and the image noise level is lower, the adaptive filtering value is larger.

[0081] Specifically, the pixel value of a pixel point may be an amplitude.

[0082] Specifically, the set logarithmic parameter may be 10.

[0083] Preferably, filtering p based on the adaptive filtering value to obtain the pixel value of p after filtering includes:

[0084] Using apxl p Represents the pixel value of p after filtering p, apxl p The calculation formula is:

[0085]

[0086] pxl ave represents the average value of the pixel values ​​in sp, σ represents the noise standard deviation, u represents the noise mean, and sp represents the side length L centered on p p The set of pixels within the range, pxl i is the pixel value of pixel i, pxl p The pixel value nsp representing the pixel point p represents the number of pixels in sp;

[0087]

[0088] adpwm represents the maximum value of the adaptive filtering value of the pixel points in the SAR image, and Lstd represents the preset length.

[0089] When filtering p, the size of the filter window is calculated based on the adaptive filter value. Therefore, when the adaptive filter value is larger, the length of the filter window is smaller, and when the adaptive filter value is smaller, the length of the filter window is larger, thereby realizing the adaptive change of the length of the filter window. When the noise level is higher, a larger window can be used to effectively reduce noise; and when the noise level is lower, the filtering efficiency can be improved by reducing the filter window.

[0090] Specifically, the preset length is 10.

[0091] Preferably, S2 includes:

[0092] In the registration process, one SAR image is selected from multiple SAR images as a reference image, and the remaining images are registered to the grid where the reference image is located.

[0093] Preferably, S3 includes:

[0094] Performing interpolation processing on the registered SAR image to obtain an interpolated image;

[0095] Get the position of each corner reflector in the interpolated image.

[0096] Preferably, S4 includes:

[0097] Generate an interferogram based on two SAR images taken adjacently in time;

[0098] The correlation formula between the phase difference ΔΦ and the elevation in the interference pattern is obtained:

[0099]

[0100] h represents the elevation, λ is the wavelength of the SAR signal, and θ is the angle between the SAR signal and the ground;

[0101] Calculate h:

[0102]

[0103] InSAR technology uses interference phase information to invert the deformation characteristics of the target. However, according to the InSAR principle, the interference phase not only includes the phase introduced by the target deformation, but also includes the phase introduced by various other factors. In order to achieve high-precision deformation measurement, the key is to accurately compensate for the phase introduced by various other factors.

[0104] Preferably, S5 includes:

[0105] The flat ground phase difference is calculated using the following formula:

[0106] ΔΦflat=Φactual-Φflat

[0107] ΔΦflat is the flat ground phase difference; Φactual is the actual phase of the corner reflector; Φflat is the theoretical flat ground phase of the corner reflector when the ground is completely flat;

[0108] Calculate the elevation phase difference:

[0109]

[0110] ΔΦheight is the elevation phase difference; hcorner is the actual elevation of the corner reflector; hflat is the assumed elevation of an ideal flat ground; λ is the wavelength of the SAR signal;

[0111] When there is a lack of InSAR data processing experience or the amount of SAR data is small, it is difficult to accurately estimate the atmospheric delay phase using the conventional linear atmospheric phase simulation method. Although this method can reduce the standard deviation of the interferogram, it cannot significantly improve the deformation time series. The deformation measurement method based on the corner reflector SAR phase signal often does not have a large amount of SAR data. Directly using the GACOS external atmosphere model can better weaken the influence of atmospheric delay to a certain extent, thereby obtaining the surface deformation phase.

[0112] Calculate the atmospheric delay phase:

[0113]

[0114] represents the atmospheric delay phase, Δd represents the atmospheric delay from the zenith;

[0115] Since GACOS is the result of zenith distance measurement, we first need to geocode the GACOS data based on the satellite orbit and the intensity image simulated by DEM. Then, we perform temporal differential processing on the GACOS data according to the baseline interferometric set generated earlier to obtain the atmospheric delay differential pair at the same time as the interferometric image pair. During spatial differential, the set reference point is regarded as the atmospheric delay spatial differential point.

[0116] The initial deformation phase is obtained by subtracting the ground phase difference, elevation phase difference and atmospheric delay phase from the interference phase of the corner reflector.

[0117] In the present invention, an interference pair can be first constructed based on temporally adjacent SAR images. In adjacent time periods, the deformation of the corner reflector is relatively small, resulting in a relatively small deformation phase, and the elevation can be estimated more accurately based on the interference phase. Subsequently, an interference pattern is reconstructed with a certain image as the main image and other images as secondary images. Based on the high-precision row and column position data and high-precision elevation data h of the corner reflector in the main image, the spatial position coordinates of each corner reflector under WGS84 are calculated in combination with the Doppler equation, the slant range equation and the ellipsoid equation. Since the position data and elevation data of the corner reflector are very accurate, the accuracy of the calculated reference phase is high, which can well compensate for the two phase terms in the interference phase.

[0118] Preferably, S6 includes:

[0119] After the above steps, after compensating the flat ground phase and elevation phase, the deformation phase is obtained. Wrapped in the interval [-π, π). In order to restore the true deformation phase of the target point The phase data needs to be unwrapped:

[0120]

[0121] n is an integer used to represent the integer ambiguity, represents the true deformation phase, represents the initial deformation phase.

[0122] In the actual data error optimization process, the reference corner reflector is connected to the other corner reflectors through the Delaunay triangulation. However, the distance between individual corner reflectors is far, and the unwrapped deformation phase may contain abnormal values.

[0123] The use of least squares estimation in the error optimization process is easily affected by outliers and may produce a certain degree of deviation from the actual deformation phase. In order to improve the reliability of the relative deformation estimation of the corner reflector, the present invention uses a robust M-estimator to optimize the deformation phase change between each arc segment. By using the idea of ​​weighted least squares method, a lower weight is given to the points with larger residuals during the iteration process, and the weights of the points with normal residuals remain unchanged. This can suppress the influence of outliers on the optimization results to a certain extent, making the final estimated value more stable and reliable.

[0124] The optimization process is as follows:

[0125] a) First, the deformation phase between the corner reflectors can be expressed as follows:

[0126] Δφ=Vα+e (4)

[0127] Where V is the correction factor, α is the unknown parameter, and Δφ is the actual observation value.

[0128] b) In robust estimation, formula (4) can be converted into a least squares solution problem, where l represents the number of observations, and the expression is shown in formula (5). As a traditional robust estimation method, when the observations follow a normal distribution, the least squares estimation is the optimal linear unbiased estimation with the best statistical properties. However, when the observations contain outliers, the assumption that the observations come from a specific normal distribution does not hold, and the least squares estimation naturally no longer has the properties of the optimal estimation, and the parameter estimation is often destroyed by a few outliers.

[0129]

[0130] In actual situations, the observed data does not follow a specific distribution. Least squares is only an ideal estimation method. Facing unknown observed data distribution patterns, reliable results cannot be obtained. In view of the above problems, the robust estimation theory can make the overall results as free from interference as possible when outliers are inevitable, and obtain the optimal or near-optimal deformation observation results in actual situations.

[0131] c) There are many estimation methods for M estimation. The easiest one to implement is based on the weighted iteration method, that is, the iterative reweighted least squares algorithm is used to solve the position parameters. In this case, equation (5) can be expressed as where e l =Δφ l -v l α, the specific formula is as follows:

[0132] α i+1 =(V T W i V) (-1) V T W i Δφ (6)

[0133] α i+1 It represents the iterative solution, i represents the number of iterations. It is not difficult to find that when the number of iterations is 0, formula (6) is the same as the most original unweighted least squares algorithm. Set the initial weight matrix to the unit matrix:

[0134] W i＝0 =I M (7)

[0135] Then, the unknown parameter α is obtained by solving i+1 , calculate the residual of each observation value, and use the HuBer loss function to update the weights of different observation results during the iteration process, as shown in formula (8)

[0136]

[0137] In the formula,

[0138] e i is the residual value after the i-th iteration, and k is a pre-given threshold value used to limit the normal range of residual variation. Residuals greater than the threshold are considered outliers, and the absolute value loss is used to reduce the weight of outliers during estimation; otherwise, they are considered normal values, and the square loss is used to obtain an unbiased estimate of the parameters.

[0139] Finally, the iterative reweighted least squares process is performed until the residual value meets the set requirements to obtain the final optimization result.

[0140] Usually k can be set to 2σ. In this case, taking the first-order derivative of equation (8), we get:

[0141]

[0142] The corresponding weight function can be solved:

[0143]

[0144] It can be observed from the Huber weight function that when all residual values ​​are between -2σ and 2σ, the Huber estimate is the same as the classical least squares estimate. However, when the absolute value of the residual value of some observations is greater than 2σ, the corresponding weight will decrease, which means that these observations have less influence on the parameter estimate. When the absolute value of the residual value of the observation is much greater than 2σ, the corresponding weight is close to zero, and the influence of these observations on the parameter estimate is almost negligible.

[0145] In summary, the use of Huber loss function when updating weights has a certain weakening effect on outliers, which enables it to better adapt to the overall distribution of data and provide a more robust model for the optimization of corner reflector shape.

[0146] Preferably, S7 includes:

[0147]

[0148] wd represents the deformation amount.

[0149] For X-band data, when the wavelength λ is 3.1 cm, in order to achieve a deformation measurement accuracy of 1 mm, The accuracy must reach 0.4rad.

[0150] Based on high-resolution time-series SAR data, the present invention firstly accurately extracts the position data and phase data of the corner reflector in each image. Subsequently, the SAR images adjacent in time are used to form an interference pair, and the elevation data of the corner reflector is more accurately estimated based on the known rough elevation of the corner reflector. Then, a SAR image is selected as the main image to construct the interference pair, and the accurate position and elevation data of the corner reflector are used to accurately compensate the flat ground phase and the elevation phase. On this basis, the atmospheric delay phase is compensated by the external atmospheric delay model GACOS, and after obtaining the high-precision deformation phase, the deformation phase of each corner reflector after unwrapping is obtained by phase unwrapping. On this basis, the M-estimator proposed by Huber is used to optimize the error of the observation value, so as to improve the stability of the observation result. Finally, deformation inversion processing is performed to extract the deformation information of each corner reflector. This method uses remote sensing non-contact measurement technology, and the deformation measurement accuracy reaches 1mm. It can accurately obtain the deformation data of key points, and provide efficient technical support for the key point deformation monitoring of infrastructure such as power towers, reservoir dams, and bridges. The present invention improves the applicability of corner reflector deformation monitoring technology in long-interval monitoring equipment scenarios and provides technical support for the engineering application of the technology.

[0151] A more specific embodiment is described below:

[0152] In this implementation case, 12 corner reflectors are installed on three power towers. In addition, a corner reflector is set around as a reference point to measure the relative deformation of the power tower. The distribution of the corner reflectors is as follows: Figure 2 In the experiment, 3m resolution COSMO-SkyMed data was used to monitor the deformation of the corner reflector. The SAR data parameters are shown in Table 1.

[0153] Table 1 Basic information of COSMO-SkyMed data used in the experiment

[0154] parameter Numeric Satellite Type COSMO-SkyMed track Lowering orbit Spatial resolution 3m×3m Polarization HH Number of images 10 Monitoring start time 20230814 Monitoring end time 20231227

[0155] Step 1:

[0156] First, Sinc interpolation is performed on the signals around the corner reflector, and the interpolation multiples in the range and azimuth directions are set to 100. The interpolation of the GCP point in the SAR image is shown in the figure below: Figure 3 As shown, the precise position coordinates of the corner reflector are then found based on each interpolated image. Figure 4 is the distribution map of corner reflector SAR images, Figure 5The precise position distribution of the corner reflector installed on tower 1 in the SAR image. Since the satellite adopts the descending right-looking shooting method, P1 is blocked and the signal is weak. There are three corner reflectors monitored this time, namely P2, P3, and P4. P10 and P11 are monitored on tower 2, and P5, P6, and P7 are monitored on tower 3. The corner reflectors not monitored above do not participate in the deformation solution.

[0157] Step 2:

[0158] Based on the 10 images in the case, the images of two adjacent dates were combined into interferometric pairs to obtain the interferometric phase and make a high-precision estimation of the elevation phase of the corner reflector. Then, the 20231102 image was selected as the main image to construct the interferometric pair, and the interferometric phase of the corner reflector was obtained. The elevation and flat ground phases were compensated using high-precision elevation data and corner reflector position data.

[0159] GACOS data based on weather models has a unique advantage in tropospheric delay correction due to its high temporal and spatial resolution and fusion of GNSS observation data, which can directly estimate the delay phase. First, the external atmospheric correction data GACOS is obtained according to the date and accurate time of image shooting. Figure 6 It is the GACOS atmospheric data on August 14, 2023. The GACOS data is processed in time according to the baseline interferometric set to obtain the atmospheric delay differential pair at the same time as the interferometric image pair. The set reference point is regarded as the atmospheric delay spatial differential point during the spatial difference, and the zenith distance atmospheric delay data is converted into the atmospheric delay phase in the InSAR processing process. Finally, the deformation phase of all corner reflector points is obtained, and the minimum cost flow method is used based on the reference point for phase unwrapping to obtain the deformation phase of each corner reflector relative to the reference corner reflector.

[0160] Step 3:

[0161] Connect each corner reflector through the Delaunay triangulation, such as Figure 7 Then, based on the unwrapped deformation phase, the M-estimator proposed by Huber was used to optimize the error of the observation value to improve the stability of the observation results. Finally, the deformation results of the corner reflectors that can be monitored on the three towers were obtained, as shown in Figure 8. From the 8 corner reflectors monitored, the maximum cumulative settlement of Tower 1, Tower 2 and Tower 3 from August 2023 to December 2023 were 13.2mm, 20.9mm and 6.8mm respectively. The overall situation is relatively stable. In the future, it is still necessary to pay attention to whether the bases of the towers have uneven settlement. Figure 8a For the electric tower, Figure 8b For the second tower, Figure 8c It is the third electric tower.

[0162] The above is only a specific embodiment of the present invention, but the protection scope of the present invention is not limited thereto. Any person skilled in the art can easily think of changes or substitutions within the technical scope disclosed by the present invention, which should be included in the protection scope of the present invention. Therefore, the protection scope of the present invention should be based on the protection scope of the claims.

Claims

1. A long-distance high-precision deformation measurement method based on corner reflector SAR phase information, characterized in that: include: S1, filtering the SAR image using an improved Lee filtering algorithm to obtain a filtered SAR image; S2, performing registration processing on the filtered SAR image; S3, obtaining the position of each corner reflector in the registered SAR image; S4, estimate the elevation of the corner reflector; S5, compensating the flat ground phase and the elevation phase to obtain the initial deformation phase; S6, unwrapping the initial deformation phase to obtain the real deformation phase of each corner reflector; S7, obtaining the deformation amount based on the real deformation phase; Among them, S1 includes: Each pixel in the SAR image is filtered in the following way to obtain a filtered SAR image: For pixel point p in the SAR image, calculate the adaptive filtering value of p; Filter p based on the adaptive filter value to obtain the pixel value of p after filtering.

2. The long-distance high-precision deformation measurement method based on corner reflector SAR phase information according to claim 1 is characterized in that: The calculation formula of adaptive filtering value is: adpw p represents the adaptive filter value, pres represents the preset logarithmic parameter, wp represents the set of pixels in a 5×5 window centered on p, nwp represents the total number of pixels in wp, pxl i Represents the pixel value of pixel i, pxl p Represents the pixel value of pixel point p.

3. The long-distance high-precision deformation measurement method based on corner reflector SAR phase information according to claim 1 is characterized in that S2 include: In the registration process, one SAR image is selected from multiple SAR images as a reference image, and the remaining images are registered to the grid where the reference image is located.

4. The long-distance high-precision deformation measurement method based on corner reflector SAR phase information according to claim 1 is characterized in that: S3 includes: Performing interpolation processing on the registered SAR image to obtain an interpolated image; Get the position of each corner reflector in the interpolated image.

5. The long-distance high-precision deformation measurement method based on corner reflector SAR phase information according to claim 1 is characterized in that S4 include: Generate an interferogram based on two SAR images taken adjacently in time; The correlation formula between the phase difference ΔΦ and the elevation in the interference pattern is obtained: h represents the elevation, λ is the wavelength of the SAR signal, and θ is the angle between the SAR signal and the ground; Calculate h:

6. The long-distance high-precision deformation measurement method based on corner reflector SAR phase information according to claim 1 is characterized in that S5 include: The flat ground phase difference is calculated using the following formula: ΔΦflat=Φactual-Φflat ΔΦflat is the flat ground phase difference; Φactual is the actual phase of the corner reflector; Φflat is the theoretical flat ground phase of the corner reflector when the ground is completely flat; Calculate the elevation phase difference: ΔΦheight is the elevation phase difference; hcorner is the actual elevation of the corner reflector; hflat is the assumed elevation of an ideal flat ground; λ is the wavelength of the SAR signal; Calculate the atmospheric delay phase: represents the atmospheric delay phase, Δd represents the atmospheric delay from the zenith; The initial deformation phase is obtained by subtracting the ground phase difference, elevation phase difference and atmospheric delay phase from the interference phase of the corner reflector.

7. The long-distance high-precision deformation measurement method based on corner reflector SAR phase information according to claim 1 is characterized in that S6 include: n is an integer used to represent the integer ambiguity, represents the true deformation phase, represents the initial deformation phase.

8. The long-distance high-precision deformation measurement method based on corner reflector SAR phase information according to claim 1 is characterized in that S7 include: wd represents the deformation amount.

Images