Blasting vibration signal inversion and surrounding rock damage intelligent identification method

Abstract

The invention provides a blasting vibration signal inversion and surrounding rock damage intelligent identification method, which comprises the following steps of: arranging a sensor at a key part of a surrounding rock to acquire a blasting vibration signal, and calibrating a surrounding rock damage index; preprocessing the blasting vibration signals to extract multi-dimensional features; after dimension reduction optimization is carried out on the multi-dimensional features, a training set is constructed with the surrounding rock damage index, and a machine learning model is trained; and inputting blasting vibration signals collected in real time into the machine learning model to identify the damage state of the surrounding rock, and adjusting blasting parameters according to an identification result. Through deep fusion of machine learning and geotechnical engineering, the key technical problems of low precision, poor adaptability, weak real-time performance and the like in the prior art are solved, experience-driven blasting construction is converted into data-driven blasting construction, reliable technical guarantee is provided for safe, efficient and intelligent construction of complex underground engineering, and the method is suitable for popularization and application. The method has wide popularization and application prospects and remarkable social and economic benefits.

Description

Technical Field

[0001] This invention belongs to the field of geotechnical engineering technology, specifically relating to a method for inverting blasting vibration signals and intelligently identifying surrounding rock damage. Background Technology

[0002] In deep underground engineering using the drill-and-blast method, the surrounding rock damage caused by blasting vibrations seriously threatens construction safety, and its accurate assessment and real-time identification are key challenges in the industry. When vibration waves generated by blasting propagate through the surrounding rock, they cause a redistribution of the internal stress field, leading to damage phenomena such as crack propagation and strength deterioration. If this damage is not identified and effectively controlled in a timely manner, it may trigger serious engineering accidents such as surrounding rock instability and collapse, endangering not only the safety of construction personnel but also causing huge economic losses and project delays. Therefore, establishing an intelligent system capable of accurately identifying the state of surrounding rock damage in real time and providing feedback to adjust blasting parameters is of significant practical importance for ensuring the safety of underground engineering construction.

[0003] Existing methods for assessing surrounding rock damage mainly include empirical formula methods, numerical simulation methods, and preliminary machine learning methods. Empirical formula methods typically rely on empirical relationships between vibration velocity, explosive charge, and distance (such as the Sadovsky formula) to predict vibration intensity and thus assess damage risk. While these methods are computationally simple, they have poor adaptability, heavily relying on specific geological conditions and blasting parameters. When the engineering environment changes, the prediction accuracy drops significantly, making it difficult to accurately reflect the nonlinear relationship between the complex and changing geological environment and the blasting dynamic response. Numerical simulation methods, by establishing finite element or discrete element models of the surrounding rock-explosive-air coupling, simulate stress wave propagation and the dynamic response of the surrounding rock during blasting. Although they can consider relatively complex geometric conditions and material properties, these methods heavily rely on the selection of empirical parameters, and are computationally expensive and time-consuming, making it difficult to meet the real-time requirements for rapid identification and dynamic feedback of surrounding rock damage at construction sites, resulting in control lag.

[0004] In recent years, machine learning technology has shown potential in geotechnical engineering monitoring, and some studies have attempted to apply it to blasting vibration analysis and surrounding rock stability evaluation. However, current applications still have significant limitations: First, most studies are limited to laboratory conditions or single geological environments, lacking adaptability verification for complex and variable engineering scenarios; second, feature extraction often relies on human experience and fails to deeply integrate multi-dimensional indicators such as time domain, frequency domain, time-frequency domain, and nonlinearity, resulting in insufficient mining of damage-sensitive information contained in vibration signals; third, most existing methods remain at the offline analysis or isolated identification stage, lacking the closed-loop control capability to automatically feed the identification results back to the construction stage to optimize blasting parameters, resulting in insufficient system intelligence and failure to form a complete closed loop from monitoring and identification to decision-making and feedback.

[0005] In summary, there is an urgent need for an innovative method that can uncover the deep characteristics of blasting vibration signals, establish a data-driven damage identification model, achieve real-time intelligent identification, and form a closed-loop control with construction feedback, so as to fundamentally solve the key technical problems of poor adaptability, low accuracy, weak real-time performance, and insufficient intelligence in existing technologies. Summary of the Invention

[0006] To address the aforementioned technical problems, this invention provides a method for inverting blasting vibration signals and intelligently identifying surrounding rock damage, comprising the following steps: Step S1: Install sensors at key locations in the surrounding rock to collect blasting vibration signals and calibrate the surrounding rock damage index; Step S2: Preprocess the blasting vibration signal to extract multidimensional features; Step S3: After dimensionality reduction optimization of the multidimensional features, construct a training set with the surrounding rock damage index and train a machine learning model; Step S4: Input the real-time collected blasting vibration signal into the machine learning model to identify the surrounding rock damage state, and adjust the blasting parameters according to the identification results.

[0007] Preferably, the method of step S1 includes: deploying at least one triaxial acceleration sensor at at least one location among the working face, the arch, the sidewalls and the bottom plate, and synchronously collecting the blasting vibration signal; The method for calibrating the surrounding rock damage index includes: quantifying the degree of surrounding rock damage and normalizing it into a damage index as label data using at least one non-destructive testing method.

[0008] Preferably, the preprocessing in step S2 includes filtering and noise reduction, amplitude normalization, and signal segmentation; The filtering and noise reduction includes: applying a bandpass filter to the blasting vibration signal, and applying wavelet threshold noise reduction when the signal quality is below a threshold. The formula for calculating the amplitude normalization is: ; in, For the original signal data points, The signal sequence to be normalized is... These are the normalized data points; The signal segmentation includes: determining the zero moment of detonation based on the trigger signal synchronized with the detonation system, and extracting the signal within a fixed window after the zero moment of detonation as the effective analysis segment.

[0009] Preferably, the multidimensional features include time-domain features, frequency-domain features, time-frequency-domain features, and nonlinear features; The time-domain features include peak vibration velocity, root mean square amplitude, total energy of the signal envelope, and energy concentration. The frequency domain features include dominant frequency, center frequency, frequency band energy percentage, and spectral entropy; The time-frequency domain features include instantaneous energy distribution, time-frequency centroid trajectory, and wavelet packet energy entropy; The nonlinear characteristics include fractal dimension, Lyapunov exponent, nonlinear energy ratio, energy decay coefficient, and frequency drift rate.

[0010] Preferably, the peak vibration velocity The calculation formula is: ; in, Indicates vibration velocity signal, Indicates time; The method for calculating the total energy of the signal envelope is as follows: The instantaneous amplitude envelope of the signal is calculated using the Hilbert transform. The expression is: ; in, This is the preprocessed vibration signal. For signal The result after Hilbert transform, This represents the instantaneous amplitude envelope of the signal.

[0011] Preferably, the formula for calculating the time-frequency centroid trajectory is: ; ; in, This represents the time-frequency representation of the signal. Indicates the amplitude in time-frequency representation; For frequency; Indicates time; For time The frequency centroid at that location; For frequency The center of mass of time; The formula for calculating the wavelet packet energy entropy is: ; in, Indicates the first i The percentage of energy of each node in the total energy; Represents the wavelet packet energy entropy; This represents the total number of wavelet packet nodes.

[0012] Preferably, the Lyapunov index The expression is: ; in, Indicates time; For time Separation vectors of adjacent orbits in the post-phase space; The separation vector of adjacent orbits in the phase space at the initial time step; The nonlinear energy ratio The calculation formula is: ; in, The original signal; For linear prediction models in time The predicted value at that location; The energy attenuation coefficient The calculation formula is: ; ; in, Indicates the initial envelope amplitude; For time The amplitude of the signal envelope at that location; The frequency drift rate The calculation formula is: ; in, The dominant frequency after damage; The dominant frequency before damage.

[0013] Preferably, the machine learning model includes at least one of the following: support vector machine model, random forest model, extreme gradient boosting model, deep neural network model, one-dimensional convolutional neural network model, and long short-term memory network model; The process of training a machine learning model includes: dividing the training set into a training subset, a validation subset, and a test subset; fitting model parameters on the training subset; optimizing hyperparameters on the validation subset using grid search, random search, or Bayesian optimization; and evaluating model performance on the test subset.

[0014] Preferably, when the machine learning model is used for classifying the damage level of surrounding rock, the performance evaluation metrics of the model include accuracy, precision, recall, F1 score, and confusion matrix. When the machine learning model is used for predicting the continuous damage index, the evaluation metrics include the root mean square error (RMSE) and the coefficient of determination (R²).

[0015] Preferably, step S4 includes: Step S41: Input the feature vector to be identified into the machine learning model, output the surrounding rock damage level or surrounding rock damage index, and calculate the prediction confidence. Step S42: When the prediction confidence level is lower than a preset threshold, a manual review process is triggered; Step S43: When the surrounding rock damage level is greater than a preset threshold, or the surrounding rock damage index exceeds a safety threshold, adjust the single-stage dosage for the next cycle according to the damage identification result using the following formula: ; in, This is the current dosage for a single dose; This is the adjusted single-stage dosage; The damage index predicted by the model; The damage threshold; The maximum permissible damage index; For adjustment coefficients; Step S44: Add the real-time collected blasting vibration signal, the extracted multidimensional features and the corresponding damage identification results to the training dataset, and perform online learning and updating of the machine learning model every few blasting cycles.

[0016] The beneficial effects of the present invention include at least the following: This invention systematically extracts time-domain, frequency-domain, time-frequency-domain, and nonlinear damage-sensitive features from blasting vibration signals, constructing a comprehensive feature set with hundreds of dimensions. Then, feature dimensionality reduction optimization techniques (correlation analysis, PCA, tree-based feature selection, etc.) are used to retain the most discriminative feature subset. Finally, machine learning models (SVM, random forest, XGBoost, deep neural networks, CNN, LSTM, etc.) are used to establish a complex mapping relationship between features and damage. This data-driven approach can adaptively learn damage patterns under different geological conditions (such as lithology, integrity, geostress, and groundwater) and blasting parameters (such as charge, detonation sequence, and borehole layout), avoiding errors caused by simplified assumptions (such as assuming the surrounding rock is a homogeneous elastic body and blasting is a point source) and empirical parameter values ​​(such as attenuation exponent and safety threshold depending on experience) in traditional methods.

[0017] This invention optimizes signal processing procedures and feature extraction algorithms, combined with an efficient machine learning inference model, enabling the entire process from vibration signal acquisition to damage status identification within seconds to tens of seconds after blasting. This real-time capability allows the system to provide timely assessment results and adjustment suggestions before the next blasting cycle, achieving dynamic feedback control and avoiding the cumulative damage and accident risks caused by analysis lag in traditional methods.

[0018] This invention, through a complete implementation plan design, forms a standardized and scalable technical system, encompassing sensor deployment, signal processing, multi-dimensional feature extraction, machine learning modeling, and final system integration and visualization. The system establishes a complete closed loop of "monitoring-identification-decision-feedback": the monitoring module collects vibration data in real time; the identification module quickly determines the damage state based on a machine learning model; the decision-making module automatically generates parameter adjustment plans based on the identification results; the feedback module applies the adjustment plans to the next construction cycle; and the online learning module continuously collects new data and updates the model. This closed-loop control mechanism realizes the transformation from "experience-driven" to "data-driven," and from passive response to proactive prevention. It can maintain optimal performance over a long period, avoiding the performance degradation caused by changes in working conditions in traditional fixed models. This greatly expands its applicability and practical value in complex underground engineering projects such as deep tunnels, underground powerhouses, and mine roadways.

[0019] This invention achieves precise control and optimized utilization of blasting energy while ensuring construction safety by accurately identifying the damage state of the surrounding rock and adjusting blasting parameters in real time. When minor damage is detected, the charge can be appropriately increased to improve construction efficiency; when damage approaches a threshold, the charge is promptly reduced or the detonation sequence is adjusted to control vibration intensity. This dynamic optimization strategy keeps the blasting vibration intensity within a reasonable range, minimizing damage to the surrounding rock and reducing disturbance to the surrounding environment. Simultaneously, by reducing over-excavation and material consumption, it indirectly reduces resource waste and carbon emissions, mitigating the impact of construction on the ecological environment.

[0020] This invention establishes a real-time intelligent identification system that can issue early warnings and take measures before damage accumulates to a dangerous level, effectively preventing major accidents such as collapses and roof falls caused by surrounding rock instability, and ensuring the safety of construction personnel and equipment. Simultaneously, by optimizing blasting parameters, the blasting effect is improved (such as uniformity of fragmentation and molding quality), reducing subsequent secondary crushing and repair work, shortening the construction period, and improving overall construction efficiency. In terms of economic benefits, it significantly reduces the total project cost.

[0021] In summary, this invention, through the deep integration of machine learning and geotechnical engineering, not only solves the technical problems of low accuracy, poor adaptability, and weak real-time performance in existing technologies, but also realizes the transformation of blasting construction from experience-driven to data-driven through an intelligent closed-loop control system. It provides reliable technical support for the safe, efficient, green, and intelligent construction of complex underground engineering projects, and has broad application prospects and significant socio-economic and environmental benefits. Attached Figure Description

[0022] Figure 1 This is a schematic diagram of the overall scheme of an embodiment of the present invention; Figure 2 This is a schematic diagram of the method flow according to an embodiment of the present invention. Detailed Implementation

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

[0024] like Figure 1 and Figure 2 As shown, this invention provides a method for inverting blasting vibration signals and intelligently identifying surrounding rock damage, aiming to solve the following technical problems: First, there are problems with poor adaptability and low accuracy. Traditional empirical formula methods rely on specific conditions and simplifying assumptions, making it difficult to accurately reflect the nonlinear relationship between complex and changing geological environments and blasting dynamic responses; numerical simulation methods rely heavily on empirical parameters and have high computational costs, making it difficult for both to achieve accurate and universal assessment of surrounding rock damage.

[0025] Second, there is the problem of weak real-time performance. Neither rough estimations based on empirical formulas nor time-consuming numerical simulations can meet the real-time requirements for rapid identification and dynamic feedback of surrounding rock damage at construction sites, resulting in delayed control.

[0026] Third, there is the problem of low intelligence and lack of closed-loop system. Most existing methods (including early machine learning applications) remain at the stage of offline analysis or isolated identification. Feature extraction relies on manual labor and is not comprehensive enough. There is a lack of closed-loop control capability to automatically feed the identification results back to the construction stage to optimize blasting parameters, resulting in insufficient system intelligence.

[0027] Therefore, this invention provides a method for inverting blasting vibration signals and intelligently identifying surrounding rock damage, including the following steps: Step S1: Data Acquisition and Monitoring Setup This step involves constructing a three-dimensional monitoring network. Based on the propagation laws of blasting dynamics and the structural characteristics of the surrounding rock, sensors are deployed at key locations such as the working face, arch, sidewalls, and floor, and even within boreholes at different depths. This principle of "key locations, three-dimensional intersection" deployment aims to create a sensing system capable of comprehensively capturing the spatial propagation and attenuation laws of blasting vibrations, laying the foundation for subsequent inversion of the spatial distribution of surrounding rock damage.

[0028] In terms of sensor selection and configuration, at least one triaxial ICP accelerometer is used as the core component. Its wide frequency range and excellent dynamic response characteristics enable it to completely record the vector information of vibration signals. The acceleration components in three orthogonal directions can comprehensively reflect the spatial characteristics of vibration, avoiding information loss caused by installation orientation deviations of uniaxial sensors. The acquisition frequency is strictly set between 1000Hz and 5000Hz according to the Nyquist sampling theorem to ensure distortion-free capture of the high-frequency components generated by blasting, while retaining sufficient frequency resolution for signal analysis. The sensor range is selected based on the estimated vibration intensity, avoiding signal clipping distortion while ensuring a high signal-to-noise ratio. Good coupling with the surrounding rock is ensured to reduce the impact of installation errors on measurement accuracy.

[0029] To achieve precise and synchronous data acquisition, the data acquisition system is linked to the on-site blasting initiation system via digital delay relays or dedicated blasting vibration recorders. At the moment of detonation, the system sends a unified trigger signal to all data acquisition channels, ensuring that all measurement points have an absolutely consistent zero-time point. This is a crucial prerequisite for accurate signal segmentation, multi-point comparison, and time-frequency analysis. The synchronous triggering mechanism also effectively distinguishes blasting vibration from environmental background noise, improving signal quality.

[0030] To construct a comprehensive database covering various operating conditions, vibration signals need to be systematically collected across multiple blasting cycles under different blasting design parameters (such as single-stage charge, total charge, detonation sequence, and borehole layout) and geological conditions (such as lithology, integrity, geostress state, and groundwater conditions). Each data record is associated with its metadata, including blasting parameters, geological conditions, measuring point location, and environmental factors, forming a structured, multi-dimensional dataset suitable for machine learning. Data collection should cover the entire process from the test blasting phase to the normal construction phase, ensuring the representativeness and diversity of the data.

[0031] To provide accurate labeling data for supervised learning models, multiple non-destructive testing methods are needed to quantify and calibrate the damage to the surrounding rock. Borehole television (in-hole imaging) is used to observe changes in parameters such as the number, width, and extension length of fractures before and after blasting, qualitatively evaluating the degree of damage. Acoustic transmission methods or cross-hole acoustic wave testing are used to measure changes in the longitudinal wave velocity of the surrounding rock, calculating the wave velocity damage factor; a greater decrease in wave velocity indicates more severe damage. Simultaneously, acoustic emission sensors are deployed to monitor micro-fracture activity in the surrounding rock; the accumulated acoustic emission energy or the number of events can serve as a quantitative indicator of damage evolution. Finally, these indicators are normalized into a comprehensive damage index ranging from 0 (no damage) to 1 (severe damage) using expert rules or weighted fusion algorithms, serving as the target variable for training the machine learning model. This multi-source information fusion calibration method can more comprehensively and accurately reflect the true damage state of the surrounding rock.

[0032] Step S2: Signal Preprocessing The core objective of signal preprocessing is to purify the raw signal, standardize the data scale, and enhance effective information. It transforms the raw vibration data collected on-site, which contains noise and non-stationary characteristics, into a high-quality feature extraction source suitable for subsequent machine learning model input. This step directly determines the reliability of feature engineering and the accuracy of model recognition.

[0033] First, adaptive filtering and noise reduction are performed. To address the low-frequency baseline drift (potentially caused by sensor zero drift or slow ground subsidence), high-frequency electromagnetic interference (from the detonator, cables, or surrounding electrical equipment), and random noise mixed in the blasting vibration signal, a multi-stage filtering strategy is employed. A Butterworth or Chebyshev bandpass filter with a passband frequency of 1-500Hz is used to preserve the dominant frequency energy of the blasting vibration while eliminating out-of-band noise. The filter order is selected based on signal characteristics and noise levels, achieving a balance between stopband attenuation and phase distortion. In cases of poor signal quality or strong interference, wavelet thresholding denoising is further introduced. By performing multi-scale wavelet decomposition on the signal, the threshold (using either a soft or hard threshold function) is adaptively selected based on noise energy in different frequency bands. The decomposed wavelet coefficients are then contracted, effectively suppressing noise while maximizing the preservation of the transient impact characteristics and abrupt changes of the blasting vibration signal. The advantage of wavelet denoising over traditional filtering lies in its time-frequency localization characteristics, which better protect the non-stationary components of the signal.

[0034] Next, data normalization and scaling are performed. To eliminate the problem of inconsistent signal amplitude dimensions caused by differences in sensor sensitivity, inconsistent installation coupling states, or differences in explosive equivalent, the preprocessed signal is normalized. A minimum-maximum normalization method is used to linearly map the amplitude of all signals to the interval [0,1] or [-1,1]. The calculation formula is as follows: ; in, For the original signal data points, The signal sequence to be normalized is... For the normalized data points, their values ​​are linearly scaled to the [0,1] interval. and Signal sequences The maximum and minimum values ​​in the data are determined. This process ensures that data collected from different measurement points and different bursts are on the same scale, avoiding the machine learning model's oversensitivity to magnitude during training, thereby improving the model's convergence speed and generalization performance.

[0035] Then, precise signal segmentation and extraction are performed. Based on a high-precision timescale triggered synchronously with the detonation system, the zero-point of the blast is determined. Using this zero point as the zero point, a fixed time window (e.g., 0-200 ms or 0-500 ms after detonation, the specific length depending on the blast scale and surrounding rock characteristics) is extracted as the effective analysis segment. This window must cover the complete process from the arrival of the stress wave at the measuring point to the significant decay of the vibration energy, ensuring that it contains the core information reflecting the dynamic response of the surrounding rock. By analyzing the energy time history curves of typical blasting events, when the vibration energy decays to less than 5% of the initial peak value and persists for a certain period of time, the effective vibration can be considered to have ended. The extraction window must be long enough to capture the complete response, but also avoid being too long to include too much invalid tail data. Signal segmentation automatically eliminates a large amount of invalid background noise before detonation and long-tail data after the vibration has completely subsided, greatly improving the computational efficiency and targeting of subsequent feature extraction.

[0036] Finally, preliminary feature enhancement calculations are performed to highlight key components of the signal that are sensitive to impairment. The instantaneous amplitude envelope of the signal is calculated using the Hilbert transform. The expression is: ; in, This is the preprocessed vibration signal. For signal The result after Hilbert transform. The Hilbert transform can be understood as delaying the phase of each frequency component of the original signal by 90 degrees. Instantaneous amplitude envelope. The algorithm describes the contour of the signal amplitude change over time, visually displaying the temporal distribution of vibrational energy and highlighting the signal's impact characteristics and periods of concentrated energy. It also calculates the instantaneous energy time history of the signal, emphasizing the concentrated periods of energy release, and estimates the signal's power spectral density, providing a preliminary observation of its dominant frequency distribution. These preliminary enhancement calculations provide directional guidance for subsequent in-depth frequency and time-frequency domain feature extraction, helping to identify key frequency components and time periods within the signal.

[0037] Step S3: Multidimensional Feature Extraction and Construction This step systematically extracts multi-dimensional and multi-scale features from the preprocessed signal to construct a high-dimensional feature space that can comprehensively characterize the blasting vibration characteristics and their correlation with surrounding rock damage. These features will serve as input to the machine learning model, and their quality and comprehensiveness directly determine the model's ability and accuracy in identifying surrounding rock damage. The feature system includes four main categories: time-domain features, frequency-domain features, time-frequency-domain features, and nonlinear features.

[0038] (1) Temporal feature extraction Temporal features are extracted directly from the amplitude variation of the signal over time, reflecting the intensity, distribution, and shape characteristics of the vibration. Peak vibration velocity is extracted. The formula for calculating the maximum impact strength is as follows: ; in, Represents vibration velocity signal (unit: m / s); Peak velocity is the most commonly used vibration intensity control index in blasting engineering. In addition, peak acceleration (PPA) and peak displacement (PPD) are extracted and calculated using the integral or differential relationship of the acceleration signal.

[0039] Obtaining the signal envelope using Hilbert transform The total energy of the signal envelope is calculated to characterize the total energy released by the vibration. Energy concentration is also calculated, defined as the proportion of the main energy to the total energy. For example, it is the ratio of the energy within a time window of ±Δt near the peak of the envelope curve to the total energy. This indicator reflects the degree of concentration of vibration energy release; a higher concentration indicates a more impactful vibration. Simultaneously, the cumulative absolute amplitude and cumulative energy are calculated, reflecting the evolution of vibration intensity over time.

[0040] (2) Frequency domain feature extraction The signal is transformed from the time domain to the frequency domain using Fast Fourier Transform (FFT) or power spectral density estimation to reveal its frequency structure characteristics. The dominant frequency and center frequency are extracted to identify the dominant vibrational frequencies; changes in these frequencies reflect alterations in the surrounding rock properties. When the surrounding rock is damaged, its natural frequencies decrease. The frequency range is divided into low-frequency, mid-frequency, and high-frequency sub-bands, and the energy proportion of each band is calculated. This helps identify frequency components related to rock mass resonance.

[0041] Simultaneously calculating the spectral width (standard deviation of the spectrum) and spectral entropy, the former measures the extent of the expansion of the spectral distribution, and the latter quantifies the complexity and uncertainty of the spectrum, these characteristics are indicative of different types of dynamic response modes of rock masses.

[0042] (3) Time-frequency domain feature extraction To address the non-stationary characteristics of blasting vibration signals (where the frequency components of the vibration change rapidly over time), relying solely on global frequency domain analysis would result in the loss of time localization information. Therefore, time-frequency analysis techniques (such as Short-Time Fourier Transform (STFT), Wavelet Transform (WT), and Wigner-Ville Distribution (WVD)) are employed to simultaneously acquire the energy distribution of the signal in both time and frequency.

[0043] Extract the instantaneous energy distribution within different time windows to capture the temporal localization characteristics of vibrational energy.

[0044] Calculate the time-frequency centroid trajectory, which reflects the path of the signal energy center as it changes over time. (Frequency centroid) and the center of time The calculation formula is: ; ; in, The time-frequency representation of a signal; Indicates frequency; Indicates time; This represents the frequency centroid at time t; Represents frequency The center of mass of time; This indicates the amplitude of the time-frequency representation.

[0045] The time-frequency centroid trajectory can reveal the migration pattern of signal energy in the time-frequency plane. Damage to the surrounding rock will cause energy to migrate to lower frequencies and subsequent time periods.

[0046] Introducing wavelet packet energy entropy This method quantizes the complexity and irregularity of signals at multiple resolutions. Wavelet packet decomposition decomposes the signal at different scales and frequency bands to obtain the energy of each node. The formula for calculating the wavelet packet energy entropy is: ; in, This represents the proportion of the energy of the i-th node in the total energy; It represents the wavelet packet energy entropy, which measures the degree of disorder in the distribution of signal energy at wavelet packet nodes; This represents the total number of wavelet packet nodes.

[0047] Wavelet packet energy entropy measures the degree of disorder in the distribution of signal energy at wavelet packet nodes and can effectively characterize the local transient characteristics of non-stationary signals.

[0048] (4) Extraction of nonlinear and damage-sensitive features The damage evolution of surrounding rock under blasting vibration is often accompanied by an enhancement of nonlinear dynamic response. Intact rock masses primarily exhibit linear elastic behavior under small strains, but as damage accumulates to a certain extent, nonlinear mechanisms such as the opening and closing of internal fissures and frictional slip begin to emerge, leading to stronger nonlinear characteristics in the vibration signal. Introducing nonlinear dynamic indices can sensitively capture this change.

[0049] The fractal dimension is calculated using box dimension or spectral dimension methods to characterize the complexity and roughness of the signal waveform; the maximum Lyapunov exponent is calculated. This reflects the system's sensitivity to initial conditions and its chaotic characteristics. Its expression is: ; in, The maximum Lyapunov exponent (in units of 1 / s) reflects the system's sensitivity to initial conditions; Indicates time; This represents the separation vector of adjacent orbits in phase space after time t; This represents the separation vector of adjacent orbits in the phase space at the initial time.

[0050] The Nonlinear Energy Ratio (NER) is defined as the proportion of nonlinear components in the total energy, and is calculated using the following formula: ; in, Indicates the nonlinear energy ratio; Represents the energy of the nonlinear component; Represents the linear component energy; Indicates the total energy of the signal; Represents the original signal; This represents the predicted value of the linear prediction model at time t.

[0051] By designing damage-sensitive features and comparing the trends of different blasts or measuring points, subtle changes in vibration signals caused by surrounding rock damage can be sensitively captured. Energy attenuation coefficient. The formula describing the decay rate of vibrational energy over time is: ; in, Indicates the energy decay coefficient; Indicates time; Indicates the initial envelope amplitude; This represents the amplitude of the signal envelope at time t.

[0052] The frequency drift rate Δf describes the relative change in the dominant frequency before and after damage, and is calculated using the following formula: ; in, Frequency drift rate; The dominant frequency after damage; The dominant frequency before damage.

[0053] By systematically extracting the above four types of features, a high-dimensional feature vector containing dozens or even hundreds of features can be constructed to comprehensively characterize the multi-dimensional characteristics of blasting vibration signals and their correlation with surrounding rock damage.

[0054] Step S4: Feature Dimensionality Reduction and Optimization When the feature dimensionality is too high, it leads to excessively long model training time, high computational resource consumption, easy overfitting, and poor model interpretability. By using feature selection and dimensionality reduction techniques, the most discriminative subset of features is selected from the hundreds of features extracted in step S3, thereby improving the generalization ability, computational efficiency, and interpretability of the machine learning model. This step achieves effective compression of the high-dimensional feature space while minimizing information loss.

[0055] First, feature correlation analysis and redundancy removal are performed. A comprehensive correlation analysis is conducted on the high-dimensional feature set constructed in step S3, and the Pearson correlation coefficients between features are calculated. This identifies and removes highly correlated redundant features. The formula for calculating the Pearson correlation coefficient is: ; in, This represents the Pearson correlation coefficient between variables X and Y; and Represents the sample mean of variables X and Y; and Let X and Y represent the values ​​of the i-th sample; Indicates the number of samples.

[0056] When the absolute value of the correlation coefficient between two features is greater than a preset threshold (e.g., 0.9), they are considered highly correlated, and only one needs to be retained. For example, peak acceleration and root mean square acceleration may be highly correlated, and peak velocity and peak displacement may also have a strong correlation; only one needs to be retained. This process can effectively reduce feature dimensionality, reduce the impact of multicollinearity on model stability, and retain the main variation information in the feature set.

[0057] Secondly, linear dimensionality reduction techniques are employed. For the remaining features, unsupervised linear dimensionality reduction methods are used, primarily Principal Component Analysis (PCA). PCA transforms potentially correlated features into a set of linearly uncorrelated principal components through orthogonal transformation, and then sorts them according to their variance contribution rate. Its mathematical expression is as follows: ; in, The new feature matrix after dimensionality reduction; Represents the original feature matrix; This represents the projection matrix.

[0058] The formula for calculating the cumulative variance contribution rate is: ; in, Let represent the variance of the i-th principal component; Indicates the original feature dimension; Indicates the number of principal components retained; Indicates the first The cumulative variance contribution rate of each principal component.

[0059] The top k principal components whose cumulative contribution rate reaches a predetermined threshold (e.g., 90%-95%) are retained as new feature inputs. This method can significantly reduce feature dimensionality while preserving the variation information of the original data to the maximum extent. In addition to PCA, other methods such as Linear Discriminant Analysis (LDA, a supervised dimensionality reduction method that finds the projection direction that maximizes inter-class divergence and minimizes intra-class divergence) or Independent Component Analysis (ICA) can also be used.

[0060] Another approach is supervised feature selection. When labeled data is sufficient, supervised feature selection methods, such as feature importance ranking based on tree models, can be used. By training a Random Forest (RF) or Extreme Gradient Boosting (XGBoost) model, the importance score of each feature in predicting damage labels can be calculated. For Random Forest, feature importance can be measured by the reduction in Gini impurity or the increase in out-of-bag (OOB) error. ; in, The importance score represents feature f; Represents a single-subject decision tree; This represents the number of decision trees in a random forest; This represents the reduction in Gini impurity caused by feature f during node n splitting.

[0061] Finally, within a deep learning framework, especially when using Convolutional Neural Networks (CNNs) or Recurrent Neural Networks (RNNs), manual feature engineering can be bypassed, allowing preprocessed vibration signals or time-frequency maps to be directly input into the network. CNNs automatically learn and extract local patterns and hierarchical features from the signal through their convolutional layers, pooling layers achieve feature dimensionality reduction and invariance, and fully connected layers perform high-level abstraction. These features are typically more discriminative and have lower dimensionality than manually designed features. The deep features extracted by CNNs (using the output of the penultimate layer or fully connected layers as feature representations) can be fused with selected key manually generated features to form the final feature representation, combining the automatic feature learning capabilities of deep learning with the interpretability advantages of manually generated features.

[0062] Step S5: Machine Learning Model Construction and Training By selecting appropriate machine learning algorithms, optimizing model parameters, and conducting rigorous performance verification, a model capable of accurately identifying the damage level of surrounding rock or predicting the continuous damage index is constructed. This embodiment transforms the features obtained from the preliminary processing into intelligent recognition capabilities with practical engineering application value.

[0063] First, the training set is prepared and the data is partitioned. The feature matrix, after feature dimensionality reduction and optimization, is combined with the corresponding surrounding rock damage labels (which can be discrete damage levels, such as "no damage," "minor damage," "moderate damage," and "severe damage," or continuous damage indices 0-1) to form a supervised learning dataset. The dataset is then cleaned to remove any missing or outlier values. Stratified sampling or time-series partitioning strategies are used to divide the dataset into a training set (typically 60%-70%), a validation set (15%-20%), and a test set (15%-20%). Stratified sampling ensures that the proportion of samples from different damage levels in each subset is consistent with the overall distribution, avoiding class distribution bias. For class imbalance problems (where the number of samples for some damage levels is much smaller than for others), oversampling techniques are used to increase minority class samples, undersampling techniques to reduce majority class samples, or class weight adjustment methods are used to ensure that the model does not favor the majority class, thereby improving the ability to identify the minority class (often the more important severe damage class).

[0064] Next, model selection and architecture design are crucial. Based on the characteristics of the problem, the size of the data, and the dimensionality of the features, a suitable machine learning model is chosen. For damage level classification problems, Support Vector Machines (SVM) can be selected (particularly suitable for small to medium-sized data, handling nonlinear problems through kernel tricks, commonly using RBF kernels or multinomial kernels), Random Forests (RF) (an ensemble learning method, robust, capable of handling high-dimensional data and nonlinear relationships, while providing feature importance assessment), Gradient Boosting Trees (such as XGBoost, LightGBM, and CatBoost, which perform excellently on structured data and have strong nonlinear fitting capabilities), and Deep Neural Networks (DNN) (multi-layer fully connected networks that learn complex nonlinear mapping relationships by stacking hidden layers). For continuous damage index prediction, corresponding regression versions of the models (SVR, RF regressors, XGBoost regressors, and DNN regression networks) can be used.

[0065] To address the sequential characteristics of vibration signals, a one-dimensional convolutional neural network (1D-CNN) can be used to automatically extract features directly from the original or preprocessed time series. 1D-CNN, by sliding convolutional kernels along the time axis, can capture hierarchical representations of local temporal patterns and features. LSTM, through its gating mechanism (input gate, forget gate, output gate), effectively solves the gradient vanishing problem of traditional RNNs, enabling it to learn the temporal evolution of vibration signals. Alternatively, a CNN-LSTM hybrid model can be constructed, first using CNN to extract local features, and then using LSTM to capture temporal dependencies.

[0066] Then, hyperparameter optimization and model training are performed. The model's hyperparameters are systematically tuned using strategies such as Grid Search, Random Search, or Bayesian Optimization.

[0067] During training, the training set is used to fit the model parameters, and performance metrics are monitored on the validation set. An early stopping strategy is employed: when the validation set performance no longer improves over several consecutive epochs (e.g., 10-20 epochs), training is stopped, and the model parameters are restored to the optimal level achieved during validation, preventing overfitting. Simultaneously, regularization techniques (L1 and L2 regularization) are applied to constrain the model parameters, and Dropout (for neural networks, randomly discarding a certain percentage of neurons during training) is used to improve the model's generalization ability, enabling it to maintain good performance even on unseen data.

[0068] Next, model validation and performance evaluation are performed. The performance of the final model is evaluated on an independent test set, which was not used at all during the entire training and parameter tuning process, and can objectively reflect the actual generalization ability of the model. For classification problems, a comprehensive evaluation is performed using multiple metrics, including accuracy (the proportion of correctly classified samples out of the total samples), precision (the proportion of true positive samples among those predicted as positive), recall (the proportion of true positive samples correctly predicted), F1 score (the harmonic mean of precision and recall), and confusion matrix. The formula for calculating the F1 score is: ; in, Indicates accuracy; Indicates recall rate; express The score is the harmonic mean of precision and recall.

[0069] A confusion matrix can visually display the classification status of each category, facilitating the analysis of which categories the model is prone to confusion in. For regression problems, indicators such as root mean square error (RMSE), mean absolute error (MAE), and coefficient of determination (R²) are used. ; ; in, Indicates the number of samples; This represents the true value of sample i; This represents the predicted value for sample i; Indicates the coefficient of determination; This represents the average of the true values.

[0070] Finally, model integration and deployment preparation are carried out. To improve the robustness and accuracy of predictions, model ensemble strategies can be considered. Voting is suitable for classification problems, where the predictions of multiple base models are voted on; stacking uses the predictions of multiple base models as new features to train a meta-model for final prediction; weighted averaging is suitable for regression problems, assigning different weights to each model based on its performance on the validation set and then calculating a weighted average. Ensemble learning can combine the advantages of multiple models, reduce the bias and variance of a single model, and improve overall prediction performance.

[0071] The finalized model and its associated preprocessing parameters are serialized and saved to form a complete inference pipeline. This standardized packaging process ensures that the model can be easily deployed to the production environment for real-time or batch processing of surrounding rock damage identification tasks.

[0072] Step S6: Real-time damage identification and intelligent feedback control The blasting vibration signals acquired in real time on-site undergo a preprocessing process identical to that used in the training phase (filtering and noise reduction, amplitude normalization, signal segmentation), feature extraction process (extracting time-domain, frequency-domain, time-frequency-domain, and nonlinear features), and dimensionality reduction process (applying saved feature scalers and dimensionality reduction matrices) to obtain feature vectors. These feature vectors are then input into a pre-trained machine learning model for forward inference. The model outputs the surrounding rock damage state caused by the current blasting cycle. This can be a damage level given in a classification format (e.g., outputting a 4-category probability distribution, taking the category corresponding to the highest probability) or a continuous damage index D∈[0,1] given in a regression format. The entire inference process should be completed within a few seconds to tens of seconds to meet real-time requirements. This process achieves end-to-end intelligent inversion from vibration signals to damage states, providing immediate and quantitative basis for subsequent decision-making.

[0073] Next, confidence assessment and uncertainty quantification are performed. For each model's prediction, its prediction confidence or uncertainty measure is output. For classification models, the predicted probability distribution can be output, and the maximum probability value is taken as the confidence level. When the maximum probability value is lower than a preset threshold, the model is considered insufficiently certain about the prediction. For regression models, Bayesian neural networks or ensemble learning methods can be used to estimate the uncertainty range of the predicted values. The prediction uncertainty of ensemble models can be quantified as follows: ; in, This indicates the number of base models in the ensemble model; This represents the predicted value of the m-th base model; This represents the average of the predictions from all base models; A quantitative indicator representing the uncertainty of prediction.

[0074] When the confidence level is lower than the preset threshold or the uncertainty is too high, the system automatically triggers the manual review process. The on-site engineer makes a comprehensive judgment by combining other monitoring data (such as displacement monitoring, stress monitoring, and video surveillance) and engineering experience to ensure the reliability of decision-making in high-risk situations and prevent engineering accidents caused by model misjudgment.

[0075] Then, intelligent decision-making and feedback control are implemented. Based on the damage level or damage index D and its confidence level output by the model, combined with the preset engineering safety threshold, the system automatically executes the corresponding decision logic. The decision rules can be designed as a hierarchical response mechanism. Based on the assessment, the decision was made to immediately halt work and initiate a detailed inspection and reinforcement process.

[0076] The core of blasting parameter adjustment is the optimization of the single-stage charge Q. Based on the damage identification results, the single-stage charge for the next cycle is adjusted using the following formula: ; in, Indicates the current blasting parameters; This indicates the adjusted blasting parameters; This represents the damage index predicted by the model; Indicates the damage threshold; Indicates the maximum permissible damage index; This represents the adjustment factor.

[0077] The larger the value, the greater the adjustment force. The formula's design logic is: when the predicted damage index D exceeds the threshold... When the dosage exceeds the limit, reduce the dosage of each segment proportionally; the greater the excess, the greater the reduction; when D is significantly lower than... hour, ≈ The formula can remain unchanged or even be appropriately increased (by adding an incremental term to the formula) to improve construction efficiency.

[0078] This decision-making process achieves an automated closed loop from identification to control, transforming passive post-event analysis into proactive real-time control, optimizing construction efficiency while ensuring the stability of the surrounding rock. The adjustment suggestions generated by the system can be pushed to blasting engineers and construction personnel in real time via mobile terminals or on-site displays, and implemented after review and confirmation by the engineers.

[0079] The system then continuously records the vibration signals, feature extraction results, model predictions, decision-making actions, and subsequent manual verification results for each blast. This newly generated labeled data is added to the training data pool, and the model is periodically subjected to incremental learning or fine-tuning. This online learning mechanism enables the model to adapt to changes in different construction stages and geological conditions, continuously improving its prediction accuracy and adaptability in specific engineering scenarios, forming a constantly evolving intelligent system.

[0080] Finally, the entire damage identification and feedback control system is integrated into the engineering management platform, providing an intuitive visual interface that displays vibration signals, damage identification results, confidence level assessments, and system decision recommendations at each monitoring point in real time. It also generates detailed analysis reports, including historical trend analysis and damage evolution patterns, providing comprehensive decision support for engineering managers and achieving transparent and digital management of the construction process.

[0081] The technical features of the above embodiments can be combined arbitrarily. For the sake of brevity, not all possible combinations of the technical features in the above embodiments are described; only preferred embodiments of the present invention are illustrated. The descriptions are relatively specific and detailed, but they should not be construed as limiting the scope of the present invention. As long as the combination of these technical features does not contradict each other, it should be considered within the scope of this specification.

[0082] It should be noted that those skilled in the art can make various modifications and improvements without departing from the inventive concept, and these all fall within the scope of protection of this invention. Therefore, the scope of protection of this invention should be determined by the appended claims.

Claims

1. A method for inverting blasting vibration signals and intelligently identifying surrounding rock damage, characterized in that: The method comprises the following steps: Step S1: arranging sensors at key positions of surrounding rock to collect blasting vibration signals and calibrating surrounding rock damage index; Step S2: pre-processing the blasting vibration signals to extract multi-dimensional features; Step S3: after dimension reduction optimization of the multi-dimensional features, constructing a training set with the surrounding rock damage index and training a machine learning model; Step S4: inputting real-time collected blasting vibration signals into the machine learning model to identify the surrounding rock damage state, and adjusting the blasting parameters according to the identification result.

2. The method according to claim 1, wherein the method is characterized in that: The method of step S1 comprises: arranging at least one three-way acceleration sensor at at least one position of the working face, the arch, the wall and the floor to synchronously collect the blasting vibration signals; The method of calibrating the surrounding rock damage index comprises: quantifying the surrounding rock damage degree by at least one non-destructive testing method and normalizing it into a damage index as label data.

3. The method according to claim 1, characterized in that: The pre-processing in step S2 comprises filtering and noise reduction, amplitude normalization and signal segmentation; The filtering and noise reduction comprises: applying a band-pass filter to the blasting vibration signals, and applying wavelet threshold denoising when the signal quality is lower than a threshold value; The calculation formula of the amplitude normalization is: ； wherein, is the original signal data point, is the signal sequence to be normalized, is the normalized data point; The signal segmentation comprises: determining the blasting zero time based on a trigger signal synchronized with the initiation system, and intercepting the signals in a fixed window after the blasting zero time as effective analysis segments.

4. The method according to claim 1, characterized in that: The multi-dimensional features comprise time domain features, frequency domain features, time-frequency domain features and nonlinear features; The time domain features comprise peak vibration velocity, root mean square amplitude, signal envelope total energy and energy concentration degree; The frequency domain features comprise main frequency, center frequency, frequency band energy proportion and spectral entropy; The time-frequency domain features comprise instantaneous energy distribution, time-frequency centroid trajectory and wavelet packet energy entropy; The nonlinear features comprise fractal dimension, Lyapunov exponent, nonlinear energy ratio, energy attenuation coefficient and frequency drift rate.

5. The method according to claim 4, characterized in that: the peak vibration velocity The formula for calculating the peak vibration velocity is: ； wherein denotes the vibration velocity signal, denotes time; The calculation method of the signal envelope total energy is: the instantaneous amplitude envelope of the signal is calculated by using Hilbert transform The expression is: ； wherein, is the pre-processed vibration signal, is the signal is the result of the Hilbert transform, is the instantaneous amplitude envelope of the signal.

6. The method according to claim 4, characterized in that: The calculation formula of the time-frequency centroid trajectory is: ； ； wherein, is a time-frequency representation of a signal; denotes the magnitude of the time-frequency representation; is a frequency; denotes a time; is a time at which the frequency centroid is located; is a frequency at which the time centroid is located; The calculation formula of the wavelet packet energy entropy is: ； wherein, represents the proportion of the energy of the i-th node in the total energy; i represents the wavelet packet energy entropy; represents the total number of wavelet packet nodes.​ 7. The method according to claim 4, characterized in that: The Lyapunov exponent The expression for the Lyapunov exponent is: ； wherein denotes time; is the time separation vector of neighboring orbits in the phase space after time is the separation vector of neighboring orbits in the phase space at initial time The non-linear energy ratio The formula for calculating the non-linear energy ratio is: ； wherein is the original signal; is the predicted value of the linear prediction model at time t; The energy attenuation coefficient The calculation formula is: ； ； wherein represents the initial envelope amplitude; is the signal envelope amplitude at time t; The frequency drift rate The formula for calculating the frequency drift rate is: ； wherein is the dominant frequency after the lesion; is the dominant frequency before the lesion.

8. The method according to claim 1, characterized in that: The machine learning model comprises at least one of a support vector machine model, a random forest model, an extreme gradient boosting model, a deep neural network model, a one-dimensional convolutional neural network model and a long short-term memory network model; The process of training the machine learning model comprises: dividing the training set into a training subset, a validation subset and a test subset, fitting model parameters on the training subset, optimizing hyperparameters on the validation subset through grid search, random search or Bayesian optimization, and evaluating model performance on the test subset.

9. The method according to claim 8, characterized in that: When the machine learning model is used for surrounding rock damage grade classification, the evaluation indicators of the model performance comprise accuracy, precision, recall, F1 score and confusion matrix; When the machine learning model is used for continuous damage index prediction, the evaluation indicators comprise root mean square error (RMSE) and determination coefficient (R²).

10. The method according to claim 1, characterized in that: Step S4 comprises: Step S41: inputting the feature vector to be identified into the machine learning model, outputting the surrounding rock damage grade or the surrounding rock damage index, and calculating the prediction confidence; Step S42: when the prediction confidence is lower than a preset threshold, triggering an artificial review process; Step S43: When the surrounding rock damage level is greater than the preset threshold, or the surrounding rock damage index exceeds the safety threshold, the single-section charge amount of the next cycle is adjusted according to the damage identification result by the following formula: ； wherein, is the current single segment dose; is the adjusted single segment dose; is the model predicted damage index; is the damage threshold; is the maximum allowed damage index; is the adjustment factor; Step S44: The real-time collected blasting vibration signal, the extracted multi-dimensional features and the corresponding damage identification result are added to the training data set, and the machine learning model is updated online every several blasting cycles.

Images