Ground deformation and fissure activity in Datong basin, China 2007–2010 revealed by multi-track InSAR

abstarct The Datong basin is one of the complex geologic environments in China. Several faults are distributed in the basin, and a number of hidden faults have been discovered within the basin. Ground fissures (GFs) and land subsidence occur in this region. In the past, ground subsidence at Datong basin has been studied using groundwater level records, intermittent small-scale levelling measurements and field investigations. However, due to the low spatial resolution and temporal sampling of these records, detailed ground deformation within the basin remains poorly understood. In this study, we exploited Envisat ASAR and ALOS PALSAR images to map the ground deformation between 2007 and 2010. Most of the deformation within the basin occurred at rates of <20 mm/year. In addition to the subsidence due to coal mining in the western mountainous areas, at least four large subsidence zones could be detected. Our analysis indicates that the ground subsidence is affected by groundwater extraction and the distribution of faults. Additionally, the vertical deformation is much larger than the deformation in the east–west direction. Obvious deformation gradients occurred on individual GFs and the activity of the GFs varied spatially. Finally, Interferometric synthetic aperture radar (InSAR) observations related to groundwater exploitation and GFs were modelled and discussed. Abbreviations: GFs: ground fissures; SAR: Synthetic Aperture Radar; ASAR: Advanced Synthetic Aperture Radar; ALOS: Advanced Land Observing Satellite; PALSAR: Phased Array type L-band Synthetic Aperture Radar; InSAR: Interferometric Synthetic Aperture Radar; ESA: European Space Agency; JAXA: Japan Aerospace Exploration Agency; SBAS: Small Baseline Subsets; SRTM: Shuttle Radar Topography Mission; DEM: Digital Elevation Model; DORIS: Dopper Orbitography and Radiopositioning Integrated by Satellite instrument; MCF: minimum cost flow; LOS: line-of-sight2D/3D: two-dimensional/three-dimensional


Introduction
Datong basin, which contains the largest deposit of coal in China, is located in the northern part of Shanxi Province (Figure 1). The urban area of Datong is surrounded on three sides by mountains, has passed only to the east and southwest, and has elevations that generally increase from southeast to northwest. The well-known Datong Volcanic Arc lies in the east of the Datong basin, and the Kouquan fault and Liulengshan piedmont fault are located in the western and eastern parts of Datong, respectively ( Figure 1). Several moderately strong earthquakes have occurred in the Datong area in its history (Wang et al. 2001;Xie et al. 2003). Since the Datong-Yanggao Earthquake (Ms 6.1) on 18 October 1989, a number of low-magnitude earthquakes have occurred in this region, indicating potentially damaging active seismicity in the area (Figure 1). In addition, ground fissures (GFs) and ground subsidence occur in this region. Since the 1980s, the groundwater has been greatly over- exploited due to the need for industrial development and domestic water, which has led to the occurrence of ground subsidence and GFs (Liu et al. 1999;Zan 2006). At present, there are more than 10 GFs in Datong with a total length of more than 34.5 km (Li and Yuan 2002;Ren et al. 2004).
As Datong is an important city in northwest China and many geological phenomena occur here, a number of studies on ground subsidence and GF activities have been conducted based on underground water level records, levelling measurement, and field investigation (Liu et al. 1999;Zan 2006). However, because of the low spatial resolution of the underground water level records and levelling measurement, detailed ground deformation remains poorly understood. The Interferometric synthetic aperture radar (InSAR) technique, which has been proven to be a useful tool for measuring many geophysical phenomena (Massonnet et al. 1993;Zebker et al. 1994;Lu et al. 2003;Zhao et al. 2012;Ji et al. 2015;Yang et al. 2016), have been used to study the ground deformation of Datong basin with Envisat ASAR images in 2014 . However, what is the horizontal displacement field in the basin? What is the relationship between the GFs and hidden faults in the region? These problems are still not well understood.
The aim of this study was to contribute to a better understanding of the ground deformation and GF activities in the Datong basin. In this work, we applied the timeseries InSAR technique to Envisat ASAR images (C-band) from the descending path and ALOS-1 PALSAR images (L-band) from the ascending path to map the ground deformation in Datong between 2007 and 2010. Ground deformation characteristics in two dimensions were obtained, and cross-correlations among the regional ground subsidence, fault activity, GFs, and underground water level were evaluated.
This article is organized as follows. First, the data collection and processing are described in Section 2. Then, the InSAR results from C-and L-band are reported and a discussion is presented in Sections 3 and 4, respectively. Finally, modelling results are discussed in Section 5, and a conclusion is provided in Section 6.

SAR data
The SAR images covering the Datong area were acquired on the descending Envisat ASAR from the European Space Agency (ESA, Paris, France) and the ascending ALOS-1 PALSAR from the Japan Aerospace Exploration Agency (JAXA, Tokyo, Japan). A total of 25 Envisat ASAR images (period of May 2007 to September 2010) and 13 ALOS-1 PALSAR images (period of June 2007 to May 2010) were collected. The coverage of the SAR images is shown in Figure 1. The radar line-of-sight look angles at scene centre were 22.8 and 38.7 degrees off nadir for ASAR and PALSAR,  13-2010.9.19 2007.6.14-2010.5.7 respectively, which means that there was a different sensitivity towards the vertical deformation. The azimuth angles were 167.9 and 10.0 degrees counter clockwise from north for ASAR and PALSAR, respectively. All of the information pertaining to the SAR images is listed in Table 1.

InSAR processing
The interferometric processing was performed with GAMMA software (ERS-ENVISAT Symposium, Gothenburg, Sweden) (Werner et al. 2000). The small baseline subsets (SBAS) differential-InSAR approach was utilized to determine the time series deformation (Berardino et al. 2002;Lanari et al. 2004). A Shuttle Radar Topography Mission (SRTM) C-band digital elevation model (DEM) with a resolution of 1 arcsecond (30 m) was used as an external DEM to remove the topographic phase from the differential interferograms. ESA DORIS precision orbits data for the Envisat satellite were employed to remove the reference phase from the differential interferograms. To eliminate the impact of temporal and spatial decorrelation (Rosen et al. 2000;Lu and Dzurisin 2014), we set spatial baselines of <650 and <2500 m and temporal baselines of <400 and <700 d for the C and L bands, respectively. Finally, 81 interferograms from C-band and 35 interferograms from L-band were generated ( Figure 2). The co-registration procedure was conducted with the aid of a DEM, which improved the coherence for the image pairs with a larger spatial baseline (Lu and Dzurisin 2010). To improve the signal-to-noise ratio, interferograms were multilooked by 2 looks in range and 10 looks in azimuth for Envisat ASAR (40 Â 40 m) and 4 looks in range and 10 looks in azimuth for ALOS-1 PALSAR (30 Â 30 m). The adaptive filter based on the local fringe spectrum was used twice for each differential interferogram (Goldstein and Werner 1998) with window dimensions of 128 Â 128 and 32 Â 32 pixels. The minimum cost flow (MCF) method was carried out to retrieve the unwrapped phase. Low-pass filtering, followed by temporal high-pass filtering, was implemented in the spatial domain to identify undesired atmospheric artefacts (Berardino et al. 2002;Lu and Dzurisin 2014). This operation also allowed us to detect possible orbital ramps. The atmospheric artefacts and the orbital ramps were removed following their identification. Most other processing steps followed the conventional SBAS InSAR approach.

Results from Envisat ASAR and ALOS PALSAR
According to the aforementioned processing procedure, time-series ground deformation maps in line-of-sight (LOS) over the Datong basin were obtained during the period of 2007-2010. Figure 3 shows the annual deformation rate calculated from the Envisat ASAR data (Figure 3(a)) and from the ALOS-1 PALSAR data ( Figure 3(b)). The positive and negative signs represent displacement towards and away from the satellite in LOS, respectively. The reference region of the deformation is the ancient town of Datong, which was considered to be stable based on a previous study (Chang'an University 2014). We averaged the deformation values on the coherent points within a 2 Â 2 km region (black box in Figure 3) to obtain the reference value. Then, the reference value was subtracted from the InSAR results. Based on the LOS deformation results from the descending Envisat (Figure 3(a)) and the results from the ascending ALOS (Figure 3(b)), the deformation characteristics are consistent between the two measurements. The deformation rate within the basin was generally small, and deformation occurred mainly in the internal region of the basin, except for the western mountain coalfield area. Most of the deformation occurred at rates of less than 20 mm/year, and the deformation was bounded by faults. At least four large subsidence zones could be detected, including Baima town, Shizhuang-Jichechang, the high-tech development zone, and Yulin County. Several local settlement locations in the lower-left corner of Figure 3(a,b) correspond to coal production plants.

InSAR accuracy assessment
The previous levelling measurement was conducted prior to the start of SAR acquisitions used in this study, and no measurement data covering our research timeframe was collected. As a result, an internal validation strategy was adopted to verify the accuracy of our InSAR results (Qu et al. 2014). After the LOS measurements were converted to vertical displacements with their corresponding incidence angles, the common part of the displacement maps between C-and L-band was selected for internal precision evaluation. Vertical deformation from the C band was subtracted from the L band measurement. Figure 4 shows a histogram of the differences. The distribution of the differences is basically a normal distribution. The standard deviation is about 5.5 mm, and the mean of the differences is À1.6 mm; large differences occurred mainly in the coal mine area in the western mountainous region and in the farmland within the basin. Basically, the L-band and C-band InSAR results are consistent.

Time series of deformation
Time series of the deformation was produced from the C-and L-band interferograms based on the SBAS InSAR approach. To analyse the deformation characteristics of the time series, four points, mountain (A), town of Baima (B), Jichechang (C), and high-tech development zone (D), were extracted (shown as triangles in Figure 3). We produced the final deformation result for each of the selected regions by averaging the deformation values of all coherence points within a box of 200 Â 200 m. The results are shown in Figure 5. It is worth noting that the two measurements were projected into the vertical direction before the comparison. We used the nearest two Envisat results in time to interpolate the Envisat deformation values corresponding to the ALOS acquisition dates. After that, we calculated the differences between the interpolated Envisat deformation estimates and the corresponding deformation values from ALOS-1 PALSAR. The standard deviations of the differences at A, B, C, and D are 1.1, 1.8, 2.9, and 3.4 mm, respectively. They have good consistency. The difference between the two results may be due to the difference in the number of coherence points within the selected region from Envisat and ALOS datasets. The interpolation of Envisat measurements to match the ALOS-1 acquisition dates might also contribute the difference. Based on Figure 5, the cumulative subsidence values around the town of Baima, Jichechang, and the high-tech development zone within nearly threeand one-half years were 60, 45, and 55 mm, respectively. The point in the mountains was basically stable. It is worth noting that the presence of horizontal displacement also has an impact on the accuracy of the converted vertical displacement.

Estimation of horizontal and vertical displacements
Due to the different perspectives on ground movement, SAR images from descending and ascending tracks can be used to recover 2D and/or 3D surface displacement fields (Wright et al. 2004;Bechor and Zebker 2006;Samsonov and Tiampo 2006;Jung et al. 2011). Because the deformation in LOS was not sensitive to displacements in the north-south direction, with the near-polar orbits of the SAR satellites, 2D deformation velocity fields in the east-west and vertical directions could be investigated with the Envisat ASAR (descending) and ALOS-1 PALSAR (ascending) data. In fact, this is reasonable for the Datong basin. Previous studies have shown that the displacement in the north-south direction is very small (Qu et al. 2013;Chang'an University 2014). The measurements equation can be written as where ½V east ; V up T is the 2D deformation vector in the east and up directions; ½V envi los ; V alos los T is the InSAR displacements in LOS form Envisat ASAR and ALOS-1 PALSAR, h 1 and h 2 are the radar incidence angles, and a 1 a 2 are the satellite heading angles. As a result, the 2D deformation velocity fields in the east-west and vertical directions can be recovered. The 2D deformation rate field from the common coverage between the Envisat ASAR and ALOS-1 PALSAR images was obtained ( Figure 6). The displacement in the east-west direction was very small compared to that in the vertical direction, and the horizontal displacement occurred mainly in the mountainous coal mine area where the vertical deformation is large (Figure 6(a)). Different directions of horizontal displacement occurred between two sides of the mining subsidence centre. The vertical displacement (Figure 6(b)) presented a deformation pattern similar with those of Figure 3(a,b), indicating that vertical motion dominates the deformation field in the Datong basin. This result is consistent with a previous investigation (Chang'an University 2014).

Mining subsidence
The Datong coalfield is one of the biggest coal production mines in China. The coalfield consists of the Jurassic coalfield and the Carboniferous-Permian coalfield with a total area of 2500 km 2 . In 2005, coal production exceeded 100 million tons. The large-scale and long-term underground extraction has caused serious environmental problems, such as mine-induced subsidence, landslides, earthquakes, buildings damage, etc. Previous studies have focused on coal mining subsidence monitoring in Datong Coalfield Yang et al. 2017). According to our InSAR results, the mining-induced subsidence rate is much larger than 20 mm/year. In addition, there are incoherent areas caused by large deformation gradients in the coal mining area, where deformation rates are likely much larger. According to our InSAR results, the mine-induced subsidence area reached 208 km 2 during the monitoring period.

Ground deformation and faults
The basement of the Datong basin has different faults (including hidden faults) in different directions, and the main fault directions are north-northeast and northeast (Liu 1980;Cheng and Yang 1996). According to the InSAR results (Figure 3), the distribution of ground subsidence is bounded by the faults. For example the subsidence funnels at the western side of the town of Baima, Shizhuang-Jichechang, and Yulin County are clearly affected by the Kouquan fault. The subsidence funnel of the hightech development zone is bounded by the Datong-Yanggao fault. Affected by the buried faults, the Shizhuang-Jichechang settlement centre appears to have a meteorlike shape (Figures 4 and 6). The profile along P1-P1' (Figure 6)

Ground deformation and groundwater
Many studies have shown that groundwater exploitation can cause an aquifer system to have reduced pore water pressure (Holzer and Johnson 1985;Woldai et al. 2009;Calderhead et al. 2010). As a result, the aquifer becomes compressed and ground subsidence occurs. The Datong basin is an inland basin, and the main water resource is groundwater. Since the 1980s, groundwater extraction has been greater than recharge, causing the groundwater level to decrease year by year (Liu et al. 1999;Zan 2006). According to a report in Datong Daily (Jin and Shuang 2007), the water level was less than 10 m below the surface in the 1970s but is now more than 60 m below. The groundwater level was particularly and obviously decreased in the area where industrial and agricultural productions are concentrated. The rate of descent was 0.52-1.24 m/year, and the well density was 11-19 wells per km 2 . The rate of descent of the groundwater level in the centre of the funnel was 3.67-12.4 m/year, resulting in the local aquifer being in a basically dry state (Jin and Shuang 2007).
At present, three groundwater fall funnels have formed in Datong. They are the Gudian-Baima, Shizhuang, and Zhijiabao funnels ( Figure 6) (Liu et al. 1999;Zan 2006). According to the InSAR results, the Gudian-Baima and Shizhuang groundwater fall funnels have good correspondence with ground settlement. However, this is not the case for the Zhijiabao groundwater fall funnel. Due to the stratigraphic structure, compressed layer thickness and other factors also affect the occurrence of ground subsidence at this location (Sun et al. 2007). We infer that the mechanical properties of the Zhijiabao area may differ from those of the other two locations. According to the underground water level records at the position of the circle with a cross in Figure 6  Due to the good correspondence between the ground deformation in the town of Baima and the groundwater fall funnel, we consider that the surface deformation in this area is caused mainly by the exploitation of groundwater. A uniformly opening sill embedded in an elastic half-space was used to explain the vertical deformation (Okada 1985;Sun et al. 2017;Kim and Lu 2018). Eight parameters were defined in the sill model: length, width, depth, strike, dip, opening and location (two parameters). The optimal parameters and associated uncertainties of the model were estimated by the downhill simplex method and Monte Carlo simulations (Press et al. 1992). The root mean square error (RMSE) between the InSAR measurements and modelled results was employed to measure the goodness of the prediction.
Because the main depth of groundwater extraction in Datong is within 200 m of the surface (Chang 2004;Zan 2006), we fixed the depth in our model at 180 m. The remaining parameters were estimated freely. After multiple iterations, the best-fit values for the model parameters (shown in Table 2) were well constrained. Figure 8 presents the InSAR measurements (a), the sill model results (b), and the residual (c) between the InSAR observations and the model deformation. Although some local residuals exist, the model fits the InSAR-derived deformation reasonably well. The best-fit model sill for the deformation is 6.3 km long, 2.7 km wide, has a strike of N39 E, and has an opening rate of 14.4 mm/year. There are still some local settlement areas that cannot be modelled. The residual settlement may be caused by localized extraction of ground water for domestic or industrial applications.

Ground fissure activities
GFs are one of the major geological hazards in Datong, and several studies on the GFs in Datong have been published (Li and Yuan 2002;Ren et al. 2004;Zan 2006;Deng 2007). According to these investigations, the occurrence of the fissures in Datong was followed by a lowering of the groundwater level. The earliest groundwater drawdown funnel formed in the Jichechang area in the late 1980s. At the same time, several GFs appeared in this region. To 1997, more than 10 GFs had formed in Datong ( Figure 6) (Li and Yuan 2002;Ren et al. 2004). According to the locations and the formation times of GFs and groundwater fall funnels, it is apparent that groundwater overexploitation has impacted GF activity.
Several published studies have reported that the formation of GFs is related to hidden faults (Peirce 1979;Lu et al. 2013;Peng et al. 2013;Brunori et al. 2015). Upward expansion of hidden faults can be promoted due to a combination of tectonic factors and human factors. If the absorption effect of the overlying strata cannot resist the  expansion of the hidden faults, GFs will be generated. As a result, the direction of a GF is generally consistent with the direction of its associated hidden fault. We investigated whether this was true for the case of Datong. According to the Geological Environment Monitoring Center of Shanxi Province (2006), there are 18 hidden faults in the Datong urban area (Figure 9). The blue line in Figure 9 outlines the investigation area of hidden faults. The directions of the faults are in the range of N30-N80 E, and the dominant direction is roughly N50 E, which is the same as that of the GFs. We superimposed the positions of the hidden faults onto the GFs, and most of the GFs corresponded to the hidden faults ( Figure 9). This result indicates that the formation of the GFs may have been controlled by the hidden faults. The profile along P2-P2' (Figure 6), which is perpendicular to the GFs, was extracted, and the result is shown in Figure 10. Five GFs are crossed by the profile. No useful InSAR results were collected near GF No. 4. Obvious deformation gradients occur near GFs No. 1 and No. 3,and gentle deformation gradients occur near GFs No. 2 and No. 5. In other words, the activity of the GFs is different at different locations. The activity of GFs is affected by tectonic factors, such as underlying concealed faults and surrounding seismic activity, and non-tectonic factors, such as groundwater exploitation, stratigraphic lithology, etc. (Lu et al. 2013). The difference in the activity of GFs at Datong also reflects the difference of tectonic and non-tectonic factors between the GFs. We collected the monitoring data of GF at Jichechang. The positions of GF monitoring stations are shown as two opposing triangles in Figure 6. The monitoring data were obtained quarterly from 2007 to 2009 by levelling measurement, which was used to monitor the vertical relative displacement between two sides of the GF. Two measuring devices are 15 m apart and are located on either side of the GF. To compare these measurements with InSAR results, time series deformation on two InSAR measurement points (S1 and S2 in Figure 6) that are close to the levelling stations are plotted ( Figure 11). Time series vertical deformation from Envisat ASAR images was used in the comparison. The InSAR-derived relative deformation difference between the two sides of GF has a similar increasing trend as that from levelling ( Figure 11). The correlation coefficient of the deformation difference derived from InSAR and that from levelling is 0.58, and the RMSE between InSAR and levelling measurements is 3.3 mm. The reason for the moderate correlation is that the two InSAR measurement points are about 500 m apart while it is 15 m between the levelling stations. As a result, the relative height difference from InSAR is greater than that from the levelling.
A creeping dislocation buried in an elastic homogeneous half-space can be used to explain the geodetic observations near a vertical strike-slip fault (Vergne et al. 2001;Biggs et al. 2008). The model that is usually used to describe the ground deformation across an active fault can be expressed as the following equation (Savage and Burford 1973): where V x ð Þ is the fault-parallel velocity, x is the distance from the fault, S is the fault slip rate, D is the locking depth, and a is a static offset.
The GFs in Datong are tectonic ones, which have the characteristics of associated hidden faults. The above model can be used to study the deformation across these GFs. As an example, we extracted a profile along F-F' (Figure 6), where there is an obvious deformation gradient. The deformation profile is shown in Figure 10. A gradient of about 6 mm/year exists between the two sides of the crack. The width of influence of the crack along this cross section is approximately 200 m, judging from the distance of the crack deformation gradient. Figure 10. Profile of the deformation rate along P2-P2' (shown in Figure 6). The blue lines show the locations of GFs. The numbers 1-5 represent the numbers of the GFs crossed by profile P2-P2', as shown in Figure 6.  . Comparison of GF activity between levelling and INSAR measurements. The red and blue circles represent InSAR-derived displacements on point S1 and point S2, respectively. The black circles show the cumulative deformation between S2 and S1. The purple squares represent the levelling measurements.
To model the ground deformation caused by the GF, we took the central position of the deformation gradient as the zero-value position of x in Equation (2). Because the deformation beyond the deformation gradient is not considered to be associated with the tectonic activity of the GF, only the deformation values close to the deformation gradient were used in the model. We found the best-fit values for the model parameters using a grid-search method. A locking depth of 44 m and a far field velocity of $7.8 mm/year are the best-fit model parameters. The best fitting model is shown as the red line in Figure 12. According to the model results, the cracks exist in the shallow soil layer within 50 m of the surface.

Conclusion
Several geological hazards occur in the Datong basin, including fault activity, GFs, and land subsidence. In this article, Envisat ASAR images (C-band) and ALOS PALSAR images (L-band) were exploited to map the ground deformation of Datong between 2007 and 2010. In addition to the mining deformation in the western mountains, four major deformation zones were found in the Datong basin. The ground deformations within the basin were considered to be related to groundwater exploitation, fault activities, and urban construction. The ground deformation within the basin is bounded by faults. The formation of GFs, one of the major hazards, is controlled mainly by hidden faults in the basin, which determine the strike and distribution of the GFs. Two-dimensional deformation velocity fields in the east-west and vertical directions were recovered from Cand L-band SAR images. The results showed that vertical displacements dominate the deformation field in the Datong basin. A simple model in a homogeneous elastic half-space was used to model the GF deformation. A locking depth of 44 m and a far field fault slip velocity of À7.8 mm/year were the best-fit model parameters. We also used a uniformly opening sill embedded in an elastic halfspace to model the vertical deformation caused by groundwater exploitation and discussed the results.

Disclosure statement
No potential conflict of interest was reported by the authors.