Subsidence Monitoring over the Southern Coalfield, Australia Using Both L-band and C-band Sar Time Series Analysis

Land subsidence is a global issue and researchers from all over the world are keen to know the causes of deformation and its further influences. This paper reports the findings from time series InSAR (TS-InSAR) results over the Southern Coalfield, Australia using both ALOS-1 PALSAR (Phased Array type L-band Synthetic Aperture Radar) and ENVISAT ASAR (Advanced Synthetic Aperture Radar) datasets. TS-InSAR has been applied to both rural and urban areas with great success, but very few of them have been applied to regions affected by underground mining activities. The TS-InSAR analysis exploited in this paper is based on GEOS-ATSA, and Measurement Point (MP) pixels are selected according to different geophysical features. Three experiment sites with different geological settings within the study zone are analysed: (1) Wollongong city, which is a relatively stable area; (2) Tahmoor town, a small town affected by underground mining activities; and (3) the Appin underground mining site, a region containing multiple underground mining activities. The TS-InSAR results show that the performance of both C-band and L-band is equally good over Wollongong, where the subsidence gradient is not significant and most subsidence rates are between´10 mm¨yr´1 to 10 mm¨yr´1. However, over the Tahmoor and Appin sites, difference in performances has been observed. Since the maximum displacement gradients that can be detected are different for L-band and C-band-based TS-InSAR methods, some rapid changes could cause the TS-InSAR to fail to estimate the correct displacements. It is well known that L-band can perform better than C-band, especially in underground mining regions and mining-affected regions where the deformation rate is much higher than city areas because of its wavelength. Statistical analyses are also conducted to further prove the above statement.


Introduction
Land subsidence is a global issue and researchers from all over the world are keen to know the causes of deformation and its further influences.Land subsidence can be divided into two categories: (1) anthropogenic-involved subsidence, the main reasons could be underground mining activities; groundwater, oil and gas extraction; and underground tunnel construction; (2) non-human-related subsidence, e.g., earthquake deformation and volcanic eruption [1].It is worth noting that ground subsidence can cause a great deal of serious environmental hazards, such as collapse of infrastructures, damage of bridges and changes of river and groundwater flow path, to name only a few.These hazards not only threaten individuals' everyday lives and properties, but also lead to massive impacts on the national economy.In Bandung Basin, Indonesia, for example, acceleration of ground subsidence has been reported by many researchers due to excessive groundwater extraction and the largest subsidence could reach to as large as ´24 cm¨yr ´1 [2].In addition, during the processing of hard rock and underground coal mining, thousands of miners die every year all over the world as a direct result of the long-term accumulated subsidence [3].However, these catastrophic events to some extent are often predictable since they often follow a natural cycle.Timely displacement information, that could provide strong indicators of these forthcoming disasters, can often be used to identify the severely affected area.One recent example is that UC Berkeley's ShakeAlert System detected the South Napa earthquake ten seconds earlier by using p-wave to monitor displacement of earth surface [4].Therefore, having the ability to measure earth surface deformation well in advance could give us a better understanding of the natural hazards.Eventually, this can help the local government or associated councils to minimize the impacts and make better decisions.
Conventionally, mapping of earth surface deformation can be achieved by using field-survey techniques, such as Global Positioning Systems (static and real-time kinematic GPS), total stations and digital levels [5].However, these methods have their limitations.They could be labor-intensive and time-consuming once the measurement regions become large.Moreover, these techniques are based on a point-by-point approach for taking measurements, which means it is very difficult to obtain a reasonable interpretable surface deformation over a large region [6].Over the last two decades, DInSAR has proven to be an effective way to monitor various earth surface subsidence activities with great successes because of its spatial coverage and precision [5,7].But there are some limitations that constrain its performance to a large extent [8], especially when monitoring low-velocity subsidence over a long period of time.The decorrelation issues caused by the large temporal and spatial baseline, often degrade the interferogram phase [9].Furthermore, atmospheric artefacts could also degrade the quality of deformation estimation and the changes in troposphere and ionosphere from one day to the next could also form different time delays similar to the deformation signals [10,11].
To overcome these problems, time series InSAR (TS-InSAR), in which jointly analysed multi-SAR images acquired on different dates, was invented in the late 1990s.Over the past 15 years, many TS-InSAR techniques were developed to minimise the error sources of traditional DInSAR technique.In general, it can be divided into three categories: (1) In the first category are techniques that make use of only one single master image to generate stacks of interferograms.These approaches estimate the ground deformation at the so-called Persistent Scatterer (PS) pixels, whose scattering characteristics remain stable even over long time intervals and large baseline separation, e.g., Persistent Scatterer InSAR (PSInSAR TM ) [8], Stanford method for Persistent Scatterers (StaMPS) [12,13], the Spatio-Temporal Unwrapping Network (STUN) [14], Stable Points Network (SPN) [15,16] and the GEOS-PSI method [17].These techniques have the advantage of associating the deformation with a specific scatter, rather than a multi-looked resolution cell.This allows displacement maps to be generated with high resolution.The achievable accuracy could be 1 mm¨yr ´1 or better where the subsidence over the study region is linear in time [18]; (2) The second category involves the usage of multi-master interferograms, and only the so-called Distributed Scatterer (DS) pixels were selected for further analysis.Unlike PS pixels whose scattering characteristics are insensitive to spatial and temporal baselines, DS pixels can only maintain their scattering characteristics with limited spatial and temporal baselines, e.g., the Small Baseline Subset (SBAS) method [19,20], the Coherent Pixel technique [21], Poly-Interferogram Rate And Time-series Estimator (π-RATE) [22,23], the Temporally Coherent Point InSAR (TCPInSAR) [6], and the Intermittent SBAS (ISBAS) [24]; (3) In the third category, single-and multi-master interferograms are combined to form the time series analysis [25], SqueeSAR method [26] and the GEOS-ATSA [2].Both DS and PS pixels are combined as MP pixels to solve the ground deformation.These different TS-InSAR techniques have been applied to measure mean velocities and cumulative deformation in various applications that include groundwater extraction-induced land subsidence, landslide and volcanic deformation, and city urbanization and expansion-induced deformation.However, it is worth noting that none of the error sources can be completely eliminated even with the TS-InSAR method, and researchers working on these techniques are trying to minimize the effect according to different geological situations.
Australia is a country with many land subsidence issues and most of them are related to anthropogenic activities including the extraction of natural resources, such as coal, natural gas and iron ore [27,28].Sydney-Gunnedah Basin, which is 500 km long and 150 km wide, has the largest coal resource storage in New South Wales (NSW).There are five major coalfields inside the basin: Newcastle, Western, Hunter, Southern and Gunnedah.It covers an area south of Berrima and Sutton Forest to north of Campbelltown, west of Tahmoor town and east of Wollongong, and has the only source of hard coking coals in NSW, which is favourable for steel production [27].Many previous InSAR studies were conducted over the Southern Coalfield [27,29,30].More specifically, both DInSAR and small stacks of InSAR methods were exploited to study the local subsidence mainly near Appin, Tahmoor and West Cliff collieries in the Southern Coalfield of NSW.Single-master-based TS-InSAR technique was not used for these analyses mainly due to two major reasons: (1) The amount of ALOS PALSAR images was limited (10 acquisitions) while normally much larger image stacks (more than 20) are required by these TS-InSAR methods to achieve precise deformation measurement; and (2) Most of these collieries in the Southern Coalfield are located within rural areas where very few PS pixels were expected by using these methods.In this paper, more advanced TS-InSAR analyses will be conducted to study the time series evolution of the earth surface subsidence, and the three study areas are the Wollongong city area, Appin mining site and Tahmoor town.
It is generally accepted that PSInSAR TM , STUN and GEOS-PSI techniques are very powerful in city areas with many manmade structures where the density of PS pixels can be high enough.Indeed, they have been successfully applied in wide urbanized areas such as Las Vegas [14] and Beijing [17].SqueeSAR and GEOS-ATSA can be applied to non-urban regions with good results because they select not only PS pixels, but also DS pixels [2,26].The TS-InSAR method utilised in this study is based on SqueeSAR and GEOS-ATSA [2].The Wollongong city area, Appin underground mining site and Tahmoor town region are our areas of interests.Appin and Tahmoor are mostly in rural areas where low density of PS pixel is expected.MP pixels are thus selected according to the varying geophysical information.More specifically, since the DS pixels selection process is time consuming, PS pixels are selected over the entire study region while DS pixels are selected around the Appin underground mining site and Tahmoor town region in order to achieve the best details.

Study Regions Within Southern Coalfield
All three regions of interest are within the Southern Coalfield (Figure 1).Wollongong city area is in the Illawarra region of New South Wales, Australia, and it is situated adjacent to the Tasman Sea and the South Coast railway line (orange rectangular box in Figure 2).Appin is a town in the Macarthur Region of New South Wales, Australia and it is situated about 35 kilometres Northwest of Wollongong.The underground coal-mining site Appin Colliery is located about 25 kilometres Northwest of Wollongong (blue rectangular box in Figure 2).Tahmoor town is located in the Macarthur Region of New South Wales and to the southwest of Appin coal-mining site.The Tahmoor Colliery, located to the south of North Bargo, is the major source of its regional economic growth (purple rectangular box in Figure 2).The topographical slopes over these three sites range between 0 and 3, thereby indicating the study regions are relatively flat (Figure 2c).

Datasets
In order to map the land displacement over the Southern Coalfield, Australia, 23 L-band ALOS-1 PALSAR scenes acquired between 29 June 2007 and 7 January 2011 and 24 C-band ASAR images acquired between 8 July 2007 and 5 September 2010, are analysed (Table 1).All the ALOS-1 images (Track 370, Frame 649) were captured in ascending orbit with the average incidence angle of 38.7 ˝while all the ASAR images (Track 381, Frame 6492) were acquired in ascending orbit with the mean incidence angle of 28.8 ˝, which means both ALOS and ENVISAT sensors are only sensitive to Line-of-Sight (LoS) displacement (Figure 2a,b).Eleven ALOS-1 images were acquired in HH and HV dual-polarisation while the other twelve were acquired in HH single polarisation.The dual-polarisation data (HH) were oversampled by a factor of two in range direction [17] and the final pixel spacing in azimuth and range were 4.82 m and 5.55 m, respectively.All the ASAR images were acquired in VV single polarisation with the azimuth and range pixel spacing of 7.80 m and 4.06 m, respectively.

Image Co-Registration and Interferogram Generation
Firstly, all images are co-registered to a common master image.An empirical SAR calibration method is then carried out for each single image within the stacks [2,31].Conventional two-pass DInSAR method is applied to form 22 PALSAR differential interferograms and 23 ASAR differential interferograms with respect to their respective master images.The topographic phase is derived using the three arc-second DEM acquired from the Shuttle Radar Topography Mission (SRTM) [32] and then removed from differential interferograms.SRTM DEM is also exploited later to geocode the TS-InSAR result from the range-Doppler radar coordinate system into the Universal Transverse Mercator (UTM) coordinate system.

Measurement Point Pixels Selection
Because of the lack of PS pixels around the Appin and Tahmoor sites, DS pixels are selected around the Appin underground mining site and Tahmoor town with 2000 pixels ˆ3000 pixels (ALOS) and 3000 pixels ˆ2000 pixels (ENVISAT).The method used to identify DS pixels and obtain the optimal phase is the same procedure used in GEOS-ATSA [2].More specifically, the two-sample Kolmogorov-Smirno (KS) test is carried out on a pixel-by-pixel basis to search for surrounding statistically homogeneous pixels (SHP) within a 20 ˆ10 window for ALOS-1 and 31 ˆ6 window for ENVISAT; image pixels with SHP higher than 20 are selected as DS candidates.In order to reduce the speckle and decorrelation noise of each DS pixel, a spatial adaptive filter is applied [26].Then the "Phase Linking" method is exploited to estimate the optimal phase value with respect to each DS pixel by using all possible interferogram combinations [33].Finally, the goodness-of-fit γ P estimated from Equation ( 1) is utilised to decide whether the image pixel should be considered as the final DS pixel (the thresholds γ P are 0.75 for ALOS-1 PALSAR and 0.55 for ENVISAT ASAR) [2,26].After that, the amplitude dispersion index (D A ) [8] is used to estimate the phase stability of each pixel.Pixels with D A < 0.25 are selected as reliable PS pixels along with DS pixels to form the total MP pixels.It is worth noting that since no MP pixels are considered as the noise source, increasing the threshold would increase the possibility of false-alarm and hence degrade the accuracy of the results.
where N is the number of SAR image stacks; ∅ nk is the spatially filtered phase; θ n and θ k are optimal phases of image n and image k, respectively; γ P is the good-of-fitness value to evaluate the quality of the optimal phase with respect to each DS pixel.

Orbit Refinement
A first-order equation is utilised to remove the linear phase "ramp" induced by residual orbital error because it is assumed that there is no deformation with a linear spatial trend in the study region [17,34].Equation ( 2) is then applied to each interferogram.
where x and y are line and pixel coordinates in the UTM system of each MP pixel, respectively, while a, b and c are equation coefficients.d is a constant value.

Network Construction and Integration
All the MP pixels are exploited to construct the Delaunay triangulation reference network and the maximum distance for any connected MP pixels is limited at 2.0 km.The Least-squares AMBiguityDecorrelation Adjustment (LAMBDA) [14] is exploited to estimate the DEM error and relative velocity for each MP pixel.During the estimation, arcs with "ensemble phase coherence" [8] less than 0.75 and 0.5 are assumed to be unreliable and removed from further analysis for both ALOS-1 PALSAR dataset and ENVISAT ASAR datasets, respectively, because these MP pixels are not likely to have high phase stability, while isolated MP pixels related to these arcs are also deleted [17].To estimate the absolute linear deformation velocity and DEM error, a maximum likelihood regression-based M-estimator [35] is exploited in this study, which belongs to the robust-fitting family [36,37].

Tropospheric Phase Delay Estimation
Traditional temporal-spatial filtering is carried out to estimate the unfavoured phase components.It is generally accepted that both tropospheric turbulence phase ϕ turb and tropospheric stratification phase ϕ stra are correlated in space and non-correlated in time; noise is non-correlated in both time and space while non-linear deformation is correlated in both time and space.A filtering operation is carried out in both time and space based on [8].After removing the modelled absolute DEM error and linear velocity, a sparse Minimum Cost Flow (MCF) method is applied to the residual phase to derive the unwrapped residual phase ϕ residual with respect to each MP pixel [2,35,38].First of all, the temporal mean residual phase ϕ residual needs to be removed from ϕ residual : φresidual " ϕ residual ´ϕresidual .A temporal high-pass filter with a defined triangular window of 360 days is applied to remove the temporally correlated component, and the resulted phase is denoted as φresidual_HP .Then a spatial low-pass filter with 2 km ˆ2 km window is exploited to remove the non-correlated component for ϕ residual and φresidual_HP , and resulting in ϕ residual_LP and φresidual_HP_LP .The tropospheric turbulence phase delay combined with tropospheric stratification phase delay is derived as: ϕ trop " ϕ turb `ϕstra " ϕ residual_LP `φ residual_HP_LP .

Non-Linear Deformation Estimation
Non-linear deformation and noise are estimated by subtracting ϕ trop from ϕ residual .Since the noise is non-correlated in the space domain while non-linear deformation is correlated in the space domain, a low-pass filter with a 1 km ˆ1 km window is applied to remove the noise component.Finally, the non-linear deformation is obtained by multiplying the filtered phase with ´λ 4π .

Experimental Results and Analyses
Images acquired on 1 January 2009 (ALOS-1 PALSAR images) and 18 January 2009 (ENVISAT ASAR scenes) are picked as master images for the two image stacks to minimize the temporal and perpendicular baselines.Figure 2a,b presents the TS-InSAR results using the ALOS-1 PALSAR and ENVISAT ASAR datasets, respectively, and both of them are resampled into a grid with a resolution of 100 m ˆ100 m for further comparison.Regions marked with red, blue and purple rectangular boxes are the Wollongong city area, Appin underground mining site and Tahmoor town region.The measured displacements are with respect to the reference point (pink point in Figure 2) centred at 34 ˝22 1 04"S and 150 ˝55 1 25"E.The total number of MP pixels obtained from ALOS-1 PALSAR dataset is 1,652,180, of which 88,003 are PS pixels.There are 403,857 MP pixels estimated from ENVISAT ASAR dataset, and 83,304 of them are PS pixels.The reason for this is because the wavelengths of the two sensors are different.L-band PALSAR sensor with longer wavelength has less decorrelation than C-band under the same baseline condition, therefore, offers a much higher density of DS pixels (the density of PS pixels are similar).The total number of refined arcs is 4,921,872 and 1,184,997 for PALSAR and ASAR datasets, respectively.Therefore, the density of MP pixels at the three sites are 355, 8004, 6686 MP/km 2 and 398, 423, 1557 MP/km 2 for ALOS-1 and ENVISAT, respectively (only PS pixels are selected in Wollongong city area while both PS pixels and DS pixels are selected in the Appin underground mining site and Tahmoor town region).To compare the time series performances between C-band and L-band products, both time series displacements in LoS direction (Figure 2a,b) are converted into vertical direction by dividing the cosine of its incidence angle.

Errors Respect to TS-InSAR Measurement
Error analysis is a very significant step of SAR interferometry application, and it can be performed in either a qualitative or quantitative manner.

Systematic Bias
It is widely accepted that the L-band sensor is less sensitive to the ground deformation as compared to C-band, and the accuracy of the TS-InSAR-estimated LoS deformation rate depends on three main components: the phase stability of MP pixel, temporal baseline distribution, and sensor wavelength.The equation to estimate the standard deviation of error for TS-InSAR-derived LoS deformation rate can be expressed as [38][39][40]: where N is the number of images in the stacks; σ v is the standard deviation of the estimated LoS displacement rate; λ is the wavelength of the sensor; σ ϕ is the standard deviation of phase residuals; and σ T is the standard deviation of the temporal baseline.Under the assumption that σ 2 T and σ ϕ are constant, the precision of the measured deformation of the C-band dataset is about four times better than equivalence for the L-band dataset.In other words, the TS-InSAR measurement derived from the ENVISAT ASAR dataset was expected to achieve more precise results than the counterpart from the ALOS PALSAR dataset.

Projection Bias
InSAR can only be utilised to measure the ground subsidence in the LoS direction because of its side-looking system, which consists of vertical, easting and northing displacement components.The equation can be expressed as follows: Within this equation, θ is the incidence angle while α is the heading angle (clockwise from north), D V is the vertical displacement, D E is the displacement in the eastern direction, while D N is the displacement in northern direction.D LoS is the displacement in the LoS direction.
Since the LoS displacement vector is insensitive to the deformation in the north-south direction with current satellite viewing geometries [41], the northing displacements were assumed negligible.For the ALOS-1 ascending pairs, the satellite heading angle is ´10 ˝and the incidence angle is 38.7 ˝.
Therefore, the coefficients are 0.78 and ´0.62 for the vertical and easting deformation, respectively.For ENVISAT ASAR ascending pairs, the heading angle is ´16 ˝, while the incidence angle is about 28.8 ˝.Thus the coefficients are 0.88 and ´0.46 in vertical and easting directions, respectively [38].In general, the 3D displacement can only be solved if at least three InSAR pairs from both ascending and descending directions are available.In this study, the vertical displacement remains the dominant displacement signal in most study regions and the easting displacement is also neglected, which means an easting displacement of 5 mm¨yr ´1 can cause a displacement error of 4 mm¨yr ´1 and 2.5 mm¨yr ´1, respectively, for ALOS-1 and ENVISAT if the LoS displacement is directly projected into the vertical direction [38].Furthermore, it is worth noting that a LoS displacement of 10 mm¨yr ´1 is equal to 12.8 mm¨yr ´1 and 11.4 mm¨yr ´1 in the vertical direction if both easting and northing displacements are neglected.
Since no GPS measurements or any field survey data are used in this paper, the LoS deformation rate derived from both datasets with respect to each MP pixel in Wollongong City are compared.

Deformation over Wollongong City
The ground cover in Wollongong City is mainly made up of low vegetation and manmade structures, where a large number of PS pixels are expected because high coherence in the interferograms can be preserved over these regions.Figures 3 and 4 show the displacement rate distribution between ALOS PALSAR and ENVISAT ASAR time series measurements in Wollongong City.It is worth noting that the resolution of both displacement maps is resampled to 100 m ˆ100 m and only common MP pixels have been taken into consideration.Figure 3 shows the time series measurement over one MP pixel marked by a red cross where a similar subsidence trend has been observed from both C-band and L-band results.Figure 4 illustrates good agreement between these two datasets.A root-mean-square (RMSE) error (Equation ( 5)) of 3.11 mm¨yr ´1 has been obtained, showing that the L-band and C-band time series subsidence is comparable over Wollongong.During the estimation, the unmatched MP pixels are not included in the calculation of RMSE.In addition, the majority of subsidence in the city is between ´10 mm¨yr ´1 to 10 mm¨yr ´1, which suggests that the local subsidence rate is quite stable.The difference between observations could be due to three reasons: (1) The different measurement precision is due to different sensors; (2) The horizontal movement is neglected; or (3) There is a mis-match in geolocation between individual MP pixels.RMSE " where P is the total number of common MP pixels within the region for RMSE estimation, D M p and D S p are the pth pixels within master and slave images, respectively.
has been observed from both C-band and L-band results.Figure 4 illustrates good agreement between these two datasets.A root-mean-square (RMSE) error (Equation ( 5)) of 3.11 mm•yr −1 has been obtained, showing that the L-band and C-band time series subsidence is comparable over Wollongong.During the estimation, the unmatched MP pixels are not included in the calculation of RMSE.In addition, the majority of subsidence in the city is between −10 mm•yr −1 to 10 mm•yr −1 , which suggests that the local subsidence rate is quite stable.The difference between observations could be due to three reasons: (1) The different measurement precision is due to different sensors; (2) The horizontal movement is neglected; or (3) There is a mis-match in geolocation between individual MP pixels.
where P is the total number of common MP pixels within the region for RMSE estimation,    and    are the pth pixels within master and slave images, respectively.

Deformation over Tahmoor
Figure 5a,b illustrates the vertical time series measurements in the Tahmoor town region with ALOS-1 PALSAR and ENVISAT ASAR, respectively.It can be observed in both figures that large areas around Tahmoor's centre are suffering from significant fall in ground elevation.In addition, L-band result indicates the western parts of Tahmoor town has experienced large subsidence with the maximum subsidence higher than 10 cm•yr −1 in vertical direction, whereas it is not found in the

Deformation over Tahmoor
Figure 5a,b illustrates the vertical time series measurements in the Tahmoor town region with ALOS-1 PALSAR and ENVISAT ASAR, respectively.It can be observed in both figures that large areas around Tahmoor's centre are suffering from significant fall in ground elevation.In addition, L-band result indicates the western parts of Tahmoor town has experienced large subsidence with the maximum subsidence higher than 10 cm¨yr ´1 in vertical direction, whereas it is not found in the C-band result.A 2 km profile line is drawn through the entire Tahmoor town (Figure 5c). Figure 5d shows the time series subsidence on this profile line for both datasets.The largest subsidence of 8.5 cm¨yr ´1 and 5.8 cm¨yr ´1 has been detected, respectively.It is clear that the left and right parts agree with each other well.However, there is a huge difference in the middle of the profile from 600 m to 1300 m, with the largest difference reaching 4.5 cm¨yr ´1 at 1000 m.The difference is mainly caused by several reasons: (1) A searching window with the same size of 100 m × 100 m is applied to find the point targets closet to the profile line for comparison, which means each value within the profile line does not reflect the subsidence of one single MP pixel, but a group of MP pixels.In addition, the numbers of MP pixels are different for both sensors within a searching window due to the threshold difference for "ensemble phase coherence" and we hold the assumption that L-band could preserve higher coherence compare to C-band; (2) The middle part may suffer from some rapid changes during the 4-year period.However, compare to the L-band sensor, pixels with large displacement rate are less likely to be considered as MP pixels for C-band sensor.In other words, within the searching window, L band TS-InSAR result contains MP pixels with larger deformation that C-band TS-InSAR cannot be captured; (3) According to the eastern part of Figure 6a-f, the loss of coherence in the high-gradient region in C-band dataset which resulted in under-estimation of the mean velocity.The difference is mainly caused by several reasons: (1) A searching window with the same size of 100 m ˆ100 m is applied to find the point targets closet to the profile line for comparison, which means each value within the profile line does not reflect the subsidence of one single MP pixel, but a group of MP pixels.In addition, the numbers of MP pixels are different for both sensors within a searching window due to the threshold difference for "ensemble phase coherence" and we hold the assumption that L-band could preserve higher coherence compare to C-band; (2) The middle part may suffer from some rapid changes during the 4-year period.However, compare to the L-band sensor, pixels with large displacement rate are less likely to be considered as MP pixels for C-band sensor.In other words, within the searching window, L band TS-InSAR result contains MP pixels with larger deformation that C-band TS-InSAR cannot be captured; (3) According to the eastern part of Figure 6a-f, the loss of coherence in the high-gradient region in C-band dataset which resulted in under-estimation of the mean velocity.

Deformation over Appin Underground Mining Site
One of the limitations of TS-InSAR techniques is related to resolving the phase ambiguity in regions suffering from rapid subsidence [42], e.g., underground mining-induced subsidence.In other words, the deformation signal may not be correctly estimated if the subsidence of the MP pixel between any two image acquisitions within the image stack is larger than one half cycle (5.9cm/46 days for L-band satellite while 1.4 cm/35 days for C-band satellite).In this circumstance, the TS-InSAR technique can only be used to identify the occurrence time and location with respect to each MP pixel, but it may not be able to calculate the correct subsidence value.The above problem can be solved if a priori knowledge is available [40], which is generally significantly limited or even unavailable in most TS-InSAR analysis cases.
Land subsidence in the Appin mine site is mostly related to underground mining activities.Figure 7a,b shows the time series displacement of L-band and C-band measurements at the Appin site and the resolution is resampled into 30 m × 30 m from original datasets, respectively.In general, the deformation patterns observed from the two measurements are not the same.In Figure 7a, several clear oval deformation patterns were found in the middle and right parts of the image, which is due to underground mining activity.However, there is no such subsidence pattern observed in Figure 7b, where the MP pixels are identified at the roads and residential areas and are reasonably stable.A possible explanation is that the deformation gradient is too large and highly non-linear over the subsidence area.According to [27], the mining-induced subsidence in the Southern Coalfield could reach up to 20 to 60 cm within 1-2 months and up to 80-100 cm in a one-year period.Therefore, no MP pixels can be obtained over the areas that were experiencing significant subsidence from C-band TS-InSAR product.

Deformation over Appin Underground Mining Site
One of the limitations of TS-InSAR techniques is related to resolving the phase ambiguity in regions suffering from rapid subsidence [42], e.g., underground mining-induced subsidence.In other words, the deformation signal may not be correctly estimated if the subsidence of the MP pixel between any two image acquisitions within the image stack is larger than one half cycle (5.9cm/46 days for L-band satellite while 1.4 cm/35 days for C-band satellite).In this circumstance, the TS-InSAR technique can only be used to identify the occurrence time and location with respect to each MP pixel, but it may not be able to calculate the correct subsidence value.The above problem can be solved if a priori knowledge is available [40], which is generally significantly limited or even unavailable in most TS-InSAR analysis cases.
Land subsidence in the Appin mine site is mostly related to underground mining activities.Figure 7a,b shows the time series displacement of L-band and C-band measurements at the Appin site and the resolution is resampled into 30 m ˆ30 m from original datasets, respectively.In general, the deformation patterns observed from the two measurements are not the same.In Figure 7a, several clear oval deformation patterns were found in the middle and right parts of the image, which is due to underground mining activity.However, there is no such subsidence pattern observed in Figure 7b, where the MP pixels are identified at the roads and residential areas and are reasonably stable.A possible explanation is that the deformation gradient is too large and highly non-linear over the subsidence area.According to [27], the mining-induced subsidence in the Southern Coalfield could reach up to 20 to 60 cm within 1-2 months and up to 80-100 cm in a one-year period.Therefore, no MP pixels can be obtained over the areas that were experiencing significant subsidence from C-band TS-InSAR product.Figure 8a shows the typical fringe patterns caused by mine subsidence derived from DInSAR interferograms with two obvious characteristics: (1) The earth surface where underground mining occurs will sink as the colour of the deformation pattern changes from yellow to blue (from the centre to the edge); and (2) The subsidence magnitude increases from the edge to the centre, therefore forming an oval or round shape (Hu et al., 2013).Figure 8b illustrates the underground mining subsidence pattern from the TS-InSAR velocity map.Regions suffering from rapid changes within a short period of time will de-correlate and therefore form gaps in the subsidence zones.Two characters are described in this paper to identify underground mining-induced subsidence regions: (1) The subsidence zones are normally rectangularly or quadrilaterally shaped (related to the shape of longwall mining sites); and (2) The magnitude of subsidence rate over the subsidence zone edges is relatively large.By using the methodology describe above, four mine subsidence bowls have been detected within Figure 7a with a maximum subsidence rate in excess of −10 cm•yr −1 .

Conclusions and Future Work
In this paper, the land surface stability over three test sites with different geological settings within the Southern Coalfield region was investigated using TS-InSAR technique.These three Figure 8a shows the typical fringe patterns caused by mine subsidence derived from DInSAR interferograms with two obvious characteristics: (1) The earth surface where underground mining occurs will sink as the colour of the deformation pattern changes from yellow to blue (from the centre to the edge); and (2) The subsidence magnitude increases from the edge to the centre, therefore forming an oval or round shape (Hu et al., 2013).Figure 8b illustrates the underground mining subsidence pattern from the TS-InSAR velocity map.Regions suffering from rapid changes within a short period of time will de-correlate and therefore form gaps in the subsidence zones.Two characters are described in this paper to identify underground mining-induced subsidence regions: (1) The subsidence zones are normally rectangularly or quadrilaterally shaped (related to the shape of longwall mining sites); and (2) The magnitude of subsidence rate over the subsidence zone edges is relatively large.By using the methodology describe above, four mine subsidence bowls have been detected within Figure 7a with a maximum subsidence rate in excess of ´10 cm¨yr ´1.
Remote Sens. 2016, 8, 543 14 of 17 Figure 8a shows the typical fringe patterns caused by mine subsidence derived from DInSAR interferograms with two obvious characteristics: (1) The earth surface where underground mining occurs will sink as the colour of the deformation pattern changes from yellow to blue (from the centre to the edge); and (2) The subsidence magnitude increases from the edge to the centre, therefore forming an oval or round shape (Hu et al., 2013).Figure 8b illustrates the underground mining subsidence pattern from the TS-InSAR velocity map.Regions suffering from rapid changes within a short period of time will de-correlate and therefore form gaps in the subsidence zones.Two characters are described in this paper to identify underground mining-induced subsidence regions: (1) The subsidence zones are normally rectangularly or quadrilaterally shaped (related to the shape of longwall mining sites); and (2) The magnitude of subsidence rate over the subsidence zone edges is relatively large.By using the methodology describe above, four mine subsidence bowls have been detected within Figure 7a with a maximum subsidence rate in excess of −10 cm•yr −1 .

Conclusions and Future Work
In this paper, the land surface stability over three test sites with different geological settings within the Southern Coalfield region was investigated using TS-InSAR technique.These three

Conclusions and Future Work
In this paper, the land surface stability over three test sites with different geological settings within the Southern Coalfield region was investigated using TS-InSAR technique.These three regions are: (1) Wollongong City, which is a relatively stable area; (2) Tahmoor Town, a small town affected by underground mining activities; and (3) Appin underground mining site, a region with many underground mining activities.Both C-band and L-band-derived measurements are compared over these three sites.Within the process, a first-order linear function is utilised to remove the linear phase "ramp" induced by residual orbital error.From the experiment result, we could see that the performance of both C-band and L-band is equally good over Wollongong City, where the dynamic range of subsidence is not significant and the subsidence rate is mostly between ´10 mm¨yr ´1 to 10 mm¨yr ´1.However, over Tahmoor and the Appin mine, their performances differ.Since the maximum displacement gradients that can be detected are different for L-band and C-band-based TS-InSAR, some rapid changes in land surface could cause the TS-InSAR to fail to estimate the correct displacements.It is well known that the L-band can perform better especially in underground mining regions and mining-affected regions where the deformation rate is much higher than city areas.L-band datasets are able to provide better results and achieve detailed deformation signals about the spatial distribution of the ground surface subsidence phenomena.The traditional DInSAR subsidence pattern is then compared with the TS-InSAR subsidence pattern induced by underground mining activities.The comparison shows a promising way to detect illegal mining sites using TS-InSAR techniques.For example, DIMDS proposed by [7] could be exploited to detect still-active illegal underground mining sites between two image acquisitions, while both still-active and ceased-illegal underground mining sites could be identified with the TS-InSAR subsidence pattern described in this paper.
TS-InSAR results derived from an L-band dataset can indeed be used to measure the subsidence around underground mining zones; possible future research is therefore to form the three-dimensional analyses and observe the time series horizontal subsidence around underground mine sites with both descending and ascending L-band image stacks.ALOS-2, the successor of ALOS-1 launched on 24 May 2014, could be utilised to conduct TS-InSAR analysis over underground mining areas.In addition, it is possible that both L-band and C-band datasets can be used together to detect the horizontal movement at certain locations, where good InSAR coherence can be preserved for C-band datasets, e.g., Sentinel-1A satellite, which has short revisit time.

Figure 2 .
Figure 2. (a) ALOS-1 PALSAR mean velocity result in LoS direction; (b) ENVISAT ASAR mean velocity result in LoS direction.The pink star represents the reference point while regions marked with red, blue and purple rectangular boxes are Wollongong city area, Appin mining sites and Tahmoor town region, respectively.The displacement results are resampled to a grid of 100 m × 100 m resolution; (c) Slope Map generated from SRTM DEM.

Figure 2 .
Figure 2. (a) ALOS-1 PALSAR mean velocity result in LoS direction; (b) ENVISAT ASAR mean velocity result in LoS direction.The pink star represents the reference point while regions marked with red, blue and purple rectangular boxes are Wollongong city area, Appin mining sites and Tahmoor town region, respectively.The displacement results are resampled to a grid of 100 m ˆ100 m resolution; (c) Slope Map generated from SRTM DEM.

Figure 3 .
Figure 3. L-band (a) and C-band (b) measurement over Wollongong city region in vertical direction (c) time series evolution over one MP pixel marked with red cross.

Figure 3 .
Figure 3. L-band (a) and C-band (b) measurement over Wollongong city region in vertical direction (c) time series evolution over one MP pixel marked with red cross.Remote Sens. 2016, 8, 543 11 of 17

Figure 4 .
Figure 4. Comparison over Wollongong between L-band and C-band TS-InSAR results.

Figure 4 .
Figure 4. Comparison over Wollongong between L-band and C-band TS-InSAR results.

Figure 5 .
Figure 5. (a) L-band and (b) C-band measurements over Tahmoor region, the resolution of each pixel is 30 m × 30 m (c) Google map over Tahmoor town (© Google Earth) (d) the profile line through Tahmoor town.

Figure 5 .
Figure 5. (a) L-band and (b) C-band measurements over Tahmoor region, the resolution of each pixel is 30 m ˆ30 m (c) Google map over Tahmoor town (© Google Earth) (d) the profile line through Tahmoor town.

Figure 7 .
Figure 7.Comparison over Appin mining site between (a) L-band and (b) C-band, the resolution of each pixel is 30 m × 30 m.

Figure 7 .
Figure 7.Comparison over Appin mining site between (a) L-band and (b) C-band, the resolution of each pixel is 30 m ˆ30 m.

Figure 7 .
Figure 7.Comparison over Appin mining site between (a) L-band and (b) C-band, the resolution of each pixel is 30 m × 30 m.