## Abstract

In this chapter, the recent advancements of differential synthetic aperture radar interferometry (DInSAR) technique are presented, with the focus on the DInSAR-based approaches leading to the generation of three-dimensional time-series of Earth’s surface deformation, based on the combination of multi-platform line-of-sight (LOS)-projected time-series of deformation. Use of pixel-offset (PO) measurements for the retrieval of North-South deformation components, which are difficult to be extracted from DInSAR data, only, is also discussed. A review of the principal techniques based on the exploitation of amplitude and phase signatures of sequences of SAR images will be first provided, by emphasizing the limitations and strength of each single approach. Then, the interest will be concentrated on the recently proposed multi-track InSAR combination algorithm, referred as minimum acceleration InSAR combination (MinA) approach. The algorithm assumes the availability of two (or more) sets of SAR images acquired from complementary tracks. SAR data are pre-processed through one of currently available multi-temporal DInSAR toolboxes, and the LOS-projected surface deformation time-series are computed. An under-determined system of linear equations is then solved, based on imposing that the 3-D displacement time-series have minimum acceleration (MA). The presented results demonstrate the validity of the MinA algorithm.

### Keywords

- ground displacement
- SAR interferometry
- pixel-offset
- minimum acceleration
- geodesy

## 1. Introduction

Over almost last two decades, the differential synthetic aperture radar interferometry (DInSAR) technique [1, 2] has evolved to become nowadays a common practice for the detection and monitoring of Earth’s crust modifications over time, both in academic and operative frameworks. DInSAR is mainly used to detect the temporal evolution of surface deformation through the generation of long-lasting displacement time-series. Several multi-temporal advanced DInSAR algorithms have been proposed in the literature [3–9]. At the present days, the availability of large archives of SAR images collected by several radar instruments operating with different wavelengths, and with complementary side-looking angle geometries has posed the problem to effectively combine the information associated with the different SAR datasets. In particular, the combination of multiple-platform line-of-sight (LOS) displacement time series can improve the ability to retrieve the three-dimensional (East-West, Up-Down, and North-South) components of the on-going surface displacement phenomena. Thus, it allows overcoming the main limitation of DInSAR, which is able only to measure the radar LOS projection of the displacement. This research field is of particular interest; and in the recent years, a few solutions have been proposed [10–22] based on the effective combination of multiple-orbit/multiple-angle DInSAR-based measurements, as well as on merging of DInSAR data products with external measurements (e.g., derived from processing GPS data).

In this chapter, first, the basic rationale of the multi-temporal DInSAR techniques for the generation of Earth’s surface displacements maps (see Section 2) is summarized; and then, the characteristics of the principal combination techniques for multi-track/multi-angle/multi-sensor SAR data recently proposed in the literature are discussed. In Section 3, the focus will be on the algorithm referred to as minimum acceleration combination technique (MinA) [23], which does not require the simultaneous process of very large sequence of differential SAR interferograms. The algorithm consists of a straightforward post-processing stage that involves the analysis of sequences of independently processed (potentially, also with different DInSAR toolboxes) multiple-platform LOS displacement time-series. Noteworthy, the adopted InSAR-combination scheme can be used in a wide context. Real SAR datasets are exploited to demonstrate the validity of the presented algorithm. Experimental results will be shown in Section 4. Conclusions and further perspectives will be provided in Section 5.

## 2. Retrieval of surface displacement components through DInSAR and pixel-offset-based techniques

In this section, the basic principles of differential SAR interferometry (DInSAR) for the detection of Earth’s surface displacement are introduced. Moreover, the potential integration of DInSAR and pixel-offset-based measurement for the measurement of large deformation signals occurring in the case of large ruptures of Earth’s crust is discussed.

### 2.1. Basics of InSAR technique

One of the major applications of the SAR technology is represented by the SAR interferometry (InSAR) technique [1, 2], which relies on the measurements of the phase difference between two complex-valued SAR images gathered from a satellite/airborne platform at different times and from different orbital positions, so as to measure the geomorphological characteristics of the ground (such as the topography height and the modifications of the surface over time, due to earthquakes, volcano eruptions, or other geophysical phenomena). Historically, InSAR has been used for the measurement of ground topography. To understand how InSAR works, let us consider the imaging geometry shown in Figure 1, and let us suppose the first SAR image (referred to as master image) is acquired from the orbital position labeled to as M, and the second image (i.e., the slave image) is taken from the orbital position labeled to as S, located at a distance *b* (usually referred to as baseline) from M. By applying simple geometrical rules, it is possible to uniquely identify the location of each ground targets on the image as well as to get an estimate of their heights (namely, *z*) relative to the reference plane. As evident by inspection of Figure 1, if a same target (namely, T) is observed form two orbital positions (master and slave), with two corresponding ranges its distances from the first (namely, *r*) and the second position (namely, *r* + *δr*) can be correctly measured and the target height can be unambiguously determined. This is obtained by finding out the solution of the following system of two equations (see Figure 1):

where *δr* and *δr* + *r* represent the radar ranges from each antenna to the target point, ϑ represents the radar side-looking angle, *α* the angle of the baseline relative to the horizontal, and *z* is the scatterer height above the flat-earth reference. *H* is the height of the sensor above the reference surface and b is the distance between the antennas, which is referred to as baseline. The ability in successfully reconstructing the unknown topography (z) is strictly dependent on the capability to precisely measure the slant-range difference *δr*, which represents one of the known terms of the system of Eqs. (1) and (2).

To understand how InSAR works, let us consider again the imaging geometry, let us assume the radar system has infinite bandwidth; under this condition the master and slave complex-valued SAR images (pixel-by-pixel) can be mathematically represented as follows:

where γ_{1} and γ_{2} are the complex reflectivity functions of the master and slave scene, respectively, and λ denotes the operative radar wavelength. An interferogram is formed on a pixel-by-pixel basis by starting from two complex-valued co-registered SAR images, as outlined in the following. For each pixel, the phase difference between the two SAR images is computed, by multiplying the first image (master) by the complex conjugate of the second image (slave) and, then, by extracting the relevant phase.

From Eqs. (3) and (4), the interferometric phase is obtained, as follows:

where the asterisk denotes the complex conjugate operation, and the symbol arg[∙] refers to the operator that extracts the phase of a complex number, which is restricted to the [− π, π] interval. However, by assuming the scattering mechanism on the ground is the same (arg[γ_{1}] = arg[γ_{2}]) between the two passages of the sensor over the illuminated area (mutually coherent observations), the measured phase difference

The observed interferometic phase ^{k}, given by:

By considering the standard interferometric configuration depicted in Figure 1 and a few mathematical calculations detailed in [26], it is possible to relate the computed interferometric phase difference to the (unknown) height topography as:

where ∆z is the surface height variation above the flat-earth reference plane,

SAR interferometry nowadays is mostly used for the detection of surface changes occurring over time. In such a case, when a slight change of the surface across the two SAR acquisition times occurs in the imaged scene, an additive term arises in the interferometric phase, which is associated with the radar line-of-sight (LOS) component of the surface displacement, in addition to a phase term that depends on topography. By the inspection of the imaging geometry depicted in Figure 2, we get:

where

In order to measure the interferometric phase term related to the surface displacement, it is thus essential to remove the interferometric phase contribution pertinent to the topography in Eq. (9). Specifically, a differential SAR interferogram is formed by synthesizing the topographic phase from an available digital elevation model (DEM) of the area (using the so so-called back-geocoding process) and by subtracting, pixel-by-pixel, these synthetic fringes to the corresponding InSAR phase, thus leaving only the terms associated with the displacement. Accordingly, computed differential SAR interferograms can be expressed as:

where k ∈ {1, .., M} specifies the considered interferometric pair (master/slave);

_{LOS} is the projection of the relevant displacement vector on the line of sight;

### 2.2. Combination of ascending/descending displacement maps

Availability of InSAR results computed from SAR data obtained from ascending and descending orbits allows also for the separation of the East-West (E-W) and the vertical components of the detected deformation. In particular, for all the pixels that are common to both radar geometries, the sum and the difference of LOS-projected deformations computed for the ascending and the descending orbits are calculated. In particular, the sum of the ascending/descending LOS-projected displacement measurements is related to the vertical component of the ground deformation, whereas the difference of the ascending/descending components gives an estimate for the E-W component of the deformation. Also, because modern spaceborne radar systems are mounted on-board satellites that fly nearly polar orbits, the North-South (N-S) component of the deformation cannot be reliably measured. To explain the rationale for the retrieval of the E-W and the vertical components of the deformation, the following assumptions are made: (i) ascending and descending radar LOS directions (

Geometric scheme to interpret the deformation component is portrayed in Figure 3. Finally, it is worth emphasizing that a fundamental advantage of InSAR technology, with respect to o GPS networks, resides in its dense sampling grid of the displacement field.

### 2.3. The advanced multi-temporal small baseline subset (SBAS) technique

In the following, the small baseline subset (SBAS) algorithm is presented. SBAS was developed in 2002 by a team of researchers from National Council Research (CNR) of Italy [19]. To introduce the rationale of the algorithm, let us consider a set of *Q* single-look-complex (SLC) SAR images acquired by a radar instrument over a certain area of interest. One of these images is selected and assumed as the reference (master) image, with respect to which all available SAR images are properly co-registered. The set is characterized by the corresponding acquisition times {t_{1}, …, t_{Q}} and the inherent perpendicular baselines vector {b_{⊥1}, …, b_{⊥Q}} estimated with respect to the reference image. Application of the standard SBAS technique starts with the generation of a set of *M* of small baseline multi-look (differential) interferograms. On these interferograms, the retrieval of the original (unwrapped) phase signals from the modulo-2π measured (wrapped) phases is carried out. The expression of the *k*th interferometric phase is as follows:

The system of Eq. (13) can be re-organized in a matrix form as:

wherein A is the incidence-like matrix directly related to the selected set of small baseline (SB) interferometric data-pairs. Let us now manipulate the system of Eq. (14) to replace the unknown phase vector *Ψ* with the mean phase velocities between adjacent time acquisitions (see [6]). Accordingly, the new unknowns become:

and the system (14) can be re-formulated as follows:

wherein **B** is the incidence matrix of the linear transformation in Eq. (15), whose detailed expression can be found in [6]. It is worth remarking that, depending upon the selection of small baseline interferometric data pairs, it is possible that SAR data could be separated into some independent subsets separated one another by large baselines. Mathematically, this leads the rank deficiency of matrix **B** and, accordingly, the system (16) can have infinite possible solutions. In order to figure it out a solution among the infinite ones, the singular value decomposition (SVD) method is applied. This allows us to evaluate the pseudo-inverse of the matrix **B**, which gives the minimum norm least-squares (LS) solution of the system (16). In this context, the minimum-norm constraint on the velocity vector **v** allows mitigating the presence of temporal discontinuities in the final result, so as to guarantee a physically sound solution. Finally, an additional integration step is necessary to compute the solution from the estimated vector **v**. After solving the system (16) an estimate of the spurious terms due to the presence of some residual topographic artifacts in the generated interferograms is usually performed [6]. Atmospheric phase screen (APS) is also estimated and filtered out [6]. The quality of retrieved LOS time-series is finally evaluated pixel-by-pixel by calculating the values of the temporal coherence factor, defined in [26].

### 2.4. Pixel-offset (PO) technique for the estimation of large rupture deformations

In areas where large and/or rapid deformation phenomena occur, the exploitation of the differential interferograms, thus the generation of displacement time series, can be however strongly limited by the presence of very high fringe rates, which in turn introduce additional difficulties in the phase unwrapping step and may lead to significant misregistration errors. A pictorial representation of how the interferometric phase degrades as the displacement amount increases is given in Figure 4. Nevertheless, the information on the occurred displacements might be preserved in the amplitude of the investigated data pair, considering the offset of the image’s pixels in range and azimuth directions.

Pixel-offset (PO) is a technique that attempts to find the same distinctive features within sub-scenes of two images relevant to the same target area. In remote sensing, this is usually performed by considering either the Fourier shift theorem [28], or normalized cross-correlation (NCC) algorithms [29]. In SAR applications, PO allows co-registering image pairs or, while already co-registered, identifying the residual shifts related to the motion of distinctive features with accuracies in the order of 1/20th of pixel.

In this scenario, the availability of a sequence of full resolution SAR data pairs, already co-registered, was assumed. For this purpose, the NCC algorithm, which is widely used for SAR images, is applied. This step, carried out on across-track and along-track directions, provides for each pixel an estimation of range and azimuth shifts, finally leading to two offset maps for each data pair. The performed NCC analysis might be evaluated through an estimator of the “goodness” of the retrieved offsets.

At this stage, in order to obtain the corresponding time series, the SBAS strategy is applied (see Section 2.3) by substituting the amplitude-driven information to the phase-driven displacement measurements. This operation was performed to the sequences of smoothed range and azimuth pixel-offset maps. Hence, the small baseline constraint in the data pair’s selection, implicit in the SBAS strategy, is convenient in order to maximize the amount of the exploitable pixels. Notice that this algorithm, which is known in the literature to as pixel offset SBAS (POSBAS) [27], is particularly attractive in the case of large deformation, because it allows us to have an estimate of the North-South deformation components, with accuracy in the order of some centimeter, whereas (as said before) the information related to North-South displacement is almost absent in the conventional DInSAR interferometric phase, being the sensors’ flight trajectories almost parallel to the North-South direction. The results of the performed experiments are shown in Section 4.

## 3. Minimum acceleration (MinA) algorithm

In this section, the minimum acceleration (MinA) combination algorithm [23] used for the extraction of displacement time-series of the Up-Down, East-West, and North-South components is detailed. To introduce the right mathematical framework, let us assume the availability of K independent sets of multiple-platform SAR data collected at ordered times *Qj* distinctive time epochs, respectively. The MinA algorithm [23] requires the preliminary generation from each single SAR data-track of the inherent LOS-projected deformation time-series. This task can be achieved by independently applying either the SBAS technique [6] or other alternative multi-temporal DInSAR approaches [4–9] to the available K SAR datasets. During this preliminary stage, the residual topography as well as the atmospheric phase screen (APS) artifacts corrupting differential SAR interferograms might be estimated and successfully filtered out from the generated LOS-projected displacement time-series [7–16]. The so-obtained LOS-projected time-series of deformation along with other ancillary information, such as the maps of temporal coherence (quantifying the goodness of obtained time-series) as well as the LOS mean deformation velocity maps are geocoded to a common spatial grid of points where to apply the subsequent combination stage. During this preliminary stage the location of high coherent targets is also identified. Henceforth, let

LOS-projected time-series of deformation are expressed with respect to the instants

Let us start by observing that a generic LOS-projected deformation measurement, namely *dLOS*, can be related to its inherent 3-D components, namely [*d*_{E − W}, *d*_{U − D}, *d*_{N − S}]*T*, as [18, 22]:

where *d*^{(j)}(*P*), *j* = 1, …, *K*(for each high coherent pixel) to their unknown 3-D components. This leads writing a system of Q-K independent linear equations with respect to the M = 3(Q-1) unknowns representing the East-West (E-W), Up-Down (U-D), and North-South (N-S) deformation velocities components between adjacent time-acquisitions, namely

where *Γj j* = 1, …, *K* are the values of temporal coherence associated to the *K* different datasets (representing a quality factor of the obtained LOS displacement time-series) and **B** is the incidence-like matrix of the linear transformation that converts LOS-projected measurements into their inherent 3-D components. It is defined by taking into account (1) as:

wherein **B**^{(j)}, *j* = 1, …, *K* is the jth incidence-like matrix of the linear transformation that relates displacement time-series with velocity deformation rates between consecutive time intervals for the jth SAR data set. Derivation of that incidence-like matrix is detailed in the paper [6]. Note that this matrix is the same as the one used in SBAS and discussed in the previous section.

The system (18) and (19) has fewer linear independent equations (*Q-K*) than unknowns (*M*), thus it is an under-determined system that does not admit a unique solution. The matrix **B** of the system has singular values that gradually decay to zero, thus rendering any solution much sensitive to noise level corrupting the vector **B** with SVD and truncating (putting to zeros) the small singular values, in such a way that they do not dominate the solution leading to spurious oscillations. In turn, Tikhonov regularization consists in replacing the solution of the system (2) and (3) by the following minimization problem:

for a suitable value of the regularization parameter α, which can effectively be found using (for instance) L-curve method [31]. The goal of L-curve is to search for a regularization parameter α for which the solution has an optimal balance between the minimization of residual norm **V**‖_{2}.

A similar regularized problem was proposed within the MSBAS algorithm [22], even though in that case the relevant system of linear equations was derived from a sequence of unwrapped multiple-tracks differential SAR interferograms. In the case of MinA algorithm, the regularization problem is achieved differently, by adding to the original system (18) and (19) other equations imposing the condition that the (unknown) 3-D (E-W, U-D, N-S) displacement time-series are with minimum curvature, that is to say the velocity deformation differences (for all the 3-D components) between consecutive time intervals is minimal. Such conditions can formally be expressed by adding to Eq. (2) the following set of 3(*Q*-2) additional equations:

Accordingly, the regularized system of linear equations becomes:

where **C** is the incidence-like matrix related to the minimum-acceleration-regularization linear transformation. The solution of Eq. (6) is finally obtained in the LS sense by applying Truncated SVD to the matrix

As earlier said, a similar system of equations was also derived within the MSBAS algorithm [22] but, at variance with MinA, it relies on searching for a minimum-velocity-norm (MN) solution considering the Tikhonov regularization. Moreover, MSBAS requires the simultaneous inversion of several (a few hundreds or more) unwrapped interferograms for the retrieval of the 3-D components of deformation, and the achieved time-series were however still affected by possible topographic and atmospheric artifacts (although considering several hundreds of interferograms various sources of noise and related artifacts are averaged and only partly filtered out) that need to be subsequently filtered out in a post-processing phase. Accordingly, even though the adopted combination strategy shares some similarities with the MSBAS algorithm, MinA does not require simultaneous processing of multi-platform/multi-angle SAR datasets and can be applied with no restrictions at all on the method (permanent scatterers and/or small baseline-oriented) used for the retrieval of LOS DInSAR time-series, making its field of applicability extremely wide.

Noteworthy, the MinA algorithm can also be extended to include azimuth- and range-pixel-offset (AZPO and RGPO) time-series, as computed using the PO-SBAS method (or alternative solutions). This case is very suitable when large deformation phenomena have to be monitored, being the accuracy of these methods is on the order of 10 cm (or larger). In this case, the strategy here adopted can be extended using AZPO and RGPO time-series of deformation, instead of the LOS deformation measurements, and applying the minimum-acceleration (MA) regularization.

The diagram block of the MinA algorithm is shown in Figure 6.

## 4. Experimental results

This section shows some experimental results obtained by applying the PO-SBAS and the MinA techniques for the estimation of three-dimensional components of terrain displacements.

### 4.1. Sierra Negra PO-SBAS results

The experiments carried out in this case were related to the area of Sierra Negra caldera, which is the most active among shield volcanoes located on Isabela Island, Galapagos Archipelago. Following an Mw 5.4 earthquake, in October 2005 Sierra Negra caldera erupted, interrupting a period of quiescence that lasted almost 30 years. The investigation of several analyses of the Sierra Negra caldera geodetic signals revealed Sierra Negra is almost continuously in an uplift phase, which started in 1992, and accelerated so as to reach about 5 m of cumulative ground displacement before the 2005 eruption. On the other hand, the October 2005 catastrophic event induced a subsidence of the inner caldera of more than 5 m [33].

Due to the large deformation dynamics affecting Sierra Negra caldera, the retrieval of ground displacements using DInSAR is a challenging task. Indeed, the application of conventional SBAS-DInSAR time series analysis on the 2003–2007 Galapagos dataset provides only a partial picture of the deformation field. In particular, a set of 25 ENVISAT SAR images were processed (see [27] for further details). Figure 7(a) and (b) shows the retrieved mean ground velocity maps relevant to the 2003–2007 period. The behavior of the northern flanks of the volcano, being the displacements still in the order of centimeters, is clearly imaged by the SBAS-DInSAR analysis, and it is in agreement with previous studies. However, due to the lack of coherence caused by the large deformation dynamics, the interferometric phase analysis is not able to measure displacements around the crater and inner caldera due to the lack of coherence caused by the large deformation dynamics.

In order to image the spatial and temporal evolution of the deformation in these areas, the SAR amplitude information was exploited. Thus, the PO-SBAS approach was applied to the same data pairs considered for the generation of the SBAS-DInSAR time series. Following the PO-SBAS steps explained, the offsets for each data pair were calculated, and a common mask of “good” pixels was selected by considering only those having a high QI value that were present at least in 70% of the whole dataset. At this stage, the PO-SBAS time-series were generated for each of the selected pixels. The accuracies for the PO-SBAS measurements, relevant to the herein analyzed test-case, were obtained by calculating the standard deviation of the measurements in an area that is known to be stable. Estimated accuracy values are in the order of 1/20th of pixel, in agreement with those expected. However, since the aim of this analysis is to emphasize the areas characterized by large deformations, the pixels whose dynamics are smaller than 1/10th of pixel were masked out. Figure 8 shows the PO-SBAS time-series and the comparison with external GPS measurements available in the area.

### 4.2. Piton de La Fournaise MinA results

To further demonstrate the capabilities of the DInSAR-driven minimum-acceleration combination algorithm, it has been applied for studying the settlements of the area of Piton de La Fournaise (Reunion Islands), which is characterized by the presence of a large volcanic system that erupted on April 3, 2007 and lead to large fractures on the ground. Such volcanic system has extensively been studied [34], however new data can provide additional information on the state of volcanism of the island. The presented experiments are based on processing three independent sets of SAR images collected by the ENVISAT/ASAR (C-band) radar instrument along ascending (48 images) and descending passes (35 images) as well as by the ALOS-1/PALSAR (L-band) sensor (11 images), spanning the 2003–2010 time interval (see Table III in [23]). These three SAR datasets were independently processed by the SBAS-DInSAR processing chain, and corresponding LOS-projected displacement time-series (and mean deformation maps) were generated. The MinA combination method is applied (as a post-processing stage) only to those pixels that remain coherent in all three independent SBAS-DInSAR processing analyses, and this permitted discriminating from the LOS-projected deformations the time-series of the 3D deformation components.

Figure 9(a)–(c) shows the maps of retrieved E-W, U-D, and N-S mean deformation velocity maps, superimposed to a gray-scale SAR amplitude image of the zone common to all the three SAR data-tracks. Also, one point, labeled to as P in Figure 9 and located in the summit area of the crater, was selected. The inherent (combined) E-W and Up-Down deformation time-series relevant to this point are shown in the plots of Figure 10. They make it evident the large cumulative E-W displacement, moving mostly eastward, affects the upper part of the Eastern flank with velocity of about 10 cm/year. This trend is abruptly interrupted by a jump of about 40 cm in correspondence of the April 2, 2007 eruption, which induced a widespread flank movement starting at the time of dike injection to feed an initial eruption, a few days before the main eruptive event; also a significant U-D signal was active even with a more moderate deformation value (around 8 cm).

## 5. Conclusions

In this chapter, a review of some existing DInSAR methods for the retrieval of the 3-D (2-D) deformation time-series is first provided. In particular, we review some recently published methods and then we focus on the MinA method. With respect to previous works, this method has the advantage to be a post-processing algorithm, thus it does not require the simultaneous processing of hundreds of differential SAR interferograms. Information on the quality of LOS-projected deformation time-series (e.g., the temporal coherence maps) as well as the *a priori* identification of very coherent targets is very proficient for the discrimination of the 3-D deformation components. One strength of the algorithm is represented by the opportunity to complement LOS measurements with other external sources of information (such as GPS/leveling data). This technique has primarily been developed as an ultimate extension of the SBAS processing chain; however, it can be used, without any further modification, to work with other general-purpose DInSAR toolboxes. Several examples are provided, thus also clarifying how this method can be easily integrated in the currently available DInSAR toolboxes.

## Acknowledgments

The author would like to thank Giuseppe Solaro, Fabiana Calò, Claudio Dema, Francesco Casu, and Riccardo Lanari, who was the co-authors of the MinA and POSBAS algorithms, which have been summarized in this Chapter. The author is also grateful to Simone Guarino, Fernando Parisi, and Maria Consiglia Rasulo for their technical support. ENVISAT data were provided by the European Space Agency, ALOS-1 data we provided by JAXA within the project entitled “Advanced Interferometric SAR Techniques for Earth Observation at L-band” of the four Research Agreement for the Advanced Land Observing Satellite-2. The DEM of the investigated zones were acquired through the SRTM archive.