Purified micro-motion H / V curve and change rate calculation method and application

Abstract

The invention relates to a method for calculating a purified micro-motion H / V curve and a change rate. The calculation method comprises the following steps: step 1, collecting surface micro-vibration signals of a two-dimensional array; 2, respectively calculating a micro-motion H / V curve and a change rate of data acquired by each three-component detector for Rayleigh surface wave purification based on the random decrement technology; 3, drawing a two-dimensional H / V curve and change rate curve contour map or a three-dimensional slice map for the H / V curve and the change rate curve of each measuring point by using an MATLAB program; and 4, according to the abnormal characteristics of the contour map of the two-dimensional H / V curve and the change rate curve or the three-dimensional slice map, deducing and explaining an abnormal target body. The method is applied to the fields of urban geological survey, engineering investigation, underground space development and the like, and a rapid, accurate and environment-friendly solution is provided. The method has the advantages of convenient field acquisition, high detection precision, strong resolution capability, shortened exploration period, reduced engineering cost, reduced exploration risk and the like.

Description

Technical Field

[0001] The present invention belongs to the technical field of engineering geophysical exploration, and particularly relates to a method and application for calculating a micro-motion H / V curve and a rate of change of Rayleigh surface wave purification based on a random reduction technology using weak vibration signals on the ground. Background Art

[0002] There is a natural weak vibration on the surface of the earth at any time and anywhere, which comes from various activities of nature and human beings. Wind, tides, changes in air pressure, volcanic activities, etc. in nature will produce vibrations; while vibrations caused by human activities include vehicle movement, factory machinery operation, and even people walking, etc. The former has a frequency of less than 1Hz, and the latter has a frequency of greater than 1Hz.

[0003] All these vibration energies will propagate far away in the form of waves. Micromotion is a complex vibration composed of body waves (P waves and S waves) and surface waves (Rayleigh waves and Love waves), and the energy of surface waves in the vertical direction accounts for more than 2 / 3 of the total signal energy. Although the amplitude and shape of micromotion signals change with time and space, they are statistically stable within a certain range of time and space, and can be described by a stationary random process in time and space.

[0004] The current microtremor detection method is based on the theory of stationary random processes to extract the surface wave dispersion curve from the microtremor signal. The mainstream method is to obtain the shear wave velocity structure of the underground medium by joint inversion of the dispersion curve and H / V (for details, please refer to the following documents: Arai H, Tokimatsu K. 2005. S-wave velocity profiling by joint inversion of microtremor dispersion curve and Horizontal-to-Vertical (H / V) spectrum [J]. Bulletin of the Seismological Society of America, 95 (5): 1766-1778. D'Amico V, Picozzi M, Baliva F, Albarello D. 2008. Ambient noise measurements for preliminary site-effects characterization in the urban area of ​​Florence, Italy [J]. Bulletin of the Seismological Society of America, 98 (3): 1373-1388.).

[0005] The essence of the micro-vibration H / V curve is that the energy distribution ratio in the vertical and horizontal directions will change when the vibration propagates in different rock and soil layer structures. The morphological values ​​of the H / V curves observed in different soil layer structures are different, that is, the energy distribution ratio can directly reflect the structural changes of the underground rock and soil layers. This is also the basis for the H / V curve to be used to detect abnormal bodies in the underground rock and soil layers.

[0006] Domestic and foreign scholars have also applied the micromotion H / V curve to the field of engineering seismic. Wang Weijun et al. (2011) used the micromotion H / V to detect the near-surface structure in Baoding area. By calculating the change in the depth of the soil-rock interface, they inferred that there were two normal faults with obvious vertical dislocations in the Baoding fault. Liang Donghui et al. (2020) also used the micromotion H / V curve to detect underground river pipelines and caves. The micromotion H / V curve is also used to obtain the thickness of the sedimentary layer. For example, Zhang Ruohan et al. (2020) used data to obtain the resonance frequency and sedimentary layer thickness distribution in Guangzhou and parts of Foshan. It is also suitable for determining the thickness of the ice and snow cover, and has been successfully applied in the detection of the internal structure of glaciers and Antarctic ice sheets (Lévêque et al., 2010; Picotti et al., 2017; Yan et al., 2018).

[0007] Traditional micro-motion detection technology is mainly based on the theory of stationary random processes. It infers the shear wave velocity structure of the underground medium by analyzing the surface wave dispersion curve and the horizontal to vertical (H / V) spectrum ratio in the micro-motion signal. These technologies have played an important role in engineering geological exploration, seismic activity assessment and underground structure research. However, when processing the H / V spectrum ratio, the existing technology fails to effectively separate and eliminate the influence of body waves and Love waves, resulting in certain limitations in the identification of underground structural changes. The conventional H / V curve analysis method is to calculate the collected three-component data according to the following formula, without Rayleigh wave purification, and is the result of all waveform calculations:

[0008]

[0009] Where: is the conventional H / V value, H body , H Love , H Rayleigh is the horizontal component amplitude spectrum of body waves, Love waves, and Rayleigh waves, V body , V Rayleigh is the vertical component amplitude spectrum of body waves and Rayleigh waves. As can be seen from the above formula, the horizontal component and the vertical component contain the energy of all types of wave fields, but in fact, it is difficult to measure the energy proportion of body waves, Love waves, and Rayleigh waves, so the joint inversion of the micromotion dispersion curve and the H / V curve is not accurate. Summary of the invention

[0010] The object of the present invention is to provide a method for calculating a purified micro-motion H / V curve and a rate of change and an application thereof.

[0011] The object of the present invention is achieved by the following technical solution: a method for calculating a purified micro-motion H / V curve and a rate of change and an application thereof, which comprises the following steps:

[0012] Step 1: Arrange three-component geophones in a two-dimensional and small-spacing manner and collect surface microvibration signals for a recording time of no less than 5 minutes;

[0013] Step 2: Calculate the micro-motion H / V curve and change rate of the data collected by each three-component geophone based on the Rayleigh surface wave purification of the random reduction technology; including the following process:

[0014] Process (1) Preprocessing the collected raw data for baseline drift correction;

[0015] Process (2). Remove short-term strong interference signals in the near field and obtain relatively stable micro-motion signals;

[0016] Process (3). Using random reduction technology to purify the three-component signal, Rayleigh wave is obtained to obtain the H / V curve; it includes the following process:

[0017] Process (1). Window the signal;

[0018] Process (2). According to the set frequency range and sampling rate, the windowed data is filtered using a Chebyshev narrowband filter to obtain a narrowband signal V f (t);

[0019] Process (3). Find the sign change points of the vertical component of the narrowband signal, and use these sign change points as the starting points of the new time window of the three components of the narrowband signal. i The following conditions need to be met:

[0020] V f (τ i )≤0, and V f (τ i +dt)>0

[0021] Where V f is the narrowband vertical component signal after the Chebyshev filter, dt is the time sampling rate,

[0022] In this way, a series of new three-component data with a time window length of Δ is formed, and the formula is as follows:

[0023] v f,i (t) = v f (τi +t)

[0024]

[0025] Where, 0≤t≤Δ, f is the center frequency of the Chebyshev filter, v f,i ,e f,i ,n f,i are the vertical component, EW component, and NS component signals after Chebyshev filtering, i is the new time window number, t is the time series, For the time shift of the two horizontal components, the vertical and horizontal components of the Rayleigh wave have a 90° phase shift, so it is equivalent to a 90° phase shift;

[0026] Process (4). Project the two horizontal component data within the time window Δ to obtain the synthesized horizontal component. The formula is as follows:

[0027] h f,i (t) = sin(θ i )·e f,i (t)+cos(θ i )·n f,i (t)

[0028] Where θ i ∈[0,2π], and let the vertical component v f,i (t) and the horizontal component h f,i (t) has the largest correlation, and the calculation formula is as follows:

[0029]

[0030] The correlation coefficient c between the normalized vertical and horizontal components 2 f,i The formula is as follows:

[0031]

[0032] Process (5). In the time window Δ, the horizontal component h f,i (t) and the vertical component v f,i (t) Using the correlation coefficient c 2 f,i Weighted, all time windows Δ are summed to obtain a new horizontal component h f,w (t) and the vertical component v f,w (t), the formula is as follows:

[0033]

[0034] Process (6). Use the following formula to get the H / V value of the current frequency;

[0035]

[0036] Process (7). Perform process (2) to process (6) on all frequencies set in process (2) to obtain an H / V curve;

[0037] Process (8): geometrically average the H / V curves of all windowed data to obtain the final H / V curve;

[0038] Process (4). Perform differential calculation on the H / V curve to obtain the H / V change rate curve;

[0039] The process (1) of performing baseline drift correction on the raw data is as follows: process 1. dividing the raw data into sections of 2000 sampling points each; process 2. performing polynomial fitting on each section of data to obtain fitting results; process 3. subtracting the fitting results of each section from the raw data of the section to obtain the result data after baseline drift correction;

[0040] Step 3: Use MATLAB program to draw two-dimensional H / V curve and change rate curve contour map or three-dimensional slice map of H / V curve and change rate curve of each measuring point;

[0041] Step 4: Infer and interpret the abnormal target body based on the abnormal characteristics of the two-dimensional H / V curve and the change rate curve contour map or the three-dimensional slice map.

[0042] The application of the calculation method based on the micro-motion H / V change rate of the present invention provides a fast, accurate and environmentally friendly solution in the fields of urban geological survey, engineering investigation and underground space development.

[0043] Compared with the prior art, the advantages of the present invention are:

[0044] 1. The present invention uses vibrations generated by human activities as the seismic source. The urban environment with dense human activities will not only not cause interference, but on the contrary provide a rich source of high-frequency signals for shallow micro-seismic exploration. Micro-seismic exploration does not require artificially excited seismic sources, will not affect the surrounding environment, has low requirements for surface conditions, and will not hinder traffic. It only needs to carry out flexible temporary traffic control during the test. After the test is completed, the control facilities can be immediately withdrawn;

[0045] 2. The random reduction technology is introduced and optimized to purify Rayleigh surface waves. The calculated H / V curve is closer to the ellipticity of the fundamental mode of Rayleigh surface waves, with clear physical meaning, and is suitable for joint inversion. At the same time, the change rate of the purified H / V curve is calculated to understand the law of the sequential change of formation wave impedance from large to small or from small to large, which can amplify the subtle signal anomalies in the soil layer and more effectively reflect the location and detailed changes of obstacles or poor geology in the soil layer;

[0046] 3. By arranging the detectors in a two-dimensional small-pitch manner, it is possible to draw three-dimensional equal-frequency slice diagrams in the YZ, XZ and YZ directions, which can intuitively present the three-dimensional distribution range of obstacles. BRIEF DESCRIPTION OF THE DRAWINGS

[0047] Figure 1 A diagram showing the steps of the calculation method of the present invention.

[0048] Figure 2 It is a flow chart of obtaining the H / V curve and the rate of change curve in step 2 of the present invention.

[0049] Figure 3 This is a schematic diagram of the detector arrangement of implementation case one.

[0050] Figure 4 This is the H / V change rate profile diagram of implementation case 1.

[0051] Figure 5 Parameters and layer diagrams of the traditional model and the model of the present invention.

[0052] Figure 6 This figure compares the H / V results of the traditional model and the model of the present invention. DETAILED DESCRIPTION

[0053] The present invention is described in detail below with reference to the accompanying drawings and implementation examples:

[0054] like Figure 1 As shown: A method for calculating a purified micro-motion H / V curve and a rate of change and its application, which includes the following steps:

[0055] Step 1: Arrange three-component geophones in a two-dimensional and small-spacing manner and collect surface microvibration signals for a recording time of no less than 5 minutes;

[0056] The two-dimensional layout can adopt regular rectangles with equal spacing or rectangles with unequal spacing, and the placement position of each detector should be conveniently recorded.

[0057] Specifically, the two-dimensional small spacing method in step 1 refers to using a three-component geophone spacing of 50 cm to 5 m according to the size of the detection target, ensuring that the three-component geophone covers the underground target, and arranging the three-component geophones in a two-dimensional plane for conveniently recording the placement positions of each three-component geophone, with a collection time of 5-15 min. The three-component geophone uses a three-component high-precision integrated digital seismograph.

[0058] Step 2: Calculate the micro-motion H / V curve and change rate of the data collected by each three-component geophone based on the random reduction technology of Rayleigh surface wave purification (such as Figure 2 It includes the following processes:

[0059] Process (1) Perform baseline drift correction preprocessing on the collected raw data, as follows:

[0060] 1) Divide the original data into sections of 2,000 sampling points each; 2) Perform polynomial fitting on each section of data to obtain fitting results; 3) Subtract the fitting results of each section from the original data of that section to obtain the result data after baseline drift correction. This process can improve the consistency of instrument response and improve the accuracy and stability of collected data.

[0061] Process (2). Remove short-term strong interference signals in the near field and obtain relatively stable micro-motion signals;

[0062] Process (3). Using random reduction technology to purify the three-component signal, Rayleigh wave, and obtain the H / V curve; it includes the following process:

[0063] Process (1). Window the signal, as follows:

[0064] The time length T of the window is set to an integer multiple of the period (the inverse of the frequency). Therefore, the time lengths of the windows of different frequencies are different. For signals of the same time length, the number of windows of different frequencies is also different:

[0065]

[0066] In the formula, n is an integer (usually 10 to 30 times), and f is the frequency.

[0067] Process (2). According to the set frequency range and sampling rate, the windowed data is filtered using a Chebyshev narrowband filter to obtain a narrowband signal V f (t); It is recommended to set the bandwidth to 0.2 times the center frequency.

[0068] Process (3). Find the sign change points of the vertical component of the narrowband signal, and use these sign change points as the starting points of the new time window of the three components of the narrowband signal. These sign change points τ i The following conditions need to be met:

[0069] V f (τ i )≤0, and V f (τ i +dt)>0

[0070] Where V f is the narrowband vertical component signal after the Chebyshev filter, dt is the time sampling rate,

[0071] In this way, a series of three-component data with a new time window length of Δ is formed, and the formula is as follows. It should be noted that the recommended value of the new time window length is: Δ=10 / f.

[0072] v f,i (t) = v f (τ i +t)

[0073]

[0074] Where, 0≤t≤Δ, f is the center frequency of the Chebyshev filter, v f,i ,e f,i ,n f,i They are the vertical component, EW component, and NS component signals after Chebyshev filtering, i is the new time window number, t is the time series, and the time shift of the two horizontal components This is equivalent to a 90° phase shift because the vertical and horizontal components of the Rayleigh wave are 90° phase shifted.

[0075] Process (4). Project the two horizontal component data within the time window Δ to obtain the synthesized horizontal component. The formula is as follows:

[0076] h f,i (t) = sin(θ i )·e f,i (t)+cos(θ i )·n f,i (t)

[0077] Where θ i ∈[0,2π], and let the vertical component v f,i (t) and the horizontal component h f,i (t) has the largest correlation, and the calculation formula is as follows:

[0078]

[0079] The correlation coefficient c between the normalized vertical and horizontal components 2 f,i The formula is as follows:

[0080]

[0081] Process (5). In the time window Δ, the horizontal component h f,i (t) and the vertical component v f,i (t) Using the correlation coefficient c 2 f,i Weighted, all time windows Δ are summed to obtain a new horizontal component h f,w (t) and the vertical component v f,w (t), the formula is as follows:

[0082]

[0083] Process (6). Use the following formula to get the H / V value of the current frequency;

[0084]

[0085] Process (7). Perform process (2) to process (6) on all frequencies set in process (2) to obtain an H / V curve;

[0086] Process (8): geometrically average the H / V curves of all windowed data to obtain the final H / V curve;

[0087] Here, in process (3), the signal is purified by Rayleigh wave, and the random reduction technology is used to suppress the interference of body waves and other surface waves, enhance the Rayleigh wave signal, suppress the influence of body waves and Love waves, etc., so that the obtained H / V curve is more similar to the ellipticity curve of Rayleigh wave.

[0088] Process (4). Perform differential calculation on the H / V curve to obtain the H / V change rate curve;

[0089] In process (4), the MATLAB program is used to convert the H / V curve of each measuring point into a change rate curve; when drawing a two-dimensional H / V change rate curve contour map or a three-dimensional slice map, the H / V curve value is selected for differential calculation. The H / V change rate map highlights the abnormal range and makes the contour boundary of the abnormality clearer, which helps to accurately analyze the location of the obstacle.

[0090] Step 3: Use MATLAB program to draw two-dimensional H / V curve and change rate curve contour map or three-dimensional slice map of H / V curve and change rate curve of each measuring point;

[0091] Step 4: Infer and interpret the abnormal target body based on the abnormal characteristics of the two-dimensional H / V curve and the change rate curve contour map or the three-dimensional slice map.

[0092] The present invention introduces and optimizes the random reduction technology, which was originally used in the field of earthquake engineering to identify the dynamic characteristics of the structure. By analyzing the response of the structure under random excitation, modal parameters such as natural frequency, damping ratio and vibration mode are identified. In the seismic design of buildings and earthquake damage assessment, these parameters are crucial to understanding the dynamic behavior of the structure. The present invention applies the random reduction technology to the micro-motion signal processing, based on the fact that the vertical signal Rayleigh wave energy is dominant and the horizontal two directions have a 900 phase difference with the vertical direction, the three-component signal is Rayleigh wave purified to suppress the influence of body waves and Love waves. The H / V calculation formula after purification is analogous to the Rayleigh wave ellipticity, and the obtained H / V is closer to the Rayleigh wave fundamental mode ellipticity, which can improve the accuracy of subsequent inversion. The purified and optimized H / V formula is as follows:

[0093]

[0094] The shape and amplitude of the H / V curve can reflect the change in the wave impedance of the underground strata. Abnormal changes in the H / V curve can indicate the existence of underground obstacles or bad geological bodies. When the underground medium contains isolated rocks or weak layers, the H / V curve will have a significant abnormal peak.

[0095] The purified H / V curve and the rate of change of the H / V curve are more sensitive to shallow structures, especially transversely inhomogeneous anomalies, such as boulders in the overburden, soft and hard interlayers in caves, or if there are obstacles such as pile foundations and box culverts, the horizontal and vertical energy will change during vibration propagation, and the H / V curve and the rate of change curve will produce obvious abnormal changes in the frequency band corresponding to the burial depth. It is precisely based on these characteristics that the present invention uses the H / V curve and the rate of change curve to amplify the subtle signal anomalies in the soil layer, thereby achieving the purpose of distinguishing obstacles. The actual measured cases have proved that this method has a good effect in underground obstacle and poor geological detection. At the same time, the purified H / V curve is closer to the fundamental mode ellipticity of the theoretical Rayleigh surface wave, with a clear physical meaning and a higher accuracy of the combined reflection.

[0096] The method of the present invention can analyze strata (soil-rock interface) with large changes in stratum wave impedance ratio within the shallow surface or cover layer depth range in urban engineering construction, or underground obstacles or unfavorable geological bodies that have an impact on engineering construction.

[0097] Calculation by building a model

[0098] The effectiveness of the present invention is illustrated below by calculating the model data. Figure 5 The two-layer model is used to obtain the simulated two-dimensional micromotion data (see Table 1). The H / V calculation results after the traditional and purified Rayleigh waves are compared (see Table 1). Figure 6 As shown in the figure, the frequency is uniformly cut to 20 Hz), it can be seen that the purified H / V value is closer to the simulation value, which can improve the accuracy of subsequent inversion results.

[0099] Table 1

[0100] Thickness(m) Vp(m / s) Vs(m / s) Density(kg / m3) Poisson 1 25.00 500.00 200.00 1900.00 0.4048 2 0.00 2000.00 1000.00 2500.00 0.3333

[0101] The present invention is described below in conjunction with implementation cases.

[0102] Example 1: Detection of H / V change rate at a construction site

[0103] Project Overview: The investigation data of the target site show that there is an uneven distribution of rounded gravel and crushed pebbles underground, as well as interlayers of crushed pebbles, rounded gravel, silty clay and silty soil. The bearing layer of the pile foundation designed for the target site is a crushed gravel layer, but the drilling information cannot effectively obtain the plane distribution of crushed pebbles and rounded gravel, and the plane distribution of the thickness of the soft soil layer in the interlayer cannot be determined. Therefore, the H / V change rate method is used for detection and identification. The nearby detailed exploration drilling holes revealed that the main soil layers from top to bottom are miscellaneous fill (stone filling), silty clay, silt (quality soil), silty clay, silty soil, rounded gravel, silty clay, (muddy) medium sand, crushed pebbles, silty clay, rounded gravel, (muddy) medium sand, crushed pebbles, silty soil, crushed pebbles, rounded gravel, silty clay, (mud-containing) medium sand, silty soil, sandy soil strongly weathered granite (syenite porphyry), and fragmented blocky strongly differentiated granite (syenite porphyry).

[0104] Step 1: This H / V rate of change method detection uses a three-component high-precision integrated digital seismograph, with a two-dimensional detector array arranged at a spacing of 0.5 to 1.5 m (see Figure 3 ), collect surface microvibration signals of 10 to 15 m;

[0105] Step 2: Calculate the H / V curve of the data collected by each three-component detector;

[0106] Step 3: Use the self-developed MATLAB program to convert the H / V change rate curve of each measuring point into a two-dimensional H / V change rate curve contour map (see Figure 4 );

[0107] Step 4: Detection results analysis and interpretation:

[0108] like Figure 3 As shown in the figure: the horizontal axis is distance; the vertical axis is frequency (Hz), high frequency reflects shallow information, low frequency reflects deep information; different grayscales represent the magnitude of the H / V change rate. From the profile, it can be seen that the two-dimensional distribution of crushed pebbles and rounded gravels is basically consistent with the existing drilling holes. It can be seen that the interlayer two-dimensional distribution of crushed pebbles, rounded gravels and silty soil and silty clay is unevenly distributed at a depth of about -33 to -41m, and there are obvious differences in weak interlayers at different positions. The detection results can guide the design and construction of pile foundations.

[0109] like Figure 4 As shown: the horizontal axis is frequency (Hz), the vertical axis is the value of H / V, and the H / V calculation results of the traditional and purified Rayleigh waves are compared.

[0110] Finally, it should be noted that the above is only a preferred embodiment of the present invention and is not intended to limit the present invention. Although the present invention has been described in detail with reference to the aforementioned embodiments, it is still possible for those skilled in the art to modify the technical solutions described in the aforementioned embodiments or to make equivalent substitutions for some of the technical features therein. Any modifications, equivalent substitutions, improvements, etc. made within the spirit and principles of the present invention should be included in the protection scope of the present invention.

Claims

1. A method for calculating a purified micro-motion H / V curve and a rate of change, characterized in that: It includes the following steps: Step 1: Arrange three-component geophones in a two-dimensional and small-spacing manner and collect surface microvibration signals for a recording time of no less than 5 minutes; Step 2: Calculate the micro-motion H / V curve and change rate of the data collected by each three-component geophone based on the Rayleigh surface wave purification of the random reduction technology; The following processes are included: Process (1) Preprocessing the collected raw data for baseline drift correction; Process (2). Remove short-term strong interference signals in the near field and obtain relatively stable micro-motion signals; Process (3). Using random reduction technology to purify the three-component signal, Rayleigh wave is obtained to obtain the H / V curve; it includes the following process: Process (1). Window the signal; Process (2). According to the set frequency range and sampling rate, the windowed data is filtered using a Chebyshev narrowband filter to obtain a narrowband signal V f (t); Process (3). Find the sign change points of the vertical component values ​​of the narrowband signal, and use these sign change points as the starting points of the new time window of the three components of the narrowband signal. These sign change points τ i The following conditions need to be met: V f (τ i ) ≤ 0, and V f (τ i + dt) > 0 Where V f is the narrowband vertical component signal after the Chebyshev filter, dt is the time sampling rate, In this way, a series of new three-component data with a time window length of Δ is formed, and the formula is as follows: v f,i (t)=v f (τ i +t) Where, 0≤t≤Δ, f is the center frequency of the Chebyshev filter, v f,i ,e f,i ,n f,i are the vertical component, EW component, and NS component signals after Chebyshev filtering, i is the new time window number, t is the time series, For the time shift of the two horizontal components, the vertical and horizontal components of the Rayleigh wave have a 90° phase shift, so it is equivalent to a 90° phase shift; Process (4). Project the two horizontal component data within the time window Δ to obtain the synthesized horizontal component. The formula is as follows: h f,i (t)=sin(θ i )·e f,i (t)+cos(θ i )·n f,i (t) Where θ i ∈[0,2π], and let the vertical component v f,i (t) and the horizontal component h f,i (t) has the largest correlation, and the calculation formula is as follows: The correlation coefficient c between the normalized vertical and horizontal components 2 f,i The formula is as follows: Process (5). In the time window Δ, the horizontal component h f,i (t) and the vertical component v f,i (t) Using the correlation coefficient c 2 f,i Weighted, all time windows Δ are summed to obtain a new horizontal component h f,w (t) and the vertical component v f,w (t), the formula is as follows: Process (6). Use the following formula to get the H / V value of the current frequency; Process (7). Perform process (2) to process (6) on all frequencies set in process (2) to obtain an H / V curve; Process (8): geometrically average the H / V curves of all windowed data to obtain the final H / V curve; Process (4). Perform differential calculation on the H / V curve to obtain the H / V change rate curve; Step 3: Use MATLAB program to draw two-dimensional H / V curve and change rate curve contour map or three-dimensional slice map of H / V curve and change rate curve of each measuring point; Step 4: Infer and interpret the abnormal target body based on the abnormal characteristics of the two-dimensional H / V curve and the change rate curve contour map or the three-dimensional slice map.

2. The method for calculating a purified micro-motion H / V curve and a rate of change according to claim 1, characterized in that: The two-dimensional small spacing method in step 1 refers to using a three-component detector spacing of 50 cm to 5 meters according to the size of the detection target, ensuring that the three-component detector covers the underground target, and arranging the three-component detectors in a two-dimensional plane for convenient recording of the placement position of each three-component detector. The acquisition time is 5-15 minutes.

3. The method for calculating the purified micro-motion H / V curve and the rate of change according to claim 1, characterized in that: Process (1) The process of performing baseline drift correction on the original data is as follows: Process 1. Divide the original data into sections of 2000 sampling points each; Process 2. Perform polynomial fitting on each section of data to obtain fitting results; Process 3. Subtract the fitting result of each section from the original data of the section to obtain the result data after baseline drift correction.

4. The method for calculating the purified micro-motion H / V curve and the rate of change according to claim 1, characterized in that: In process (4), the MATLAB program is used to convert the H / V curve of each measuring point into a change rate curve; when drawing a two-dimensional H / V change rate curve contour map or a three-dimensional slice map, the H / V curve value is selected for differential calculation. The H / V change rate map highlights the abnormal range and makes the contour boundary of the abnormality clearer, which helps to accurately analyze the location of the obstacle.

5. The method for calculating the purified micro-motion H / V curve and the rate of change according to claim 1, characterized in that: The three-component geophone adopts a three-component high-precision integrated digital seismograph.

6. Application of a method for calculating a purified micro-motion H / V curve and a rate of change according to any one of claims 1 to 5, characterized in that: It is used in urban geological surveys, engineering investigations and underground space development. It is used to analyze strata with large changes in stratum wave impedance ratio within the shallow surface or cover layer depth range in urban engineering construction, or underground obstacles or unfavorable geological bodies that have an impact on engineering construction.

Images