Comparison of PSI and DInSAR approach for the subsidence monitoring caused by coal mining exploitation

Underground coal exploitation has its reflection in ground movements such as subsidence, sinking or shaking. These cause buildings and infrastructure damage, therefore it is important to measure the magnitude of deformation. Last decades, Differential Interferometric Synthetic Aperture Radar (DInSAR) captured considerable attention as a tool for deformation monitoring. The results of conventional DInSAR, which utilizes two SAR images, are degraded due to atmospheric, topographic and orbital errors. To overcome these limitations, various stacking-based methods have been introduced. Therefore, the goal of presented study is to compare Persistent Scatterer Interferometry (PSI) as stacking-based method with classical DInSAR for monitoring of subsidence caused by underground coal exploitation. Deformations in the areas of active mining exploitation are characterised typically by rapid non-linear movement. The comparison has been performed for the area of active exploitation in Rydułtowy mine located in Upper Silesian Coal Basin (USCB) in Poland. Results from two separate PSI and DInSAR processing portray similar deformation pattern over the study area. Unfortunately, due to the temporal decorrelation, PSI clearly demonstrate smaller information coverage in respect to DInSAR results. Additionally, due to the applied linear deformation model, PSI failed in displacement estimation with magnitude higher than 12 cm. In contrast, DInSAR thanks 6-day temporal baseline and no assumption for a deformation model, was able to capture the maximum magnitude of subsidence reaching 86 cm/year. However, these results are affected by atmospheric artefacts which in presented case study can reach even 14 cm/year. To achieve few cm level of accuracy and to estimate high deformation magnitude such as in presented study case (1m/year), integrated use of both InSAR techniques seems to be the reasonable solution. * Corresponding author


INTRODUCTION
Poland is one of the biggest coal producer in Europe. Upper Silesian Coal Basin (USCB) which is located in Southern Poland covers an area of almost 6000km 2 (Przyłucka et al., 2016). Coal exploitation in the USCB dated back the nineteenth century (Ilieva et al., 2019). Coal exploitation performed underground has its consequence in ground movements such as subsidence, sinking or shaking (Mutke et al., 2019). This ground deformation often results in buildings and infrastructure destruction (Przyłucka et al., 2019). Because the USCB is relatively densely populated and urbanised, monitoring of the ground deformation in this area is especially crucial for safety reason. Since the conventional geodetic survey is point-based and requires time and resources to perform measurements, alternative techniques, mainly Differential Interferometric Synthetic Aperture Radar (DInSAR) captured considerable attention in the last two decades. However, the most SAR sensors (e.g., ERS-1,2, Envisat, TerraSAR-X, and Sentinel-1) are C-band or X-band. Interferograms generated by SAR acquisitions with C-band or X-band are very sensitive to the temporal decorrelation. This significantly reduces the possibility of InSAR applications in rural areas. Since Sentinel-1 is C-band sensor, temporal decorrelation is important issue, however, temporal baseline of 6-days for Sentinel-1 A\B satellites increases coherence and offer the possibility to monitor ground deformation also in rural areas be means of DInSAR. Seeing that, the conventional DInSAR utilizes two SAR images, the results are affected by atmospheric, topographic and orbital errors (Osmanoğlu et al., 2016). To reduce the above mentioned errors, various stacking-based methods have been introduced. Among them, Persistent Scatterer Interferometry (PSI) is based on points or stable natural reflectors which remain stable over the long time interval (Crosetto et al., 2016). In the PSI processing chain, Atmospheric Phase Screen (APS) is estimated and removed (Ferretti et al., 2001). Various approaches and implementations of PSI have been proposed in last years (see Feretti et al., 2000;Ferretti et al., 2001;Hooper et al., 2004;Warner et al., 2003;Kampes and Hanssen (2004), Crosetto et al. (2005) Kampes and Adam, 2005). The PSI approach introduced by Ferretti et al., 2000 has an advantage that persistent scatterers are not affected by the baseline decorrelation. All SAR images can be utilized for the interferograms generation even if the baseline between these two images is longer than a critical baseline. This approach makes it possible to achieve sub-meter digital elevation model (DEM) precision and millimetres level of deformation measurements (Ferretti et al., 2001). Moreover, the analysis is performed at the full spatial resolution and is focused on pixels containing a single dominant scatterer (Pepe and Calò, 2017). Therefore, PS are smaller than the resolution cell of SAR image (Ferretti et al., 2001). However, the PSI approach utilizes the linear deformation model. As a consequence, in the areas where displacements are characterized by non-linear motion, there will be no PSs. Having considered the advantages and limitation of DInSAR and PSI, the goal of this study is to compare results of both methods for the monitoring of subsidence caused by the underground coal exploitation. The comparison is performed for the area of active exploitation in the Rydułtowy mine located in USCB.

Study area
Study area selected for this research is the oldest mine in the USCB, Rydułtowy mine, which covers 32 km 2 . The coal exploitation over the study area reaches 1000 m and sometimes even 1200 m under the ground level (Pawluszek-Filipiak and Borkowski, 2020). Although, the mining exploitation is carried out very deep, considerable ground deformations are observed often in the form of building and infrastructure damages (Przyłucka et al., 2015). According to Pawluszek-Filipiak and Borkowski (2020), the subsidence can reach even one m/year in the selected study area. USCB is also very densely populated area, therefore, deformation monitoring is very important to mitigate negative effects of the subsidence. The location of the study area with Sentinel-2 image is presented in figure 1.

Data used
In this study, we utilize in total 62 ascending Sentinel-1 images. Data has been acquired from the Scientific data Hub platform available at the Copernicus webpage: https://scihub.copernicus.eu/dhus/#/home. More detailed specification of the used data is presented in Table 1.

METHODOLOGY
The methodology of the presented research is based on the InSAR processing of Sentinel-1 data. The InSAR processing was performed by conventional DInSAR and PSI techniques, separately. Results of both techniques have been compared with each other. The overall view of applied methodology is shown in figure 2 while more detailed description of intermediate steps is presented in the following subsections.

Consecutive DInSAR processing
DInSAR processing relies on calculation of differential interferograms between two neighbouring SAR acquisitions (e.g. φ1-2, φ2-3, φ3-4,… φn-1,n). In this study 62 SAR images have been used to calculate 61 consecutive differential interferograms. ALOS-3D digital elevation model has been used to simulate phase corresponding to the topography, which has been afterwards removed during interferometric processing. Phase unwrapping procedure has been carried out using the Minimum Cost Flow function and then converted to displacements. The displacements estimated from each SAR pair, have been accumulated to provide complete time-series of interferometric results (e.g.,φ1-2, φ2-3, φ3-4, …, φn-1, n). The relationship between SAR images, baselines and time are presented in figure 3. It can been observed that maximal spatial baseline between two SAR images in presented stack is almost 200m (image 39 and 40).

PSI processing
In this study, we utilized classical PSI approach introduced by Ferretti et al., (2000Ferretti et al., ( , 2001. The PSI technique was developed to employ the persistent scatterers (PS), whose sizes are smaller than the resolution cell of SAR image (Ferretti et al., 2001). For 61 slaves images, 62 full-resolution interferograms have been created in reference to the same master image. Figure 4 presents the spatial and temporal relation between the master image and slave images. The maximum temporal and spatial baselines are equal to 228 days and 130m, respectively. Afterwards, pixels with a relatively constant amplitude in time are selected as potential scatterers (PSC). PSC with amplitude dispersion higher than 0.25 are considered for the first iteration. In the next step, the topographical component of the phase, APS and baseline error contributions are estimated as a linear phase term both in range and azimuth direction. All these components are simultaneously estimated by minimizing the temporal phase residuals for pixels that are coherent enough (Ferretti et al., 2001). In presented approach, constant linear velocity model is adopted for target motion. Continuously, 1D+2D unwrapping approach is applied for sparse grid points previously selected based on their coherence. After APS removal, finally motion of each pixel in the image was estimated and more PSs were identified. Afterwards, both topographic components and PS velocity were estimated on a PS-by-PS basis. Finally, PS with coherence higher than 0.6 have been identified. More detailed description of the applied PSI technique can be found in Ferretti et al.( 2000Ferretti et al.( , 2001.

PSI and DInSAR comparison
Achieved results from PSI and DInSAR processing have been compared with each other. The comparison was performed for the PS points, which exhibit maximum deformation gradient with corresponding DInSAR pixel. The root mean square error (RMSE) has been calculated for the differences between DInSAR and PSI for the selected points in time-series results. Additionally, for the last observed date, we transformed our PS results into the raster format and we calculated differences between these two results in order to assess reliability of the achieved results.

RESULTS
Results of the consecutive DInSAR and PSI processing have been presented in figures 5a and 5b, respectively. As expected, the DInSAR result provides better information coverage in respect to PSI results. The reason for that is that in DInSAR processing the temporal baseline is always 6 days (with some exception of 12 days when there is lack of SAR image). On the contrary, in PSI processing, all images are coregistered to the one master image. The temporal baseline between master and the first slave image reaches even 228 days. Therefore, temporal decorrelation affected the information coverage in PSI results. Additionally, it is clearly visible that PSI provides empty "holes" in the centre of the maximum subsidence. DInSAR was able to retrieve the maximum deformation in LOS direction reaching -0.86 m while PS was only able to retrieve -0.12 m. For the PS point with maximum detected deformation magnitude, deformation time series have been extracted ( figure  6). Having observed deformation behaviour in time, it can be noticed that PSI and DInSAR results correspond with each other very fairly in the beginning and in the end of the investigated period. Differences reaching 2/3 cm can be clearly observed between 14.10. 2018 and 13.12.2018. Calculated RMSE between DInSAR and PSI in time-series for this point is 15mm. The explanation for this can be atmospheric artefact, which were not modelled and removed from DInSAR processing. Also non-linear deformation trend can be observed. Having considered, the APS estimation in PSI processing, PSI results are almost "atmospheric free". Therefore, to evaluate the error rate which can corresponds to atmospheric error, we simply calculate the difference between PSI and DInSAR results. In order to performed it, we transform PS points into the raster format. Then, DInSAR raster has been subtracted from PSI raster and achieved differences are presented in figure 7. The maximum and minimum values of 0.12 and -0.14 can be observed for the PSI-DInSAR differences. These differences clearly show some clusters of colours especially violet cluster on the left-bottom part of figure 7. These clusters have an extent reaching 2 km width and 3 km height, therefore it clearly shows that the differences mainly corresponds to accumulated atmospheric artefacts, which were not removed in the DInSAR processing. According to Gomba et al.,2016, in such a small study area ionospheric error can be neglected, thus these error mainly corresponds to tropospheric component. Moreover, displacement error related to the topography can be calculated according to Long et al. (2018). Taking into account the slant range from ground point do the satellite platform about 700 km, incidence angle about 36°, the normal baseline equal 200 m (maximal value in our study) and DEM error equal 10 m, the error due to topography is equal 2 mm and can be also neglected. In figure 8 the histogram of the differences between PSI and DInSAR is shown. It is clear that the distribution of the differences is almost perfect symmetric and close to the normal distribution. However the mean value indicates that the DInSAR result is in average 22 mm above the PSI result. It can be interpreted as an error due to atmospheric artefacts. The standard deviation is equal 23 mm (Figure 8). Therefore RMSE for the whole scene can be calculated to be equal of 32 mm.

DISCUSSION AND CONCLUSSION
The results of conventional DInSAR, which utilizes two SAR images, are affected by atmospheric, topographic and orbital errors. To reduce these adverse effects various stacking-based methods have been introduced such as Persistent Scatterer Interferometry (PSI). PSI is based on points or stable natural reflectors which remain stable over long temporal series. Since PSI utilized stable reflectors, its application is very limited by temporal decorrelation. Additionally, PSI approach uses linear deformation model. Therefore, in areas where displacements are characterized by non-linear motion, there are none PSs. From another point of view, temporal baseline of 6-days for Sentinel-1 A\B satellites increases coherence and offer the possibility to monitor ground deformation also in rural areas.
Having considered the advantages and limitations of DInSAR and PSI, the goal of presented study was to compare PSI and DInSAR results for monitoring subsidence cause by underground coal exploitation. The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XLIII-B3-2020, 2020 XXIV ISPRS Congress (2020 edition) The comparison has been performed for the area of active exploitation in Rydułtowy mine located in USCB. Results from two separate processing of PSI and DInSAR portray similar deformation patters over the study area. Unfortunately, due to the temporal decoration, PSI results clearly demonstrate small information coverage in respect to DInSAR results. Additionally, due to the applied linear deformation model, PSI failed in the displacement estimation with magnitude higher than 12 cm/year. For the point with the maximum PS-velocity, extracted time series from DInSAR and PSI corresponds with each other quite fairly. In contrast, DInSAR thanks to 6-day temporal baseline and none deformation model assumption, was able to capture the maximum magnitude of subsidence reaching 86 cm/year. Besides the possibility to catch big deformation gradient, DInSAR results are affected by atmospheric influence. In presented case when 61 LOS deformation maps were accumulated, cumulative error in some places can reach even 14 cm/year. We interpret this error as due to atmosphere. This clearly demonstrate that in order to achieve cm level of accuracy of DInSAR, atmospheric modelling needs to be carried out. Nevertheless, when the displacement is much higher than atmospheric artefact, DInSAR can be applied. However in many mining application few cm level accuracy is sufficient. To achieve this level of the accuracy and in the same time estimate high magnitude of the deformation (1m/year) such as in presented study case, the reasonable solution seems to be using of both InSAR techniques and then to integrate results as presented in studies of Przyłucka et al. (2015); Pawluszek-Filipiak and Borkowski (2020).