Distributed spaceborne D-InSAR time sequence deformation inversion method

Abstract

The invention discloses a distributed spaceborne D-I nSAR time sequence deformation inversion method, and the method can solve the influence of atmospheric space-time non-stationary disturbance through the correlation of different components in the space-time dimension, employing a space-time dimension iterative algorithm, and introducing particle filtering, thereby achieving the inversion of more accurate topographic information. According to the method, based on particle filtering and a trend-signal-noise model, a deformation phase and an atmospheric phase are roughly separated in a distributed dimension, and then a final result of the deformation phase is obtained through alternate iteration of a time dimension and an interferogram space dimension. Experimental results prove that a deformation inversion result of the method can obtain higher precision than that of a traditional window-based filtering method.

Description

Technical Field

[0001] The present invention belongs to the technical field of synthetic aperture radar, and in particular relates to a distributed spaceborne SAR time series deformation inversion method. Background Art

[0002] Differential Interferometry (D-InSAR) is a method for detecting surface deformation using phase differences between two or more Synthetic Aperture Radar (SAR) images. Distributed spaceborne SAR is a new radar system composed of multiple SAR satellites working together. Applying D-InSAR technology to distributed spaceborne SAR systems offers advantages over traditional single-satellite D-InSAR, such as improved deformation measurement data acquisition efficiency and increased observation dimensionality. Consequently, many countries have launched their own distributed spaceborne D-InSAR systems, such as Germany's TerraSAR-X and TanDEM-X dual-satellite system and China's Lutan-1.

[0003] While distributed spaceborne D-InSAR systems offer the aforementioned advantages, they are severely impacted by non-stationary, non-Gaussian atmospheric phase disturbances. This is particularly true when subjected to severe atmospheric turbulence. Because atmospheric conditions vary with geography and time, the atmospheric delay phases in the differential interferograms of each satellite observed at different times and orbital positions. These non-Gaussian, non-stationary characteristics introduce additional interferometric phase differences, disrupting the accuracy of the distributed spaceborne D-InSAR collaborative deformation inversion.

[0004] Existing methods that use external data to compensate for atmospheric phases cannot fully meet the compensation requirements for high-precision deformation measurements from distributed spaceborne D-InSAR under non-stationary, non-Gaussian conditions with severe atmospheric turbulence. Furthermore, methods that use dense formations to eliminate atmospheric effects impose strict constraints on distributed configurations, limiting the observational capabilities of deformation inversion.

[0005] First, external data (such as numerical meteorological models or GPS measurements) are often limited by spatial resolution, coverage, and accuracy, and cannot meet the atmospheric compensation requirements of distributed spaceborne D-InSAR. Second, the dense formation method of eliminating atmospheric effects requires the distributed spaceborne SAR system to have a dense formation configuration, with each satellite very close in space. This allows the influence of atmospheric phase to be eliminated through inter-satellite phase difference. However, this method has strict requirements on the configuration of the distributed system and can only be applied to terrain mapping, not deformation monitoring. Therefore, other more sophisticated atmospheric phase correction methods are needed to achieve high-precision deformation inversion of distributed spaceborne D-InSAR. Summary of the Invention

[0006] To solve the above problems, the present invention provides a distributed spaceborne D-InSAR time series deformation inversion method. Based on particle filtering and the trend-signal-noise model, the method coarsely separates the deformation phase and the atmospheric phase in the distributed dimension. Then, the final result of the deformation phase is obtained by alternating iterations in the time dimension and the interferogram space dimension. Experimental results show that the deformation inversion results of the proposed method can achieve higher accuracy than the traditional window-based filtering method.

[0007] The technical solution of the present invention is:

[0008] A distributed spaceborne D-InSAR time series deformation inversion method, comprising:

[0009] Step 1: Coarse separation of distributed deformation. If there is no time-dimensional deformation model, first use the spatial K-means clustering algorithm to obtain some cluster points. In the distributed dimension, use the least squares method to roughly separate the deformation and atmospheric phase. Then, use polynomial fitting of the time series deformation variables to obtain the time-dimensional deformation model. In addition, the time-dimensional variance-covariance matrix is ​​initialized using data obtained from some discrete ground observation stations in the target area:

[0010] Finally, in this section, the correlation of parameters in different dimensions is first analyzed, as shown in Table 1.

[0011] Table 1 Correlation of parameters in the STI-PF algorithm in different dimensions

[0012]

[0013] Combining the relevant models and assumptions for each of these dimensions, the distributed dimension can be used to separate deformation, the associated atmospheric phase, and the unrelated components of the atmospheric phase, as well as other noise phases. In the temporal dimension, the associated deformation components and the temporally unrelated atmospheric phases can be precisely separated. In the interferogram dimension, the atmospheric phase and the noise phase can be further separated. The final deformation inversion result is obtained by iterating the temporal and spatial dimensions of the interferogram.

[0014] Step 2, time-dimensional TSN-PF: Combine the time-dimensional deformation state equation obtained in the previous step with the time-dimensional TSN model (as shown in formula (10)) to achieve the separation of deformation and atmospheric phase; further use the spatial distribution of the separated atmospheric phase to estimate the variance-covariance matrix of the spatial-dimensional atmospheric phase:

[0015] For differential interferometry phase data, decomposing it into long-term trends, abrupt signals, and additive noise is a common modeling approach, which helps to better understand and analyze the characteristics of interferometry data. This decomposition can be achieved using the Trend-Signal-Noise (TSN) model:

[0016] y =t+ s + n (1)

[0017] in, y The observation vector t is a deterministic but unknown trend (space or time), s is a zero-mean random signal vector, n is a zero-mean noise vector; the underlined variables are random variables, and the rest are deterministic variables. n It has non-stationary characteristics, so we introduce particle filtering and design the TSN-PF high-precision three-dimensional deformation inversion method.

[0018] The unknown trend t in Equation (1) is generally expressed as t = Ax using the unknown parameter vector x, where A is a design matrix. Assuming that the signal s and the noise n are uncorrelated, their variance-covariance matrices (VCMs) are expressed as Q ss and Q nn The complete form of the trend-signal-noise model can be written as:

[0019]

[0020] Among them, Q yy yes y The VCM can be obtained through the error propagation law.

[0021] The trend-signal-noise model based on particle filtering can be divided into the following three steps:

[0022] 1. Calculate the estimated value of x:

[0023] In this step, s and n Merge processing, that is SN = s + n . So the original trend-signal-noise model can be transformed into:

[0024] y =Ax+ SN (3)

[0025] Using this as the state equation, combined with the measurement equation, the particle filter process is performed to obtain the estimated value of x The main process of particle filtering is as follows:

[0026] A. Particle set initialization: k = 0; for i = 1, 2, ..., N s, the sampling particles are generated by the prior p(x0)

[0027] B. For k = 1, 2, ..., loop through the following steps:

[0028] 1) Importance sampling: for i=1,2,...,N s , generate sampling particles from the importance probability density Calculate particle weights And normalize it;

[0029] 2) Resampling: Particle set Resample to obtain the particle set

[0030] 3) Output: Calculate the state estimate at time k:

[0031] 2. Conduct s and n Prediction:

[0032] When the estimated value of x OK, then you can get s and n The optimal linear unbiased prediction of :

[0033]

[0034] 3. Estimate the atmospheric spatial dimension variance-covariance matrix:

[0035] Assuming that the atmosphere is stationary in the spatial dimension, the spatial variance-covariance matrix of the atmosphere can be estimated by sampling the spatial distribution of the atmosphere at different times.

[0036] Step 3, spatial-dimensional TSN-PF: The atmospheric phase estimated in the previous step and its spatial-dimensional variance-covariance matrix are combined with the spatial-dimensional TSN model to accurately separate the atmospheric phase and the noise phase. After the separation is completed, the atmospheric phase of different time series at the same location is used to estimate the time-dimensional variance-covariance matrix of the atmospheric phase:

[0037] The distribution of the atmosphere in the spatial dimension can be written as a fixed spatial trend, a sudden change signal, and additive noise. Therefore, it can be modeled using TSN and solved using particle filtering. The TSN model is the same as in step 2.

[0038] 1. Calculate the estimated value of x:

[0039] The particle filter-based solution method is the same as that in Part 2.

[0040] 2. Conduct s and nPrediction:

[0041] According to the spatial dimension variance-covariance matrix of the atmosphere obtained in Part 2, the phase estimation of random atmosphere and noise can be obtained.

[0042] 3. Estimate the variance-covariance matrix of the deformation time dimension:

[0043] Assuming that the spatial dimension of the deformation is stationary, the temporal dimension variance-covariance matrix of the deformation can be estimated based on the temporal dimension distribution of the deformation at different vacancy positions.

[0044] Step 4: Compare the atmospheric phase time dimension variance-covariance matrix obtained above with the input prior atmospheric phase time dimension variance-covariance matrix; if the difference between the two is greater than the set threshold, return to step 2 and re-perform the time dimension TSN-PF processing; if the difference is less than or equal to the threshold, directly output the final deformation estimation result:

[0045] If the difference between the variance-covariance matrix calculated in step 2 and the variance-covariance matrix initialized in step 1 is higher than the threshold, the variance-covariance matrix in step 3 is substituted into step 2 and iterated; otherwise, the iteration is terminated to obtain the final deformation inversion result.

[0046] Beneficial effects:

[0047] 1. To address the difficulty in decoupling deformation signals from atmospheric disturbances in existing technologies, a trend-signal-noise model was designed. This model combines deformation, atmosphere, and noise to accurately characterize these two factors.

[0048] 2. Innovatively introduce distributed dimensions and design a distributed dimension-space dimension-time dimension information iterative processing algorithm to significantly improve the utilization of multi-dimensional information;

[0049] 3. To address the failure of traditional Kalman filtering under non-stationary and non-Gaussian assumptions, a particle filter is introduced to dynamically track the non-stationary characteristics of the atmospheric phase screen, and maintain high-precision deformation inversion under severe convective weather conditions. BRIEF DESCRIPTION OF THE DRAWINGS

[0050] Figure 1 A schematic flow chart of the method proposed in the present invention;

[0051] Figure 2 It is the topographic map of the experimental scene;

[0052] Figure 3 Schematic diagram of slope calculation for Ritter algorithm;

[0053] Figure 4 is the scene gradient field map;

[0054] Figure 5 is the scene slope field map;

[0055] Figure 6 This is the time sequence diagram of landslide simulation;

[0056] Figure 7 The comparison results of STI-PF inversion with a priori time-dimensional deformation model (left is the true value of scene deformation, right is the inversion result);

[0057] Figure 8 Error distribution histograms for the inverted deformation using the STI-PF method (left) and the standard window-based filtering method (right). DETAILED DESCRIPTION

[0058] The present invention will be further described below with reference to the accompanying drawings and embodiments.

[0059] The overall processing flow of the present invention is as follows Figure 1 As shown in the figure, through the joint processing of distributed dimension, time dimension and space dimension, the deformation and atmospheric phase are first roughly separated by using the distributed dimension, and then the deformation and atmospheric phase and other phase components are separated by particle filtering in time and space dimensions.

[0060] The input of the method of the present invention is the center position of the formation satellite imaging, satellite imaging parameters (bandwidth, synthetic aperture time, wavelength, etc.), the target scene DEM, the center position of the target scene, the target scene size, the observation time interval, and the atmospheric phase time-dimensional variance-covariance matrix; the output is the three-dimensional deformation result. The specific steps are as follows:

[0061] Step 1: Coarse separation of distributed deformation. If there is no time-dimensional deformation model, first use the spatial K-means clustering algorithm to obtain some cluster points. In the distributed dimension, use the least squares method to roughly separate the deformation and atmospheric phase. Then, use polynomial fitting of the time series deformation variables to obtain the time-dimensional deformation model. In addition, the time-dimensional variance-covariance matrix is ​​initialized using data obtained from some discrete ground observation stations in the target area:

[0062] In the distributed dimension, the three stars are close to each other and the observation time is close, so the deformation variables they observe can be assumed to be constant and a deterministic function; the atmospheric phase part can be assumed to be a random process with partial correlation. It can be expressed as:

[0063]

[0064] in, is the deformation phase of the winding, are the relevant components in the atmospheric phase; both are considered to be constant parts in the distributed dimension, that is, is the uncorrelated part of the atmospheric phase, For other noise phases, the two parts are considered to be unrelated random processes, that is, Based on this assumption, the least squares method can be used to extract The model corresponding to formula (5) can be written as:

[0065]

[0066] Among them, x cons is the parameter to be estimated, including deformation and related components in the atmosphere; C is the design matrix:

[0067]

[0068] in, and They represent the unit slant range vector of the i-th satellite in the system to the target scene and the three-dimensional deformation vector of the scene respectively.

[0069] According to the LS algorithm, the solution of formula (6) can be obtained as:

[0070]

[0071] From formula (8), the deformation and the related parts in the atmosphere can be inverted.

[0072] Step 2: Time-dimensional TSN-PF. Combine the time-dimensional deformation state equation obtained in the previous step with the time-dimensional TSN model (as shown in formula (10)) to separate the deformation from the atmospheric phase. Further, use the spatial distribution of the separated atmospheric phase to estimate the variance-covariance matrix of the spatial atmospheric phase:

[0073] For D-InSAR measurements, according to the geometric relationship, the differential interferometry phase model can be expressed as:

[0074]

[0075] The subscript _ represents a random variable, and the others are deterministic components; φ defo To represent the surface deformation, it is assumed to be a random variable because it can be regarded as a deterministic function in time (linear function, quadratic function or periodic function trend) Plus a random part with zero mean and second order stationary in time It is the phase caused by the reference DEM error; and They represent the atmospheric phase of the main and auxiliary images respectively. Since all time series images share the same main image, the atmospheric phase of the main image is constant. The orbital errors of the main and auxiliary images are and The same applies; is the noise phase.

[0076] The state equation is obtained by fitting the coarse separation results obtained in step 3. Because the surface deformation is decomposed into deterministic and random components, the state equation only needs to model the deterministic deformation and can generally be set as a linear function, quadratic function, sine and cosine function, etc. At the same time, due to the resampling step in PF, it does not require the noise to be a temporally stationary random process.

[0077] The corresponding temporal trend-signal-noise model is:

[0078]

[0079] At the same time, this equation is the measurement equation in time. The deformation part is dispersed in x p,r and s p,r The two parts respectively contain the deterministic part of the deformation (that is, the part that can be represented by the model) and the random part (the part that cannot be represented by the model); at the same time, x p,r Also contains a and The part ξ that is composed of the sum of s p,r that is This part is a second-order stationary and zero-mean normal distribution in time; where A is the design matrix, which is related to the spatiotemporal baseline of the observed time series and is expressed as:

[0080]

[0081] Among them, B T is the time baseline; [·] 2π is the winding phase.

[0082] Noise component n p,r It can be expressed as:

[0083]

[0084] Finally, using the TSN-PF algorithm in step 2, we can get the estimated value It contains the deterministic component of deformation; it also obtains and Among them, the prior information of the deformation part that cannot be represented by the model can be easily obtained from some discrete ground observation stations in the target area (such as GNSS measurement, laser equipment, etc.).

[0085] The output includes deterministic deformations relative to a global reference point, random deformations, DEM phase residuals, and estimates of the phase sum of the APS, satellite orbits, and noise.

[0086] Step 3, spatial-dimensional TSN-PF: The atmospheric phase and its spatial-dimensional variance-covariance matrix estimated in the previous step are combined with the spatial-dimensional TSN model (as shown in formula (14)) to accurately separate the atmospheric phase and the noise phase. After the separation is completed, the atmospheric phase of different time series at the same location is used to estimate the time-dimensional variance-covariance matrix of the atmospheric phase:

[0087] Through the TSN-PF algorithm in the time dimension, we have obtained the time series However, these results are all relative to the same global reference point. Next, we further separate the APS from the orbit error and phase noise using the spatial TSN-PF algorithm.

[0088] Assume that k represents the kth auxiliary image, then n k It can be expressed as:

[0089]

[0090] Based on the above components, we can once again construct a trend-signal-noise model:

[0091] n k =Ry k + ν k + μ k (14)

[0092] Where R is the spatial design matrix with M rows and 4 columns, y k is a 4×1 parameter vector that includes three parameters that determine the linear surface trend caused by the troposphere and orbit errors, as well as a coefficient related to vertical stratification; ν k is the phase including the turbulent component; finally, μ k is the noise component. According to the TSN-PF algorithm, we can solve and

[0093] Step 4: Compare the obtained atmospheric phase time-dimension variance-covariance matrix with the input prior atmospheric phase time-dimension variance-covariance matrix. If the difference between the two is greater than the set threshold, return to step 2 and re-perform the time-dimension TSN-PF processing. If the difference is less than or equal to the threshold, directly output the final deformation estimation result.

[0094] Next, an implementation example is given with specific parameters.

[0095] Step 1: Coarse separation of distributed deformation. If there is no time-dimensional deformation model, first use the spatial K-means clustering algorithm to obtain some cluster points. In the distributed dimension, use the least squares method to roughly separate the deformation and atmospheric phase. Then, use polynomial fitting of the time series deformation variables to obtain the time-dimensional deformation model. In addition, the time-dimensional variance-covariance matrix is ​​initialized using data obtained from some discrete ground observation stations in the target area:

[0096] First, a simulation deformation scene is generated. According to the current mainstream research on landslide-prone conditions, factors such as slope angle, height range, rock properties of slope surface, and slope direction or oblique direction of slope surface in mountainous areas may all lead to landslide formation. For the sake of simplicity, this invention only uses height range and slope angle as the decisive factors for spatial landslide generation. The terrain of the target scene and the corresponding contour terrain are as follows: Figure 2 shown.

[0097] The principle of slope calculation by Ritter algorithm is as follows Figure 3 As shown, if the slope of the center point e is calculated, the elevation of the four pixels directly adjacent to the center point e is considered, then the slope S e It can be expressed as:

[0098]

[0099] Where d represents the pixel size. For example, if SRTM-30m is used, d is 30. (e1-e3) and (e2-e4) are the elevation differences in two directions. If SRTM-30m is used, the two correspond to the elevation differences in longitude and latitude, respectively.

[0100] For a surface f, the gradient is defined as:

[0101]

[0102] The gradient field and slope field of the study area have been established. Figure 4 、 Figure 5 shown

[0103] According to the time history curve relationship of the corresponding surface deformation variable in the landslide disaster, combined with the spatial relationship between height and slope angle, the landslide deformation can be simulated as follows: Figure 6 shown.

[0104] Next, we generate time series data. In this example, we simulate imaging using a distributed spaceborne SAR system consisting of three satellites at daily intervals. The system parameters are shown in Table 2. The simulation process is based on the input satellite imaging center position, satellite imaging parameters (including bandwidth, synthetic aperture time, wavelength, etc.), and the target scene center position. Three SAR images are generated every day for seven days, resulting in a total of 21 SAR images.

[0105] Table 2 Distributed spaceborne SAR system parameters

[0106]

[0107] Subsequently, the 21 SAR images were precisely registered using either cross-correlation or interferometric data registration methods assisted by an external digital elevation model (DEM). After registration, differential interferometry processing was performed using the three SAR images acquired on the first day as primary images and the 18 images acquired on the remaining six days as auxiliary images, resulting in 18 differential interferograms.

[0108] The 18 differential interferograms were divided into six groups based on acquisition time, with each group containing three images collected each day. Within each group, a rough separation of deformation and atmospheric phase was performed using Equation (8). After separation, the deformation results were fitted to obtain the time-dimensional deformation equation of state.

[0109] Step 2, time-dimensional TSN-PF; combine the time-dimensional deformation state equation obtained in the previous step with the time-dimensional TSN model (as shown in formula (10)) to achieve the separation of deformation and atmospheric phase; further use the spatial distribution of the separated atmospheric phase to estimate the variance-covariance matrix of the spatial-dimensional atmospheric phase.

[0110] Step 3, spatial-dimensional TSN-PF: The atmospheric phase and its spatial-dimensional variance-covariance matrix estimated in the previous step are combined with the spatial-dimensional TSN model (as shown in formula (14)) to accurately separate the atmospheric phase and the noise phase. After the separation is completed, the atmospheric phase of different time series at the same location is used to estimate the time-dimensional variance-covariance matrix of the atmospheric phase.

[0111] Step 4: Compare the atmospheric phase time dimension variance-covariance matrix obtained above with the input prior atmospheric phase time dimension variance-covariance matrix; if the difference between the two is greater than the set threshold, return to step 4 and re-perform the time dimension TSN-PF processing; if the difference is less than or equal to the threshold, directly output the final deformation estimation result.

[0112] The experimental results are as follows Figure 7As shown in the figure, the deformation inverted by STI-PF is very close to the reference deformation in the simulation experiment. In fact, on the seventh day, the RMSE between the deformation inverted by the STI-PF method and the simulated reference deformation was approximately 4.5 mm. This demonstrates that the STI-PF method can accurately recover deformation information and has excellent deformation monitoring capabilities. Secondly, the STI-PF method does not require prior unwrapping and can directly process using the wrapped phase. This feature is very important when dealing with complex nonlinear perturbations. Typically, in InSAR processing, unwrapping algorithms are required to convert complex phase information into deformation information. However, unwrapping algorithms themselves can also introduce errors, which can affect the accuracy of deformation monitoring. The advantage of the STI-PF method is that it can directly process using the wrapped phase, which not only saves steps but also reduces the impact of unwrapping errors on deformation monitoring results. Third, the STI-PF method can adaptively adjust to non-stationary and non-Gaussian atmospheric phases. In InSAR processing, varying atmospheric conditions can cause nonlinear perturbations in the phase, which can affect the accuracy of deformation monitoring results. Moreover, the atmospheric phase is non-stationary in both space and time. Traditional processing methods often require the direct assumption that the atmospheric phase is stationary and Gaussian. However, the STI-PF method can be adaptively adjusted according to the resampling steps, so it also works well for non-stationary and non-Gaussian atmospheric phases. Figure 8 The error distribution histogram obtained by inverting deformation using the STI-PF method shows that its error range is approximately -6mm to 6mm, with an RMSE of approximately 4.5mm. This demonstrates that the STI-PF method can effectively improve the accuracy and reliability of deformation monitoring. In contrast, the error distribution histogram obtained by inverting deformation using the traditional window-based filtering method shows a relatively dispersed error distribution, with an error range of approximately -10mm to 10mm and an RMSE of approximately 8.8mm. Compared with traditional methods, the STI-PF method can reduce inversion error by nearly 50%, achieving higher accuracy and broader application prospects.

[0113] Of course, the present invention may have many other embodiments. Without departing from the spirit and essence of the present invention, those skilled in the art may make various corresponding changes and modifications based on the present invention, but these corresponding changes and modifications should all fall within the scope of protection of the claims attached to the present invention.

Claims

1. A distributed spaceborne D-InSAR time series deformation inversion method, characterized by: The following steps are involved: The first step is to roughly separate the distributed deformation. If there is no time-dimensional deformation model, a spatial K-means clustering algorithm is first used to obtain some cluster points. In the distributed dimension, the least squares method is used to roughly separate the deformation and atmospheric phase. Then, a polynomial fitting of the time series deformation variables is used to obtain the time-dimensional deformation model. In addition, data obtained from some discrete ground observation stations in the target area are used to initialize the time-dimensional variance-covariance matrix. The second step is the time-dimensional TSN-PF. The time-dimensional deformation state equation obtained in the previous step is combined with the time-dimensional TSN model (as shown in formula (10)) to achieve the separation of deformation and atmospheric phase. The spatial distribution of the separated atmospheric phase is further used to estimate the variance-covariance matrix of the spatial-dimensional atmospheric phase. Step 3: Spatial TSN-PF. The atmospheric phase and its spatial variance-covariance matrix estimated in the previous step are combined with the spatial TSN model (as shown in formula (14)) to accurately separate the atmospheric phase and the noise phase. After the separation is completed, the atmospheric phase of different time series at the same location is used to estimate the temporal variance-covariance matrix of the atmospheric phase. In the fourth step, the obtained atmospheric phase time dimension variance-covariance matrix is ​​compared with the input prior atmospheric phase time dimension variance-covariance matrix; if the difference between the two is greater than the set threshold, return to step 4 and re-perform the time dimension TSN-PF processing; if the difference is less than or equal to the threshold, the final deformation estimation result is directly output.

2. The distributed spaceborne D-InSAR time series deformation inversion method according to claim 1, characterized in that: In step 1, if there is no deformation model in the time dimension, the spatial K-means clustering algorithm is first used to obtain some cluster points. In the distributed dimension, the least squares method is used to roughly separate the deformation and the atmospheric phase. Then, the deformation model in the time dimension is obtained by fitting the time series deformation variables with polynomials.

3. The distributed spaceborne D-InSAR time series deformation inversion method according to claim 1, characterized in that: In step 2, the complete form of the trend-signal-noise model is: y =Ax+ s + n Q yy =Q ss +Q nn in, y The observation vector t is a deterministic but unknown trend (space or time), s is a zero-mean random signal vector, n is a zero-mean noise vector, Q yy yes y The VCM can be obtained through the error propagation law.

4. The distributed spaceborne D-InSAR time series deformation inversion method according to claim 1, characterized in that: In step 2, the trend-signal-noise model based on particle filtering can be divided into three steps:

1. Calculate the estimated value of x; 2. s and n prediction; 3. Estimate the atmospheric spatial dimensional variance-covariance matrix.

5. The distributed spaceborne D-InSAR time series deformation inversion method according to claim 1, characterized in that: In step 2, the main process of particle filtering is as follows: A. Particle set initialization: k = 0; for i = 1, 2, ..., N s , the sampling particles are generated by the prior p(x0) B. For k = 1, 2, ..., loop through the following steps: 1) Importance sampling: for i=1,2,...,N s , generate sampling particles from the importance probability density Calculate particle weights And normalize it; 2) Resampling: Particle set Resample to obtain the particle set 3) Output: Calculate the state estimate at time k:

6. The distributed spaceborne D-InSAR time series deformation inversion method according to claim 1, characterized in that: In step 2, s and n The prediction method is: When the estimated value of x Confirm, then you will get s and n The optimal linear unbiased prediction of : in and They are s and n The optimal linear unbiased prediction, Q ss and Q yy They are y and VCM of x, A is the design matrix.

7. The distributed spaceborne D-InSAR time series deformation inversion method according to claim 1, characterized in that: In step three, the trend-signal-noise model is the same as that in step two.

8. The distributed spaceborne D-InSAR time series deformation inversion method according to claim 1, characterized in that: In step 4, the atmospheric phase time dimension variance-covariance matrix is ​​compared with the input prior atmospheric phase time dimension variance-covariance matrix; if the difference between the two is greater than the set threshold, the process returns to step 4 and re-performs the time dimension TSN-PF processing; if the difference is less than or equal to the threshold, the final deformation estimation result is directly output.

Images