Unmanned aerial vehicle surface subsidence monitoring method and system based on GM (1, 1)-AR model

Abstract

The invention provides an unmanned aerial vehicle surface subsidence monitoring method and system based on a GM (1, 1)-AR model, and belongs to the field of coal mining subsidence monitoring, and the method comprises the steps: S1, collecting multi-stage image data, carrying out the post-processing, generating a digital elevation model, extracting the point coordinates of ground feature points in the digital elevation model, and obtaining the elevation values of the feature points; s2, inputting the feature point elevation value into a GM (1, 1) model to obtain a feature point elevation fitting value, calculating to obtain a residual value based on the feature point elevation value and the feature point elevation fitting value, inputting the residual value into an AR model to obtain an optimized residual value, and calculating to obtain a corrected feature point elevation based on the feature point elevation fitting value and the optimized residual value; s3, performing elevation calculation based on the corrected feature points to obtain surface subsidence data; and correcting the elevation data of the feature points obtained by the unmanned aerial vehicle, thereby solving the problem that the unmanned aerial vehicle photogrammetry technology precision cannot meet the requirement of obtaining the surface subsidence parameter and cannot be suitable for the area needing accurate measurement.

Description

Technical Field

[0001] The present invention relates to the technical field of coal mining subsidence monitoring, and particularly to an unmanned aerial vehicle (UAV) surface subsidence monitoring method and system based on the GM(1,1)-AR model. Background Technique

[0002] As one of the main energy sources in China at present, coal resources play a crucial role in the national economic development and people's livelihood. The large-scale exploitation of coal resources has caused changes in the original underground rock formation structure, resulting in surface subsidence above, damaging surface farmland, buildings (structures), water system structures, etc., and destroying the ecological environment of the mining area. Accurately monitoring the coal mining subsidence area is a difficult problem that urgently needs to be solved in the current industry.

[0003] As a new emerging surface monitoring means, UAV photogrammetry technology plays an increasingly important role in surface subsidence monitoring with its characteristics of flexibility, high efficiency, precision and comprehensive coverage. The UAV photogrammetry technology collects photo data by carrying a high-resolution camera and flying over the area to be measured according to the specified flight route. These data can generate results such as Digital Elevation Model (DEM) and Digital Orthophoto Map (DOM) after post-processing. The DEM can accurately describe the spatial distribution of the geomorphic form, and by subtracting multiple-phase DEMs, the change of surface subsidence can be intuitively shown. In surface subsidence monitoring, UAV photogrammetry technology has significant advantages. First, its non-contact measurement method ensures the safety of the operation and avoids the risk of personnel entering dangerous areas in traditional measurement methods. Second, the UAV operation process is simple, and it can quickly obtain data of a large area, improving the monitoring efficiency. However, the accuracy of UAV photogrammetry technology cannot meet the requirements for obtaining surface subsidence parameters and is not applicable in areas that require precise measurement. Summary of the Invention

[0004] The technical problem to be solved by the present invention is how to solve the problem that the accuracy of UAV photogrammetry technology cannot meet the requirements for obtaining surface subsidence parameters and is not applicable in areas that require precise measurement.

[0005] The present invention solves the above technical problem through the following technical solutions: An unmanned aerial vehicle surface subsidence monitoring method based on the GM(1,1)-AR model, the method includes:

[0006] S1. Collect multi-phase image data and perform post-processing to generate a digital elevation model, extract the point coordinates of ground feature points in the digital elevation model, and obtain the elevation values of the feature points;

[0007] S2, input the feature point elevation value into the GM (1,1) model to obtain the feature point elevation fitting value, calculate the residual value based on the feature point elevation value and the feature point elevation fitting value, input the residual value into the AR model to obtain the optimized residual value, and calculate the corrected feature point elevation based on the feature point elevation fitting value and the optimized residual value;

[0008] S3. Calculate the surface subsidence data based on the corrected feature point elevations.

[0009] The present invention uses a feature point elevation data sequence to establish a GM (1, 1) model for data processing, obtains a feature point elevation data sequence and an internal residual sequence, uses an AR model to model the internal residual sequence, obtains an optimized residual value, and calculates a corrected feature point elevation based on the feature point elevation fitting value and the optimized residual value, thereby correcting the feature point elevation data obtained by a drone, and obtaining accurate surface subsidence data based on the corrected feature point elevation calculation, thereby solving the problem that the accuracy of drone photogrammetry technology cannot meet the requirements for obtaining surface subsidence parameters and cannot be applied to areas that require precise measurement.

[0010] Preferably, the process of collecting multiple phases of image data includes:

[0011] S111. Collect mining area data to obtain the scope of the surface area to be measured above the coal mining face;

[0012] S112. According to the scope of the area to be measured, plan the positions of the image control points for drone photogrammetry on professional mapping software, deploy the image control points on site and measure the coordinates of the image control points;

[0013] S113. Determine the flight parameters of the drone according to the range of the area to be measured, and collect multiple periods of image data over the area to be measured by the drone.

[0014] Preferably, the process of post-processing the multi-period image data includes:

[0015] S121, importing multi-period image data, and using the Brown-Conrady model to correct image distortion to obtain preprocessed image data;

[0016] S122, under the same spatial reference, hierarchically storing and displaying the preprocessed image data using a bottom-up sampling method according to different resolutions;

[0017] S123, performing matching of image points of the same name on the preprocessed image data to obtain relative coordinates of relative orientation of the image data;

[0018] S124, performing aerial triangulation on the relative coordinates of the relative orientation of the image data to obtain the absolute coordinates of the image data.

[0019] Preferably, the feature point elevation fitting value Z (1) The calculation method is:

[0020] Z (1) ={z (1) (2),z (1) (3),…,z (1) (n)}, where n represents the number of original data involved in modeling, n = m + 1, Z (0) is the numerical value of the characteristic point height sequence, a represents the development coefficient, and u represents the gray action amount.

[0021] Preferably, the residual value δ is calculated as follows:

[0022] δ=Z (1) -Z (0) ={δ(2),δ(3),...,δ(n)}.

[0023] Preferably, the optimized residual value δ (1) The calculation process includes:

[0024] S21. Input the residual value δ into the AR model, and obtain the mathematical model as follows:

[0025] δ (1) (n) = φ 1 δ(n-1)+φ 2 δ(n-2)+…+φ p δ(np)+e n

[0026] Among them, φ is the AR model parameter, φ=[φ 1 φ 2 … φ p ] T , e n is a white noise with a mean of 0 and a variance of σ, p is the model order; n ≥ 3, when n = 2, δ (1) (2) = δ(2);

[0027] The optimized residual value δ (1) ={δ (1) (2) δ (1) (3), ..., δ (1) (n)}:

[0028]

[0029] Preferably, the corrected feature point elevation Z (2) ={z (2) (2),z (2) (3),...,z (2) (n)}:

[0030]

[0031] The present invention also provides an unmanned aerial vehicle surface subsidence monitoring system based on the GM (1, 1)-AR model, the system comprising:

[0032] The data acquisition module is used to collect and post-process multiple image data to generate a digital elevation model, extract the point coordinates of ground feature points in the digital elevation model, and obtain the elevation values ​​of the feature points;

[0033] The data correction module is used to input the feature point elevation value into the GM (1,1) model to obtain the feature point elevation fitting value, calculate the residual value based on the feature point elevation value and the feature point elevation fitting value, input the residual value into the AR model to obtain the optimized residual value, and calculate the corrected feature point elevation based on the feature point elevation fitting value and the optimized residual value;

[0034] The subsidence value calculation module is used to calculate the surface subsidence data based on the modified feature point elevation.

[0035] Preferably, the process of collecting multiple phases of image data includes:

[0036] S111. Collect mining area data to obtain the scope of the surface area to be measured above the coal mining face;

[0037] S112. According to the scope of the area to be measured, plan the positions of the image control points for drone photogrammetry on professional mapping software, deploy the image control points on site and measure the coordinates of the image control points;

[0038] S113. Determine the flight parameters of the drone according to the range of the area to be measured, and collect multiple periods of image data over the area to be measured by the drone.

[0039] Preferably, the process of post-processing the multi-period image data includes:

[0040] S121, importing multi-period image data, and using the Brown-Conrady model to correct image distortion to obtain preprocessed image data;

[0041] S122, under the same spatial reference, hierarchically storing and displaying the preprocessed image data using a bottom-up sampling method according to different resolutions;

[0042] S123, performing matching of image points of the same name on the preprocessed image data to obtain relative coordinates of relative orientation of the image data;

[0043] S124, performing aerial triangulation on the relative coordinates of the relative orientation of the image data to obtain the absolute coordinates of the image data.

[0044] Preferably, the feature point elevation fitting value Z (1) The calculation method is:

[0045] Z (1) ={z (1) (2),z (1) (3),…,z (1) (n)}, where n represents the number of original data involved in modeling, n = m + 1, Z (0) is the numerical value of the characteristic point height sequence, a represents the development coefficient, and u represents the gray action amount.

[0046] Preferably, the residual value δ is calculated as follows:

[0047] δ=Z (1) -Z (0) ={δ(2),δ(3),...,δ(n)}.

[0048] Preferably, the optimized residual value δ (1) The calculation process includes:

[0049] S21. Input the residual value δ into the AR model, and obtain the mathematical model as follows:

[0050] δ (1) (n) = φ 1 δ(n-1)+φ 2 δ(n-2)+…+φ p δ(np)+e n

[0051] Among them, φ is the AR model parameter, φ=[φ 1 φ 2 …φp]T,e n is a white noise with a mean of 0 and a variance of σ, p is the model order; n ≥ 3, when n = 2, δ (1) (2) = δ(2);

[0052] The optimized residual value δ (1) ={δ (1) (2) δ (1) (3), ..., δ (1) (n)}:

[0053]

[0054] Preferably, the corrected feature point elevation Z (2) ={z (2) (2),z (2) (3),...,z (2) (n)}:

[0055]

[0056] The advantages provided by the present invention are:

[0057] (1) The present invention uses the feature point elevation data sequence to establish the GM (1, 1) model for data processing, obtains the feature point elevation data sequence and the internal residual sequence, uses the AR model to model the internal residual sequence, obtains the optimized residual value, calculates the corrected feature point elevation based on the feature point elevation fitting value and the optimized residual value, realizes the correction of the feature point elevation data obtained by the drone, and obtains accurate surface subsidence data based on the corrected feature point elevation calculation, thereby solving the problem that the accuracy of drone photogrammetry technology cannot meet the requirements for obtaining surface subsidence parameters and cannot be applied to areas that require precise measurement.

[0058] (2) The present invention can obtain data with accuracy that meets the requirements of aerial triangulation solution by post-processing the collected image data, including correcting image distortion, hierarchical storage and display of images according to resolution, and matching image points with the same name, thereby improving the efficiency of data processing. BRIEF DESCRIPTION OF THE DRAWINGS

[0059] Figure 1 A flow chart of a method for monitoring ground subsidence using an unmanned aerial vehicle based on a GM (1, 1)-AR model provided in an embodiment of the present invention;

[0060] Figure 2 A schematic diagram of the on-site deployment of image control points in the UAV surface subsidence monitoring method based on the GM (1, 1)-AR model provided in an embodiment of the present invention;

[0061] Figure 3 A schematic diagram of correcting image distortion in a UAV surface subsidence monitoring method based on a GM (1, 1)-AR model provided in an embodiment of the present invention;

[0062] Figure 4 A schematic diagram of feature point selection in a UAV surface subsidence monitoring method based on a GM (1, 1)-AR model provided in an embodiment of the present invention;

[0063] Figure 5 A schematic diagram of the modified characteristic point elevation calculation in the UAV surface subsidence monitoring method based on the GM (1, 1)-AR model provided in an embodiment of the present invention. DETAILED DESCRIPTION

[0064] In order to make the purpose, technical solution and advantages of the present invention more clearly understood, the technical solution of the present invention is clearly and completely described below in combination with specific embodiments and with reference to the accompanying drawings. Obviously, the described embodiments are part of the embodiments of the present invention, not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without creative work are within the scope of protection of the present invention.

[0065] like Figure 1 As shown, this embodiment provides a UAV surface subsidence monitoring method based on the GM (1, 1)-AR model, comprising the following steps:

[0066] Step S1, collecting multiple image data above the area to be measured, post-processing the multiple image data to generate a digital elevation model (DEM), extracting the point coordinates of ground feature points in the digital elevation model, and obtaining a ground feature point elevation data sequence;

[0067] The process of collecting multi-phase image data specifically includes:

[0068] S111. Collect mining area data to obtain the scope of the surface area to be measured above the coal mining face;

[0069] S112. According to the scope of the area to be measured, plan the image control point positions of the drone photogrammetry on professional mapping software, such as Figure 2 As shown, the image control points are arranged on site and the coordinates (x, y, z) of the image control points are measured;

[0070] S113. Determine the flight parameters of the UAV according to the scope of the area to be measured. The flight parameters of the UAV include the route, altitude, heading overlap rate, lateral overlap rate and other parameters. Use UAV photogrammetry to carry out surface subsidence monitoring. The UAV collects multiple periods of image data over the area to be measured.

[0071] The process of post-processing multi-phase image data specifically includes:

[0072] S121, importing multiple image data taken by the drone during on-site operation into a computer, using the Brown-Conrady model to correct image distortion, and obtaining pre-processed image data; refer to Figure 3 , image distortion includes radial distortion and tangential distortion. Radial distortion includes barrel distortion and pincushion distortion. The way to correct radial distortion of the image is:

[0073]

[0074] The way to correct the tangential distortion of the image is:

[0075]

[0076] In formula (1) and formula (2), (x, y) is the coordinate of the distorted image, (x', y') is the coordinate of the image after distortion correction, r is the distance from the pixel to the center of the image, and k is 1 , k 2 , k 3 is the radial distortion parameter, p 1 , p 2 is the tangential distortion parameter.

[0077] The cameras carried by drone photogrammetry are all non-metric cameras. Since the lens group used in non-metric cameras has not been strictly calibrated, there are often effects such as lens distortion and focal length changes, so image distortion needs to be corrected.

[0078] S122. Under the same spatial reference, the pre-processed image data is hierarchically stored and displayed by adopting a bottom-up sampling method according to different resolutions. The present invention establishes an influence pyramid and hierarchically stores and displays the data by adopting a bottom-up sampling method according to different resolutions. When operating the data, the user can select a certain level to operate according to the needs, which greatly improves the overall efficiency.

[0079] S123, the preprocessed image data may have a certain degree of overlap. The present invention performs same-name image point matching on the preprocessed image data, searches for similar feature points on two photos with a certain degree of overlap, and uses a grayscale-based least squares algorithm or a feature-based image matching strategy or a feature-based SIFT matching algorithm to perform same-name image point matching, and finally obtains relative coordinates of relative orientation of the image data;

[0080] S124, performing aerial triangulation on the relative coordinates of the relative orientation of the image data to obtain the absolute coordinates of the image data.

[0081] In DJI’s mapping data processing software, aerial triangulation is performed. During free network adjustment, the effects of camera focal length, lens distortion, etc. can be constrained. After manually piercing control points, aerial triangulation is performed and control network adjustment is performed to make the data meet the aerial triangulation accuracy requirements. After completing the aerial triangulation, the absolute coordinates of the photo points are calculated, and the photos are stitched to generate a digital elevation model DEM.

[0082] Step S2, inputting the feature point elevation value into the GM (1, 1) model to obtain the feature point elevation fitting value, calculating the residual value based on the feature point elevation value and the feature point elevation fitting value, inputting the residual value into the AR model to obtain the optimized residual value, and calculating the corrected feature point elevation based on the feature point elevation fitting value and the optimized residual value;

[0083] like Figure 4As shown, according to the layout of the working face, a strike line and a dip line are selected, and a feature point is selected every 5 meters on the strike line and every 5 meters on the dip line.

[0084] The GM(1,1) model is a grey system model, and its basic form is:

[0085]

[0086] Where m = 1, 2, ..., n-1, n is the number of original data involved in modeling, Z (1) is the fitting value of the feature point elevation, Z (0) is the numerical value of the feature point height sequence, a is the development coefficient, and u is the gray action.

[0087] Z (0) ={z (0) (1),z (0) (2),…,z (0) (n)}(4)

[0088] The development coefficient a and the ash action u can be obtained by the following formula:

[0089] Y GM =B GM X GM (5)

[0090] X GM =(B GM T B GM ) -1 B GM T Y GM (6)

[0091] Among them, Y GM is a constant term, Y GM =[z (0) (2),z (0) (3),...,z (0) (n)] T , X GM is the parameter matrix to be solved, X GM =[a,u] T , B GM is the coefficient matrix of the parameters to be solved, B GM The expression is:

[0092]

[0093] After determining the development coefficient a and gray action u of formula (3), the characteristic point elevation fitting value Z is calculated according to formula (3): (1) :

[0094] Z (1) ={z (1) (2),z (1) (3),…,z (1) (n)} (8)

[0095] The residual value δ is calculated as: δ = Z (1) -Z (0) ={δ(2),δ(3),...,δ(n)}

[0096] The AR model is a process that uses itself as a regression variable, that is, a linear regression model that uses a linear combination of data from several previous periods to describe variables at a certain moment in the future. It is a form of time series.

[0097] After the residual value δ obtained by the GM(1,1) model is expressed using the AR model, its specific mathematical model is as follows:

[0098] δ (1) (n) = φ 1 δ(n-1)+φ 2 δ(n-2)+…+φ p δ(np)+e n (9)

[0099] Among them, δ (1) is the optimized residual sequence, φ is the AR model parameter, φ=[φ 1 φ 2 … φ p ] T , e n is a white noise with a mean of 0 and a variance of σ. p is the model order, which indicates how long the model considers in the past. n ≥ 3. When n = 2, δ (1) (2) = δ(2);

[0100] The present invention selects the least square method to solve the AR model parameter φ, and the specific mathematical model is as follows:

[0101] Y=Hφ+ε (10)

[0102] in:

[0103] ε=[e p+1 e p+2 … e N ] T (11)

[0104] φ=[φ 1 φ 2 … φ p ] T (12)

[0105] Y=[y p+1 y p+2 … y N ] T (13)

[0106] In formula (11), ε represents the noise matrix; Y represents the true value of the AR model estimate, and H represents the coefficient matrix of the AR model parameters. The specific structure is as follows:

[0107]

[0108] According to the least squares principle, the calculation formula of AR(P) model parameters is as follows:

[0109] φ=(H T H)- 1 H T Y(15)

[0110] The optimized residual value δ (1) ={δ (1) (2) δ (1) (3), ..., δ (1) (n)}:

[0111]

[0112] The optimized residual value δ (1) Fitting value Z with feature point elevation (1) The sum is the corrected feature point elevation Z (2) :Z (2) ={z (2) (2),z (2) (3),...,z (2) (n)}:

[0113]

[0114] The present invention combines the GM (1, 1) model and the AR model to form a GM (1, 1) + AR model, uses the feature point elevation data sequence to establish the GM (1, 1) model for data processing, obtains the feature point elevation data sequence and the internal residual sequence, the internal residual is the difference between the known feature point elevation value and the feature point elevation fitting value, uses the AR model to model the internal residual sequence, obtains the optimized residual sequence value, the sum of the optimized residual sequence value and the feature point elevation fitting value is the final feature point elevation value, uses the GM (1, 1) + AR model to correct the feature point elevation data obtained by the unmanned aerial vehicle, and finally obtains accurate working face strike line and inclination line settlement data that can meet the parameter requirements, thereby solving the problem that the accuracy of the unmanned aerial vehicle photogrammetry technology cannot meet the requirements for obtaining surface settlement parameters and cannot be applied to areas that require precise measurement.

[0115] Step S3: Calculate the surface subsidence data based on the corrected feature point elevations.

[0116] The obtained corrected multi-period trend lines and characteristic point elevations of the inclination lines are subjected to difference processing.

[0117] Assume that the trend line and the trend line are observed for N periods respectively. After optimization processing using the GM(1,1)+AR model, the N-period characteristic point elevation values ​​of the trend line and the N-period characteristic point elevation values ​​of the trend line are obtained. The corrected characteristic point elevation value of the trend line is D, and the corrected characteristic point elevation value of the trend line is E. Then the final settlement value ΔD of the surface trend line and the final settlement value ΔE of the trend line are:

[0118] ΔD=D N -D 1 (18)

[0119] ΔE=E N -E 1 (19)

[0120] Among them, D 1 is the elevation value of the characteristic point corrected in the first phase of the trend line, D 2 is the elevation value of the characteristic point corrected at the Nth stage of the trend line, E 1 is the elevation value of the characteristic point corrected in the first phase of the trend line, E N It is the elevation value of the characteristic point corrected for the Nth period of the trend line.

[0121] The present invention also provides an unmanned aerial vehicle surface subsidence monitoring system based on the GM (1, 1)-AR model, comprising:

[0122] The data acquisition module is used to collect multi-phase image data and perform post-processing to generate a digital elevation model, extract the point coordinates of ground feature points in the digital elevation model, and obtain the elevation values ​​of the feature points; the process of collecting multi-phase image data includes:

[0123] S111. Collect mining area data to obtain the scope of the surface area to be measured above the coal mining face;

[0124] S112. According to the scope of the area to be measured, plan the positions of the image control points for drone photogrammetry on professional mapping software, deploy the image control points on site and measure the coordinates of the image control points;

[0125] S113. Determine the flight parameters of the drone according to the range of the area to be measured, and collect multiple periods of image data over the area to be measured by the drone.

[0126] The process of post-processing multi-period image data includes:

[0127] S121, importing multi-period image data, and using the Brown-Conrady model to correct image distortion to obtain preprocessed image data;

[0128] S122, under the same spatial reference, hierarchically storing and displaying the preprocessed image data using a bottom-up sampling method according to different resolutions;

[0129] S123, performing matching of image points of the same name on the preprocessed image data to obtain relative coordinates of relative orientation of the image data;

[0130] S124, performing aerial triangulation on the relative coordinates of the relative orientation of the image data to obtain the absolute coordinates of the image data.

[0131] The data correction module is used to input the feature point elevation value into the GM (1,1) model to obtain the feature point elevation fitting value, calculate the residual value based on the feature point elevation value and the feature point elevation fitting value, input the residual value into the AR model to obtain the optimized residual value, and calculate the corrected feature point elevation based on the feature point elevation fitting value and the optimized residual value;

[0132] Among them, the feature point elevation fitting value Z (1) The calculation method is:

[0133] Z (1) ={z (1) (2),z (1) (3),…,z (1) (n)}, where n represents the number of original data involved in modeling, n = m + 1, Z (0) is the numerical value of the characteristic point height sequence, a represents the development coefficient, and u represents the gray action amount.

[0134] The residual value δ is calculated as:

[0135] δ=Z (1) -Z (0) ={δ(2),δ(3),...,δ(n)}.

[0136] The optimized residual value δ (1) The calculation process includes:

[0137] S21. Input the residual value δ into the AR model, and obtain the mathematical model as follows:

[0138] δ (1) (n) = φ 1 δ(n-1)+φ 2 δ(n-2)+…+φ p δ(np)+e n

[0139] Where φ represents the AR model parameter, φ=[φ 1 φ 2 … φ p ] T , e n represents white noise with a mean of 0 and a variance of σ, p is the model order, n ≥ 3, when n = 2, δ (1) (2) = δ(2);

[0140] The optimized residual value δ (1) ={δ (1) (2) δ (1) (3), ..., δ (1) (n)}:

[0141]

[0142] Corrected feature point elevation Z (2) ={z (2) (2),z (2) (3),...,z (2) (n)}:

[0143]

[0144] The subsidence value calculation module is used to calculate the surface subsidence data based on the modified feature point elevation.

[0145] The above embodiments are only used to illustrate the technical solutions of the present invention, rather than to limit the same. Although the present invention has been described in detail with reference to the aforementioned embodiments, those skilled in the art should understand that the technical solutions described in the aforementioned embodiments may still be modified, or some of the technical features thereof may be replaced by equivalents. However, these modifications or replacements do not deviate the essence of the corresponding technical solutions from the spirit and scope of the technical solutions of the embodiments of the present invention.

Claims

1. A UAV surface subsidence monitoring method based on the GM (1, 1)-AR model is characterized by: Methods include: S1, collecting multiple image data and performing post-processing to generate a digital elevation model, extracting the point coordinates of ground feature points in the digital elevation model, and obtaining the elevation values ​​of the feature points; S2, input the feature point elevation value into the GM (1,1) model to obtain the feature point elevation fitting value, calculate the residual value based on the feature point elevation value and the feature point elevation fitting value, input the residual value into the AR model to obtain the optimized residual value, and calculate the corrected feature point elevation based on the feature point elevation fitting value and the optimized residual value; S3. Calculate the surface subsidence data based on the corrected feature point elevations.

2. The UAV surface subsidence monitoring method based on the GM (1, 1)-AR model according to claim 1 is characterized in that: The process of acquiring multi-phase image data includes: S111. Collect mining area data to obtain the scope of the surface area to be measured above the coal mining face; S112. According to the scope of the area to be measured, plan the positions of the image control points for drone photogrammetry on professional mapping software, deploy the image control points on site and measure the coordinates of the image control points; S113. Determine the flight parameters of the drone according to the range of the area to be measured, and collect multiple periods of image data over the area to be measured by the drone.

3. The UAV surface subsidence monitoring method based on the GM (1, 1)-AR model according to claim 1 is characterized in that: The process of post-processing multi-period image data includes: S121, importing multi-period image data, and using the Brown-Conrady model to correct image distortion to obtain preprocessed image data; S122, under the same spatial reference, hierarchically storing and displaying the preprocessed image data using a bottom-up sampling method according to different resolutions; S123, performing matching of image points of the same name on the preprocessed image data to obtain relative coordinates of relative orientation of the image data; S124, performing aerial triangulation on the relative coordinates of the relative orientation of the image data to obtain the absolute coordinates of the image data.

4. The UAV surface subsidence monitoring method based on the GM (1, 1)-AR model according to claim 1 is characterized in that: Feature point elevation fitting value Z (1) The calculation method is: Z (1) ={z (1) (2),z (1) (3),…,z (1) (n)}, where n represents the number of original data involved in modeling, n = m + 1, Z (0) is the numerical value of the characteristic point height sequence, a represents the development coefficient, and u represents the gray action amount.

5. The method for monitoring ground subsidence using an unmanned aerial vehicle based on the GM (1, 1)-AR model according to claim 4 is characterized in that: The residual value δ is calculated as: δ=Z (1) -Z (0) = {δ(2),δ(3),...,δ(n)}.

6. The method for monitoring ground subsidence using an unmanned aerial vehicle based on the GM (1, 1)-AR model according to claim 5 is characterized in that: The optimized residual value δ (1) The calculation process includes: S21. Input the residual value δ into the AR model, and obtain the mathematical model as follows: d (1) (n)=φ1δ(n-1)+φ2δ(n-2)+…+φ p δ(np)+e n Where φ is the AR model parameter, φ=[φ1φ2…φp]T, e n is a white noise with a mean of 0 and a variance of σ, p is the model order; n ≥ 3, when n = 2, δ (1) (2) = δ(2); The optimized residual value δ (1) ={δ (1) (2) δ (1) (3), ..., δ (1) (n)}:

7. The UAV surface subsidence monitoring method based on the GM (1, 1)-AR model according to claim 6 is characterized by: Corrected feature point elevation Z (2) ={z (2) (2),z (2) (3),...,z (2) (n)}:

8. The UAV surface subsidence monitoring system based on GM (1, 1)-AR model is characterized by: The system includes: The data acquisition module is used to collect and post-process multiple image data to generate a digital elevation model, extract the point coordinates of ground feature points in the digital elevation model, and obtain the elevation values ​​of the feature points; The data correction module is used to input the feature point elevation value into the GM (1,1) model to obtain the feature point elevation fitting value, calculate the residual value based on the feature point elevation value and the feature point elevation fitting value, input the residual value into the AR model to obtain the optimized residual value, and calculate the corrected feature point elevation based on the feature point elevation fitting value and the optimized residual value; The subsidence value calculation module is used to calculate the surface subsidence data based on the modified feature point elevation.

9. The UAV surface subsidence monitoring system based on the GM (1, 1)-AR model according to claim 8 is characterized by: The process of acquiring multi-phase image data includes: S111. Collect mining area data to obtain the scope of the surface area to be measured above the coal mining face; S112. According to the scope of the area to be measured, plan the positions of the image control points for drone photogrammetry on professional mapping software, deploy the image control points on site and measure the coordinates of the image control points; S113. Determine the flight parameters of the drone according to the range of the area to be measured, and collect multiple periods of image data over the area to be measured by the drone.

10. The UAV surface subsidence monitoring system based on the GM (1, 1)-AR model according to claim 8 is characterized in that: The process of post-processing multi-period image data includes: S121, importing multi-period image data, and using the Brown-Conrady model to correct image distortion to obtain preprocessed image data; S122, under the same spatial reference, hierarchically storing and displaying the preprocessed image data using a bottom-up sampling method according to different resolutions; S123, performing matching of image points of the same name on the preprocessed image data to obtain relative coordinates of relative orientation of the image data; S124, performing aerial triangulation on the relative coordinates of the relative orientation of the image data to obtain the absolute coordinates of the image data.

Images