Open access peer-reviewed chapter

Generation of Earth’s Surface Three-Dimensional (3-D) Displacement Time-Series by Multiple-Platform SAR Data

Written By

Antonio Pepe

Submitted: 26 April 2017 Reviewed: 28 September 2017 Published: 20 December 2017

DOI: 10.5772/intechopen.71329

From the Edited Volume

Time Series Analysis and Applications

Edited by Nawaz Mohamudally

Chapter metrics overview

1,930 Chapter Downloads

View Full Metrics

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 [39]. 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 [1022] 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.

Advertisement

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):

r + δr 2 = r 2 + b 2 2 r b sin ϑ α E1
z = H r cosϑ E2

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).

Figure 1.

SAR interferometric configuration. The dashed lines show that radar signal paths for the first interferogram pair formed by antennas at M and S.

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:

γ ̂ 1 = γ 1 exp j 4 π λ r E3
γ ̂ 2 = γ 2 exp j 4 π λ r + δr E4

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:

arg γ ̂ 1 γ ̂ 2 = arg j 4 π λ δr + arg γ 1 arg γ 2 E5

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 ψ ~ k (where k identifies a specific interferometric pair of a multiple baseline configuration) depends upon on the range difference δr, only:

ψ ~ k = arg exp j 4 π λ δr E6

The observed interferometic phase ψ ~ k is 2π ambiguous, and the obtained image is called an interferogram. Since the ambiguity of the phase, which is measured modulo 2π, the information on range difference δr is retrieved from the interferogram by applying the phase unwrapping operation [24, 25], thus estimating the inherent absolute interferometric phase ψk, given by:

ψ k = 4 π λ δr E7

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:

ψ k ψ 0 k + ψ k z z = 4 π λ b k sin ϑ 0 k α k 4 π λ b k r sin ϑ 0 k E8

where ∆z is the surface height variation above the flat-earth reference plane, ϑ 0 k is the side-looking angle of each point in the image, assuming zero local height, b k = b k cos ϑ 0 k α k represents the projection of the baseline in the direction perpendicular to the line-of-sight from the radar to the target. The first term ψ 0 k = 4 π λ b k sin ϑ 0 k α k in Eq. (8) accounts for the phase contribution corresponding to the flat-Earth (z = 0), which is the term that is present in the absence of any height elevation on the ground. From Eq. (8), it is clear that the sensitivity of the InSAR measurement can be improved by increasing the baseline. However, the perpendicular baseline, cannot exceed the limiting value (critical baseline) for which the variation of the phase difference across one single ground range resolution element is 2π.

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:

ψ k = ψ k z z + ψ k d Los k d Los k = 4 π λ b k r sin ϑ 0 k z + 4 π λ d Los k E9

Figure 2.

Differential SAR interferometry imaging geometry.

where d Los k represents the projection of the surface displacement vector onto LOS (range) direction pertinent to the kth interferometric pair. Note that the presence of the flat earth phase contribution was neglected, for the sake of convenience.

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:

ψ k = ψ disp k + ψ topo k + ψ orb k + ψ prop k + ψ scat k E10

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

ψ disp k = 4 π λ d Los k represents possible displacement of the scatterer between observations, where dLOS is the projection of the relevant displacement vector on the line of sight;

ψ topo k = 4 π λ b k r sin ϑ 0 k z ~ represents the residual-topography induced phase due to a non-perfect knowledge of the actual height profile (i.e., the DEM errors z ~ );

ψ orb k accounts for possible inaccurate orbital information;

ψ prop k denotes the phase components due to the variation of propagation conditions (pertinent to the change in the atmospheric and ionospheric dielectric constant) between the two master/slave acquisitions;

ψ scat k accounts for change in scattering behavior [13].

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 ( d LOS Asc and d LOS Desc , respectively) lay wholly in the east–z plane and (ii) the sensor side-looking angle is approximately the same along the ascending and descending orbits. Both these assumptions are acceptable. If we refer to the same homologous pixel imaged by the ascending and descending orbits, the E-W and Up-Down components of the measured surface deformation can be estimated from the ascending/descending LOS measurement (e.g., the LOS-projected rates of deformation) as follows:

d LOS East d LOS Desc d LOS Asc 2 sin ϑ E11
d LOS Up d LOS Desc + d LOS Asc 2 cos ϑ E12

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.

Figure 3.

Combination scheme of ascending and descending displacement fields.

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 {t1, …, tQ} 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 kth interferometric phase is as follows:

ψ k = 4 π λ d LOS k 4 π λ b k r sin ϑ 0 k z ~ + ψ orb k + ψ prop k + ψ scat k E13

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

A Ψ = Ψ E14

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:

v = v 1 = Ψ 2 Ψ 1 t 2 t 1 v Q 2 = Ψ Q 1 Ψ Q 2 t Q 1 t Q 2 T E15

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

B v = Ψ E16

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.

Figure 4.

Pictorial representation of the effects due to different amounts of deformations. (a) Plot of the deformation time-series showing the temporal evolution of displacement in the Sierra Negra Caldera, as computed from a sequence of ENVISAT DInSAR interferograms. (b)–(d) Differential interferograms relevant to different deformation regimes and time epochs: when the deformation rate is low the information conveyed in the phase is fully preserved (b). As deformation rate increases the fringe spatial frequency increases (c), and in the occurrence of the large rupture of terrain the fringe rate is so high that the corresponding interferometric phase (d) is completely corrupted by decorrelation noise, thus making the use of phase not effective. This figure is a re-adaptation of Figure 1 originally presented in [27].

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.

Advertisement

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 t j = t 0 j t 1 j t Q j 1 j T j = 1 , , K over the same area on the ground, consisting of 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 [49] 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 [716]. 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 d j = d 0 j d 1 j d Q j 1 j T j = 1 , , K be the geocoded LOS-projected deformation time-series relevant to a generic pixel P that belongs to the group of high-coherent pixels common to all the available K SAR datasets.

LOS-projected time-series of deformation are expressed with respect to the instants t 0 j , j = 1 , , K , which are singularly taken as reference for each dataset, that is to say d 0 j = 0 , j = 1 , , K . Let us now describe the algorithm works, and Q = j = 1 K Q j be the total number of the available SAR images collected at the “whole” ordered times T = j = 1 K t j = T 0 T 1 T Q 1 T .

Let us start by observing that a generic LOS-projected deformation measurement, namely dLOS, can be related to its inherent 3-D components, namely [dE − W, dU − D, dN − S]T, as [18, 22]:

d L O S = d i ^ L O S = sin ϑ cos ϕ d E a s t W e s t cos ϑ d U p D o w n + sin ϑ sin ϕ d N o r t h S o u t h E17

where i ̂ LOS is the LOS-direction versor. Note that θ and ϕ represent the radar side-looking and the satellite heading angles, respectively; the imaging geometries for ascending/descending data-tracks are shown in Figure 5. Extending to our case what originally proposed in [6] and subsequently adapted in [22], let us relate the available LOS deformations 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 V E = [ V E 1 , V E 2 , , V E Q 1 ] T , V U = [ V U 1 , V U 2 , , V U Q 1 ] T and V N = [ V N 1 , V N 2 , , V N Q 1 ] T . This system of linear equations can be expressed using matrix formalism as follows

B V E V U V N = Γ 1 d 1 d 0 1 Γ 2 d 2 d 0 2 Γ K d K d 0 K = d ˜ E18

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:

B = B 1 Γ 1 sin ϑ 1 cos ϕ 1 B 1 Γ 1 cos ϑ 1 B 1 Γ 1 sin ϑ 1 sin ϕ 1 B 2 Γ 2 sin ϑ 2 cos ϕ 2 B 2 Γ 2 cos ϑ 2 B 2 Γ 2 sin ϑ 2 sin ϕ 2 B K Γ K sin ϑ K cos ϕ K B K Γ K cos ϑ K B K Γ K sin ϑ K sin ϕ K E19

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.

Figure 5.

SAR data acquisition geometries for descending (a) and ascending (b) orbits, respectively.

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 d ˜ . It represents a canonical example of a linear discrete ill-posed problem whose meaningful solution can be obtained by replacing the “original” linear system (18) and (19) by a nearby system that is less sensitive to perturbations of the right-hand side of the system, and considers the solution of this new system as a good approximation of the original one. This operation is known as regularization and can be performed using truncated singular value decomposition (TSVD) [30], Tikhonov regularization [31], maximum entropy principle [32]. TSVD practically consists in decomposing the matrix 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:

min V B V d ˜ 2 + α 2 V 2 E20

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 B V d ˜ 2 and the “weight” of the minimum-norm velocity regularization ‖V2.

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:

V E i + 1 V E i i = 1 , , Q 2 V U i + 1 V U i i = 1 , , Q 2 V N i + 1 V N i i = 1 , , Q 2 E21

Accordingly, the regularized system of linear equations becomes:

B C V E V U V N = d ˜ 0 E22

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 Ω = B C . Once the problem in Eqs. (21) and (22) is solved, the 3-D velocity deformation components are independently time-integrated to recover the relevant 3-D displacement time-series.

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.

Figure 6.

MinA diagram block. K independent sets of LOS measurements of the surface ground displacement, as obtained by processing data from K multi-angle/multi-sensor SAR systems also using independent DInSAR toolboxes, are combined (taking into account the associated quality reconstruction maps). Combination is based on the application of minimum-acceleration constraints on the achievable Up-Down, East-West, and North-South time-series of deformation. Use of external data is helpful for better constraining the solutions.

Advertisement

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.

Figure 7.

SBAS-DInSAR results. (a) Mean deformation velocity map of the Galapagos Islands retrieved by applying the SBAS technique and (b) zoom of the study area, (c) Comparison between PO-SBAS and GPS measurements corresponding to the GV01 station.

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.

Figure 8.

Examples of PO-SBAS time-series in azimuth (right) and range (left) directions, respectively. (a)-(f) Comparison between PO-SBAS and GPS measurements in the proximity of selected GPS stations. (g)-(h) AZO and RGO displacement mean velocity maps. The figure is adapted from [27].

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).

Figure 9.

Geocoded maps of the Up-Down (a), east–west (b), and north–south (c) mean velocity deformation.

Figure 10.

DInSAR results retrieved for the Piton de la Fournaise study area. Zoom view of the Up-Down (a) and East-West (b) mean deformation displacement maps. (c) and (d) are the MinA-driven time-series obtained by combining the LOS time-series for the Up (c) and East-West (d) components, respectively.

Advertisement

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.

Advertisement

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.

References

  1. 1. Massonnet D, Feigl KL. Radar interferometry and its application to changes in the earth’s surface. Reviews of Geophysics. 1998;36:441-500
  2. 2. Massonnet D, Rossi M, Carmona C, Adragna F, Peltzer G, Feigl K, Rabaute T. The displacement field of the landers earthquake mapped by radar interferometry. Nature. 1993;364:138-142
  3. 3. Ferretti A, Prati C, Rocca F. Permanent scatterers in SAR interferometry. IEEE Transactions on Geoscience and Remote Sensing. 2001;39(1):8-20
  4. 4. Werner C, Wegm¨uller U, Strozzi T, Wiesmann A. Interferometric point target analysis for deformation mapping, In: Proceedings of the Geoscience and Remote Sensing Symposium; 21-25 July 2003; Toulouse, France; 2003;7:4362-4364
  5. 5. Hooper A, Zebker H, Segall P, Kampes BM. A new method for measuring deformation on volcanoes and other natural terrains using InSAR persistent scatterers. Geophysical Research Letters. 2004;31(23):L23611. DOI: 10.1029/2004GL021737
  6. 6. Berardino P, Fornaro G, Lanari R, Sansosti E. A new algorithm for surface deformation monitoring based on small baseline differential SAR interferograms. IEEE Transactions on Geoscience and Remote Sensing. 2002;40(11):2375-2383
  7. 7. Usai S. A least squares database approach for SAR interferometric data. IEEE Transactions on Geoscience and Remote Sensing. 2003;41(4):753-760
  8. 8. Mora O, Mallorqui JJ, Broquetas A. Linear and nonlinear terrain deformation maps from a reduced set of interferometric SAR images. IEEE Transactions on Geoscience and Remote Sensing. 2003;41(10):2243-2253
  9. 9. Crosetto M, Crippa B, Biescas E. Early detection and in-depth analysis of deformation phenomena by radar interferometry. Engineering Geology. 2005;79(1/2):81-91
  10. 10. Wright TJ, Parsons BE, Lu Z. Toward mapping surface deformation in three dimensions using InSAR. Geophysical Research Letters. 2004;31. DOI: L01607, 10.1029/2003GL018827
  11. 11. Gray L. Using multiple RADARSAT InSAR pairs to estimate a full three-dimensional solution for glacial ice movement. Geophysical Research Letters. 2011;38(5):L05502
  12. 12. Gudmundsson S, Sigmundsson F, Carstensen J. Three-dimensional surface motion maps estimated from combined interferometric synthetic aperture radar and GPS data. Journal of Geophysical Research. 2002;107(B10):2250-2264
  13. 13. Spata A, Guglielmino F, Nunnari G, Puglisi G. SISTEM: A new approach to obtain three-dimensional displacement maps by integrating GPS and DInSAR data. Presented at the FringeWorkshop; November 30–December 4 2009; Frascati, Italy; 2009
  14. 14. Fialko Y, Simons M, Agnew D. The complete (3-D) surface displacement field in the epicentral area of the 1999 M(w)7.1 Hector mine earthquake, California, from space geodetic observations. Geophysical Research Letters. 2001;28(16):3063-3066
  15. 15. Fialko Y, Sandwell D, Simons M, Rosen P. Three-dimensional deformation caused by the Bam, Iran, earthquake and the origin of shallow slip deficit. Nature. 2005;435(7040):295-299. DOI: 10.1038/nature03425
  16. 16. Jun H, Wei LZ, Jun ZJ, Chong RX, XiaoLi D. Inferring three-dimensional surface displacement field by combining SAR interferometric phase and amplitude information of ascending and descending orbits. Science China Earth Sciences. 2010;53(4):550-560. DOI: 10.1007/s11430-010-0023-1
  17. 17. Shirzaei M. A seamless multitrack multitemporal InSAR algorithm. Geochemistry, Geophysics, Geosystems. 2015;16:1656-1669. DOI: 10.1002/2015GC005759
  18. 18. Hu J, Ding X, Li Z, Zhu J, Sun Q, Zhang L. Kalman-filter based approach for multisensor, multitrack, and multitemporal InSAR. IEEE Transactions on Geoscience and Remote Sensing. 2013;51(7):4226-4239
  19. 19. Manzo M et al. Surface deformation analysis in the Ischia Island (Italy) based on spaceborne radar interferometry. Journal of Volcanology and Geothermal Research. 2006;151:399-416
  20. 20. Gourmelen N, Amelung F, Casu F, Manzo M, Lanari R. Mining related ground deformation in Crescent Valley, Nevada: Implications for sparse GPS networks. Geophysical Research Letters;34:L09309. DOI: 10.1029/2007GL029427
  21. 21. Ozawa T, Ueda H. Advanced interferometric synthetic aperture radar (InSAR) time series analysis using interferograms of multiple-orbit tracks: A case study on Miyake-jima. Journal of Geophysical Research. 2011;116:B12407. DOI: 10.1029/2011JB008489
  22. 22. Samsonov S, d’Oreye N. Multidimensional time-series analysis of ground deformation from multiple InSAR data sets applied to Virunga Volcanic Province. Geophysical Journal International. 2012;191:1095-1108
  23. 23. Pepe A, Solaro G, Calò F, Dema C. A minimum acceleration approach for the retrieval of multiplatform InSAR deformation time-series. IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing, August 2016;9(8)
  24. 24. Gens R. Two-dimensional phase unwrapping for radar interferometry: Developments and new challenges. International Journal of Remote Sensing. 2003;24(4):703-710
  25. 25. Chen CW, Zebker HA. Phase unwrapping for large SAR interferograms: Statistical segmentation and generalized network models. IEEE Transactions on Geoscience and Remote Sensing. 2002;40(8):1709-1719
  26. 26. Pepe A, Lanari R. On the extension of the minimum cost flow algorithm for phase unwrapping of multi-temporal differential SAR interferograms. IEEE Transactions on Geoscience and Remote Sensing. 2006;44(9):2374-2383
  27. 27. Casu F, Manconi A, Pepe A, Lanari R. Deformation time-series generation in areas characterized by large displacement dynamics: The SAR amplitude pixel-offset SBAS technique. IEEE Transactions on Geoscience and Remote Sensing. 2011;49(7):2752-2763
  28. 28. Michel R, Avouac JP, Tabouri J. Measuring ground displacement from SAR amplitude images: Application to the Landers earthquake. Geophysical Research Letters. 1999;26(7):875-878
  29. 29. Scambos T, Dutkiewicz M, Wilson J, Bindschadler R. Application of image cross-correlation to the measurement of glacier velocity using satellite image data. Remote Sensing of Environment. 1992;42:177-186
  30. 30. Hansen PC. The truncated SVD as a method for regularization. BIT. 1987;27:354-553
  31. 31. Rezghi M, Hosseini SM. A new variant of the L-curve for Tikhonov regularization. Journal of Computational and Applied Mathematics. 2009;231:914-924
  32. 32. Skilling J, Gull SF. Algorithms and applications, in maximum-entropy and bayesian methods in inverse problems. In: Fundamental Theories of Physics. Vol. 14. Netherlands: Springer; ISBN: 978-94-017-2221-6
  33. 33. Geist D, Harpp K, Naumann T, Poland M, Chadwick W, Hall M, Rader E. The 2005 eruption of Sierra Negra volcano, Galapagos. Bulletin of Volcanology. 2008;70(6):655-673
  34. 34. Peltier A, Bianchi M, Kaminski E, Komorowski J-C, Rucci A, Staudacher T. PSInSAR as a new tool to monitor pre-eruptive volcano ground deformation: Validation using GPS measurements on piton de la Fournaise. Geophysical Research Letters. 2010;37. DOI: L12301, 10.1029/2010GL043846

Written By

Antonio Pepe

Submitted: 26 April 2017 Reviewed: 28 September 2017 Published: 20 December 2017