Phase unwrapping method based on Bayesian-Sade integral Kalman filtering
By employing the Bayesian-Sadder integral Kalman filter algorithm, combined with phase gradient estimation and quality map, the problem of independent step error propagation in phase unwrapping and denoising in interferometric synthetic aperture radar is solved, achieving efficient and accurate phase unwrapping.
Patent Information
- Application Number
- CN202411515164.8
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- Filing Date
- 2024-10-29
- Publication Date
- 2025-02-07
AI Technical Summary
In existing interferometric synthetic aperture radar data processing, phase unwrapping and denoising are independent steps, which leads to error propagation. The integrated method has a long execution time and poor performance.
A phase unwrapping algorithm based on Bayes-Sade integral Kalman filtering is adopted. By calculating the phase gradient estimate and quality map, the algorithm unwrappes high-quality pixels to low-quality pixels using the Bayes-Sade integral Kalman filtering framework, establishes the state equation and observation equation, and calculates the unwrapping estimate.
It improves the accuracy and efficiency of phase unwrapping, reduces global noise propagation, and shortens processing time.
Smart Images
Figure CN119395698A_ABST
Abstract
Description
Technical Field
[0001] The invention relates to the technical field of interferometric synthetic aperture radar signal processing, and in particular to a phase unwrapping method based on Bayesian-Saad integral Kalman filtering. Background Art
[0002] Phase unwrapping is a key step in Interferometric Synthetic Aperture Radar (InSAR) data processing and is related to the application and development of InSAR technology. In recent years, many practical phase unwrapping methods have been proposed. Among them, the path tracking algorithm can achieve phase unwrapping by integrating along the tracking path. When unwrapping the wrapped phase with a high signal-to-noise ratio, the path tracking algorithm can obtain better results. However, when processing high-noise interferograms, errors will propagate along the selected path. The minimum norm method does not need to identify residuals or adopt a path tracking strategy, so its running time is acceptable. The minimum norm method is more efficient, but phase smoothing occurs during the unwrapping process.
[0003] In the traditional interferometric synthetic aperture radar data processing flow, phase denoising and unwrapping are regarded as two independent steps: denoising first and then unwrapping. However, this separate processing inevitably introduces additional errors because each step relies on different algorithms and assumptions. If only the denoising or unwrapping algorithm is improved, the unwrapping step will also be affected by the algorithmic errors in the denoising step. To address this problem, some integrated denoising and unwrapping methods have emerged in recent years. These methods perform phase unwrapping and phase noise filtering simultaneously, and therefore have received increasing attention. These integrated methods transform the phase unwrapping problem into a state estimation problem, and also filter the noise during the unwrapping process to reduce the global transmission of noise. However, the current integrated denoising and unwrapping methods cannot completely replace the traditional cascade processing framework of filtering first and then unwrapping, and their execution time is relatively long. Summary of the invention
[0004] The invention provides a phase unwrapping algorithm based on Bayesian-Saad integral Kalman filtering, which belongs to a state estimation method and is used to solve the phase unwrapping of InSAR.
[0005] To achieve the above object, the present invention provides a phase unwrapping algorithm based on Bayesian-Saad integral Kalman filtering, comprising the following steps:
[0006] S1. Calculate the phase gradient estimate of each pixel in the distance / azimuth direction and the quality map corresponding to the interference pattern.
[0007] S2. Based on the quality map, guide the Bayesian-Saad integral Kalman filter framework to disentangle from high-quality pixels to low-quality pixels.
[0008] S3. Based on the phase gradient estimate, the state equation and observation equation of Bayesian-Saad integral Kalman filter phase unwrapping are established. The state prediction mean and its output covariance are first calculated, and then the observation prediction mean and its output covariance, input-output covariance are calculated. Finally, the unwrapping estimate is calculated. Based on the quality map obtained in S1, high-quality pixels are unwrapped to low-quality pixels, and the final unwrapping phase is obtained.
[0009] Optionally, phase unwrapping is the process of recovering the true phase from a known two-dimensional phase value matrix modulo 2π, so phase unwrapping can be expressed as:
[0010]
[0011] In the formula, and φ(i,j) are the true phase and winding phase of pixel (i,j) respectively, and k(i,j) is an integer representing the fuzzy number of pixel (i,j).
[0012] Optionally, the S1 obtains the phase gradient estimation value of each pixel in the distance / azimuth direction through a phase gradient estimator based on an enhanced matrix bundle model.
[0013] Optionally, in S1, the quality map is phase derivative variance. When the phase derivative variance is used to represent the quality, the quality of the pixel point (m, n) is represented by q(m, N) = PDV m,n , PDV m,n is the phase derivative variance of pixel (m,n); therefore, the smaller the value of q(m,n), the higher the quality of pixel (m,n).
[0014] Optionally, the S2, according to the quality map: PDV obtained by S1, unwraps high-quality pixels to low-quality pixels, that is, q(m,n) in PDV, unwraps pixels with low q values to pixels with high q values.
[0015] Optionally, the S3 includes:
[0016] S3.1. Based on the phase gradient estimation value in S1, the state equation and observation equation of the Bayesian-Saad Kalman filter phase unwrapping are established to obtain the observed values of the unwrapped phase and observation coefficient of the pixel point;
[0017] Equation of state:
[0018]
[0019] Observation equation:
[0020]
[0021] In the state equation, is the true phase gradient estimate of the interference pattern at pixel k, ω(k) is the process noise, which is a Gaussian white noise with zero mean and a variance of Q(k); in the observation equation, v(k) is the observation noise, which is a Gaussian white noise with zero mean and a covariance of R(k).
[0022] S3.2. Based on S3.1, update the state equation to obtain the predicted value of the sampling point, the first step predicted value of the unwrapped phase and the corresponding error covariance matrix;
[0023] First set the sigma point:
[0024]
[0025] Among them, c = τ + D, where τ is the scaling parameter, D is the dimension, and here D = 1; the weights of these sigma points are:
[0026]
[0027] Get the state sampling point according to the sigma point:
[0028] X i (k-1)=x(k-1)+L x ξ (i) ,i=0,1,2
[0029] Where L x It is obtained by Cholesky decomposition of the covariance matrix P(k-1), that is, P(k-1) = L x L x T .
[0030] To obtain the state prediction value, input the sampling point into the state equation:
[0031]
[0032] Use the state prediction value to obtain the state prediction mean and its output covariance:
[0033]
[0034] State prediction mean is also the predicted value at pixel k, where:
[0035]
[0036] α is a positive integer
[0037] In the present invention, α1=α2=1, α3=2.
[0038] A=E ξ [φ(ξ)φ(ξ) T ]
[0039]
[0040] [D] n =E ξ [k(ξ,ξ (n) ,θ)φ(ξ) T ]
[0041] [K] nm = k(ξ (n) ,ξ (m) ,θ)
[0042] Among them, function k is the kernel function, and its formula is:
[0043]
[0044] In the above formula, θ=[α, l]. In the present invention, the parameter of state prediction is θ x , the parameter of observation prediction is θ y .
[0045] S3.3, based on S3.1-S3.2, obtaining the observed value and the corresponding covariance matrix;
[0046] Get the observation sampling points according to the sigma point:
[0047]
[0048] Where L y From the covariance matrix P x The Cholesky decomposition of (k) is obtained, that is, P x (k) = L y L y T .
[0049] To obtain the observed predicted value, input the sampling point into the observation equation:
[0050]
[0051] Use the observed predicted values to obtain the observed predicted mean and its output covariance:
[0052]
[0053] The input-output covariance is:
[0054]
[0055] in:
[0056] W c =BΦ -1
[0057] B=E ξ [ξφ(ξ) T ]
[0058] Get the Kalman gain matrix:
[0059] G k =P xy (k) / P y (k)
[0060] Finally determine the updated estimates and covariance matrix:
[0061]
[0062] in is the unwrapping estimate at pixel k.
[0063] S3.4, re-assign the pixel points, and repeat S3.1-S3.3 to obtain the final unwrapped phase.
[0064] The technical effects of the present invention are:
[0065] Simulated data with different coherence coefficients and real data from the Three Gorges region of China are tested to evaluate the effect of the phase unwrapping method based on the Bayesian-Saade integral Kalman filter. BRIEF DESCRIPTION OF THE DRAWINGS
[0066] Figure 1 is a flow chart of a phase unwrapping method based on Bayesian-Saad integral Kalman filtering in an embodiment of the present invention;
[0067] Figure 2 Schematic diagram of simulated interference phase data in an embodiment of the present invention, wherein (a) is the interference phase of simulated data with a coherence coefficient of 0.8, with a size of 256*256, and (b) is the real phase of the simulated data;
[0068] Figure 3 for Figure 2 (a) is a schematic diagram of the unwrapping of the phase using the unscented Kalman filter method. Figure 2 (a) is a schematic diagram of the unwrapping result, and (b) is a phase unwrapping method based on Bayesian-Saad integral Kalman filtering of the present invention. Figure 2 (a) Schematic diagram of the disentanglement results;
[0069] Figure 4 This is a schematic diagram of the phase of the measured data in the embodiment of the present invention, and the phase data of the Three Gorges region in China is selected;
[0070] Figure 5 For Figure 4 Schematic diagram of the unwrapping results, where (a) is the phase unwrapping method of the unscented Kalman filter. Figure 4 (b) is a schematic diagram of the unwrapping result of the phase unwrapping method based on Bayesian-Saad integral Kalman filtering of the present invention. Figure 4 (c) is a schematic diagram of the unwinding result of Figure (a), (d) is a schematic diagram of the rewinding result of Figure (b). DETAILED DESCRIPTION
[0071] In order to make the above-mentioned purpose of the present invention more obvious and understandable, the present invention is further described in detail below with reference to the accompanying drawings and specific implementation methods.
[0072] Reference Figure 1 As shown, the present invention proposes a phase unwrapping method based on Bayesian-Saad integral Kalman filtering, which specifically includes the following steps:
[0073] S1. Calculate the phase gradient estimate of each pixel in the distance / azimuth direction and the quality map corresponding to the interference pattern.
[0074] The quality graph that can represent the quality of the interference pattern is calculated. It should be noted that:
[0075] The commonly used quality map is the phase derivative variance (Phase Derivative Variance), which is defined as:
[0076]
[0077] PDV m,n represents the phase derivative variance value at the (m,n) position calculated in the k×k window, is the phase gradient of the winding phase in the vertical direction, is the phase gradient of the winding phase in the horizontal direction, and Displayed in a k×k window and The mean of .
[0078] S2. Based on the quality map, guide the Bayesian-Saad integral Kalman filter framework to disentangle from high-quality pixels to low-quality pixels.
[0079] S3. Based on the phase gradient estimate, the state equation and observation equation of Bayesian-Saad integral Kalman filter phase unwrapping are established. The state prediction mean and its output covariance are first calculated, and then the observation prediction mean and its output covariance, input-output covariance are calculated. Finally, the unwrapping estimate is calculated. Based on the quality map obtained in S1, high-quality pixels are unwrapped to low-quality pixels, and the final unwrapping phase is obtained.
[0080] S3.1. Based on the phase gradient estimation value in S1, the state equation and observation equation of the Bayesian-Saad Kalman filter phase unwrapping are established to obtain the observed values of the unwrapped phase and observation coefficient of the pixel point;
[0081] Equation of state:
[0082]
[0083] Observation equation:
[0084]
[0085] In the state equation, is the true phase gradient estimate of the interference pattern at pixel k, ω(k) is the process noise, which is a Gaussian white noise with zero mean and a variance of Q(k); in the observation equation, v(k) is the observation noise, which is a Gaussian white noise with zero mean and a covariance of R(k).
[0086] S3.2. Based on S3.1, update the state equation to obtain the predicted value of the sampling point, the first step predicted value of the unwrapped phase and the corresponding error covariance matrix;
[0087] First set the sigma point:
[0088]
[0089] Among them, c = τ + D, where τ is the scaling parameter, D is the dimension, and here D = 1; the weights of these sigma points are:
[0090]
[0091] Get the state sampling point according to the sigma point:
[0092] X i (k-1)=x(k-1)+L x ξ (i) ,i=0,1,2
[0093] Where L x It is obtained by Cholesky decomposition of the covariance matrix P(k-1), that is, P(k-1) = L x Lx T .
[0094] To obtain the state prediction value, input the sampling point into the state equation:
[0095]
[0096] Use the state prediction value to obtain the state prediction mean and its output covariance:
[0097]
[0098] State prediction mean is also the predicted value at pixel k, where:
[0099]
[0100] W=Φ -T AΦ -1
[0101] [Φ] nm =φ n (ξ (m) )
[0102]
[0103] φ(ξ)=[ξ α1 ,..,ξ αn ], α is a positive integer
[0104] In the present invention, α1=α2=1, α3=2.
[0105] A=E ξ [φ(ξ)φ(ξ) T ]
[0106]
[0107] [D] n =E ξ [k(ξ,ξ (n) ,θ)φ(ξ) T ]
[0108] [K] nm = k(ξ (n) ,ξ (m) ,θ)
[0109] Among them, function k is the kernel function, and its formula is:
[0110]
[0111] In the above formula, θ=[α, l]. In the present invention, the parameter of state prediction is θ x, the parameter of observation prediction is θ y S3.3, based on S3.1-S3.2, obtaining the observed value and the corresponding covariance matrix;
[0112] Get the observation sampling points according to the sigma point:
[0113]
[0114] Where L y From the covariance matrix P x The Cholesky decomposition of (k) is obtained, that is, P x (k) = L y L y T .
[0115] To obtain the observed predicted value, input the sampling point into the observation equation:
[0116]
[0117] Use the observed predicted values to obtain the observed predicted mean and its output covariance:
[0118]
[0119] The input-output covariance is:
[0120]
[0121] in:
[0122] W c =BΦ -1
[0123] B=E ξ [ξφ(ξ) T ]
[0124] Get the Kalman gain matrix:
[0125] G k =P xy (k) / P y (k)
[0126] Finally determine the updated estimates and covariance matrix:
[0127]
[0128] in is the unwrapping estimate at pixel k.
[0129] The pixel points are re-assigned, and the steps S3.1-S3.3 are repeatedly executed to obtain the final unwrapped phase of the interference pattern.
[0130] In order to verify the technical effect of the present invention, the unscented Kalman filter unwrapping method and the phase unwrapping method of the present invention are used to perform unwrapping experiments on the same interference pattern; the simulation data used are as follows: Figure 2 As shown, Figure 2 (a) is the interference pattern of simulated data with a coherence coefficient of 0.8. Figure 2 (b) is the real phase diagram of the simulated data; the real data used are Figure 4 Shown are the phase data for the Three Gorges region of China; Figure 3 The phase unwrapping method based on the Bayesian-Saad integral Kalman filter of the present invention is respectively used to unwrap the phase of the phase unwrapping method based on the Bayesian-Saad integral Kalman filter. Figure 2 (a) The untangling result; Figure 5 The phase unwrapping method based on the Bayesian-Saad integral Kalman filter of the present invention is respectively used to unwrap the phase of the phase unwrapping method based on the Bayesian-Saad integral Kalman filter. Figure 4 In the simulation data experiment, the interference patterns of simulated data with coherence coefficients of 0.2 to 0.9 were tested. Figure 2 Only the simulated data interference diagram with a coherence coefficient of 0.8 is listed in the figure. In the simulated data experiment, the present invention uses the root mean square error (RMSE) between the unwrapped phase and the true phase to represent the accuracy of phase unwrapping. Table 1 lists the simulated data interference diagram with a coherence coefficient of 0.2 to 0.9, and the root mean square errors of the unscented Kalman filter unwrapping method and the phase unwrapping method based on the Bayesian-Saade integral Kalman filter of the present invention are used respectively. In Table 1, the unscented Kalman filter unwrapping method is represented by UKF, and the phase unwrapping method based on the Bayesian-Saade integral Kalman filter is represented by BSQKF.
[0131] Table 1 Root mean square error (RMSE) between the unwrapped phase of simulated data and the true phase
[0132]
[0133] According to the results in Table 1, it can be seen that when the coherence coefficient is 0.3-0.6, the root mean square error of the phase unwrapping method based on the Bayesian-Saad integral Kalman filter of the present invention is smaller than the root mean square error of the unscented Kalman filter unwrapping method.
[0134] The present invention calculates the residual count of the result after re-wrapping of the real data unwrapping result to judge the effect of unwrapping the real data. The residual count of the original wrapped phase of the real data is 10163, the residual count of the re-wrapping of the unscented Kalman filter phase unwrapping method is 483, and the residual count of the re-wrapping of the phase unwrapping method based on the Bayesian-Saade integral Kalman filter of the present invention is 3893.
[0135] In summary, the phase unwrapping method based on Bayesian-Saad integral Kalman filtering of the present invention can obtain relatively reliable phase unwrapping results.
[0136] Finally, it should be noted that the embodiments used above are only used to illustrate the technical solutions of the present invention. Obviously, they cannot be used to limit the scope of rights of the present invention. Ordinary technicians in this field can understand that all or part of the processes of the above embodiments and equivalent replacements made according to the claims of the present invention still fall within the scope covered by the present invention.
Claims
1. A phase unwrapping method based on Bayesian-Saad integral Kalman filtering, characterized in that: The following steps are involved: S1. Calculate the phase gradient estimate of each pixel in the distance / azimuth direction and the quality map corresponding to the interference pattern. S2. Based on the quality map, guide the Bayesian-Saad integral Kalman filter framework to disentangle from high-quality pixels to low-quality pixels. S3. Based on the phase gradient estimate, the state equation and observation equation of Bayesian-Saad integral Kalman filter phase unwrapping are established. The state prediction mean and its output covariance are first calculated, and then the observation prediction mean and its output covariance, input-output covariance are calculated. Finally, the unwrapping estimate is calculated. Based on the quality map obtained in S1, high-quality pixels are unwrapped to low-quality pixels, and the final unwrapping phase is obtained.
2. The phase unwrapping method based on Bayesian-Saad integral Kalman filtering according to claim 1, characterized in that: Phase unwrapping is the process of recovering the true phase from the known two-dimensional phase value matrix modulo 2π. Therefore, phase unwrapping can be expressed as: In the formula, and φ(i,j) are the true phase and winding phase of pixel (i,j) respectively, and k(i,j) is an integer representing the fuzzy number of pixel (i,j).
3. The phase unwrapping method based on Bayesian-Saad integral Kalman filtering according to claim 1, characterized in that: In S1, the phase gradient estimation value of each pixel in the distance / azimuth direction is obtained by a phase gradient estimator based on an enhanced matrix bundle model.
4. The phase unwrapping method based on Bayesian-Saad integral Kalman filtering according to claim 1, characterized in that: In S1, the quality map is phase derivative variance. When phase derivative variance is used to represent quality, the quality of pixel point (m, n) is represented by q(m, n)=PDV. m,n , PDV m,n is the phase derivative variance of pixel (m,n); therefore, the smaller the value of q(m,n), the higher the quality of pixel (m,n).
5. The phase unwrapping method based on Bayesian-Saad integral Kalman filtering according to claim 1, characterized in that: In S2, according to the quality map: PDV obtained in S1, high-quality pixels are unwrapped to low-quality pixels, that is, q(m,n) in PDV, pixels with low q values are unwrapped to pixels with high q values.
6. The phase unwrapping method based on Bayesian-Saad integral Kalman filtering according to claim 1, characterized in that: The S3 includes: S3.
1. Based on the phase gradient estimation value in S1, the state equation and observation equation of the Bayesian-Saad Kalman filter phase unwrapping are established to obtain the observed values of the unwrapped phase and observation coefficient of the pixel point; Equation of state: Observation equation: In the state equation, is the true phase gradient estimate of the interference pattern at pixel k, ω(k) is the process noise, which is a Gaussian white noise with zero mean and a variance of Q(k); in the observation equation, v(k) is the observation noise, which is a Gaussian white noise with zero mean and a covariance of R(k); S3.
2. Based on S3.1, update the state equation to obtain the predicted value of the sampling point, the first step predicted value of the unwrapped phase and the corresponding error covariance matrix; First set the sigma point: Among them, c = τ + D, where τ is the scaling parameter, D is the dimension, and here D = 1; the weights of these sigma points are: Get the state sampling point according to the sigma point: X i (k-1)=x(k-1)+L x ξ (i) ,i=0,1,2 Where L x It is obtained by Cholesky decomposition of the covariance matrix P(k-1), that is, P(k-1) = L x L x T . To obtain the state prediction value, input the sampling point into the state equation: Use the state prediction value to obtain the state prediction mean and its output covariance: State prediction mean is also the predicted value at pixel k, where: W=Φ -T AF -1 [F] nm =φ n (x) (m) ) α is a positive integer In the present invention, α1=α2=1, α3=2. A=E ξ [φ(ξ)φ(ξ) T ] [D] n =E ξ [k(ξ,ξ (n) ,θ)φ(ξ) T ] [K] nm =k(ξ (n) ,x (m) ,i) Among them, function k is the kernel function, and its formula is: The above In the present invention, the parameter of state prediction is θ x , the parameter of observation prediction is θ y ; S3.3, based on S3.1-S3.2, obtaining the observed value and the corresponding covariance matrix; Get the observation sampling points according to the sigma point: Where L y From the covariance matrix P x The Cholesky decomposition of (k) is obtained, that is, P x (k) = L y L y T . To obtain the observed predicted value, input the sampling point into the observation equation: Use the observed predicted values to obtain the observed predicted mean and its output covariance: The input-output covariance is: in: IN c =BΦ -1 B=E ξ [ξφ(ξ) T ] Get the Kalman gain matrix: G k =P xy (k) / P y (k) Finally determine the updated estimates and covariance matrix: in is the unwrapping estimate at pixel k; S3.4, re-assign the pixel points, and repeat S3.1-S3.3 to obtain the final unwrapped phase.