Method for detecting and analyzing rock mass damage of mine blasting surrounding rock based on ultrasonic detection

Abstract

This invention discloses a method for detecting and analyzing rock mass damage in mine blasting based on ultrasonic testing, belonging to the field of mine blasting surrounding rock damage detection technology. The method involves deploying an ultrasonic sensor array within a monitoring borehole preserving the surrounding rock, acquiring background noise signals before blasting, calculating the background cross-correlation function between each sensor pair, and extracting the background surface wave dispersion curve. After triggering the blast, the method continuously acquires the blasting response signal, calculates the response cross-correlation function, and extracts the response surface wave dispersion curve. The difference between the two dispersion curves is processed to obtain the dispersion curve offset, and the change in the shear wave velocity profile of the surrounding rock around the monitoring borehole is inverted to determine the depth boundary and spatial distribution range of the damage zone within the surrounding rock. By analyzing the difference in dispersion characteristics before and after blasting, the wave velocity change caused by the surrounding rock damage is obtained, and the distribution of the damage zone is determined.

Description

Technical Field

[0001] This invention relates to the field of damage detection technology for surrounding rock in mine blasting, specifically a method for detecting and analyzing damage to surrounding rock mass in mine blasting based on ultrasonic testing. Background Technology

[0002] Damage detection in surrounding rock masses during mine blasting typically relies on methods such as core sampling, borehole video recording, sonic logging, cross-hole sonic tomography, and microseismic monitoring. Core sampling and borehole video recording can directly obtain information about the rock mass structure near the borehole wall, but they mainly reflect local characteristics of the borehole surface or core samples, making it difficult to provide continuous information on the changes in the mechanical properties of the rock mass with depth within a certain radial range around the borehole. Sonic logging and cross-hole sonic tomography require active seismic sources or transmitting transducers, placing high demands on borehole conditions and on-site construction coordination. Repeated observations before and after blasting in boreholes that retain surrounding rock are difficult to implement on-site, and the results are mostly single wave velocity values ​​or local profiles, insufficiently describing the continuity of the spatial boundaries of the damaged area. Microseismic monitoring can record the location of rupture events and the released energy, but it is not sensitive to changes in the wave velocity field within the already formed damaged area and cannot directly characterize the amplitude and distribution range of shear wave velocity changes caused by surrounding rock damage. When performing non-destructive and repeatable wave velocity change detection in the same surrounding rock area before and after blasting, conventional methods are affected by changes in observation conditions, coupling state, or excitation mode. Single wave velocity test results are difficult to directly compare, and wave velocity changes caused by surrounding rock damage are easily mixed with system biases. Background noise cross-correlation and surface wave dispersion extraction techniques can recover medium wave velocity information using environmental vibrations. However, in mine blasting surrounding rock damage detection, the difference in surface wave dispersion in the background noise before blasting and the response signal after blasting has not yet been systematically used for damage identification. There is a lack of methods to simultaneously acquire both types of signals using an ultrasonic sensor array within the borehole containing the surrounding rock monitoring data and extract the dispersion curve offset. Furthermore, there is a lack of systematic methods for inverting the change in shear wave velocity profile of the surrounding rock around the borehole from the dispersion curve offset, and determining the depth boundary and spatial distribution range of the damaged area accordingly. This results in the determination of damage boundaries within the surrounding rock after blasting relying heavily on experience. Summary of the Invention

[0003] To address the challenges of continuously acquiring changes in shear wave velocity within the surrounding rock before and after blasting, and the influence of subjective experience on determining the depth boundary and spatial extent of the damaged zone in mine blasting damage detection, this paper proposes an ultrasonic-based method for detecting and analyzing damage in the surrounding rock mass during mine blasting. This method involves deploying an ultrasonic sensor array within a monitoring borehole that preserves the surrounding rock. Background noise signals before blasting and blasting response signals after blasting are collected. The cross-correlation function between each sensor pair is calculated to extract the background surface wave dispersion curve and the response surface wave dispersion curve. The difference between the two dispersion curves is processed to obtain the dispersion curve offset. The change in the shear wave velocity profile of the surrounding rock around the monitoring borehole is then retrieved, determining the depth boundary and spatial extent of the damaged zone within the surrounding rock.

[0004] To achieve the above objectives, the present invention provides the following technical solution: The present invention provides a method for detecting and analyzing rock mass damage in mine blasting based on ultrasonic detection, comprising: deploying an ultrasonic sensor array in a monitoring borehole with preserved surrounding rock, and collecting background noise signals received by the ultrasonic sensor array before blasting; calculating the background cross-correlation function between each sensor pair based on the background noise signal, and extracting the background surface wave dispersion curve from the background cross-correlation function; triggering blasting and continuously collecting blasting response signals received by the ultrasonic sensor array; calculating the response cross-correlation function between each sensor pair based on the blasting response signal, and extracting the response surface wave dispersion curve from the response cross-correlation function; performing difference processing on the background surface wave dispersion curve and the response surface wave dispersion curve to obtain the dispersion curve offset; inverting the change in shear wave velocity profile of the surrounding rock around the monitoring borehole based on the dispersion curve offset; and determining the depth boundary and spatial distribution range of the damage zone inside the surrounding rock based on the change in shear wave velocity profile. This method utilizes the difference in surface wave dispersion characteristics between background noise and blasting response to highlight the changes in rock wave velocity caused by blasting, thereby enabling quantitative identification and spatial localization of surrounding rock damage.

[0005] As a technical solution of this invention, when deploying an ultrasonic sensor array within a monitoring borehole that preserves the surrounding rock and collects background noise signals, multiple ultrasonic sensors are arranged at preset equal intervals along the borehole axis and fixedly coupled using a grouting coupling medium. During a quiet period before blasting and without external vibration sources, each ultrasonic sensor is activated to synchronously collect environmental vibration signals, obtaining a time-series signal of a preset duration as the background noise signal. The equal-interval deployment and grouting coupling ensure consistent coupling between each sensor and the borehole wall, improving the spatial comparability of the background noise signals.

[0006] As a technical solution of this invention, when extracting the background surface wave dispersion curve, the background noise signal is divided into multiple background noise segments according to a preset time window length. Each background noise segment undergoes spectral whitening processing to obtain a whitened background noise segment. Cross-correlation is performed on the whitened background noise segments corresponding to two sensors within the same sensor pair to obtain multiple cross-correlation waveforms. These multiple cross-correlation waveforms are then superimposed and averaged to obtain a background cross-correlation function. Time-frequency analysis is performed on the background cross-correlation function to extract the curve showing the change in phase velocity with frequency as the background surface wave dispersion curve. Through segmented cross-correlation and superposition averaging, random noise and non-steady-state interference can be suppressed, improving the stability of the background dispersion curve.

[0007] As a technical solution of this invention, when triggering a blast and continuously acquiring the blast response signal received by the ultrasonic sensor array, the blast initiation time is taken as the trigger zero point. Each ultrasonic sensor is activated to synchronously acquire the blast vibration signal according to a preset sampling frequency and preset acquisition duration. The blast vibration signals acquired by each ultrasonic sensor within the preset acquisition duration are arranged in chronological order to obtain the blast response time sequence as the blast response signal. The blast response signal is divided into multiple blast response segments according to a preset time window length. Each blast response segment undergoes spectral whitening processing to obtain a whitened blast response segment. Cross-correlation calculations are performed on the whitened blast response segments corresponding to two sensors within the same sensor pair to obtain multiple blast cross-correlation waveforms. The multiple blast cross-correlation waveforms are superimposed and averaged to obtain the response cross-correlation function. Time-frequency analysis is performed on the response cross-correlation function, and the curve of phase velocity changing with frequency is extracted as the response surface wave dispersion curve. The background noise and the blast response use the same time-frequency processing parameters, ensuring that the two dispersion curves are comparable at the same frequency coordinate.

[0008] As a technical solution of this invention, when obtaining the dispersion curve offset through differential processing, under the same frequency coordinates, the phase velocity values ​​in the response surface wave dispersion curve and the phase velocity values ​​in the background surface wave dispersion curve are subtracted from each other at each frequency point to obtain the phase velocity difference at each frequency point. These phase velocity differences are then arranged in ascending order of frequency to form the dispersion curve offset. Differential processing can eliminate systemic biases caused by differences in initial geological conditions of the monitoring borehole and sensor coupling, highlighting the changes in dispersion characteristics caused by blasting damage.

[0009] As one technical solution of this invention, when inverting the change in shear wave velocity profile, a radial layered model of the surrounding rock around the borehole is established, and an initial shear wave velocity is assigned to each layer in the radial layered model. Using the dispersion curve offset as the target data, based on ray theory, the shear wave velocity of each layer in the radial layered model is iteratively adjusted until the residual between the theoretical dispersion offset calculated in the forward model and the dispersion curve offset is less than a preset threshold. The difference between the shear wave velocity of each layer and the corresponding initial shear wave velocity at the end of the iteration is output as the change in shear wave velocity profile. As a further implementation scheme, the disturbance of the shear wave velocity of each layer relative to the initial shear wave velocity is used as the parameter to be inverted. A Jacobian matrix between the disturbance and the dispersion curve offset is established. The least squares algorithm is used to calculate the adjustment step size of the parameter to be inverted and update the shear wave velocity of each layer. The updated shear wave velocity of each layer is substituted into the surface wave dispersion forward modeling equation to calculate the theoretical dispersion offset. The residual between the theoretical dispersion offset and the dispersion curve offset is repeatedly calculated until the residual is less than a preset threshold, at which point the iteration is terminated. By using ray-theoretic forward modeling and least-squares iterative inversion, the spatial resolution and inversion stability of the shear wave velocity profile variation can be improved.

[0010] As a technical solution of this invention, when determining the depth boundary and spatial distribution range of the damage zone within the surrounding rock, the upper and lower boundaries of continuous segments in the shear wave velocity profile variation exceeding a preset threshold are marked as the depth boundary of the damage zone along the depth direction of the monitoring borehole. The depth boundaries of multiple monitoring boreholes located at different spatial positions are obtained, and these depth boundaries are fitted onto a surface in the same three-dimensional geological model to obtain the spatial distribution range of the damage zone. This extends the damage assessment from a single borehole to a three-dimensional description of the damage range over a spatial region, improving the engineering applicability of damage detection.

[0011] As a technical solution of this invention, based on the acquisition of the blasting response signal, the blasting response signal is further processed by first-arrival wave pickup and amplitude normalization. Specifically, the ratio of the short-term average value to the long-term average value of each sampling point in the blasting response signal is calculated, and the moment when the ratio first exceeds a preset trigger threshold is marked as the arrival time of the first-arrival wave. Based on the arrival time of the first-arrival wave, an effective signal segment is extracted from the blasting response signal, the maximum absolute value of the effective signal segment is calculated, and the amplitude of each sampling point within the effective signal segment is divided by the maximum absolute value to obtain the amplitude-normalized signal. Through first-arrival wave pickup and amplitude normalization, the influence of differences in blasting energy and propagation path on the surface wave dispersion curve of the response can be reduced, improving the comparability of the background and response signals and the reliability of the dispersion curve offset.

[0012] The technical effects and advantages provided by the present invention in the above technical solution are as follows: An array of ultrasonic sensors was deployed in the monitoring borehole while preserving the surrounding rock to collect background noise signals before blasting and response signals after blasting. Cross-correlation calculations were performed on the records of each sensor pair, and surface wave dispersion curves were extracted. Dispersion information was obtained using the same array and the same frequency coordinates before and after blasting. The dispersion curve offset obtained by difference processing can weaken the relatively fixed influence of sensor coupling conditions and local borehole structure, highlighting the changes in surface wave phase velocity caused by blasting disturbance. The background cross-correlation function and the response cross-correlation function were obtained by fragment whitening, cross-correlation superposition averaging, and time-frequency analysis. The dispersion curve offset was obtained by subtracting frequency points at the same frequency coordinates. Compared with damage identification based solely on wave velocity testing at a single moment after blasting or borehole wall images, this processing method obtains a continuous dispersion offset distribution along the borehole depth direction, which is more sensitive to wave velocity changes caused by surrounding rock damage. This allows damage identification to extend beyond the borehole wall surface to the rock mass within a certain radial range around the borehole.

[0013] Using the dispersion curve offset as the target data for inversion, and combining the radial stratification model of the surrounding rock around the borehole with ray theory, the perturbation of the shear wave velocity in each layer of the model relative to the initial shear wave velocity is iteratively adjusted. By establishing the Jacobian matrix between the perturbation and the dispersion curve offset, the adjustment step size is calculated using the least squares algorithm. The updated shear wave velocities of each layer are then substituted into the surface wave dispersion forward modeling equation to obtain the theoretical dispersion offset, until the residual between the theoretical dispersion offset and the dispersion curve offset is less than a preset threshold. This inversion process utilizes the sensitivity of different frequency surface waves to different depth ranges, transforming the dispersion offset into a quantitative distribution of wave velocity changes in the surrounding rock mass. Boundary marking is performed on continuous segments along the borehole depth where the rate of change in the shear wave velocity profile exceeds a preset threshold, thus obtaining the depth boundary of the damaged zone. The depth boundaries of multiple monitoring boreholes are fitted to a three-dimensional geological model surface to form the spatial distribution range of the surrounding rock damage zone. Compared to methods such as core drilling and borehole television imaging that rely on local or surface information, the distribution of the damage zone obtained is quantitatively constrained by the wave velocity change field, which can reflect the continuous morphology of the damaged body inside the surrounding rock after blasting and reduce the influence of human experience judgment. Attached Figure Description

[0014] To more clearly illustrate the technical solutions in the embodiments of this application or the prior art, the drawings used in the embodiments will be briefly introduced below. Obviously, the drawings described below are only some embodiments recorded in this invention. For those skilled in the art, other drawings can be obtained based on these drawings.

[0015] Figure 1 This is a flowchart of a method for detecting and analyzing damage to surrounding rock mass in mine blasting based on ultrasonic testing. Figure 2This is a flowchart for extracting the surface wave dispersion curve of blasting vibration signal response; Figure 3 This is a flowchart for calculating the offset of the surface wave dispersion curve; Figure 4 It is a profile and variation diagram of radial transverse wave velocity of the surrounding rock. Detailed Implementation

[0016] To make the objectives, technical solutions, and advantages of the embodiments of the present invention clearer, 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, 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.

[0017] See Figure 1 This invention provides a method for detecting and analyzing rock mass damage in mine blasting based on ultrasonic detection, comprising: deploying an ultrasonic sensor array in a monitoring borehole containing the remaining surrounding rock, and collecting background noise signals received by the ultrasonic sensor array before blasting; calculating the background cross-correlation function between each sensor pair based on the background noise signal, and extracting the background surface wave dispersion curve from the background cross-correlation function; continuously collecting blasting response signals received by the ultrasonic sensor array after blasting is triggered; calculating the response cross-correlation function between each sensor pair based on the blasting response signal, and extracting the response surface wave dispersion curve from the response cross-correlation function; performing difference processing on the background surface wave dispersion curve and the response surface wave dispersion curve to obtain the dispersion curve offset; inverting the change in shear wave velocity profile of the surrounding rock around the monitoring borehole based on the dispersion curve offset; and determining the depth boundary and spatial distribution range of the damage zone inside the surrounding rock based on the change in shear wave velocity profile.

[0018] In practice, multiple ultrasonic sensors from the ultrasonic sensor array are sequentially arranged at predetermined equal intervals on the inner wall of the monitoring borehole along its axial direction. The predetermined equal interval is set to 1 meter, determined based on the depth resolution requirements of the axial damage zone in the monitoring borehole and the total length of the ultrasonic sensor array. The monitoring borehole is a borehole preserved within the surrounding rock, and its diameter is determined based on the outer diameter of the ultrasonic sensors and the thickness of the grouting coupling medium. Each ultrasonic sensor is connected to a ground acquisition station via a cable. Cement grout is used as the grouting coupling medium, injected from the bottom of the monitoring borehole towards the borehole opening. This grout fills the gaps between the ultrasonic sensors and the inner wall of the monitoring borehole. After the cement grout solidifies, each ultrasonic sensor is fixedly coupled to the inner wall of the monitoring borehole.

[0019] In practice, during the quiet period before blasting and without external vibration sources, each ultrasonic sensor synchronously acquires environmental vibration signals to obtain a time-series signal of a preset duration, which serves as the background noise signal. Each ultrasonic sensor begins sampling simultaneously upon receiving a synchronous trigger command from the acquisition station. The preset sampling frequency is set to 10kHz, determined based on the upper limit of the ultrasonic sensor's frequency response range of 5kHz and the Nyquist sampling theorem. The preset duration is set to 30 minutes, determined based on the time required for coherence accumulation in the surface wave frequency band of the background noise signal. During the acquisition process, each ultrasonic sensor continuously records the environmental vibration signal, forming a background noise time series corresponding to each ultrasonic sensor. The background noise time series corresponding to all ultrasonic sensors constitute the background noise signal.

[0020] In the specific implementation, the background noise signal is divided into multiple background noise segments according to a preset time window length. The preset time window length is set to 60s, determined based on the inverse relationship between frequency resolution and time window length in surface wave dispersion analysis. During the segmentation process, the overlap length between adjacent background noise segments is set to half the preset time window length to improve the average number and stability of the background cross-correlation function. Spectral whitening processing is performed on each background noise segment to obtain whitened background noise segments. The spectral whitening processing is implemented as follows: a Fast Fourier Transform is performed on the background noise segment to obtain its spectrum; the amplitude of each frequency component in the background noise segment spectrum is normalized, while retaining the phase information of each frequency component to obtain the whitened background noise segment spectrum; an Inverse Fast Fourier Transform is performed on the whitened background noise segment spectrum to obtain the whitened background noise segment. The spectral whitening processing uses the following formula: in, Indicates the first Discrete frequency points, Indicates discrete frequency index. Take a positive integer; This indicates the frequency of the background noise segment after Fourier transform. Complex spectral values ​​at; Represents the complex spectrum value The model; This indicates the background noise segment at a frequency after spectral whitening. The complex spectral values ​​at the specified location. Spectral whitening is performed on each background noise segment to obtain the corresponding whitened background noise segment.

[0021] In practice, cross-correlation is performed on the whitened background noise segments corresponding to the two sensors within the same sensor pair to obtain multiple cross-correlation waveforms. The cross-correlation operation is as follows: for the whitened background noise segments corresponding to the first and second sensors within the same sensor pair, the whitened background noise segment corresponding to the first sensor is fixed, while the whitened background noise segment corresponding to the second sensor is time-shifted point by point. At each time-shifted position, the values ​​of the two whitened background noise segments at the overlapping sampling points are multiplied point by point and then accumulated to obtain the cross-correlation value at the corresponding time-shifted position. The cross-correlation values ​​corresponding to all time-shifted positions constitute the cross-correlation waveform. The cross-correlation operation is repeated on different background noise segments for the two sensors within the same sensor pair to obtain multiple cross-correlation waveforms. These multiple cross-correlation waveforms are then superimposed and averaged to obtain the background cross-correlation function. The superposition and averaging method is as follows: the cross-correlation values ​​of the multiple cross-correlation waveforms at the same time-shifted position are added, and the sum is divided by the total number of cross-correlation waveforms to obtain the value of the background cross-correlation function at the corresponding time-shifted position.

[0022] In practice, time-frequency analysis is performed on the background cross-correlation function to extract the curve of phase velocity versus frequency, which serves as the background surface wave dispersion curve. The time-frequency analysis employs the phase-shifting method. A Fourier transform is performed on the background cross-correlation function to obtain the background cross-correlation spectrum. Within this spectrum, for each discrete frequency component, the phase difference between the two sensors within the same sensor pair is extracted. Based on the phase difference, the sensor spacing between the two sensors within the same sensor pair, and the angular frequency corresponding to the discrete frequency component, the phase velocity value corresponding to each discrete frequency component is calculated. The phase velocity value is calculated as: sensor spacing multiplied by angular frequency, then divided by the phase difference. The phase velocity values ​​corresponding to each discrete frequency component are arranged from low to high frequency to form the background surface wave dispersion curve.

[0023] In specific implementation, please refer to Figure 2 Using the blast initiation moment as the trigger zero point, each ultrasonic sensor synchronously acquires blast vibration signals according to a preset sampling frequency and preset acquisition duration. The preset sampling frequency is set to 10kHz, determined based on the upper limit of the ultrasonic sensor's frequency response range of 5kHz and the Nyquist sampling theorem. The preset acquisition duration is set to 30 minutes, determined based on a preset time window length of 60 seconds, the same as that used for background noise signal processing, allowing the blast vibration signal to be divided into multiple blast response segments. After receiving the blast initiation trigger signal, the ground acquisition station sends a synchronous acquisition command to each ultrasonic sensor, and each ultrasonic sensor simultaneously begins sampling, continuously recording the blast vibration signal. The blast vibration signals acquired by each ultrasonic sensor within the preset acquisition duration are arranged in chronological order, forming the blast response time sequence corresponding to each sensor. The blast response time sequences corresponding to all sensors constitute the blast response signal.

[0024] In practical implementation, the ratio of the short-term average to the long-term average of each sampling point in the blasting response signal is calculated. The short-term average is calculated using a short-term time window, with a length of 0.01 s, determined based on the time scale corresponding to the upper limit of the ultrasonic sensor's frequency response range of 5 kHz. The long-term average is calculated using a long-term time window, with a length of 0.5 s, determined based on the dominant frequency range and noise fluctuation period of the mine blasting vibration signal. For the blasting response time series... The time corresponding to each sampling point The ratio of the short-term average to the long-term average is calculated using the following formula: in, Indicates the blast response time series. The time corresponding to each sampling point It is a positive integer; Indicates The average absolute value of the amplitude of the sampling points of the blast response time series within a short time window of 0.01s, centered at the blast point; Indicates The average absolute value of the amplitude of the sampling points of the blast response time series within a long time window of 0.5s, centered at the blast point; Indicates sampling point time The ratio of the short-term average to the long-term average at a given location. The moment when the ratio of the short-term average to the long-term average first exceeds a preset trigger threshold is marked as the arrival time of the first arrival wave. The preset trigger threshold is set to 3, determined based on the upper limit of the statistical distribution of the ratio of the short-term average to the long-term average in the background noise signal before blasting, so that the background noise before blasting will not trigger the marking of the arrival time of the first arrival wave.

[0025] In practice, a valid signal segment is extracted from the blasting response signal based on the arrival time of the first arrival wave. The valid signal segment is extracted as follows: starting from the arrival time of the first arrival wave, the blasting response signal 2 seconds after the arrival time of the first arrival wave is taken as the valid signal segment. The 2-second extraction length is determined based on the period corresponding to the lower limit of the ultrasonic sensor's frequency response range (1kHz) and the time required for the surface wave to propagate in the medium surrounding the monitoring borehole. The maximum absolute value of the valid signal segment is calculated, and the amplitude of each sampling point within the valid signal segment is divided by the maximum absolute value to obtain the amplitude-normalized signal. The amplitude range of each sampling point in the amplitude-normalized signal is between -1 and 1.

[0026] In the specific implementation, the blast response signal is divided into multiple blast response segments according to a preset time window length. The preset time window length is set to 60s, the same as the preset time window length used in background noise processing. The overlap length between adjacent blast response segments is set to half the preset time window length to improve the average number and stability of the response cross-correlation function. Spectral whitening processing is performed on each blast response segment to obtain whitened blast response segments. The spectral whitening processing method is as follows: a fast Fourier transform is performed on the blast response segment to obtain the spectrum of the blast response segment; the amplitude of each frequency component in the spectrum of the blast response segment is normalized, while retaining the phase information of each frequency component to obtain the spectrum of the whitened blast response segment; an inverse fast Fourier transform is performed on the spectrum of the whitened blast response segment to obtain the whitened blast response segment.

[0027] In practical implementation, cross-correlation calculations are performed on the whitening blasting response segments corresponding to the two sensors within the same sensor pair to obtain multiple blasting cross-correlation waveforms. The cross-correlation calculation method is as follows: for the whitening blasting response segments corresponding to the first sensor and the second sensor within the same sensor pair, the whitening blasting response segment corresponding to the first sensor is fixed, and the whitening blasting response segment corresponding to the second sensor is time-shifted point by point. At each time-shifted position, the values ​​of the two whitening blasting response segments at the overlapping sampling points are multiplied point by point and then accumulated to obtain the cross-correlation value at the corresponding time-shifted position. The cross-correlation values ​​corresponding to all time-shifted positions constitute the blasting cross-correlation waveform. The cross-correlation calculation is repeated on different blasting response segments for the two sensors within the same sensor pair to obtain multiple blasting cross-correlation waveforms. The multiple blasting cross-correlation waveforms are superimposed and averaged to obtain the response cross-correlation function. The superposition and averaging method is as follows: the cross-correlation values ​​of the multiple blasting cross-correlation waveforms at the same time-shifted position are added, and the sum is divided by the total number of blasting cross-correlation waveforms to obtain the value of the response cross-correlation function at the corresponding time-shifted position.

[0028] In practical implementation, time-frequency analysis is performed on the response cross-correlation function to extract the curve of phase velocity as a function of frequency, which serves as the response surface wave dispersion curve. The time-frequency analysis employs the phase-shift method. A Fourier transform is performed on the response cross-correlation function to obtain the response cross-correlation spectrum. In the response cross-correlation spectrum, for each discrete frequency component, the phase difference between the two sensors within the same sensor pair is extracted. Based on the phase difference, the sensor spacing between the two sensors within the same sensor pair, and the angular frequency corresponding to the discrete frequency component, the phase velocity value corresponding to each discrete frequency component is calculated. The phase velocity value is calculated as follows: sensor spacing multiplied by angular frequency, then divided by phase difference. The phase velocity values ​​corresponding to each discrete frequency component are arranged in ascending order of frequency to form the response surface wave dispersion curve.

[0029] In specific implementation, please refer to Figure 3The background surface wave dispersion curve and the response surface wave dispersion curve are obtained. The background surface wave dispersion curve and the response surface wave dispersion curve are obtained from the above embodiment. The frequency range intersection of the background surface wave dispersion curve and the response surface wave dispersion curve is determined. The discrete frequency points of the background surface wave dispersion curve within the frequency range intersection are used as reference frequency coordinates, and the frequency values ​​in the reference frequency coordinates are arranged in ascending order. For the response surface wave dispersion curve, at each frequency value in the reference frequency coordinates, if the response surface wave dispersion curve has corresponding phase velocity data at the corresponding frequency value, it is directly used; if the response surface wave dispersion curve does not have corresponding phase velocity data at the corresponding frequency value, the phase velocity value of the response surface wave dispersion curve at the corresponding frequency value in the reference frequency coordinates is obtained by linear interpolation of the phase velocity values ​​of two adjacent frequency values ​​in the response surface wave dispersion curve. After frequency coordinate unification processing, the background surface wave dispersion curve and the response surface wave dispersion curve correspond to each other at the same frequency coordinates.

[0030] In practice, the phase velocity values ​​in the response surface wave dispersion curve and the phase velocity values ​​in the background surface wave dispersion curve are subtracted from each other frequency point to obtain the phase velocity difference at each frequency point. The subtraction at each frequency point is performed using the following formula: in, Represents the first frequency in the reference frequency coordinate system. One frequency value, It is a positive integer. Take from 1 , This represents the total number of frequency points in the reference frequency coordinate system. This indicates the frequency response surface wave dispersion curve at frequency The phase velocity value at that location, in meters per second; This indicates the background surface wave dispersion curve at frequency [frequency missing]. The phase velocity value at that location, in meters per second; Represents frequency The phase velocity difference at each point, expressed in meters per second. A negative value indicates that the surface response wave dispersion curve is at a frequency of The phase velocity value at that point is lower than that of the background surface wave dispersion curve at frequency. The phase velocity value at that location.

[0031] In practice, the phase velocity differences at each frequency point are arranged in ascending order of frequency to form the dispersion curve offset. Each element in the dispersion curve offset is... The order of the elements is as follows Arranged from low to high, they form the dispersion curve offset with frequency as the independent variable and phase velocity difference as the dependent variable.

[0032] In the specific implementation, a radially layered model of the surrounding rock around the monitoring borehole is established. Using the borehole axis as the central axis, the surrounding rock is divided into multiple concentric cylindrical layers radially outwards. The number of layers in the radially layered model is set to 10, determined based on the length of the ultrasonic sensor array and the radial influence range of the damage zone around the monitoring borehole. The thickness of each layer is set to 0.5m, determined based on the 1m spacing between adjacent ultrasonic sensors in the ultrasonic sensor array and the effective detection depth of the dispersion curve. An initial shear wave velocity is assigned to each layer in the radially layered model. The initial shear wave velocity is determined based on the reference shear wave velocity profile obtained by inverting the background surface wave dispersion curve before blasting. The reference shear wave velocity profile is interpolated according to the position of each layer in the radially layered model to obtain the initial shear wave velocity value corresponding to each layer.

[0033] In practice, the perturbation of the shear wave velocity in each layer of the radially layered model relative to the initial shear wave velocity is used as the parameter to be inverted. The shear wave velocity perturbations of all 10 layers are arranged in sequence to form the parameter vector to be inverted. The initial value of the parameter vector to be inverted is a zero vector.

[0034] In practical implementation, using the dispersion curve offset as the target data, the shear wave velocity of each layer in the radially layered model is iteratively adjusted based on ray theory. The dispersion curve offset is obtained from the above embodiment, where each element represents the phase velocity difference at different frequencies. According to ray theory, the sensitivity of the surface wave phase velocity to changes in the shear wave velocity depends on the energy distribution of the surface wave field in each radial layer and the ray path. When calculating the theoretical dispersion offset using forward modeling, the current shear wave velocity of each layer in the radially layered model is substituted into the surface wave dispersion forward modeling equation to calculate the theoretical phase velocity values ​​at different frequencies. Then, the theoretical phase velocity values ​​are subtracted from the corresponding phase velocity values ​​in the background surface wave dispersion curve to obtain the theoretical dispersion offset.

[0035] In practice, a Jacobian matrix is ​​established between the parameters to be inverted and the dispersion curve offset. The first [missing information] in the Jacobian matrix... Line number Column elements represent the first Layer shear wave velocity disturbance on the first The partial derivative of the dispersion offset at the nth frequency point. The Jacobian matrix elements are calculated as follows: with respect to the nth frequency point... A small perturbation is applied to the shear wave velocity of the layer, with the perturbation magnitude set to 0.001 of the initial shear wave velocity value. The surface wave dispersion forward modeling is then performed again to obtain the theoretical dispersion offset after the perturbation. This theoretical dispersion offset is then compared with the theoretical dispersion offset without the perturbation on the [missing information - likely a specific date or time]. Dividing the difference at each frequency point by the small perturbation yields the Jacobian matrix at the nth frequency point. Line number Column elements.

[0036] In practice, the least squares algorithm is used to calculate the adjustment step size of the parameters to be inverted based on the Jacobian matrix. In each iteration, the adjustment step size of the parameters to be inverted is calculated using the following formula: in, This represents the adjustment step size vector for the parameters to be inverted. The vector length is 10, and the unit is meters per second. Let the Jacobian matrix be of dimension ? , This represents the total number of frequency points in the dispersion curve offset; Representing the Jacobian matrix The transpose of the matrix, with dimension ; The damping factor ranges from 0.01 to 1.0. In this embodiment, the damping factor is set to 0.1, which is determined based on the singular value distribution of the Jacobian matrix and the noise level of the dispersion curve offset data. Represents the identity matrix with dimension . ; This represents the dispersion curve offset residual vector, obtained by subtracting the dispersion curve offset from the theoretical dispersion offset calculated in the forward modeling at each frequency point. The vector length is... After calculating the adjustment step size vector according to the above formula, each element in the adjustment step size vector is added to the shear wave velocity of the corresponding layer of the radial layered model to obtain the updated shear wave velocity of each layer.

[0037] In practice, the updated shear wave velocities of each layer are substituted into the surface wave dispersion forward modeling equation to calculate the theoretical dispersion offset. Each element in the theoretical dispersion offset is the difference between the theoretical phase velocity value obtained from the updated shear wave velocity forward modeling and the phase velocity value at the corresponding frequency in the background surface wave dispersion curve. The residual between the theoretical dispersion offset and the dispersion curve offset is calculated. The residual is calculated using the root mean square error (RMSE). The square of the difference between the theoretical dispersion offset and the dispersion curve offset at each frequency point is taken, averaged, and then the square root of the average is obtained. The RMS error is compared with a preset threshold, which is set to 0.05 of the standard deviation of the phase velocity difference at all frequency points in the dispersion curve offset. If the root mean square error is less than the preset threshold, the iteration terminates; if the root mean square error is not less than the preset threshold, the dispersion curve offset residual vector is recalculated based on the current Jacobian matrix and the updated theoretical dispersion offset, the step size is recalculated and the shear wave velocity of each layer is updated, and the above process of calculating the theoretical dispersion offset, calculating the residual, and judging whether it is less than the preset threshold is repeated until the root mean square error is less than the preset threshold and the iteration terminates.

[0038] In practice, the output is the difference between the shear wave velocity of each layer and the corresponding initial shear wave velocity at the end of the iteration. The shear wave velocity value of each layer in the radially layered model at the end of the iteration is subtracted layer by layer from the corresponding initial shear wave velocity value to obtain the change in shear wave velocity for each layer. The changes in shear wave velocity for each layer are arranged in radial sequence and output as the change in shear wave velocity profile.

[0039] See Figure 4 In the figure, the horizontal axis represents radial distance in meters, ranging from 0 to 5 meters, and the vertical axis represents shear wave velocity in meters per second, with values ​​ranging from approximately 2800 to 3100 meters per second. The legend shows three curves: the "initial shear wave velocity profile" (solid blue line), the "iteration termination shear wave velocity profile" (dashed orange line), and the "change in shear wave velocity profile" (dashed green line).

[0040] The blue solid line represents the initial shear wave velocity profile of the surrounding rock around the borehole established in the above embodiment. Along the radial distance from the inside to the outside, the shear wave velocity is maintained at about 3000 meters per second. The wave velocity is basically stable from 0 meters to about 0.6 meters in the radial direction. The wave velocity then slowly decreases, and a more obvious wave velocity drop trend appears at about 2.5 meters, and then tends to stabilize.

[0041] The orange dashed line represents the terminated shear wave velocity profile obtained through iterative inversion calculations in the above embodiment. It largely coincides with the initial velocity profile in the radial range of 0 to approximately 0.4 meters, exhibiting a consistent velocity level. In the radial range of approximately 0.4 to 2.5 meters, the orange dashed line is significantly lower than the blue solid line, indicating a decrease in shear wave velocity within this range, with a minimum velocity of approximately 2850 meters per second. In the radial range of 2.5 to 5 meters, the orange dashed line again tends to coincide with the blue solid line, and the velocity returns to the initial level.

[0042] The green dashed line represents the change in shear wave velocity profile, i.e., the difference between the terminated profile and the initial profile, with the ordinate fluctuating around zero. In the radial range of 0 to approximately 0.4 meters, the change is essentially zero, indicating no significant change in wave velocity. In the radial range of 0.4 to 2.5 meters, the change is negative, reaching a minimum of approximately -150 meters per second, indicating a significant decrease in the shear wave velocity of the surrounding rock in this range, reflecting substantial damage to the surrounding rock. In the radial range of 2.5 to 5 meters, the change gradually returns to zero, indicating that the wave velocity has returned to normal and the surrounding rock has not suffered significant damage.

[0043] In specific implementation, the above embodiment obtains the changes in shear wave velocity profiles at different depths within the monitoring borehole. Different sensor pairs in the ultrasonic sensor array within the monitoring borehole correspond to different depth intervals, and the changes in shear wave velocity profiles of the surrounding rock within the corresponding depth interval are obtained through inversion for each sensor pair. The changes in shear wave velocity profiles corresponding to all sensor pairs within the monitoring borehole are arranged from shallowest to deepest according to the depth of the sensor pair, forming a sequence of shear wave velocity profile changes along the depth direction of the monitoring borehole. Each depth position in the sequence corresponds to a shear wave velocity change value; the depth position is the depth of the center point of the sensor pair, with the depth unit being meters and the shear wave velocity change unit being meters per second.

[0044] In practice, the rate of change of the shear wave velocity profile is calculated along the depth direction of the monitored borehole. For every two adjacent depth positions in the shear wave velocity profile change sequence, the difference in the shear wave velocity profile change is calculated, and then divided by the depth interval between the adjacent depth positions to obtain the rate of change between the corresponding adjacent depth positions. The rate of change is calculated using the following formula: in, Indicates the first [number] along the direction of the monitored borehole depth. Location of each depth sampling point It is a positive integer. Take from 1 , This represents the total number of depth sampling points in the sequence of changes in the shear wave velocity profile within the borehole. Indicates the first [number] along the direction of the monitored borehole depth. The location of each depth sampling point is shown in meters. Indicates the location of the depth sampling point The change in transverse wave velocity profile at a given location, in meters per second; Indicates the location of the depth sampling point The change in transverse wave velocity profile at a given location, in meters per second; Indicates the location of the depth sampling point The rate of change of shear wave velocity profile along the depth direction is expressed in meters per second per meter. Each calculated rate of change value is mapped to the midpoint depth of adjacent depth locations, forming a sequence of rates of change along the depth direction.

[0045] In the specific implementation, a preset change threshold is set to 0.2 m / s / m. The preset change threshold is determined based on the following: the natural gradient of shear wave velocity along the depth direction in undamaged surrounding rock is less than 0.1 m / s / m; blasting damage causes a decrease in shear wave velocity in the surrounding rock, and the gradient of the shear wave velocity profile change along the depth direction at the boundary of the damaged zone exceeds 0.2 m / s / m. Therefore, the preset change threshold is set to 0.2 m / s / m. The change rate sequence is traversed, and continuous depth segments with absolute change rates exceeding the preset change threshold are identified as candidate damaged zone segments. Within a continuous depth segment, all absolute change rate values ​​exceed the preset change threshold, and the absolute change rate values ​​of adjacent segments on both sides of the continuous depth segment are less than or equal to the preset change threshold.

[0046] In practice, depth boundaries are marked for continuous depth segments. The upper boundary of a continuous depth segment is the depth value corresponding to the shallowest depth position within that segment, and the lower boundary is the depth value corresponding to the deepest depth position within that segment. If multiple continuous depth segments exist, the upper and lower boundaries of each segment are marked separately, and these two boundaries are used together as the depth boundary of the damaged area. Along the borehole depth direction, the depth boundary marking result of the damaged area includes the upper and lower boundary depth values ​​of one or more damaged segments.

[0047] In practice, the depth boundaries of multiple monitoring boreholes located at different spatial positions are obtained. These boreholes are arranged in a spatial grid around the mine blasting area. The spatial coordinates of each borehole are obtained through mine surveying, including east, north, and elevation coordinates. Within each borehole, the depth boundary of its respective damaged zone is determined using the aforementioned depth boundary marking method, including the upper and lower boundary depth values ​​of each damaged section. The upper boundary depth values ​​of all monitoring boreholes are converted to upper boundary elevation coordinates in three-dimensional space, and the lower boundary depth values ​​of all monitoring boreholes are converted to lower boundary elevation coordinates in three-dimensional space.

[0048] In practice, the boundaries of each depth are fitted to surfaces within the same three-dimensional geological model. The three-dimensional geological model uses a three-dimensional block model of a mining area, with its spatial extent extending 10m beyond the monitoring borehole distribution area. The block size in the three-dimensional block model is 1m × 1m × 1m. Surface fitting employs radial basis function interpolation, using a Gaussian function with a shape parameter set to 1.0. The shape parameter value is based on the ratio of the monitoring borehole spacing to the surface undulation scale of the damaged area boundary. With a monitoring borehole spacing of 5m and a damaged area boundary surface undulation scale of 5m, a shape parameter of 1.0 ensures the smoothness of the interpolated surface matches the actual undulation of the damaged area boundary. Surface fitting is divided into upper boundary surface fitting and lower boundary surface fitting: the three-dimensional spatial points of the upper boundary of each monitoring borehole are used as data points for upper boundary surface fitting, and radial basis function interpolation is performed to obtain the upper boundary surface; similarly, the three-dimensional spatial points of the lower boundary of each monitoring borehole are used as data points for lower boundary surface fitting, and radial basis function interpolation is performed to obtain the lower boundary surface. The three-dimensional spatial region between the upper boundary surface and the lower boundary surface is the spatial distribution range of the damage zone.

[0049] The above description is merely a specific embodiment of this application, but the scope of protection of this application is not limited thereto. Any changes or substitutions that can be easily conceived by those skilled in the art within the scope of the technology disclosed in this application should be included within the scope of protection of this application.

Claims

1. A method for detecting and analyzing damage to surrounding rock mass in mine blasting based on ultrasonic testing, characterized in that, The method includes: An ultrasonic sensor array was installed in the monitoring borehole where the surrounding rock was preserved, and the background noise signal received by the ultrasonic sensor array before blasting was collected. Based on the background noise signal, calculate the background cross-correlation function between each sensor pair, and extract the background surface wave dispersion curve from the background cross-correlation function; After the blast is triggered, the blast response signal received by the ultrasonic sensor array is continuously collected; Based on the blasting response signal, calculate the response cross-correlation function between each sensor pair, and extract the response surface wave dispersion curve from the response cross-correlation function; The difference between the background surface wave dispersion curve and the response surface wave dispersion curve is processed to obtain the dispersion curve offset. Based on the dispersion curve offset, the change in the shear wave velocity profile of the surrounding rock around the monitoring borehole is inverted. Based on the change in shear wave velocity profile, the depth boundary and spatial distribution range of the damage zone inside the surrounding rock are determined.

2. The method for detecting and analyzing rock mass damage in mine blasting based on ultrasonic testing according to claim 1, characterized in that, An ultrasonic sensor array was deployed in the monitoring borehole where the surrounding rock was preserved, and the background noise signal received by the ultrasonic sensor array before blasting was collected, specifically: Along the axial direction of the monitoring borehole, multiple ultrasonic sensors in the ultrasonic sensor array are sequentially arranged at a preset equal interval on the inner wall of the monitoring borehole, and each ultrasonic sensor is fixedly coupled to the inner wall of the monitoring borehole through a grouting coupling medium. During the quiet period before the blast and in the absence of external vibration sources, each of the ultrasonic sensors is activated to synchronously collect environmental vibration signals and obtain a time-series signal of a preset duration, which serves as the background noise signal.

3. The method for detecting and analyzing rock mass damage in mine blasting based on ultrasonic testing according to claim 1, characterized in that, Based on the background noise signal, the background cross-correlation function between each sensor pair is calculated, and the background surface wave dispersion curve in the background cross-correlation function is extracted, specifically as follows: The background noise signal is divided into multiple background noise segments according to a preset time window length, and each background noise segment is subjected to spectral whitening processing to obtain whitened background noise segments. The whitened background noise segments corresponding to two sensors within the same sensor pair are cross-correlated to obtain multiple cross-correlation waveforms. The multiple cross-correlation waveforms are then superimposed and averaged to obtain the background cross-correlation function. Time-frequency analysis is performed on the background cross-correlation function to extract the curve of phase velocity versus frequency, which is used as the background surface wave dispersion curve.

4. The method for detecting and analyzing rock mass damage in mine blasting based on ultrasonic testing according to claim 1, characterized in that, After the blast is triggered, the blast response signal received by the ultrasonic sensor array is continuously collected, specifically: Using the moment of blasting initiation as the trigger zero point, each of the ultrasonic sensors is activated to synchronously collect blasting vibration signals according to the preset sampling frequency and preset acquisition duration; The blasting vibration signals collected by each ultrasonic sensor within the preset acquisition time are arranged in chronological order to obtain the blasting response time sequence corresponding to each sensor, which is used as the blasting response signal.

5. The method for detecting and analyzing rock mass damage in mine blasting based on ultrasonic testing according to claim 3, characterized in that, Based on the blast response signal, the response cross-correlation function between each sensor pair is calculated, and the response surface wave dispersion curve in the response cross-correlation function is extracted, specifically as follows: The blasting response signal is divided into multiple blasting response segments according to the preset time window length, and each blasting response segment is subjected to spectral whitening processing to obtain whitened blasting response segments. The whitening blasting response segments corresponding to two sensors within the same sensor pair are cross-correlated to obtain multiple blasting cross-correlation waveforms. The multiple blasting cross-correlation waveforms are then superimposed and averaged to obtain the response cross-correlation function. Time-frequency analysis is performed on the response cross-correlation function to extract the curve of phase velocity as a function of frequency, which is used as the response surface wave dispersion curve.

6. The method for detecting and analyzing rock mass damage in mine blasting based on ultrasonic testing according to claim 1, characterized in that, The difference between the background surface wave dispersion curve and the response surface wave dispersion curve is processed to obtain the dispersion curve offset, specifically: Under the same frequency coordinate, the phase velocity value in the response surface wave dispersion curve is subtracted from the phase velocity value in the background surface wave dispersion curve at each frequency point to obtain the phase velocity difference at each frequency point. The phase velocity differences at each frequency point are arranged in order of increasing frequency to form the dispersion curve offset.

7. The method for detecting and analyzing rock mass damage in mine blasting based on ultrasonic testing according to claim 1, characterized in that, Based on the dispersion curve offset, the inversion of the shear wave velocity profile change of the surrounding rock around the monitoring borehole specifically includes: A radial stratification model of the surrounding rock around the monitoring borehole is established, and an initial shear wave velocity is assigned to each layer in the radial stratification model; Using the dispersion curve offset as the target data, the transverse wave velocity of each layer in the radial layered model is iteratively adjusted based on ray theory until the residual between the theoretical dispersion offset calculated by forward modeling and the dispersion curve offset is less than a preset threshold. The difference between the shear wave velocity of each layer and the corresponding initial shear wave velocity at the end of the output iteration is used as the change in the shear wave velocity profile.

8. The method for detecting and analyzing rock mass damage in mine blasting based on ultrasonic testing according to claim 1, characterized in that, Based on the changes in the shear wave velocity profile, the depth boundary and spatial distribution range of the damage zone within the surrounding rock are determined, specifically as follows: Along the depth direction of the monitoring borehole, the upper and lower boundaries of the continuous segments in the change rate of the shear wave velocity profile exceeding a preset change threshold are marked as the depth boundaries of the damage zone. The depth boundaries of each of the multiple monitoring boreholes arranged in different spatial locations are obtained, and the depth boundaries are fitted to the same three-dimensional geological model to obtain the spatial distribution range of the damage zone.

9. The method for detecting and analyzing rock mass damage in mine blasting based on ultrasonic testing according to claim 1, characterized in that, After triggering the blast, following the continuous acquisition of the blast response signal received by the ultrasonic sensor array, the process further includes first-arrival wave pickup and amplitude normalization of the blast response signal, specifically: Calculate the ratio of the short-term average value to the long-term average value of each sampling point in the blasting response signal, and mark the moment when the ratio first exceeds the preset trigger threshold as the arrival time of the first wave; Based on the arrival time of the first wave, an effective signal segment is extracted from the blast response signal, and the maximum absolute value of the effective signal segment is calculated. The amplitude of each sampling point within the effective signal segment is divided by the maximum absolute value to obtain the amplitude-normalized signal.

10. The method for detecting and analyzing rock mass damage in mine blasting based on ultrasonic testing according to claim 7, characterized in that, Using the dispersion curve offset as the target data, and based on ray theory, the shear wave velocity of each layer in the radial layering model is iteratively adjusted until the residual between the theoretical dispersion offset calculated in forward modeling and the dispersion curve offset is less than a preset threshold. Specifically: The perturbation of the shear wave velocity of each layer in the radially layered model relative to the initial shear wave velocity is used as the parameter to be inverted, and the Jacobian matrix between the perturbation and the dispersion curve offset is established. The least squares algorithm is used to calculate the adjustment step size of the parameters to be inverted based on the Jacobian matrix, and the shear wave velocity of each layer is updated according to the adjustment step size. Substitute the updated shear wave velocities of each layer into the surface wave dispersion forward modeling equation to calculate the theoretical dispersion offset. Repeat the calculation of the residual between the theoretical dispersion offset and the dispersion curve offset until the residual is less than the preset threshold and the iteration is terminated.

Images