Performance Evaluation of Cubature Kalman Filtering and Extended Kalman Filtering Based Phase Unwrapping For Insar

Phase unwrapping is one of the most critical steps in InSAR data processing. Filter-based phase unwrapping methods can perform phase unwrapping and noise removal simultaneously, and show a clear advantage in dealing with high noise and density fringe interferograms. Based on the precision estimation theory of nonlinear functions in Extended Kalman Filter (EKF) and Cubature Kalman Filter (CKF) algorithms, we derive and analyze in detail the estimated precision of the EKF and the CKF on nonlinear observation equations in the phase unwrapping model. Results show that both EKF and CKF algorithms have similar precision in areas with better phase quality. However, in poorer quality areas, the precision of the two shows a certain difference and there is no consistent superiority for each algorithm. Experiments using simulated and measurement data confirm the validity of our theoretical analysis.


Introduction
Phase unwrapping is one of the most critical steps in InSAR data processing and its precision directly affects the precision of elevation measurements. Filter-based unwrapping methods transform phase unwrapping into state estimation, and have become the latest approach, attracting many scholars' attention. This method is not affected by phase residual points, and avoids the phase loss and distortion caused by traditional methods by filtering the noisy interferograms before phase unwrapping is carried out, which results in simultaneous noise filtering and phase unwrapping.
Krä mer and Loffeld proposed an InSAR phase unwrapping method based on Kalman filters for the first time in 1996 [1]. Then, they put forward a method combining Extended Kalman Filters (EKFs) with local slope estimation [2] in 1997. Since then, the filter-based phase unwrapping algorithms have been gradually paid more attention. In 1999, Kim and Griffiths came up with the idea of using EKF to achieve multi-baseline phase unwrapping, focusing on the advantages of the EKF in phase unwrapping and the high precision potential provided by the powerful information fusion capabilities of the EKF algorithm [3]. In 2008, specific state and observation models were presented and the EKF algorithm was elaborated by Loffeld et al. [4]. The algorithm transformed phase unwrapping into state estimation, integrating state and observation models, and achieving the purpose of phase unwrapping and noise elimination simultaneously. In 2009, Osmanoglu proposed an improved EKFPU algorithm, in which a simple linear operation was applied to the gradient estimation, and error detection and correction were added into it [5]. A 3-D phase unwrapping method based on EKF was also then investigated to  [6]. Chirico et al also applied EKF to Multichannel interferometric phase unwrapping [7]. Guolin Liu proposed Kalman filter phase unwrapping algorithms considering terrain factors by introducing related control variables to the state space model of the Kalman filter [8]. The model error and topographical factors were considered in this algorithm and the experimental results indicated that the algorithm can deal with steep terrain and slope situations effectively. An EKF multi-baseline phase unwrapping algorithm was proposed by Xie and Pi in 2011 [9]. In all these EKF-based phase unwrapping algorithms, a non-linear observation model performs the linearization approximation for the EKF, which always leads to a loss of high-order phase information. To address this, Xie and Pi presented a new phase unwrapping algorithm based on an Unscented Kalman Filter (UKFPU) and they improved and applied it to other algorithms [10][11][12]. Then, a Cubature Kalman Filtering-based phase unwrapping method (CKFPU) was studied by Liu [13]. In practical terms, noisy interferograms are usually complex and their statistic model is not easy to establish. For this, Gao et al. proposed an adaptive unscented Kalman filter phase unwrapping method which is shown to yield the greatest precision and the greatest robustness to noise [14]. In 2019, Gao et al. applied UKFPU to multi-baseline SAR interferograms combining a refined Two-stage Programming Approach [15]. All of these algorithms were founded on the assumption that the noise affecting both the evolution and measurement stages is Gaussian. Martinez-Espla et al. presented a particle-filter phase-unwrapping algorithm which was not subject to this constraint [16]. Then, Xie et al. made improvements on this algorithm [17,18]. Because of the complexity of the measured noise data, the determination of the statistical model is a difficult problem. Although some models have been presented [4,7,10,16,18], they are either approximate or subjective, or mainly emphasize the theoretical nature of the algorithm. So far, the theoretical and applied research of filter-based phase unwrapping algorithms has mainly focused on EKFPU, CKFPU and UKFPU, with a Gaussian model hypothesis for noise, and good results have been achieved [4,7,10,11,13,14].
It has been demonstrated that UKFs, or CKFs, are effective and efficient tools in many nonlinear fields [19][20][21][22][23][24]. Compared to EKFs, they improve the estimation precision when the system model has severe nonlinearities. Nevertheless, for different practical problems (system models), the degree of improvement in precision is different. If the improvement is small or even negligible, the EKF algorithm may be preferred since it strikes a balance between computational complexity and precision; otherwise, a UKF or a CKF will be selected. From [24], it is known that when  is 0, the nonlinear estimation mean and estimation precision of the UKF and CKF methods are the same. However, the CKF has the advantages of fewer parameters, fewer sampling particles, stricter theory, and simpler calculations. Therefore, in this paper we compare the precision of EKF and CKF methods in phase unwrapping applications. Some instructive conclusions are drawn to develop the theory and application of filter-based phase unwrapping methods. This paper is organized as follows: first, a phase unwrapping system model, EKFPU and CKFPU are introduced in Section 2. Section 3 elaborates on EKF and CKF error analysis and analyzes the errors of the nonlinear function in the phase unwrapping model. Tests and analysis are given in Sections 4 and 5. The discussion and conclusion are given in Sections 6 and 7, respectively.

Phase Unwrapping System Model
A simple and effective phase unwrapping system model can be expressed as [13]: represents the estimation error of the phase gradient, which is generally Gaussian white noise. and are the observation errors of the real and imaginary parts of the plural observation, respectively, and are considered to be zero-mean Gaussian white noise [13].

EKFPU Algorithm
From equations (1) and (2), the calculation steps of the one-dimensional EKFPU algorithm are as follows [4,13]: First, we calculate the predicted values of the state vector and its variance matrix: where ,kk x is the phase estimate of the current pixel and , kk P is an estimated variance matrix, whose initial value can be selected based on empirical values.

1,kk 
x is the one-step prediction phase value of the pixel to be unwrapped, whose state covariance matrix is is the linearized observation matrix and I is the unit matrix. The calculation is as follows: is the measured value, and () hx is a phase wrapped nonlinear measurement function.

CKFPU Algorithm
The calculation steps of the one-dimensional CKFPU algorithm are [13]: (1) First, we calculate the initial cubature sampling points and corresponding weights based on the spherical-radial cubature rule: Then, we calculate the propagated cubature points through the nonlinear state equation: (12) to obtain the prediction state and prediction variance: (3) Measurement update First, we perform factorization: (14) and then obtain the cubature points We then calculate the propagated cubature points through nonlinear measurement equations , , (16) Measurement Prediction, Innovation Variance and Covariance Estimation are obtained through: The gain matrix, updated state and covariance are: Since the state prediction equation is linear in our phase unwrapping model, it can be simplified along with the covariance calculation in the CKFPU algorithm, which is consistent with the state prediction method of EKF (Equations. (3) and (4)).
For two-dimensional phase unwrapping, the unwrapping path of EKFPU and CKFPU algorithms can be found in [13]. This paper mainly studies the precision comparison of the EKFPU and CKFPU algorithms. Except 2.2 and 2.3, the other pre-treatment process is identical.

Numerical Stability Analysis of Nonlinear Functions in Phase Unwrapping Models
Here, we analyze the numerical stability of the EKF and the CKF from the point of Taylor expansion of nonlinear functions [24].
We assume an N-dimensional vector~( , ), x N x P , where 1 11 Then, the Taylor expansion of the nonlinear function () gx near the mean is: The true mean of () gxis where 2,3, k  , and  represents the partial conductance of .
The true estimation precision of () gx is For the phase unwrapping model in this paper, the nonlinear functions are the sine and cosine functions and the state vector is one-dimensional. In the following, estimation precision was analyzed for the sine and cosine functions in the EKFPU and CKFPU algorithms.

Estimation Precision of The Sine Function
The Taylor expansion of the sine function sin( ) yx  near the mean is: The mean is:    x  , the difference between the two is no more than 0.46. When 0 P  , 0 W  , so the difference between them approaches zero.

Estimation Precision of Cosine Function
The Taylor expansion of cos( ) yx  near the mean is:  x  , the difference between the two is no more than 0.54. When 0 P  , 0 V  , so the error between them approaches zero. As can be seen from the above, the advantage of CKF over EKF is mainly reflected in the second order terms: sin( ) 2!
x P . Since the nonlinearity of the sine and cosine functions is weak and periodic, CKF appears to be superior to EKF when this part plays a major role in the overall error term. When this part tends to zero (even if P is large, the absolute value of this part may be close to zero due to the periodicity of the sine and cosine function), the superiority of CKF is diminished. Therefore, the smaller the state variance P is, the more similar the estimation precision of the two algorithms. As P increases, two phenomena appear gradually: on the one hand, when the second-order term plays a major role, the CKF is superior to EKF. On the other hand, due to the periodicity of the sine and cosine function, this superiority is not consistent.

Simulated Experiments
This experiment aims to compare the unwrapping precision of CKFPU and EKFPU for two-dimensional simulated data. To ensure a fair comparison, a coherence coefficient was selected to guide the unwrapping path and the maximum likelihood method was used to estimate the phase gradient. The statistical characteristics of the process noise and observation noise can be found in [15]. The same setting will be applied on the following real data experiment. The steep mountainous terrain scene shown in Figure 1 (a) was generated by the peaks function in Matlab, and has large terrain fluctuation and high steepness. The simulation parameters are shown in Table 1. The true phase, interferogram with noise and coherence map for the scene are shown in  Figure 1(b)-(d) (the parameters are same as the previous except the vertical baseline). The noise was generated by software provided by the Delft University of Technology in the Netherlands, reflecting the degree of geometrical incoherence; that is to say, the more serious the geometric incoherence, the greater the noise. It is clear that the data is characterized by complex and dense stripes.  In addition, the phase difference between the EKFPU and UKFPU (named EKF-CKF unwrapped difference) maps and their statistical histograms are shown in Figure 2 (g) and (h). It can be seen from Figures 2 that the unwrapped results are mostly continuous and the errors are small in most areas. The number of pixels with large errors is very small and these pixels are distributed in areas with dense fringe and high noise. Compared to the original interferograms, the fringe details of the rewrapped maps are preserved well and a lot of the noise has been removed. From Figures 2 (e)-(h), we know that in areas where the stripes are sparse and the noise is low (in sparse stripe areas), the unwrapping precision of both methods is about the same. In areas where the fringes change quickly and the noise is large (in dense stripe areas), the precision of both methods has a relatively larger deviation. This is consistent with the previous theoretical analysis. For the simulated data, the larger the phase gradient is, the higher the noise is, and the larger the variance of the state parameters is, thus a greater difference between the EKFPU and CKFPU is obtained. Moreover, although the difference of the two methods in the dense stripe areas is relatively larger, it is not significant compared to the magnitude of the real phase.

TerraSAR-X Data
Two TerraSAR-X images covering the Taiyuan Gujiao Mining Area, acquired on January 4, 2013 and January 15, 2013, were selected as the measurement data. The parameters of the TerraSAR-X are listed in Table 2.

Experiment Area
The target area was located in a coal mining working face. It was active during 4 January 2013 and 15 January 2013 by collecting and analyzing the corresponding coal mining archive and GPS measurement data. There was about 20-30 cm (80-100 radians in the line of sight direction) of subsidence during this repeated cycle. Additionally, most study area was covered with vegetative, which could lead to serious speckle noise. Therefore, the interferogram to test the performance of the EKFPU and CKFPU in this paper is 2×2 looks (Figure 3(a)), and its coherence map is shown in Figure  3 Figure 4 (a)-(d) show the unwrapped and rewrapped phase maps of the two methods. The EKF-CKF phase difference map and its statistical histogram are also shown in Figure 4(e) and (f). As can be seen from Figure 4(a)-(d), the two methods perform very similarly in terms of both unwrapping ability (the unwrapped phase range) and the shape of unwrapped maps. They were similar to the actual subsidence of the working face, although there is a large difference between the unwrapped phase of the two methods and the actual subsidence in the center of the working face. This is mainly due to the decorrelation caused by overly fast subsidence. As can be seen from Figures 4 (e) and (f) (the EKF-CKF unwrapped difference and its histogram), the results of the two methods are basically similar except for a few pixels, while most error areas are distributed within 0.5. From the coherence map, we see that the pixels with large unwrapped difference based on the two methods are all distributed in the regions with poor coherence. This is because the variance of the state parameter is large in the region with poor coherence. According to the theoretical analysis, the greater the variance, the greater the difference between the two methods and greater instability will occur.   regions with good phase quality. Additionally, the precision of the two methods shows a certain difference in the area where the phase quality is poor. It should be pointed out that the severe decorrelation area is generally full of aliased and distorted stripes. So, a satisfying unwrapped result will not be obtained regardless of the method (this situation requires additional research). Therefore, in practical applications an appropriate method can be selected according to the actual situation of the interferogram. An adaptive unwrapping model can also be designed for an interferogram with a large change in phase quality, selecting different methods based on the noise level. It should also be noted that this paper compares the nonlinear unwrapping precision problem caused by the sine and cosine functions. If there are other nonlinear functions in the system equations, additional analysis will be required.