^{1}

^{2}

^{3}

In view of the limitations of traditional InSAR technology in selecting stable target point for orbit refining and surface subsidence inversion in complicated mining area, this paper proposes a time-series InSAR mining area subsidence monitoring method based on the fusion of multi threshold targets. On the basis of the traditional technology, the deviation threshold parameters, the regional window threshold parameters and the coherence threshold parameters are set to extract the relatively stable target points on the ground. Applying this method and traditional InSAR method to practical cases, the monitoring results of surface subsidence in the study area are obtained and verified. The results show that: (1) there are three mining subsidence areas in the mining area, the maximum annual average subsidence rate is -156 mm/a, and the maximum subsidence is -376 mm. Compared with the optical image data, the location of the mining subsidence area is consistent with the mining work area of the coal mine; (2) The absolute average difference of subsidence in the mining area using the two methods shall not exceed 12 mm. It shows that the InSAR method of fusing multi threshold targets can not only effectively overcome the limitations of traditional InSAR, but also ensure high accuracy, and has more advantages in the monitoring of surface subsidence in mining areas.

Long-term underground mining of resources will lead to surface subsidence and deformation, which will lead to geological disasters such as surface collapse, cracks and landslides. It will not only bring serious economic losses, but also endanger the lives and safety of residents. In order to effectively control the occurrence of disasters in the mining area, it is necessary to continuously and accurately monitor the surface changes in the mining area through monitoring means, to comprehensively and timely obtain the status of subsidence. Although the traditional coal mine deformation monitoring method can accurately obtain the deformation information of the measuring points [1], this method that only relies on limited discrete observation data cannot truly reflect the dynamic process of surface deformation. In recent years, the D-InSAR (differential interferometric synthetic aperture radar) of synthetic aperture radar (SAR) has developed rapidly with its outstanding advantages [2], [3], [4]. In particular, the introduction of multiple time series InSAR methods [5], [6], [7] has effectively remedied the problem of D-InSAR method, such as the monitoring accuracy, and the shortcomings that its reliability is affected by time, space de-coherence and atmospheric effects [8], [9], [10].

Affected by its own special conditions, the surface subsidence in the mining area is characterized by complex terrain conditions, large gradient, high subsidence rate and small deformation range, which makes it difficult to monitor [11], [12]. Small Baseline Subsets (SBAS) technology can effectively solve the problem of slow ground deformation monitoring in the mining area with its characteristics [13], [14]. In the processing, the target points are generally selected manually for orbit refining and settlement inversion [15], [16]. Although this method can effectively remove the residual phase and leveling effect, if the research area is not examined in advance, the blind selection based on ground target point information will lead to errors [17], [18]. Due to the mining is prone to bring a large number of surface subsidence, which makes the surface pixels less coherent, thus it is not easy to precisely obtain stable pixel points [19], [20]. Therefore, the selection of target points is a key research issue in current InSAR monitoring of the mining area.

To solve the above problems, this paper proposes an InSAR subsidence monitoring method based on fusion of multi-threshold targets. By setting the deviation threshold parameter, the region window threshold parameter and the coherence threshold parameter, the most stable target point is extracted and integrated into InSAR processing. Taking Yaojie mining area in Lanzhou as the research area, the application effect of InSAR method based on fusion of multi threshold target in subsidence monitoring of mining area is explored.

In this paper, the basic principle of InSAR method that fuses multi threshold targets is described firstly; then, taking Yaojie coal mine as the research object, the surface subsidence monitoring is carried out with the method in this paper; finally, the proposed method is compared with the traditional InSAR method, which proves the effectiveness of the proposed method.

The interference phase of the ground object generally contains a lot of noise phases. The ground object with relatively pure interference phase can be selected by measuring the noise contribution. SAR image is composed of multiple pixels. Assuming that the standard deviation of Gaussian noise of any pixel is* σ*, and the corresponding amplitude is in Rice distribution. Namely:

where, *I* is Bessel function, and *g* refers to the positive real number reflected energy of ground target. If the distribution is approximately close to Gaussian distribution, *g/σ*_{n}>4 is required. At this time, the amplitude dispersion index can be approximately considered as:

where, *σ*_{∂} is the amplitude standard deviation corresponding to the ground object, and *σ*_{nI} represents the standard deviation of the imaginary part. *m*_{∂} refers to the mean time series amplitude corresponding to the ground object and *D*_{∂} means the time series amplitude dispersion index.

The mathematical expression of coherence coefficient is:

where, *S*_{1} and *S*_{2} represent a set of interference pair images.* M* and *N* are the size of the local window. *i* and *j *represent the pixel coordinates of the extracted pixel. * is the multiplication of conjugate complex numbers.

In this paper, a method is proposed to first screen the ground object by setting a reasonable deviation index threshold parameter, and then extract the ground object with relatively pure interference phase. Secondly, since the coherence performance directly reflects the quality of the interference phase generated by data processing, the coherence threshold parameter is set to screen the ground object for the second time, and most of the pixels with poor coherence are eliminated, and finally the detection results are obtained.

When calculating the coherence coefficient, the setting of the size of the local moving window is one of the important factors to ensure the accuracy of the coherence coefficient. In general, the larger the window is, the higher the reliability is. However, the resolution is reduced to a certain extent. It is easy to misjudge nearby unstable points as stable target points, and it is easy to miss individual target points that exist separately in discrete distribution. Therefore, in order to balance the two mutually influencing factors of image resolution and coherence coefficient estimation accuracy, the study area is divided into windows before the second screening. In this paper, the whole region is divided into several sub-regions by setting the threshold parameter of the local moving window, and then the ground objects in each sub-region are filtered for coherence, and the target points larger than the coherence threshold value are extracted as reference points. Finally, all sub-regions are tessellated.

Assuming that SAR image interferes at *t*_{A} and *t*_{B} (*t*_{A}>*t*_{B}) to obtain the *j*_{th} interferogram, and the coordinates in the azimuth and range directions are (*x*,*r*), then the phase at this point can be expressed as:

where, λ is the radar wavelength, d(tB,x) and d(tA,x) are the cumulative deformation of the radar line of sight (LOS) direction, and d(t0,x)=0.

At the pixel point (*x*,*r*), it is assumed that the phase M interferogram is correspondingly generated by the primary image IE and secondary image IS, which can be expressed as:

The acquisition of the main image of the *j*_{th} (*j*=1, 2, ..., M) interferogram is later than that of the auxiliary image, that is, *IE*_{j}>*IS*_{j}, then the observation equation of its phase is as follows:

The above formula represents M observation equations with N unknown interference phases, and its matrix form is:

A is an approximate incidence matrix of M×N, and the interferogram generated in the data can determine matrix A. The data are all in a small baseline subset. When M is close to or less than N, the least square solution of the equation in the sense of the minimum norm can be obtained by singular value decomposition, and finally the time series shape variable can be obtained.

Yaojie Coal Mine is located at the junction of Gansu and Qinghai, with a longitude range of $102^{\circ} 50^{\prime} \sim 102^{\circ} 56^{\prime}$ and a latitude range of $36^{\circ} 20^{\prime} \sim 36^{\circ} 26^{\prime}$. The administrative division belongs to Honggu District of Lanzhou City, which is located in the west of the Loess Plateau with an altitude of 1,800~2,400 m. The terrain is relatively complex [21].

There are three mining areas in the mining area, including Yaojie No. 3 Mine, Jinhe Coal Mine and Haishiwan Coal Mine. Due to continuous underground mining for many years, geological disasters such as surface collapse and landslide occur frequently in the mining area, which seriously affect the safety of various construction facilities and mining production in the mining area [22]. The information of the study area is shown in Figure 1.

In this study, the single view complex data of Sentinel-1 A satellite with short time and space baseline in the interference wide range mode (IW) is selected. The satellite is equipped with 5.6cm C-band synthetic aperture radar. The revisit period is 12 days, and the image resolution is 5 m × 20 m. The time span of the experimental image selected in this paper is from January 2018 to April 2020, the radar incidence angle of the image is 39.06 o, and the azimuth resolution and range resolution are 13.91 m and 3.70 m respectively. In addition, in order to improve the accuracy of image registration, precise orbital ephemeris data is added to reduce the orbital error. The reference DEM data adopts 90 m resolution SRTM DEM data.

The traditional InSAR method uses manually selected target points. However, it is likely to increase the error if the selection is made without grasping the information of ground target points. In this paper, a fused multi-threshold target point extraction method is proposed.

After multiview registration of the original single view complex data, the intensity information map of the study area is obtained. All image data are analyzed using the method proposed in this paper to extract ground target points. First, in order to extract the best point as the target point to the greatest extent, the deviation threshold parameter is set as 3.2, and the first screening is conducted. The results are shown in subgraph (a) of Figure 2. A total of 29,453 target points were screened out in the experiment. Secondly, due to the large window of the study area, when calculating the pixel coherence, it is easy to generate clustered control points in the high coherence area, and it is easy to misjudge nearby unstable points as control points in the low coherence area, and it is easy to miss individual control points that exist separately in discrete distribution. Therefore, in order to balance high image resolution and high coherence coefficient accuracy, local window threshold parameters are set to decompose the entire region into several sub-regions, and then the coherence calculation is carried out. After repeated experiment, the region window threshold is set as 16 sqkm, and the entire region is divided into 132 sub-regions. Screening is carried out within each sub-region. Finally, the coherence threshold parameter is set as 0.9, and the second screening is carried out based on the results of the first screening. The screened ground control points are shown in subgraph (b) of Figure 2. Finally, 103 ground control points are screened in the experiment in total. The final screening results will be used as the target points for subsequent orbital refining for land subsidence inversion.

The InSAR method fusing multi-threshold target was used to analyze and process the data in the study area. The spatial baseline threshold was set at 5% of the critical limit, and the time baseline threshold was set at 180 days. The obtained time and space baseline connection diagram is shown in Figure 3.

After the interference processing, all the interference image pairs obtained are checked, and the image pairs with poor coherence are removed to avoid errors in deformation extraction. Then the target point selection method, that is, the InSAR method fusing multi threshold target proposed in this paper replaces the manual point selection method, and the InSAR orbit refining and re-flattening are carried out. Finally, the singular value decomposition method is used to estimate the deformation rate and residual terrain, and the atmospheric phase and residual terrain phase are removed through the space domain low-pass filtering and time domain high-pass filtering algorithms to obtain the temporal deformation rate.

After the above processing, the annual average rate chart of radar line of sight (LOS) direction in the study area from January 2018 to April 2020 is obtained. According to the formula ∆h=∆r/cosθ (θ is the radar incidence angle, ∆r is the settlement in the direction of radar line of sight, and ∆h is the settlement in the vertical direction), and the annual average vertical settlement rate in the study area is obtained, as shown in Figure 4. The monitoring results show that the surface subsidence in most areas of the study area is relatively stable, and the average annual subsidence rate is mainly -24~13 mm/a. There are three areas, A, B and C, with obvious surface subsidence, which is shown as a funnel like subsidence gradually increasing from the edge area to the center. The subsidence rate is -25~-156 mm/a for subsidence and 14~87 mm/a for uplift, which is basically consistent with the results obtained in literature [16]. According to the comparison with optical image data, the three subsidence funnels in the figure are basically consistent with the actual mining face of the coal mine. Figure 5 shows the side view of the profile along HD direction. Based on the comparative analysis of the underground mining data in the mining area, the underground mining in the mining area is the main reason for the subsidence of the three subsidence areas.

In order to analyze the temporal evolution law of the ground subsidence in the mining area, the time series analysis and calculation of the whole study area were carried out to obtain the settlement accumulation of each time point in the study period relative to the first phase image (2018-01-02). The settlement accumulation diagram of 11 important time nodes was taken, as shown in Figure 6. The results show that there is no obvious settlement in the study area before May 2018, and the settlement amount is concentrated between -5~-15mm. Since May 2018, the study area has seen obvious settlement, and the settlement value and settlement range gradually increase with time, and the settlement value reached at -376 mm in April 2020. Among them, in September 2018, January 2019 and April 2020, the settlement degree is relatively large, and the settlement values are 58 mm, 54 mm and 45 mm higher than before.

In order to further analyze the settlement development process of mining subsidence area, three areas (A, B, C) with the most obvious settlement are selected in the study area, and the time sequence settlement trend is analyzed. The analysis results are shown in Figure 7. Compared with the underground mining data of coal mines in the study area, the three areas are located in the three mining areas such as Yaojie No. 3 Coal Mine, Jinhe Coal Mine and Haishiwan Coal Mine. The three mining areas showed a linear subsidence trend in the whole study period, and the overall subsidence rate of Area A (No. 3 Mine) and Area B (Jinhe Coal Mine) was relatively stable. Zone C (Haishiwan Coal Mine) shows a phased change trend. The subsidence trend is relatively serious before March 2018, with a subsidence value of -78 mm. The subsidence trend is relatively stable after April 2018.

In order to verify the effectiveness of InSAR method integrating multi threshold targets in this paper, the same original data is used, and the traditional time series InSAR method is used to obtain the time series surface subsidence values of study areas A, B, and C. The monitoring results using the two methods are compared and analyzed. The comparison results are shown in in subgraphs (a), (b), and (c) of Figure 8.

The results in Figure 8 show that the maximum difference between the two monitoring results in Area A is 9.3 mm, the minimum is 0.12 mm, and the average value of the difference is 0.76 mm. The dispersion degree of the difference is calculated. The standard deviation range is 0.06-4.66 mm, and the average standard deviation is 2.06 mm. The maximum difference between the two monitoring results in Area B is 18.8 mm, the minimum is 0.39 mm, and the average value of the difference is 6.8 mm. The dispersion degree of the difference is calculated, and the standard deviation range is 0~9.4 mm, and the average standard deviation is 3.96 mm. The maximum difference between the two monitoring results in area C is 42.2 mm, the minimum is 0.8 mm, and the average value of the difference is 28.7 mm. The dispersion degree of the difference is calculated. The standard deviation range is 0.4~21.1 mm, and the average standard deviation is 12.3 mm. From the above analysis, it can be concluded that the monitoring results obtained by the InSAR method fusing multi threshold targets are roughly consistent with those obtained using the traditional time series InSAR method. Based on the monitoring accuracy of previous InSAR method, it can be considered that the results obtained by the method proposed in this paper are reliable.

Due to the complex surface environment of the mining area, it is extremely difficult to effectively extract stable target points for InSAR processing. Therefore, this paper proposes an InSAR monitoring method based on the fusion of multi threshold targets. This method is used to inverse the time series settlement information such as the annual average settlement rate and the time series cumulative settlement value in the study area, and the monitoring results using this method is compared with the monitoring results of the traditional InSAR method. The following conclusions are drawn from the study: (1) Through comparative analysis of the monitoring results of the two methods, the average values of the settlement differences of the three mining areas are 0.76 mm, 6.8 mm and 28.7 mm respectively, and the average values of the standard deviation are 2.06 mm, 3.96 mm and 12.3 mm respectively. The monitoring results using the two methods are relatively consistent, which proves that the accuracy of the InSAR method based on the fusion of multi threshold targets and the traditional InSAR method in extracting the surface settlement of the mining area can reach the sub-centimeter level. (2) In order to extract more target points with stable surface phase, the traditional InSAR method is more rigorous in selecting the research area, which makes it have many limitations in practical application. The method proposed in this paper can not only ensure the accuracy, but also overcome the limitations of the traditional InSAR method in the practical selection of target points. Especially when a large number of surface subsidence occurs in the mining area, resulting in poor surface pixels, the method proposed in this paper has better applicability, and can provide a new monitoring means for the follow-up monitoring and analysis of mining area subsidence. (3) The monitoring results show that there are three mining subsidence basins in the study area, respectively located in No. 3 Coal Mine (Area A), Jinhe Coal Mine (Area B) and Haishiwan Coal Mine (Area C). The maximum subsidence value is - 376 mm, the maximum annual average subsidence rate is -156 mm/a, and the subsidence trend is gradually increasing, so it is necessary to monitor the mining area for a long time series to provide disaster early warning.

The data used to support the findings of this study are available from the corresponding author upon request.

The authors declare that they have no conflicts of interest.

*Acadlore Trans. Geosci.*, 1(1), 22-32. https://doi.org/10.56578/atg010104

*Acadlore Trans. Geosci.*, vol. 1, no. 1, pp. 22-32, 2022. https://doi.org/10.56578/atg010104

*Acadlore Transactions on Geosciences*, v 1, pp 22-32. doi: https://doi.org/10.56578/atg010104

*Acadlore Transactions on Geosciences*, 1, (2022): 22-32. doi: https://doi.org/10.56578/atg010104