A method for monitoring large deformation surface area of ground subsidence

Abstract

The present application relates to a kind of ground subsidence large deformation surface area monitoring method, comprising the following steps: according to SBAS-InSAR technology obtains the subsidence data of n-1 period in study area;Using dynamic subsidence rate interpolation method obtains the theoretical subsidence value of SBAS-InSAR each period corresponding monitoring point position;According to each period SBAS-InSAR deformation original result and the theoretical subsidence amount of corresponding period monitoring point position, it is calculated to obtain the correction coefficient of each period SBAS-InSAR deformation original result in corresponding period monitoring point position;Using Kriging interpolation method is solved to obtain the correction coefficient of each period SBAS-InSAR deformation original result in corresponding period arbitrary point position in each period SBAS-InSAR deformation original result in corresponding period monitoring point position correction coefficient, finally obtains the corrected mining area dynamic subsidence result.The present application effectively overcomes the "phase ambiguity" and other problems existing in large deformation monitoring of InSAR technology, improves the monitoring accuracy of mining area large subsidence area.

Description

Technical Field

[0001] This invention belongs to the field of mining subsidence monitoring, and particularly relates to a method for monitoring large deformation areas of surface subsidence. Background Technology

[0002] Coal mines occupy a vital position in China's energy structure and are one of the pillars of the national economic development. With rapid economic development, the demand for coal is constantly rising, and coal mining activities are becoming increasingly frequent. However, the coal mining process inevitably triggers surface subsidence, a phenomenon that has a significant impact on the surrounding ecological environment and infrastructure. Subsidence can lead to structural damage to buildings, road collapses, changes in groundwater levels, and water pollution, posing potential threats to the lives and safety of surrounding residents. Furthermore, geological disasters caused by subsidence can cause enormous economic losses and social problems, urgently requiring effective monitoring methods to prevent and control these risks.

[0003] Traditional methods for monitoring subsidence in mining areas mainly rely on ground monitoring points and manual observation. While these methods provide relatively accurate subsidence data in local areas, the sparse distribution of monitoring points and the limited frequency of manual observation often make it difficult to comprehensively and timely reflect the subsidence situation of the entire mining area. This results in a lag in the judgment of subsidence trends, affecting the effectiveness of emergency decision-making.

[0004] In recent years, with the continuous development of remote sensing technology, Synthetic Aperture Radar Interferometry (InSAR) technology has gradually become an important means of monitoring overall subsidence in mining areas. InSAR technology utilizes radar imagery acquired by satellites and analyzes the phase differences of ground reflection signals to efficiently and accurately acquire surface deformation information over a large area. Compared with traditional ground monitoring methods, InSAR has advantages such as wide-area coverage, high spatial resolution, and high temporal frequency, enabling it to capture the overall subsidence trend of mining areas in a timely manner. However, when faced with large deformations, the phase information of InSAR may exhibit a "phase ambiguity" phenomenon, meaning that the phase values ​​cannot accurately reflect the true surface deformation. This is because InSAR technology relies on the measurement of phase differences; when the deformation exceeds the phase period of the radar wave, accurate calculation is difficult, leading to distortion of the monitoring results. Summary of the Invention

[0005] To address the aforementioned problems, this invention provides a method for monitoring large deformation areas of surface subsidence. This method combines subsidence data from known ground monitoring points with SBAS-InSAR results using interpolation techniques to generate more accurate monitoring results for mining area subsidence. This method offers the following advantages: First, it enables dynamic and continuous monitoring of mining area subsidence, capturing deformation changes in a timely manner; second, it globally corrects the SBAS-InSAR results using monitoring data from ground monitoring points, improving the accuracy of the monitoring results; finally, this method reduces the time and cost of manual monitoring, improves monitoring efficiency and response speed, and provides a scientific basis and strong support for safe production and environmental protection in coal mines. The application of this innovative technology will provide important guarantees for the sustainable development and safety management of the coal mining industry.

[0006] To achieve the above objectives, this application provides the following technical solution:

[0007] A method for monitoring large deformation areas of surface subsidence, characterized by comprising the following steps:

[0008] S1) Obtain n SAR images during the mining period in the study area, and obtain the original SBAS-InSAR deformation results of the study area for n-1 periods based on SBAS-InSAR technology;

[0009] S2) Collect settlement data from known monitoring points in the study area, and use dynamic settlement rate interpolation to analyze and process the settlement data of the monitoring points to obtain the theoretical settlement values ​​of the monitoring point locations for each period of SBAS-InSAR.

[0010] S3) Based on the original deformation results of SBAS-InSAR in each period and the theoretical settlement of the monitoring point location in the corresponding period, the settlement deviation of the monitoring point location in each period of SBAS-InSAR is calculated and used as the correction coefficient of the original deformation results of SBAS-InSAR in each period at the monitoring point location in the corresponding period.

[0011] S4) The Kriging interpolation method is used to calculate the correction coefficients of the original SBAS-InSAR deformation results at the corresponding monitoring point locations in each period, and the correction coefficients of the original SBAS-InSAR deformation results at any point location in the corresponding period are obtained.

[0012] S5) Based on the original SBAS-InSAR deformation results of each period and the correction coefficient of any point in the corresponding period, the corrected SBAS-InSAR deformation result of that point is obtained, thus obtaining the corrected dynamic subsidence result of the mining area.

[0013] Further, step S2) specifically includes:

[0014] Obtain the location information of known monitoring points in the study area and the settlement monitoring value S at each time node t;

[0015] The settlement rate of the monitoring point in adjacent monitoring cycles is calculated using the following expression:

[0016]

[0017] Where i represents the monitoring round at each monitoring point, and t i This represents the time node corresponding to the i-th round of monitoring. Indicates the monitoring point at t i+1 -t i Settlement rate over time period This means that at t=t i Settlement monitoring values ​​at the monitoring points;

[0018] Based on the dynamic settlement rate of each monitoring cycle, the acquisition time T of n long-term SAR image data is calculated. x At that time, the theoretical settlement value of each monitoring point is calculated using the following expression:

[0019]

[0020] Among them 1 <x≤n,t i ≤T x <t i+1 , Representing the time node t=T x The theoretical settlement value corresponding to the monitoring point.

[0021] Further, step S3) specifically includes:

[0022] Calculate the theoretical settlement of monitoring points at each period of SBAS-InSAR The calculation expression is:

[0023]

[0024] Based on the original deformation results of SBAS-InSAR at various periods Theoretical settlement at monitoring points during the corresponding period The settlement deviation of monitoring points at various periods in SBAS-InSAR is calculated using the following expression:

[0025]

[0026] in This represents the correction coefficient for the original SBAS-InSAR deformation results at the corresponding monitoring point locations for each period.

[0027] Further, step S4) specifically includes:

[0028] The correction coefficients of the original SBAS-InSAR deformation results for each period at the corresponding monitoring point locations are analyzed. Assuming that k monitoring points are deployed in the study area, the expression for calculating the correction coefficient of any point in the study area for each period is as follows:

[0029]

[0030] in Indicates t = T x The correction coefficient at any position (x,y), θ j Indicates the weighting coefficient. Indicates t = T x The correction coefficient corresponding to the location of monitoring point j.

[0031] Furthermore, step S4) specifically includes:

[0032] Calculate the weighting coefficients, assuming that the correction coefficients at any point (x, y) in the study area have the same expected value c and variance σ. 2 That is, it holds true for any point:

[0033]

[0034] Calculate the pairwise distance D and the corresponding semivariance γ(D) between any two known monitoring points. The distances are obtained from the coordinate values.

[0035]

[0036] Where D 1,2 Δx represents the straight-line distance between monitoring point 1 and monitoring point 2, and Δy represents the difference between the horizontal and vertical coordinates of the two monitoring points, respectively.

[0037] The semivariance function is defined as:

[0038]

[0039] in This represents the semivariance between monitoring point 1 and monitoring point 2;

[0040] Fit the relationship between all obtained (D, γ) values ​​using a spherical model, with the following function:

[0041]

[0042] Where C0 represents the nugget constant, C0+C represents the sill value, C represents the arch height, and a represents the range.

[0043] For any point in the study area, first calculate its distance D and corresponding semivariance C from all known monitoring points using the above method;

[0044] The goal of the Kriging interpolation method is to find a set of weighting coefficients that minimizes the difference between the estimated and true values ​​of the correction coefficients at any point in the study area, i.e.:

[0045]

[0046] in, ε(x) represents the estimated value of the correction coefficient at any point in the study area. i ,y i () represents the true value of the correction coefficient at any point in the study area;

[0047] Simultaneously satisfying the unbiased estimation condition:

[0048]

[0049] By constructing a short covariance matrix and a covariance vector, and combining them with unbiased constraints, a system of linear equations can be established, and the optimal weight coefficient vector can be obtained by solving it.

[0050] Further, step S5) specifically includes:

[0051] The original SBAS-InSAR deformation results for n-1 periods in the study area were corrected, and the calculation expression is as follows:

[0052]

[0053] in This indicates the corrected SBAS-InSAR deformation results.

[0054] Beneficial effects:

[0055] (1) This invention uses the settlement data of known ground monitoring points and SBAS-InSAR results to construct a correction model, and corrects the SBAS-InSAR results for each period. This effectively overcomes the problem of "phase ambiguity" in InSAR technology in large deformation monitoring and improves the monitoring accuracy of large settlement areas in mining areas.

[0056] (2) The present invention adopts the dynamic settlement rate interpolation method to further study the local settlement law of the monitoring point. It can interpolate the settlement rate of the monitoring point at different time periods and accurately reflect the dynamic changes of settlement in the mining area.

[0057] (3) In this invention, the Kriging interpolation method assigns different weights to each monitoring point to minimize interpolation error. This method ensures that the interpolation results are both highly accurate and avoid systematic biases, making the calibration model and SBAS-InSAR results more reliable. It can effectively capture the spatial variation trend of subsidence data, enabling the method to exhibit different subsidence patterns in different areas within the mining area, adapting to the different impacts of geological conditions and mining activities on ground deformation, thereby improving interpolation accuracy. Attached Figure Description

[0058] Figure 1 This is a flowchart of the method for monitoring large deformation areas of surface subsidence according to the present invention;

[0059] Figure 2 This is a flowchart of the Kriging interpolation method of the present invention;

[0060] Figure 3 This is a diagram showing the experimental results of the research area of ​​this invention. Detailed Implementation

[0061] The technical solutions provided by the present invention will be described more clearly below with reference to the embodiments and accompanying drawings, but the scope of protection claimed by the present invention is not limited to the following embodiments.

[0062] To verify the effectiveness of the method of this invention, we selected the Zhenchengdi mining area as the study area and took the correction process of the last SBAS result data (i.e., total settlement) as an example. The invention will be further described in detail below with reference to the accompanying drawings. The specific implementation steps are as follows:

[0063] (1) Data preparation: The mining area started mining on March 25, 2016 and ended mining in April 2017. 35 SAR data were obtained from March 2016 to September 2017 in the study area, and monitoring data from known monitoring points in the mining area were collected.

[0064] (2) SBAS-InSAR processing of the acquired multi-temporal SAR image data: First, select the study area and crop the SAR image to generate a short baseline connection map, and select the super master image.

[0065] Interferometric processing is performed to eliminate topographic phase errors and reduce noise effects, obtaining high-precision ground deformation information. Through orbital refinement and secondary inversion calculations, atmospheric effects and systematic errors are eliminated, yielding time-series surface deformation results for the mining area, such as... Figure 3 As shown, the maximum total settlement in the SBAS results obtained by processing 35 SAR data of the study area is 169.6 mm.

[0066] (3) The settlement data of each monitoring point were analyzed. First, the settlement rate of the monitoring point in different time periods was calculated. The rate was calculated by the two periods of settlement data. The settlement value of the monitoring point at the corresponding time of the SAR image was generated by the dynamic rate interpolation method. The maximum settlement value of the monitoring point was 1410.3 mm, indicating that the SBAS monitoring results in the area with severe settlement were unreliable.

[0067] Spatially align the monitoring points and SBAS results to obtain the settlement value corresponding to the SBAS at the monitoring point location, and calculate the correction coefficient for the monitoring point location based on the correction model.

[0068] (4) The correction coefficients for the entire study area at each period were calculated using the Kriging interpolation method, such as... Figure 2 As shown: Calculate the distance D between each pair of monitoring points and the corresponding semivariance γ(D);

[0069] The relationship between distance D and corresponding semivariance γ(D) at this time node is fitted using a spherical model, and the distance and semivariance between each point in the study area and the monitoring point are calculated based on the trained spherical model.

[0070] The optimal weighting coefficient θ for that point can be obtained by solving the objective function. j And calculate the correction coefficient;

[0071] Based on the above process, the corresponding correction coefficients are calculated for each point in the study area for each period, resulting in the area correction coefficients for the study area. Figure 3 As shown;

[0072] The same process was used to build a correction model for data from different periods, and the correction coefficients for any point in the study area were calculated.

[0073] (5) The obtained correction coefficients and SBAS results are superimposed, such as... Figure 3 As shown, the corrected SBAS results can more comprehensively and accurately represent the overall settlement of the mining area. The SBAS settlement data for each period were corrected according to the correction coefficients for different periods, resulting in more accurate dynamic monitoring results of the study area.

[0074] (6) To verify the accuracy of the corrected results, we used a monitoring point that was not involved in the calculation. Figure 3 The control points in the middle circle were used for verification. The measured data and the corrected SBAS results at five time points were selected for analysis. As shown in Table 1, the maximum error between the corrected SBAS results and the measured values ​​(all of which are cumulative settlement) is 4.8 mm, and the root mean square error is 3.79 mm, which verifies the effectiveness of the method of the present invention.

[0075] Table 1. Measured Values ​​and Correction Results

[0076] time 2016-10-02 2016-11-19 2017-02-11 2017-03-19 2017-05-06 Measured value 1182.3mm 1283.2mm 1300.8mm 1303.2mm 1307.2mm Corrected SBAS results 1180.6mm 1286.4mm 1296.0mm 1301.5mm 1302.8mm

[0077] The embodiments described above are merely examples of several implementations of the present invention, and while the descriptions are relatively specific and detailed, they should not be construed as limiting the scope of the invention. It should be noted that those skilled in the art can make various modifications and improvements without departing from the concept of the present invention, and these modifications and improvements all fall within the scope of protection of the present invention.

Claims

1. A method for monitoring large deformation areas of surface subsidence, characterized in that, Includes the following steps: S1) Obtain n SAR images during the mining period in the study area, and obtain the original SBAS-InSAR deformation results of the study area for n-1 periods based on SBAS-InSAR technology; S2) Collect settlement data from known monitoring points in the study area, and use dynamic settlement rate interpolation to analyze and process the settlement data of the monitoring points to obtain the theoretical settlement values ​​of the monitoring point locations for each period of SBAS-InSAR. S3) Based on the original deformation results of SBAS-InSAR in each period and the theoretical settlement of the monitoring point location in the corresponding period, the settlement deviation of the monitoring point location in each period of SBAS-InSAR is calculated and used as the correction coefficient of the original deformation results of SBAS-InSAR in each period at the monitoring point location in the corresponding period. S4) The correction coefficients of the original SBAS-InSAR deformation results at the corresponding monitoring point positions in each period are calculated using the Kriging interpolation method to obtain the correction coefficients of the original SBAS-InSAR deformation results at any point position in the corresponding period. S5) Based on the original SBAS-InSAR deformation results of each period and the correction coefficient of any point in the corresponding period, the corrected SBAS-InSAR deformation result of that point is obtained, thus obtaining the corrected dynamic subsidence result of the mining area. Step S2) specifically includes: Obtain the location information and time nodes of known monitoring points in the study area. Settlement monitoring values ; The settlement rate of the monitoring point in adjacent monitoring cycles is calculated using the following expression: in This indicates the monitoring rounds at each monitoring point. Indicates the first The corresponding time nodes for the round monitoring Indicates that the monitoring point is at Settlement rate over time period Indicates in Settlement monitoring values ​​at the monitoring points; Based on the dynamic settlement rate of each monitoring cycle, calculate the acquisition time points of n long-term SAR image data. At that time, the theoretical settlement value of each monitoring point is calculated using the following expression: in , , Indicates time node The theoretical settlement value corresponding to the monitoring point.

2. The method for monitoring large deformation areas of surface subsidence according to claim 1, characterized in that, Step S3) specifically includes: Calculate the theoretical settlement of monitoring points at each period of SBAS-InSAR The calculation expression is: = - Based on the original deformation results of SBAS-InSAR at various periods Theoretical settlement at monitoring points during the corresponding period The settlement deviation of monitoring points at various periods in SBAS-InSAR is calculated using the following expression: in This represents the correction coefficient for the original SBAS-InSAR deformation results at the corresponding monitoring point locations for each period.

3. The method for monitoring large deformation areas of surface subsidence according to claim 2, characterized in that, Step S4) specifically includes: The correction coefficients of the original SBAS-InSAR deformation results for each period at the corresponding monitoring point locations are analyzed. Assuming that k monitoring points are deployed in the study area, the expression for calculating the correction coefficient of any point in the study area for each period is as follows: ； in express any position The correction factor, Indicates the weighting coefficient. express Time monitoring point The correction coefficient corresponding to the position.

4. The method for monitoring large deformation areas of surface subsidence according to claim 3, characterized in that, Step S4) specifically also includes: Calculate the weighting coefficients, assuming any point in the study area. The correction coefficients all have the same expectation. and variance That is, it holds true for any point: , Calculate the pairwise distance D between known monitoring points and the corresponding semivariance. The distance is calculated based on the coordinate values: in This represents the straight-line distance between monitoring point 1 and monitoring point 2. and These represent the differences between the horizontal and vertical coordinates of the two monitoring points, respectively. The semivariance function is defined as: in This represents the semivariance between monitoring point 1 and monitoring point 2; All obtained The relationship between the two is fitted using a spherical model, and the function is as follows: ； in Represents the nugget constant. Indicates the base value. Indicates the height of the arch. Indicates range change; For any point in the study area, first calculate its distance to all known monitoring points using the method described above. and the corresponding semivariance ; The goal of the Kriging interpolation method is to find a set of weighting coefficients that minimizes the difference between the estimated and true values ​​of the correction coefficients at any point in the study area, i.e.: in, This represents the estimated value of the correction coefficient at any point in the study area. This represents the true value of the correction coefficient at any point in the study area; Simultaneously satisfying the unbiased estimation condition: By constructing a short covariance matrix and a covariance vector, and combining them with unbiased constraints, a system of linear equations can be established, and the optimal weight coefficient vector can be obtained by solving it.

5. The method for monitoring large deformation areas of surface subsidence according to claim 3, characterized in that, Step S5) specifically includes: The original SBAS-InSAR deformation results for n-1 periods in the study area were corrected, and the calculation expression is as follows: in This indicates the corrected SBAS-InSAR deformation results.

Images