A positioning method for deep underground space based on InSAR technology

Abstract

The application discloses a positioning method for deep underground space based on InSAR technology, and particularly relates to obtaining ground surface time-series deformation and cumulative deformation results based on a satellite platform SAR, inverting a ground engineering azimuth angle based on ground surface time-series deformation data, and an underground space inversion algorithm based on an IPIM model; the InSAR technology and the IPIM model are combined to effectively improve the inversion precision of the deep underground space; the application proposes to determine the azimuth angle of the underground space in advance based on the texture information of the time-series LOS ground surface deformation, which can reduce the required inversion parameters and improve the inversion efficiency; a hybrid intelligent genetic particle swarm optimization algorithm is used to invert the optimal characteristic parameters of the underground goaf; for the inversion precision of the deep underground space, the average inversion precision can reach 3.62% under the condition of known geological parameters; and the application can invert the underground space parameters with high precision.

Description

Technical Field

[0001] This invention relates to the field of positioning methods for deep underground spaces, and more specifically to a positioning method for deep underground spaces based on InSAR technology. Background Technology

[0002] The implementation of underground engineering projects disrupts the original stress balance of the surrounding rock strata, causing stress redistribution until a new equilibrium is reached. This process also involves rock movement and fracturing, inevitably inducing surface deformation. This provides the possibility for underground spatial positioning based on surface deformation.

[0003] Currently, underground space exploration and location technologies primarily rely on geophysical methods, including microgravity techniques, transient electromagnetic methods, ground-penetrating radar monitoring, and temperature measurements. These methods utilize the physical and chemical changes in rock masses and the Earth's surface induced by underground engineering to invert and locate the position and spatial distribution of underground spaces. However, all of these methods require significant human, material, and financial resources, are inefficient, difficult to implement, and depend on prior information. These limitations make it difficult for geophysical methods to achieve large-scale underground space location and inversion.

[0004] Interferometric synthetic aperture radar (InSAR) is an active microwave remote sensing technology that uses phase interferometry to acquire data on Earth's surface deformation. This technology offers advantages such as a large monitoring range, high measurement accuracy, low cost, all-weather and 24 / 7 capability, and a large amount of archived data. InSAR has become a powerful tool for monitoring geological hazards, including earthquakes, volcanoes, iceberg drift, landslides, coal mine subsidence, and coal fires. Applying InSAR technology to the inversion and location of underground spaces can significantly compensate for the shortcomings of traditional geophysical methods.

[0005] InSAR-based underground space positioning methods include center positioning and spatial morphology inversion. Center positioning typically uses the point of maximum deformation as the center of the underground space. This method is only applicable to horizontally shaped underground spaces and struggles to achieve good positioning accuracy for inclined underground spaces. Spatial morphology inversion is based on the functional relationship between surface deformation and underground space; therefore, the accuracy of this functional relationship directly affects the inversion accuracy. Existing spatial morphology inversion algorithms are mostly based on probability integral models, which are well-suited for shallow, large-scale, and highly deformed underground space positioning; however, these models are not suitable for deep underground engineering. Summary of the Invention

[0006] To address the aforementioned technical shortcomings, the purpose of this invention is to provide a positioning method for deep underground spaces based on InSAR technology. This method involves setting up a sampling unit and a heating unit within a housing. A battery pack provides electrical energy to drive the heating unit and the motor inside the sampling unit to extract and heat the asphalt to be sampled, ensuring the fluidity of the asphalt and facilitating manual operation and transportation.

[0007] To solve the above-mentioned technical problems, the present invention adopts the following technical solution:

[0008] This invention provides a method for locating deep underground spaces based on InSAR technology, comprising the following steps:

[0009] S1: Obtain the temporal and cumulative deformation results of the land surface based on satellite platform SAR;

[0010] S2: Inversion of the azimuth angle of underground engineering based on surface time-series deformation data, specifically:

[0011] S21. First, by setting a deformation threshold, the time-series deformation result map is converted into a binary map;

[0012] S22. Use an edge detection algorithm to extract texture from a binary image;

[0013] S23. Use two-dimensional fast Fourier transform to convert the texture map into a spectrum map;

[0014] S24. Perform Hough transform on the spectrogram and combine it with local peak detection to obtain the main direction of each time-series deformation map;

[0015] S25. Average all the obtained main directions and combine them with the temporal deformation development trend to finally determine the azimuth of the underground space;

[0016] S3: Underground space inversion algorithm based on IPIM model:

[0017] The IPIM model, based on stochastic medium theory, can effectively describe the functional relationship between deep underground space and surface deformation. This model assumes the underground space is rectangular in shape and further subdivides it into numerous micro-units. It assumes that these micro-units are independent and move independently, and obtains the surface deformation induced by the entire underground space through an accumulation method. The specific functional expression is as follows:

[0018]

[0019] in:

[0020] W max =mqcosα

[0021] l = L1 - s1 - s2

[0022]

[0023]

[0024]

[0025] Wherein, W and U are the surface subsidence and horizontal displacement, respectively;

[0026] m represents the height of the underground space;

[0027] q and b are the settlement coefficient and the horizontal displacement coefficient, respectively;

[0028] α, β, and θ represent the dip angle of the underground space, the main influence angle, and the influence propagation angle, respectively.

[0029] Burial depths of H and H1;

[0030] L1 and L2 are the main length and perpendicular length of the underground space, respectively.

[0031] l is the length calculated in the main direction;

[0032] s1 and s2 are the inflection point offset distances in the main direction, respectively. Since the length perpendicular to the main direction is small, there is no inflection point offset distance.

[0033] S is the horizontal distance from the center of the underground space to the calculation origin;

[0034] ρ(y) is the horizontal offset function of the unit at different locations from the center of its subsidence basin;

[0035] r(y) is the influence radius function of elements at different locations;

[0036] When the underground space is horizontal, the angle of inclination should be ignored;

[0037] The random medium theory posits that there are no connections between elements, and their movements are independent. Therefore, the deformation at any point on the Earth's surface can be represented as:

[0038]

[0039] Among them, U E and U N These represent horizontal displacements in the east-west and north-south directions, respectively.

[0040] and These represent the angles between the coal seam mining direction and the east and north directions, respectively.

[0041] Preferably, the specific process of step 1 is as follows:

[0042] S11. Collect SAR image data of the monitored area acquired by the satellite platform;

[0043] S12. After performing high-precision registration, cropping, and filtering preprocessing on the acquired SAR image data, time-series D-InSAR technology is used to obtain differential interferograms of nearby time SAR images; high-precision processing requires at least one-thousandth of a pixel.

[0044] S13. Select a highly coherent and stable point as a reference point to perform phase unwrapping operation and obtain the temporal deformation results of the land surface; highly coherent means that the coherence of urban areas is greater than 0.9 and the coherence of non-urban areas is greater than 0.7.

[0045] S14. The obtained deformation results are superimposed sequentially to obtain the cumulative surface deformation result;

[0046] Preferably, step S2 consists of the following steps:

[0047] Step 21: Set a deformation threshold for each scene's surface deformation result. The deformation threshold should be set to the regional average surface deformation or 0. Convert the deformation results into a binary map according to the following rules:

[0048]

[0049] In the formula, Value is the pixel value of the transformed binary image, and D... value D represents the surface deformation value. thre The set deformation threshold;

[0050] Step 22: Edge detection algorithms can effectively remove irrelevant information, significantly reduce data volume, and improve processing speed. Binarization is also used to better achieve edge detection. Commonly used edge detection operators include: Roberts Cross operator, Prewitt operator, Sobel operator, Kirsch operator, Marr-Hildreth operator, Canny operator, Laplacian operator, etc. The Sobel operator is used to extract texture from binary images; the formulas for horizontal and vertical edge detection are as follows:

[0051]

[0052]

[0053]

[0054] In the formula, G x For horizontal edge detection images, G y A is the image for vertical edge detection, A is the image for gradient calculation, and G is the gradient calculation result, i.e., the texture of the binary image.

[0055] Step 23: The two-dimensional Fast Fourier Transform (FFT) can transform an image from the spatial domain to the frequency domain. This not only reflects the image's relative structure but also facilitates the implementation of algorithms such as convolution and target recognition. Furthermore, based on frequency information, it can remove low-frequency and noise signals from the image signal, and is also more conducive to extracting the principal azimuth angle of underground spaces. The formula for the two-dimensional FFT is as follows:

[0056]

[0057] In the formula, M and N are the row and column data of the image;

[0058] Step 24: Use Hough transform to detect texture maps and spectrograms, and combine it with local peak detection algorithm to determine the deformation center and the main deformation method;

[0059] Step 25: Take the average value of the main directions obtained in Step 4, and combine it with the development trend of temporal deformation to determine the final azimuth of the underground space, as shown in the following formula:

[0060]

[0061] In the formula, The final azimuth angle. The main direction of each image is denoted by n, where n is the number of images, and θ is determined based on the development trend of temporal deformation, taking values ​​of 0°, 90°, 180°, or 270°.

[0062] Preferably, in step S3, there are many underground space parameters to be inverted, and they have a highly nonlinear relationship with surface deformation. Therefore, an intelligent algorithm is introduced in this study to obtain the global optimal solution. Based on the intelligent genetic particle swarm optimization algorithm, the underground space inversion algorithm implementation process is as follows:

[0063] Step 1: Input geological parameters and main azimuth angle, and preset the search domain for underground space parameters, as well as the parameter settings for the intelligent genetic particle swarm optimization algorithm;

[0064] The geological parameters q, b, β, θ, s1, s2 and the principal azimuth angle finally obtained in step S2 are used. In the known information input in the inversion algorithm, the geological parameters are determined based on geological conditions, field measurement data, and engineering construction methods. The underground spatial parameters H, L1, L2, m, α, and X are set through analogy. c Y c Upper and lower limits; where X c Y c The geographic coordinates of the underground space center; the parameters of the intelligent genetic particle swarm optimization algorithm include population size, learning factor, maximum and minimum inertia factor, crossover probability, mutation probability, maximum number of iterations, and termination condition;

[0065] Step 2: Randomly generate population particles and calculate the three-dimensional coordinates of the Earth's surface based on IPIM;

[0066] Based on the search domain of underground space parameters set in step one, the number of particles of the set population size is randomly generated, and the three-dimensional deformation value of the surface is calculated based on the IPIM model. The calculation formula is as follows:

[0067]

[0068] Step 3: Constructing a model of the relationship between three-dimensional surface deformation and line-of-sight deformation;

[0069] The surface deformation information acquired by SAR satellites is the line-of-sight deformation information D. LOS Step two simulates the three-dimensional deformation W and U of the Earth's surface. N U E The information, based on the SAR satellite's incident angle and heading angle, establishes the relationship between the two as follows:

[0070] D LOS =Wcosθ h -sinθ h (U N cosα h +U E sinα h )

[0071] Where, θ h and α h Given the incident angle and heading angle of the SAR satellite, respectively, the functional relationship between the deep subsurface space and LOS deformation can be constructed by combining the IPIM model:

[0072]

[0073] Step 4: Calculate the population particle fitness function;

[0074] D obtained from SAR data inversion LOS And the D obtained in step three los Construct the particle fitness function, as shown in the following formula:

[0075]

[0076] Where fit represents the particle fitness value; The LOS deformation of surface point i is calculated using Formula 4. denoted as LOS deformation of surface point i obtained by SAR technology inversion, and n is the total number of pixels in the study area;

[0077] Step 5: Population particle update;

[0078] The population particles are sorted according to their fitness values. Half of the particles with poor fitness are discarded, while the particles with good fitness are regenerated into a new particle population using particle update algorithm and genetic algorithm respectively, so as to maintain the set population size.

[0079] Then, repeat steps three through five until the termination condition is met or the maximum number of iterations is reached, thus obtaining the optimal parameters for the underground space.

[0080] The beneficial effects of this invention are as follows:

[0081] 1. InSAR technology can acquire surface deformation information over a wide area with high precision, while IPIM can better describe the relationship between surface deformations in deep underground space. Combining the two can effectively improve the inversion accuracy of deep underground space.

[0082] 2. In order to reduce the impact of surrounding projects, this invention proposes to determine the azimuth of underground space in advance based on the texture information of temporal LOS surface deformation. This can reduce the required inversion parameters and improve the inversion efficiency.

[0083] 3. The optimal characteristic parameters of underground mining areas are inverted by using a hybrid intelligent genetic particle swarm optimization algorithm. For the inversion accuracy of deep underground space, the average inversion accuracy can reach 3.62% under known geological parameters. This invention can invert underground space parameters with high accuracy. Attached Figure Description

[0084] To more clearly illustrate the technical solutions in the embodiments of the present invention or the prior art, the drawings used in the description of the embodiments or the prior art will be briefly introduced below. Obviously, the drawings described below are only some embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0085] Figure 1 This is a flowchart of the deep subsurface space inversion method based on InSAR and IPIM in this embodiment;

[0086] Figure 2 This is a schematic diagram of the IPIM model in this embodiment;

[0087] Figure 3 This is a surface deformation map obtained based on InSAR technology in this embodiment;

[0088] Figure 4 This is a diagram showing the actual data positioning results in this embodiment. Detailed Implementation

[0089] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0090] Taking a known underground mining space in Fengfeng Mining Area, Hebei Province as an example, this study uses seven RADARSAT-2 datasets (from September 19, 2015 to March 5, 2016) for experimental verification; see [link to relevant documentation]. Figure 1 A method for locating deep underground spaces based on InSAR technology, with the following specific implementation plan:

[0091] S1: Obtain the temporal and cumulative deformation results of the land surface based on satellite platform SAR;

[0092] Based on the selected SAR data, six interferometric pairs were formed using D-InSAR technology. The specific processing flow included: precise registration of the SAR data, followed by interferometric processing to generate interferograms; filtering and removal of flat-land effects from the generated interferograms; selecting high-coherence regions as stable areas, and performing phase unwrapping on pixels with coherence higher than 0.3 to obtain LOS surface deformation information; and converting the surface deformation information from the radar coordinate system to the geographic coordinate system based on geocoding. Due to the small study area, the influence of atmospheric noise was ignored during processing; the final surface deformation results are as follows: Figure 3 As shown.

[0093] S2: Inversion of the azimuth angle of underground engineering based on surface time-series deformation data, specifically:

[0094] The specific steps of step S2 are as follows:

[0095] Step 21: Set a deformation threshold for each scene's surface deformation result. The deformation threshold should be set to the regional average surface deformation or 0. Convert the deformation results into a binary map according to the following rules:

[0096]

[0097] In the formula, Value is the pixel value of the transformed binary image, and D... value D represents the surface deformation value. thre The set deformation threshold;

[0098] Step 22: Edge detection algorithms can effectively remove irrelevant information, significantly reduce data volume, and improve processing speed. Binarization is also used to better achieve edge detection. Commonly used edge detection operators include: Roberts Cross operator, Prewitt operator, Sobel operator, Kirsch operator, Marr-Hildreth operator, Canny operator, Laplacian operator, etc. The Sobel operator is used to extract texture from binary images; the formulas for horizontal and vertical edge detection are as follows:

[0099]

[0100]

[0101]

[0102] In the formula, G x For horizontal edge detection images, G y A is the image for vertical edge detection, A is the image for gradient calculation, and G is the gradient calculation result, i.e., the texture of the binary image.

[0103] Step 23: The two-dimensional Fast Fourier Transform (FFT) can transform an image from the spatial domain to the frequency domain. This not only reflects the image's relative structure but also facilitates the implementation of algorithms such as convolution and target recognition. Furthermore, based on frequency information, it can remove low-frequency and noise signals from the image signal, and is also more conducive to extracting the principal azimuth angle of underground spaces. The formula for the two-dimensional FFT is as follows:

[0104]

[0105] In the formula, M and N are the row and column data of the image;

[0106] Step 24: Use Hough transform to detect texture maps and spectrograms, and combine it with local peak detection algorithm to determine the deformation center and the main deformation method;

[0107] Step 25: Take the average value of the main directions obtained in Step 4, and combine it with the development trend of temporal deformation to determine the final azimuth of the underground space, as shown in the following formula:

[0108]

[0109] In the formula, The final azimuth angle. The main direction of each image is denoted by n, where n is the number of images, and θ is determined based on the development trend of temporal deformation, taking values ​​of 0°, 90°, 180°, or 270°. This step yields an azimuth angle of 236.93°, while the actual value is 236.00°, a difference of only 0.93°.

[0110] S3: Underground space inversion algorithm based on IPIM model:

[0111] The IPIM model, based on stochastic medium theory, can effectively describe the functional relationship between deep underground space and surface deformation. This model assumes the underground space is rectangular in shape and further subdivides it into numerous micro-units. It assumes that these micro-units are independent and move independently, and obtains the surface deformation induced by the entire underground space through an accumulation method. See the schematic diagram of the IPIM model in this embodiment. Figure 2 The specific function expression is as follows:

[0112]

[0113] in:

[0114] W max =mqcosα

[0115] l = L1 - s1 - s2

[0116]

[0117]

[0118]

[0119] Wherein, W and U are the surface subsidence and horizontal displacement, respectively;

[0120] m represents the height of the underground space;

[0121] q and b are the settlement coefficient and the horizontal displacement coefficient, respectively;

[0122] α, β, and θ represent the dip angle of the underground space, the main influence angle, and the influence propagation angle, respectively.

[0123] Burial depths of H and H1;

[0124] L1 and L2 are the main length and perpendicular length of the underground space, respectively.

[0125] l is the length calculated in the main direction;

[0126] s1 and s2 are the inflection point offset distances in the main direction, respectively. Since the length perpendicular to the main direction is small, there is no inflection point offset distance.

[0127] S is the horizontal distance from the center of the underground space to the calculation origin;

[0128] ρ(y) is the horizontal offset function of the unit at different locations from the center of its subsidence basin;

[0129] r(y) is the influence radius function of elements at different locations;

[0130] When the underground space is horizontal, the angle of inclination should be ignored;

[0131] The random medium theory posits that there are no connections between elements, and their movements are independent. Therefore, the deformation at any point on the Earth's surface can be represented as:

[0132]

[0133] Among them, U E and U N These represent horizontal displacements in the east-west and north-south directions, respectively.

[0134] and These represent the angles between the coal seam mining direction and the east and north directions, respectively.

[0135] In step S3, there are many underground space parameters to be inverted, and they have a highly nonlinear relationship with surface deformation. Therefore, an intelligent algorithm is introduced to obtain the global optimal solution. Based on the intelligent genetic particle swarm optimization algorithm, the underground space inversion algorithm implementation process is as follows:

[0136] Step 1: Input geological parameters and main azimuth angle, and preset the search domain for underground space parameters, as well as the parameter settings for the intelligent genetic particle swarm optimization algorithm;

[0137] The geological parameters q, b, β, θ, s1, s2 and the principal azimuth angle finally obtained in step S2 are used. In the known information input in the inversion algorithm, the geological parameters are determined based on geological conditions, field measurement data, and engineering construction methods. The underground spatial parameters H, L1, L2, m, α, and X are set through analogy. c Y c Upper and lower limits; where X c Y c The geographic coordinates of the underground space center; the parameters of the intelligent genetic particle swarm optimization algorithm include population size, learning factor, maximum and minimum inertia factor, crossover probability, mutation probability, maximum number of iterations, and termination condition;

[0138] Step 2: Randomly generate population particles and calculate the three-dimensional coordinates of the Earth's surface based on IPIM;

[0139] Based on the search domain of underground space parameters set in step one, the number of particles of the set population size is randomly generated, and the three-dimensional deformation value of the surface is calculated based on the IPIM model. The calculation formula is as follows:

[0140]

[0141] Step 3: Constructing a model of the relationship between three-dimensional surface deformation and line-of-sight deformation;

[0142] The surface deformation information acquired by SAR satellites is the line-of-sight deformation information D. LOS Step two simulates the three-dimensional deformation W and U of the Earth's surface. N U E The information, based on the SAR satellite's incident angle and heading angle, establishes the relationship between the two as follows:

[0143] D LOS =Wcosθ h -sinθ h (U N cosα h +U E sinα h )

[0144] Where, θ h and α h Given the incident angle and heading angle of the SAR satellite, respectively, the functional relationship between the deep subsurface space and LOS deformation can be constructed by combining the IPIM model:

[0145]

[0146] Step 4: Calculate the population particle fitness function;

[0147] D obtained from SAR data inversion LOS And the D obtained in step three los Construct the particle fitness function, as shown in the following formula:

[0148]

[0149] Where fit represents the particle fitness value; D losi D is the LOS deformation of surface point i calculated by formula 4. LOSi denoted as LOS deformation of surface point i obtained by SAR technology inversion, and n is the total number of pixels in the study area;

[0150] Step 5: Population particle update;

[0151] The population particles are sorted according to their fitness values. Half of the particles with poor fitness are discarded, while the particles with good fitness are regenerated into a new particle population using particle update algorithm and genetic algorithm respectively, so as to maintain the set population size.

[0152] Then, repeat steps three through five until the termination condition is met or the maximum number of iterations is reached, thus obtaining the optimal parameters for the underground space.

[0153] In this embodiment, the geological parameters required for the inversion were determined based on the engineering geological conditions, on-site measured data, and construction method, as follows: q = 0.40, b = 0.25, β = 60.81°, θ = 85°, s1 = s2 = 20m. Analogous to similar engineering conditions in the surrounding area, the parameter search domain for this underground space was set as follows: H ∈ [300m, 1300m], L1 ∈ [400m, 1000m], L2 ∈ [0m, 300m], m ∈ [0m, 8m], α ∈ [0°, 45°], X c ∈[600m,2000m], Y c ∈[300m,1000m]. In the embodiment, the population size of GA-PSO is set to 500, the learning factor is 1.50, the maximum and minimum inertia factors are 0.9 and 0.5 respectively, the crossover probability is 0.5, the mutation probability is 0.1, the maximum number of iterations is 1000, and the fitness value is less than 0.1.

[0154] After the above steps, the location and inversion of the underground space were finally achieved, and the results are as follows: Figure 4 As shown in the figure. The specific results for each underground space geometric parameter are shown in Table 1:

[0155] Table 1. Results and errors of underground space parameter inversion

[0156]

[0157] Except for the dip angle of the underground space, the relative errors of other underground space geometric parameters are all below 5%, with an average value of 3.62%. This indicates that the method proposed in this invention can effectively and accurately invert the geometric parameters of underground space.

[0158] Obviously, those skilled in the art can make various modifications and variations to this invention without departing from its spirit and scope. Therefore, if these modifications and variations fall within the scope of the claims of this invention and their equivalents, this invention also intends to include these modifications and variations.

Claims

1. A positioning method for deep underground space based on InSAR technology, characterized in that: Includes the following steps: S1: Obtain the temporal and cumulative deformation results of the land surface based on satellite platform SAR; S2: Inversion of the azimuth angle of underground engineering based on surface time-series deformation data, specifically: S21. First, by setting a deformation threshold, the time-series deformation result map is converted into a binary map; S22. Use an edge detection algorithm to extract texture from a binary image; S23. Use two-dimensional fast Fourier transform to convert the texture map into a spectrum map; S24. Perform Hough transform on the spectrogram and combine it with local peak detection to obtain the main direction of each time-series deformation map; S25. Average all the obtained main directions and combine them with the temporal deformation development trend to finally determine the azimuth of the underground space; S3: Underground space inversion algorithm based on IPIM model: The PIM model assumes that the underground space is rectangular in shape and further subdivides it into numerous tiny units. It assumes that these tiny units are independent and move independently, and obtains the surface deformation induced by the entire underground space through an accumulation method. The specific function expression is as follows: in: in, and These are surface subsidence and horizontal displacement, respectively. The height of the underground space; and These are the settlement coefficient and the horizontal displacement coefficient, respectively. , and These are the inclination angle of the underground space, the main influence angle, and the influence propagation angle; and Burial depth; , These are the main directional length and the length perpendicular to the main directional length of the underground space, respectively. Calculate the length in the main direction; , These are the inflection point offsets in the main direction, respectively. Since the length perpendicular to the main direction is small, there is no inflection point offset. This represents the horizontal distance from the center of the underground space to the calculation origin. This is the horizontal offset function between the unit at different locations and the center of its subsidence basin; The influence radius function for elements at different locations; When the underground space is horizontal, the angle of inclination is ignored; Deformation at any point on the Earth's surface is represented as: in, and These represent horizontal displacements in the east-west and north-south directions, respectively. and These represent the angles between the coal seam mining direction and the east and north directions, respectively. The specific steps of step S2 are as follows: Step 21: Set a deformation threshold for each scene's surface deformation result. The deformation threshold should be set to the regional average surface deformation or 0. Convert the deformation results into a binary map according to the following rules: In the formula, The pixel values ​​of the converted binary image. This represents the surface deformation value. The set deformation threshold; Step 22: Use the Sobel operator to extract texture from the binary image. The formulas for horizontal and vertical edge detection are as follows: In the formula, For horizontal edge detection images, Here, A is the image for vertical edge detection, and A is the image for gradient calculation. This is the result of gradient calculation, i.e., the texture of the binary image; Step 23: Use the two-dimensional fast Fourier transform (FFT) method to transform the image from the spatial domain to the frequency domain. The formula for the two-dimensional fast Fourier transform is as follows: In the formula, M and N are the row and column data of the image; Step 24: Surface deformation induced by underground engineering usually manifests as circular or elliptical shapes. Hough transform can quickly and effectively identify circular or elliptical structures contained in images. Therefore, Hough transform is used to detect texture maps and spectrograms, and combined with local peak detection algorithms, the deformation center and main deformation method are determined. Step 25: The principal directions of the images obtained in Step 4 are not the same. Therefore, it is necessary to take the average of the principal directions obtained in Step 4 and combine them with the development trend of temporal deformation to determine the final azimuth of the underground space. The formula is as follows: In the formula, The final azimuth angle. The main direction of each scene's image. For the number of images, It is determined based on the developmental trend of temporal deformation. , , or .

2. The positioning method for deep underground space based on InSAR technology as described in claim 1, characterized in that: The specific process of step 1 is as follows: S11. Collect SAR image data of the monitored area acquired by the satellite platform; S12. After performing high-precision registration, cropping, and filtering preprocessing on the acquired SAR image data, time-series D-InSAR technology is used to obtain differential interferograms of nearby time SAR images; high-precision processing requires at least one-thousandth of a pixel. S13. Select a highly coherent and stable point as a reference point to perform phase unwrapping operation and obtain the temporal deformation results of the land surface; highly coherent means that the coherence of urban areas is greater than 0.9 and the coherence of non-urban areas is greater than 0.

7. S14. The obtained deformation results are superimposed sequentially to obtain the cumulative surface deformation results.

3. The positioning method for deep underground space based on InSAR technology as described in claim 1, characterized in that: In step S3, an intelligent algorithm is introduced to obtain the global optimal solution. Based on the intelligent genetic particle swarm optimization algorithm, the implementation process of the underground space inversion algorithm is as follows: Step 1: Input geological parameters and main azimuth angle, and preset the search domain for underground space parameters, as well as the parameter settings for the intelligent genetic particle swarm optimization algorithm; Geological parameters , , , , , The principal azimuth angle finally obtained in step S2 In the known information input in the inversion algorithm, the geological parameters are determined based on geological conditions, field measurement data, and engineering construction methods. Underground space parameters are set through analogy. , , , , , , Upper and lower limits; among which, , The geographic coordinates of the underground space center; the parameters of the intelligent genetic particle swarm optimization algorithm include population size, learning factor, maximum and minimum inertia factor, crossover probability, mutation probability, maximum number of iterations, and termination condition; Step 2: Randomly generate population particles and calculate the three-dimensional coordinates of the Earth's surface based on IPIM; Based on the search domain of underground space parameters set in step one, the number of particles of the set population size is randomly generated, and the three-dimensional deformation value of the surface is calculated based on the IPIM model. The calculation formula is as follows: Step 3: Constructing a model of the relationship between three-dimensional surface deformation and line-of-sight deformation; The surface deformation information acquired by SAR satellites is line-of-sight deformation information. Step two simulates the three-dimensional deformation of the Earth's surface. , , The information, based on the SAR satellite's incident angle and heading angle, establishes the relationship between the two as follows: in, and Given the incident angle and heading angle of the SAR satellite, respectively, the functional relationship between the deep subsurface space and LOS deformation can be constructed by combining the IPIM model: ； Step 4: Calculate the population particle fitness function; Surface deformation information obtained based on InSAR technology And the calculation obtained in step three Construct the particle fitness function, as shown in the following formula: in, This represents the particle fitness value; For surface points predicted based on the IPIM model LOS deformation; Surface points obtained by InSAR technology inversion LOS deformation; n is the total number of pixels in the study area; Step 5: Population particle update; The population particles are sorted according to their fitness values. Half of the particles with poor fitness are discarded, while the particles with good fitness are regenerated into a new particle population using particle update algorithm and genetic algorithm respectively, so as to maintain the set population size. Then, repeat steps three through five until the termination condition is met or the maximum number of iterations is reached, thus obtaining the optimal parameters for the underground space.

Images