Kalman Filter Application to GBSAR Interferometry for Slope Monitoring

The use of Kalman filtering techniques for landslide monitoring has proved effective as a tool for estimating and predicting land displacements. Ground-Based Synthetic Aperture Radars (GBSAR) are popular remote sensing instruments able to provide displacement maps of the investigated area, with submillimeter precision. These instruments outperform other sensors in several respects, such as all-weather and all-day monitoring. However, in some cases, for instance in vegetated scenarios, the displacement is affected by a significant uncertainty due to the decorrelation of the radar signal. In such a case, to retrieve any reliable information, noise must be filtered out using appropriate methods. Given the success of kinematic modeling of landslide movements through Kalman filtering, this technique seems to be the optimal candidate for processing the displacement measured by interferometric GBSAR data. This paper investigates this idea, by applying Kalman Filter to GBSAR measurements acquired in two different campaigns: a landslide monitoring in north Spain, and a sliding glacier monitoring in the Alpes, Italy. A proper initialization of the algorithm parameters is fundamental for a correct application of the Kalman filter. In this work, we present a strategy that exploits information from coherent pixels for tuning the filtering parameters and optimizing the filter performance on areas with low coherence.

with an accuracy of a fraction of the signal wavelength. The 31 achievable high resolution coupled with the return time of 32 few minutes, make them optimal candidates for real time 33 monitoring, which is of great importance on active landslides. 34 These sensors are now widely used for environmental mon-35 itoring, especially in areas not directly accessible with other 36 instruments [5], [6]. 37 The quality of the interferometric radar data strongly 38 depends on the characteristics of the investigated scenario. 39 For instance, vegetated areas are typically characterized by 40 low coherence [7], [8], i.e., the signal is subject to decorre-41 lation effects, which results in noisier interferometric data. 42 In these cases, standard data interpolation or filtering tech-43 niques could lead to an underestimation of the actual result, 44 thus, to a misinterpretation of the data. Therefore, it is of 45 paramount importance to find a way to correctly weight the 46 obtained measurements. To this end, the authors of this article 47 identified the Kalman Filter as a tool to perform dynamic 48 106 Compared to the study illustrated in [20], we apply Kalman 107 algorithms to GBSAR measurements acquired in different 108 scenarios, thus, characterized by a great variety of backscat-109 tered signals. The case studies reported in this paper are those 110 of a sliding Glacier in the Alpes, Monte Rosa, Italy, and a slow 111 active landslide near Formigal, Spain. These two case studies 112 provide useful insights into the application of Kalman filter 113 on interferometric GBSAR data.

114
The paper is organized as follows: section II reviews the 115 basic concepts of GBSAR interferometry, while in section III 116 the mathematical formulation of the Kalman filter is pre-117 sented. Finally, in section IV, the experimental results are 118 presented and discussed. 120 Interferometry techniques allow to determine the displace-121 ment of a target between two radar acquisitions. The displace-122 ment is obtained by manipulating the phase of the complex 123 valued GBSAR images. The phase ϕ of a complex valued 124 image I , relative to a pixel, is the sum of three terms where ϕ dist is a phase contribution which depends on the 127 relative distance between radar and the imaged target, ϕ atm is 128 the phase contribution due to the atmospheric conditions, and 129 ϕ noise is the noise term. By subtracting the phases of images 130 acquired at different times, we get 131 ϕ = ϕ dist + ϕ atm + ϕ noise .

II. GBSAR INTERFEROMETRY
(2) 132 After estimation and removal of the atmospheric contri-133 bution, the phase difference is directly related to the target 134 displacement in the radar line of sight, through the equation 135 where λ is the wavelength associated to the central frequency 137 of the radar signal. Whether the noise term is present or not 138 determine the possibility of correctly measuring the displace-139 ment of a target. In fact, interferometric analysis is usually 140 performed on pixels characterized by high signal quality, for 141 which the noise phase term ϕ noise in (1) can be considered 142 negligible. For these pixels, by cumulating the displacements 143 retrieved using consecutive GBSAR acquisitions, it is pos-144 sible to obtain the target movement over time. However, 145 if the noise term greatly affects the interferometric phase, it is 146 difficult to retrieve the correct displacement information. The 147 presence of noise can be due for instance because of resid-148 uals of the atmospheric phase screen correction or to non-149 compensated phase wrapping. The latter effect is common in 150 vegetated scenarios, where decorrelation can cause a loss of 151 phase information, which makes measurement data difficult 152 to interpret.

153
In order to single out high quality pixels in a GBSAR 154 image, the amplitude dispersion index parameter (D A ) is often 155 used [22], [23]. It can be defined for each pixel of a GBSAR image as where σ A and µ A are respectively the standard deviation 159 and mean value of the pixel signal amplitude time series.

160
It quantifies the pixel signal quality and is used in radar 161 interferometry to preliminary select candidate areas to carry 162 out the analysis. Indeed, this parameter is related to the cor-163 responding coherence value [24], but is easier to calculate.

164
The D A provides us also with an estimation of the

170
Since σ ϕ can be regarded as an estimate of the interferometric

176
In this section the mathematical formulation of the KF algo-177 rithm, relative to our case studies, is reviewed.

178
The FK algorithm allows to estimate/predict an unknown The variable x k to be estimated is called the state variable.

185
In the first phase, the state variable predictionx k , at time k, 186 is determined using the estimate obtained at the previous time Here, A is a matrix that recursively relates the state variable to the selected pixel, at time k. That is, The matrix A represents the kinematic model describing 203 the target movement. In the case studies analyzed in this 204 paper, we are interested in the ground displacement due to 205 the sliding of a glacier or a landslide.

206
In both cases, the ground motion results from the sum 207 of gravitational and frictional forces acting on the selected 208 target of the slope/glacier. Furthermore, in interferometric 209 analysis a single pixel is monitored during a certain period. 210 Therefore, the retrieved displacement is relative to a point of 211 the scenario with fixed coordinates. The point is characterized 212 by a constant inclination. Thus, in most cases, the dynamics 213 will evolve towards a stationary situation, where the ground 214 is subject to constant motion. In this case, the matrix A reads 215 where t is the value of the time interval between successive 217 acquisitions.

218
From a practical viewpoint, the prediction process of 219 the KF algorithm consists of two steps; the computation of 220 the state variable predictionx k , and the computation of the 221 covariance state variable predictionP k . The state variable 222 prediction for the k-th time step, is given by the transition 223 equation with matrix A defined in (7), and x k−1 the final estimate from 226 the previous iterative step. The prediction of the state variable 227 covariance matrix is given by Where t is the time interval between consecutive acquisi-234 tions, and σ 2 w is the covariance of the white noise acceleration 235 process. One of the main challenges to properly apply KF, 236 is the choice of σ 2 w .

238
In the estimation process the measurement z k is used to 239 correct the predictionx k , and to obtain the final estimate 240 x k . In fact, each element z k of the measurement collection, 241 is related to the variable x k , through the so-called measure-242 ment equation where H is the state-to-measurement vector, and e k is a 245 zero-mean gaussian error, called the measurement error, 246 whose variance is R = σ 2 e .

247
The estimation process involves three steps: the computa-248 tion of the so-called Kalman gain K k , the correction of the 249 prediction estimatex k to get the final estimate x k , and the 250 correction of the covariance prediction. The Kalman Gain is 251 given by the following expression Here, R is the covariance of the measurement error e k , and 254P k is the prediction of the state variable covariance matrix, 255 calculated in (9). In our case, the state-to-measurement vector 256 H is equal to The correction of the predicted valuex k is realized by 259 properly weighting the deviation of the measurement z k from 260 the predicted value, with the Kalman gain. Specifically, the 261 final estimate of the state variable is given by 263 Finally, the correction of the covariance prediction is cal-264 culated with the following expression Moreover, the large distance of targets from the sensor, leads 310 to a decrease in the intensity of the backscattered signal, 311 hence, to a decrease of the signal to noise ratio. For these 312 reasons, we estimated the measurement uncertainty directly 313 from the scenario under investigation, as the average of the 314 standard deviation of the displacement of coherent pixels of 315 the scenario. Doing so, we obtained a measurement uncer-316 tainty of the order of 1 mm. Nevertheless, we experienced to 317 substantial deviation by tuning the measurement uncertainty 318 σ e from 0.1 mm to 1 mm, as can be seen in Figure 1. 319 In this figure, we show results of KF filtering of a measured 320 interferometric displacement, with fixed Q, and σ e = 0.1 mm 321 and σ e = 1 mm. This evidences that the algorithm is robust 322 to a fine tuning of the R parameter within this range. 323 We then proceed as follows: by keeping fixed the R value, 324 determined by σ e = 1 mm, we are able to properly tune σ 2 w , 325 and hence, Q values, considering high quality pixels. Once 326 determined the proper choice for Q for the investigated sce-327 nario, we apply KF to noisier pixels in the image. This time, 328 we keep the Q parameter constant, and tune the measurement 329 error variance R, according to the noise that affects the data, 330 based on (5).

332
In order to study the applicability of KF techniques to 333 GBSAR interferometric data, measurements acquired in two 334 VOLUME 10, 2022  Our purpose is to use the KF to extrapolate displace-338 ment information from noisy pixels with low signal quality, 339 by tuning the algorithm parameters based on the information 340 retrieved from high quality neighboring pixels. As already 341 said, to identify high-and low-quality pixels we used the D A 342 parameter. Pixels whose D A is less than 0.15 are identified as 343 stable targets, whose signal has a high quality.   In what follows we show the results obtained for a stable 361 pixel, with D A = 0.11 (target A in Fig. 2), and a pixel with 362 D A = 0.35 (target B in Fig. 1), at about the same range value. 363 We calculated the displacement of the targets over time, using 364 the cumulative interferometric phase according to (3).  line represents the filtered variable, while the red points 368 represent the measured displacement before filtering. In this 369 case, the variable σ w which defines the Q covariance matrix 370 was set equal to 3.6 mm/min 2 . It can be noted that, the 371 blue line representing the filtered displacement over time, 372 partially deviates from the measured values. This is better 373 evidenced by the results shown in Fig. 4, where the deviation 374 of filtered values from the measured ones is represented. The 375 residual clearly shows systematic behaviour. This suggests 376 that this filter does not adequately reproduce real motion. 377 In this case, the model has been assigned too much weight 378 with respect to the measurements. Therefore, we repeated the 379   Fig. 5 and Fig. 6, for σ w = 36 mm/min 2 . In this 382 case, the measurement trend is better reproduced, despite 383 the noise being filtered, and the residuals appear symmetric, 384 which suggests a more reliable estimation.

385
Once optimized the KF parameters, for high quality pixels 386 of the scenario, we repeated the processing for low coherence 387 pixels, by scaling the measurement error variance R, based on 388 the corresponding D A value, according to (5). Fig. 7 shows 389 the results obtained for target B (see Fig. 2). In this case, the measured displacement was noisier, and the filter success-391 fully reduced the noise.

392
Since this scenario was subject to a rapid movement, 393 the displacement trend was already evident even before 394 the filtering operation. However, the filtered displacement 395 is smoother, as the KF is able to filter out non-physical 396 fluctuations.

398
To further test the performance of the proposed method, 399 we analyzed a second series of GBSAR measurements, 400 acquired in Spain, near Formigal, with a C-band system. 401 The system return time was of about one hour. Aim of the 402 experimental campaign was to monitor a slope landslide, just 403 above a road. In this case, the slope was subject to a slow 404 movement, with a velocity of about 1 mm per day. Moreover, 405 the slope was partially covered by grass and low vegetation. 406 These two aspects make this scenario the perfect candidate for 407 our study, as it comprises both high-and low-quality pixels, 408 subject to slow movements, hence, more difficult to measure. 409 Fig. 8 shows the cumulative displacement map, obtained 410 using interferometric techniques, after atmospheric phase 411 compensation, for a time series lasting 8 days. Only pixels 412 with D A < 0.4 are shown in the image.

413
As done for the Macugnaga case study, KF was applied 414 first to high quality pixels, in order to optimize the algo-415 rithm parameters for this scenario. Here, we report the results 416 obtained for target A (see Fig. 8), which has D A = 0.15. 417 Fig. 9 shows the result of the filtering operation once the 418 parameter Q has been optimized for this scenario. In this 419 case, the covariance noise model σ w has been set equal to 420 0.04 mm/min 2 . Fig. 10 shows the residual between the mea-421 sured and filtered displacements.  Once optimized the algorithm parameters, the KF was 423 applied to noisier targets, characterized by a higher D A .

424
As done for the Macugnaga case study analysis, the R param-425 eter was tuned according to the D A value of the pixels.

426
In Fig. 11 and Fig.12   It must be noticed that for short time intervals (for instance 438 time below 20 hours in Fig. 11), the filtered displace-439 ment trend shows a probably non-physical peak. However, 440 for times beyond 20 hours, such peaks are less and less 441 pronounced. This is because the KF algorithm relies on the 442 pixel's displacement history and, in case of noisier mea-443 surements, it takes longer to converge. Therefore, in case 444 of low coherence pixels care must be taken when starting 445 to trust the filtering result. For instance, in case of target B 446 (Fig. 11), measurements acquired in the first 20 hours should 447 be discarded for the analysis. However, this is not a problem 448