A multiscale approach for detection and mapping differential subsidence using multi-platform InSAR products

Detecting and mapping subsidence is currently supported by interferometric synthetic aperture radar (InSAR) products. However, several factors, such as band-dependent processing, noise presence, and strong subsidence limit the use of InSAR for assessing differential subsidence, which can lead to ground instability and damage to infrastructure. In this work, we propose an approach for measuring and mapping differential subsidence using InSAR products. We consider synthetic aperture radar (SAR) data availability, data coverage over time and space, and the region's subsidence rates to evaluate the need of post-processing, and only then we interpret the results. We illustrate our approach with two case-examples in Central Mexico, where we process SAR data from the Japanese ALOS (L-band), the German TerraSAR-X (X-band), the Italian COSMO-SkyMed (X-band) and the European Sentinel-1 (C-band) satellites. We find good agreement between our results on differential subsidence and field data of existing faulting and find potential to map yet-to-develop faults.

1 Introduction

Several locations around the world experience land subsidence due to groundwater extraction (e.g. Gambolati and Teatini, 2015; Semple et al., 2017). Central Mexico alone, has more than twenty urban areas reported as subsiding (Brunori et al., 2015; Cabral-Cano et al., 2008; Chaussard et al., 2017; Pacheco-Martínez et al., 2015). More importantly, spatial variation of subsidence in several of such locations has led to differential subsidence and, consequently, to ground faulting and infrastructure damage (Avila-Olivera and Garduño-Monroy, 2008; Pacheco-Martínez et al., 2013).

City-scale subsidence patterns and rates have been characterized since the 1940's using multiple geodetic techniques. In recent decades, InSAR measurements have allowed the mapping of large subsiding areas at the scale of entire cities and with frequent observations (e.g. Amelung et al., 1999; Hoffmann et al., 2001). However, several factors limit the application and interpretation of InSAR results for differential subsidence mapping, such as data availability, processing particularities, data integration, and signal interpretation in the presence of strong subsidence.

In this work, we illustrate an approach focused specifically on detecting and mapping differential subsidence based on SAR data from different platforms, elevation data and basic field information. We illustrate our approach with two case studies in Central Mexico, which have different subsidence characteristics and different data availability. The ultimate goal of our approach is to share the expertise we have acquired after studying several cases of differential subsidence for expanding its application to other areas with data availability constrains.

2 Materials and Methods

We propose using SAR data, a digital elevation model (DEM) from the Shuttle Radar Topography Mission (SRTM) mission, and surface field data when available.

We rely on four SAR datasets, of which two cover Aguascalientes and two Mexico City. The two datasets covering Aguascalientes consist on 34 ALOS PALSAR scenes acquired form August 2007 to March 2011 and six TerraSAR-X acquired from December 2009 and September 2012. The two datasets covering Mexico City consist of 144 Sentinel-1 scenes acquired from 2014 to 2017 and 21 COSMO-SkyMED scenes acquired from December 2011 to June 2012. We process the SAR using InSAR techniques, which relies on the existence of at least two satellite acquisitions to perform pixel-wise phase differences to generate a map of interferometric phase differences, so-called interferogram (Fig. 1).

Figure 1 Schematic illustration describing the concept of repeat-pass InSAR. The interferogram was calculated from a COSMO-SkyMed SAR pair acquired on 3 April and 30 June 2012 over Mexico City. Point A is a stable area known as Chimalhuache hill. Point B is located close to Peñón de los Baños. Point C is located over Benito Juárez International Airport. Shaded relief from SRTM data.

We perform single-interferogram processing using the Delft object-oriented radar interferometric software (Doris) (Kampes et al., 2004) and time series processing using either the Stanford Method for Persistent Scatterers (StaMPS) (Hooper et al., 2007) or the Miami InSAR time-series software in Python (MintPy) (Yunjun et al., 2019). Additionally, we perform post-processing of the InSAR velocity maps results using subsidence gradient (Cabral-Cano et al., 2008) or band pass filtering (Solano-Rojas, 2018).

Figure 2 Case example of differential subsidence analysis over Aguascalientes Valley. (a) Satellite image (source: © Esri, Digital Globe, GeoEye, Earthstar Geographics, CNES Airbus, USDA, USGS, ADX, Getmapping, Aerogrid, IGN, swisstopo, and the GIS community) overlaid by locations of groundwater-extraction wells in the area. Red lines represent faults mapped in the city. (b) Slopes map from SRTM topographic data. Red square shows the limits of (d)–(g). (c) Velocity map from 34 ALOS PALSAR scenes acquired form August 2007 to March 2011. (d)–(g) High resolution interferograms of Aguascalientes' urban area from TerraSAR-X data. Each fringe corresponds to 1.55 cm of subsidence.

3 Case studies

3.1 Aguascalientes

The valley of Aguascalientes experiences subsidence rates of up to 100 [mm yr⁻¹] due to overexploitation of its aquifer system (Chaussard et al., 2014). The city is located within a tectonic graben (Fig. 2a). The topographic slopes map of the area shows flat topography in the middle of the graben, which forms the valley (Fig. 2b). Mapped surface faults in the region show that the graben-delimiting faults are orientation roughly NE–SW, in agreement with the topographic slopes. The InSAR velocity map from ALOS data (Fig. 2c) shows a clear subsidence pattern within the graben, in agreement with previous observations (Pacheco-Martínez et al., 2015). However, the graben limits do not exactly coincide with the transition between subsiding and stable areas (see northern portion of the eastern main fault). There seems to be, instead, a better agreement between the topographic slope map and the faults location.

High-resolution single interferograms produced with pairs of TerraSAR-X (TSX) images (Fig. 2d–g) reveal a-few-fringes patterns in the easternmost portion of the study area, as opposed to the many-fringes pattern observed in the centre of the interferograms. Evidently, interferograms comprising a longer time span shows a more complex interferometric pattern, due to the accumulation of subsidence (compare Fig. 2d vs Fig. 2f, for instance). The discontinuities in the phase coincide well with the presence of the field-surveyed surface faulting. As a matter of fact, the interferometric pattern shows additional phase discontinuities that may correspond to not-yet mapped or not-yet-developed surface faults. We interpret that these phase discontinuities.

Figure 3 Case example of differential subsidence analysis over South Mexico City. (a) Velocity map from 144 Sentinel-1satellites acquired from 2014 to 2017. (b) Horizontal gradient calculated from (a). (c) Areas of differential subsidence calculated from high-resolution COSMO-SkyMED data. (d) 10-year worth of subsidence-related shallow faults mapped by the city's government (CENAPRED, 2017). Shaded relief from SRTM data.

3.2 Mexico City

Mexico City is one of the largest and most populated urban areas in the world. The city was built on a highly compressible lacustrine sediment sequence, which has been subjected to fast land subsidence with rates exceeding 350 mm yr⁻¹, in response to aggressive groundwater extraction,

Our Sentinel-1 InSAR velocity map from 2014 to 2017 (Fig. 3a) shows that Sierra de Santa Catarina is stable, while its surroundings subside rapidly ($\sim -$400 mm yr⁻¹). We calculate a gradient velocity map, which shows the highest values around Sierra de Santa Catarina (Fig. 3b). Additionally, we calculate a velocity map from the COSMO-SkyMed constellation data and obtain the boundaries of sharp transition between highly-subsiding and highly-uplifting (Fig. 3c). The areas of pronounced differential displacements agree well with the location of pre-existing, subsidence-related faults identified during 10 years of field surveys (CENAPRED, 2017).

4 Conclusions

The proposed approach shows the capability of multiplatform, multiresolution InSAR for detecting and mapping differential subsidence and surface faulting. The strategic advantage of this approach is its reliance on mostly available products (i.e. SAR data, SRTM DEM, basic field surveys). Medium-resolution results from ALOS and Sentinel-1 satellites depict well the larger-scale subsidence patterns. Detailed faulting cartography in the urban area coincide with phase discontinuities in this interferograms produced with higher-resolution data from TerraSAr-X, with potential for mapping of not-reported faults and improved understanding of the subsurface.

Additionally, post-processing techniques from both low resolution (such as subsidence gradient) and high resolution (such as band-pass filtering) datasets provide a tool for improving the understanding of differential subsidence, as shown in the Mexico City study case. Differential subsidence shows to be particularly large in areas of sharp geotechnical transitions denoted by topographic slope changes.