Abstract
Based on the same blasting vibration signal, the S transform and HHT method have been applied for the time and frequency analysis. Both the S transform and HHT method are analyzed and compared from the point of decomposition path and analysis function. It expressed that these two methods are complete and effective and have the characteristics of adaptive multiresolution. Especially, the time resolution of the HHT method is invariable and its accuracy is high. Moreover, its frequency resolution can be adaptively adjusted with the inherent characteristics of the signal. However, the time and frequency resolution of the S transform are limited by the uncertainty principle, and the resolution improvement is at the cost of reducing the resolution of another. Furthermore, the HHT method has much more adaptability and is superior to S transform in the analysis of nonstationary signals. In advance, the waveform curve and power spectral density of each detonator are extracted, which indicates that it is not scientific to use different segment detonators for the same circle of the auxiliary hole. Under the premise of the same charge quantity, the vibration energy peak value produced by the same segment detonator is less than that of several segments of detonators. Therefore, it is believed that the HHT method has broad application prospects in the study of wave propagation, attenuation law, structural dynamic response characteristics, and blasting vibration damage in rocks.
1. Introduction
Blasting vibration is an unavoidable hazard to the surrounding environment in the blasting excavation process. It is of significant importance to study the blasting effect caused by blasting vibration [1–5]. Sadowski formula is mainly used for regression analysis of blasting vibration signals. However, the relationship between frequency and energy distribution of blasting waves is less [6–9]. Meanwhile, through the engineering practice, it is shown that if the vibration velocity of a single particle is used as the sole index to measure the blasting vibration intensity, the actual damage of buildings cannot be well reflected in many cases. The other vital influence factors are not considered, such as the frequency, duration, and energy of blasting vibration [10–13].
Blasting vibration signals could be analyzed and processed by various digital signal processing methods, and the time and frequency can be extracted. Then the main characteristic information of the signals can be comprehensively analyzed by time-frequency processing methods. Those have been widely used now, such as the traditional Fourier transform, short-time Fourier transform, and wavelet transform [14–17].
The Fourier transform is limited to the local frequency domain and could not be analyzed from the perspective of the time domain. Short-time Fourier transform (STFT), which is based on Fourier transform, is a commonly used time-frequency analysis method. The time window is loaded onto the time signal first for the STFT, and then Fourier transform on the data signal of the time window is performed. Then, the short-time spectrum of the signal could be obtained from the angle of frequency and time.
Wavelet analysis is a method to process the local time-frequency characteristics of data signals. The window area could be kept unchanged, and the window shape can be changed gradually with the analysis situation. Both the time window and frequency window could be transformed. S transform adopts Gauss window function, and the window width is proportional to the reciprocal of the frequency, which avoids the selection of window function and improves the defect of fixed window width, and the phase spectrum of each frequency component in time-frequency representation keeps direct contact with the original signal.
A new method for analyzing and processing the time-frequency distribution of nonstationary data signal is proposed, which is HHT (Hilbert–Huang Transform) method [18, 19]. It consists of two parts, i.e., Huang transform and Hilbert spectrum analysis. For the Huang transform, the EMD (Empirical Mode Decomposition) method is used to decompose the signal into a finite number of IMFs (Intrinsic Mode Function). It is considered that any complex function is composed of some IMF components, which are uncorrelated, simple, and nonsinusoidal functions. IMF components could be extracted by EMD decomposition. The original vibration signal is composed of n-order IMF components and residuals after decomposition [20, 21]. For the Hilbert spectrum analysis, it is to transform the Hilbert spectrum with a finite number of IMF components obtained by EMD so as to obtain the distribution spectrum of energy in the time-frequency plane, namely, the Hilbert spectrum. Energy and instantaneous frequency, rather than the frequency and energy of the whole signal, are obtained by the Fourier transform. By integrating the Hilbert spectrum, the energy spectrum could be well obtained, which can be used to analyze the energy distribution of vibration signals varied with time and frequency. It breaks the limitation of the uncertainty principle and can accurately express all kinds of information on the time-frequency plane [22, 23].
The basic properties of those three different transforms are compared in Table 1.
From Table 1, it can be shown that the basic functions of the three transform methods are different in signal analysis. Correspondingly, the decomposition method of the signal is different, and the result is different. Therefore, it cannot be made a unified analysis for the signal. The sine or cosine function is adopted as the basic function in Fourier transform, while the basis function of wavelet transform is preset and could not be changed with the parameters of the wave itself, so the results are quite different when choosing various wavelet basic functions. It means that the signal analysis has to be met some difficulties. For the HHT method, it is based on the data signal itself to decompose. Moreover, there is no choice to select any predefined function so that better analysis results can be obtained by the EMD decomposition. The IMF components obtained by EMD decomposition are stationary and have stronger local characteristics than those obtained by Fourier transform and wavelet transform. HHT transform method has strong automatic adaptability, high decomposition efficiency, and better suitability to sudden change and fast attenuation of data signals.
The time-frequency processing methods of signal mainly include wavelet transform and EMD. However, wavelet transform is based on the analysis of specified wavelet basis. Once the wavelet basis is specified, the decomposed mode is fixed. For different signal analysis, it is better to use different wavelet basis to achieve a better suitable processing effect. The EMD method can just solve this problem. EMD is equivalent to an adaptive orthogonal basis time-frequency signal processing method. For the processing of an unknown nonlinear and nonstationary signal, the decomposition can be started directly without prior analysis. According to some certain fixed modes and levels, without human intervention and setting, the EMD method would be automatically divided.
Once there is an intermittent signal in the actual signal, frequency aliasing would occur by the EMD decomposition method. It is shown that there are multiple scale components in an IMF component or a scale component exists in multiple IMF components. To solve this issue, the noise assistant processing method, EEMD (Ensemble Empirical Mode Decomposition), and CEEMD (Complementary Ensemble Empirical Mode Decomposition) are proposed. EEMD is to add white noise to the original signal to change the distribution of the extreme point of the signal. In contrast, CEEMD adds a set of noise signals to the original signal to change the distribution of the extreme point of the signal.
When the same number and amplitude of white noise are added, the residual noise of EEMD will decrease with the average number of integration. The residual noise of CEEMD is maintained to a small extent, regardless of the average number of integration. In a sense, EEMD and CEEMD can save calculation time when the interference of small residual noise is guaranteed.
Compared with EMD, the EEMD algorithm can overcome the phenomenon of mode aliasing and reconstruct the signal better, but the residual white noise brings noise interference to the signal. CEEMD algorithm can not only overcome the phenomenon of mode aliasing but also could accurately reconstruct the original signal without energy loss. Therefore, compared with EMD, the white noise added in the decomposition process by those two methods can partially eliminate the noise, but it can not completely offset it.
S transform is a further optimization of the wavelet transform, which reflects that the frequency is related to the resolution in time and frequency domain, and the local phase information and the time and frequency resolution improvement could be significantly reflected.
The relevant theoretical analysis, comparison, and summary of rock blasting vibration signal are carried out here. Wavelet S transform and HHT transform methods are selected to effectively extract the time-frequency and energy distribution characteristics of the blasting vibration signal. Meanwhile, the waveform curve and power spectral density of each detonator are extracted. Especially, for the same circle of the auxiliary hole, those two transform methods would be applied to identify the blasting vibration effect, which is of great significance to guide the blasting design.
2. Blasting Vibration Signal
Based on mastering the influencing factors of blasting vibration, the blasting vibration signals in the test are analyzed. Figure 1 shows the layout of the blasting holes in this section. The excavation section area is 12.15 m2, and the cutting method is a straight and parallel hole cutting way [24]. The cutting hole depth is 2.2 m, and the depth of the other holes is 2.0 m. In order to ensure the cutting quality, six holes are evenly distributed around the empty hole. The more holes and the bigger holes diameter are, the more advantageous it is to cut blasting because much more auxiliary free surface and compensation space could be provided. Full-face blasting excavation is adopted here. The millisecond detonators are used, i.e., MS1, MS3, MS5, MS7, MS8, MS9, MS10, MS11, MS13, MS15, and MS19, respectively. The emulsion explosive is used and the weight of each roll is 300 g. The charge quantity and other blasting parameters are shown in Table 2 [24].
The intuitive analysis method of blasting vibration signal is to directly analyze the measured waveform and determine the characteristic quantity of blasting vibration from the waveform diagram. Blasting vibration waveform can be analyzed intuitively using a blasting vibration instrument system, and three parameters of vibration amplitude, main frequency, and duration could be obtained.
The frequency of vibration wave is related to the rock property. Generally, the harder the rock medium is, the higher the frequency of the seismic wave is. And the softer the medium is, the lower the frequency is. For important safety projects involving blasting vibration, suitable parameters of blasting seismic wave should be determined by actual observation through a small-scale blasting test.
Referenced values for blasting protection of buildings have been proposed in criteria for the safety of blasting vibration (GB6722-2014) [25], but some problems should be discussed. Firstly, besides the characteristics of blasting vibration motion itself, the factors determining the dynamic response of building structures to blasting vibration are closely related to the inherent dynamic response of structures. It has an obvious selective amplification effect on blasting seismic response. The blasting vibration component consistent with the natural frequency of the structure would be amplified to the greatest extent, in which the response of the structure to blasting vibration motion is directly determined.
The blasting vibration signal is selected here to analyze the typical characteristics, as shown in Figure 2.
3. Theory of S Transform and HHT Analysis
3.1. Theory of S Transform
S transform is a time-frequency reversible analysis method of time-window Fourier transform proposed by Stockwell et al. It is an extension of continuous wavelet transform based on Morlet wavelet.
The continuous wavelet transform of a function h(t) can be set as follows:
where and are translation factor and scaling factor, respectively.
If the continuous wavelet transform is multiplied by a phase term, that is,
then, its wavelet generating function is
Since the allowable wavelet function in formula (3) could not satisfy the zero-average condition, it means that formula (2) is not a continuous wavelet transform in a strict sense. Specifically, the continuous S-transform should be defined as follows:
It is worth mentioning that the continuous S transform can also be derived from the time-window Fourier transform, whose positive transform is as follows:
If the normalized Gaussian window function is used in the time-window function and the scaling and translation t are carried out, it could beThen,
Formula (7) is further simplified, order , and then the formula (4) of continuous S transform could be obtained.
As can be seen from the above, a Gaussian window function with variable width has been adopted in the S transform, and the window width is in direct proportion to the reciprocal of the frequency. Therefore, the time window in the low-frequency band is wider and the high-frequency resolution can be obtained, while the time window in the high-frequency band is narrower, which a high time resolution could be obtained. That is to say, it has self-adaptability to some extent.
3.2. Theory of HHT Analysis
HHT analysis method is mainly composed of EMD and Hilbert transform. In order to study transient and nonstationary phenomena, the frequency must be a function of time. IF (Instantaneous Frequency) is an intuitive and basic physical concept in the HHT method. Significantly, EMD is proposed by Huang N E et al. and the EMD decomposition is the core of HHT method.
The original signal can be composed of n-order IMF components and residuals :
Signals can be expressed as contours in the time-frequency plane. The amplitude distribution on the processed time-frequency plane is called the Hilbert spectrum.
Based on the time scale (time difference between adjacent peak points) characteristics of the signal itself, the EMD method decomposes the signal into a set of IMFs with different time scales, which could satisfy two kinds of definitions. Each IMF can be considered as a modal function inherent in the signal.
At the same time, IMF can be arbitrarily combined to highlight certain characteristics and achieve the analysis purpose. This filter is based on the IMF component and has self-adaptability. It fully retains the inherent nonlinear and nonstationary characteristics of the signal itself. When the noise is one or more IMF components, the filter is directly used to eliminate it.
4. S Transform and HHT Analysis
4.1. S Transform Analysis of Blasting Vibration Signal
The time-frequency distribution of each detonator after S transform is shown in Figure 3, such as the frequency range and magnitude in the time-frequency diagram.
The energy distribution pattern and frequency distribution range of each detonator through S transform could be clearly revealed. Meanwhile, the waveform curve and power spectral density of each detonator are extracted, as shown in Figure 4.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)
(j)
(k)
The following can be seen from the waveform and spectrum curve of each detonator signal:(1)The interval time of millisecond detonator has a great influence on the spectrum distribution of blasting vibration signals, and the peak power spectrum of signal with large vibration velocity is not necessarily large. Although there are differences in waveforms between kinds of detonator segments, the frequency range covered by the power spectrum curve is approximately the same, indicating that the blasting vibration signal belongs to the shock wave signal with the broad-frequency domain.(2)The frequency of full-face millisecond blasting is in the range of 0 ∼ 1000 Hz under the geological condition, and the main frequency is 0 ∼ 200 Hz, while the central frequency is above 100 Hz, which is consistent with the frequency range of the main vibration phase monitored by measuring instrument.(3)The main frequency and the duration time of the main frequency are different for kinds of signals with various interval times. The phase difference of superimposed signals with different interval times makes the signals form new peaks or valleys, and the interval time between two adjacent peaks or valleys has been changed. Then the energy distribution in different frequency bands has been inevitably changed with the main frequency and duration time of blasting vibration signal, which will be clarified in the follow-up analysis.
The vibration components response of different frequency bands can be obtained, as the vibration components intensity of the corresponding frequency band can be directly reflected by the frequency band energy. Moreover, the energy duration is related to the accumulation effect of low-cycle fatigue damage and the nonelastic failure of the structure, so the normalized frequency and the duration of frequency band are taken as indicators to study the energy characteristics of blasting vibration signals. Time-frequency distributions of energy through the S transform for each detonator are shown in Figure 5.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)
(j)
(k)
The following can be seen:(1)The delay interval has a great influence on the total energy of the signal. The energy of the blasting vibration signal with a large peak particle velocity must be large.(2)The energy duration at different frequencies is various. Under this test condition, the energy duration in the frequency band decreases gradually from low frequency to high frequency.(3)The charge of MS7, MS8, MS9, MS10, and MS11 detonator segments are 8.1 kg, 3.6 kg, 5.4 kg, 4.5 kg, and 4.5 kg, respectively, while the corresponding energy peaks are 0.35, 0.35, 0.35, 0.37, and 0.38, respectively. In contrast, the MS8 detonator charged the smallest explosive with a relatively high peak energy. Therefore, the selection of the MS8 detonator must be careful.(4)From the time-frequency spectrum of detonator vibration signals, it can be seen that the frequency distribution range is wide. The influence of detonator segments on vibration signal is more reflected in the time-frequency amplitude. The energy of the MS11 detonator is the largest, even more than the value of the MS1 detonator in the cutting hole, while the charge is only 50% of its charge, which indicates that it is not scientific to use different segment detonators for the same circle of the auxiliary hole. Moreover, the vibration reduction is also invalid. This also shows that under the premise of the same charge quantity, the vibration energy peak value produced by the same segment detonator is less than that of several segments of detonators.
It is not suitable for the analysis of nonstationary vibration signals by Fourier transform. Wavelet and wavelet packet analysis uses local basic quantities and wavelet basic functions as analysis tools. It can be more intuitive and effective for the analysis of nonstationary signals, but it has high requirements for selecting wavelet basic functions, which means the adaptability is not strong enough. Compared with the traditional blasting signal analysis theory, the HHT analysis method is a newer method. Especially, its unique time-frequency analysis ability could meet the requirement of nonstationary signals such as blasting vibration. It differs from the Fourier transform in that the complex harmonic components are not needed to reveal the information contained in the signal. Meanwhile, it differs from wavelet and wavelet packet methods in that it does not depend on the selection of wavelet bases and has better adaptability.
4.2. HHT Analysis of Blasting Vibration Signal
The following is the filtering and denoising effect of blasting vibration signal based on EMD method. The sampling frequency is 8 kHz and the sampling time is 2.2 s. Firstly, the original signal is decomposed by EMD, and then the IMF components c1 ∼ c19 are obtained, as shown in Figure 6. It can be seen that after EMD decomposition, each IMF component shows a certain scale mode in the time domain, and each component is arranged in order of its specific frequency from high to low, which reflects the multiresolution characteristics and adaptability of the EMD method.
Any signal can be decomposed into several IMF, which is called empirical mode decomposition (EMD). The specific decomposition method is as follows: firstly, for any signal to be analyzed, all the maximum points and all the minimum points are found, and then two H-Spline curves are used to fit these extreme points to get the upper and lower envelope lines. All the signal data are between the two envelope lines, and then the average value of the upper and lower envelope lines could be obtained. The difference between the signal value and the mean value is recorded as the first component. For each IMF component, the Hilbert spectrum, Hilbert transform, and corresponding analytic signal are obtained.
It can be seen from Figure 6 that the IMF component of the blasting vibration signal decomposed by EMD is very rich. After signal decomposition, to determine the main information components contained in the signal and avoid the loss of effective information, it is necessary to determine the dominant components of the signal so as to realize the signal reconstruction process. The cross-correlation coefficient is used here to objectively evaluate the correlation between the IMF component and original signal, as shown in Table 3.
It can be seen from Table 3 that the correlation coefficient between IMF1∼IMF6 components and the original signal is much larger than the other order IMF components so it can be determined as the dominant component of the signal. Then the above six order IMF components are selected to synthesize the blasting vibration signal. The reconstructed signal and the separated self-similar high-frequency noise from the signal IMF1 to IMF6 components are illustrated in Figures 7 and 8, respectively.
Most of the effective information of the signal has been contained in the filtered signal, and it can be determined as the characteristic signal, which provides the premise for further analysis.
5. Comparison and Discussion of Blasting Vibration Signal Analysis Methods
As the inheritance and development of wavelet transform and short-time Fourier transform, S transform uses Gauss window function whose frequency is proportional to the reciprocal of window width, eliminating the selection of window function and improving the defect of fixed window width. Furthermore, in the time-frequency representation, the phase spectrum of each frequency component keeps direct contact with the original signal, inducing more features that can be obtained in the analysis. Meanwhile, the feature extracted by the S transform is insensitive to noise information.
The special feature of the S transform is that it not only keeps the direct relationship with the Fourier transform but also has different resolutions at different frequencies. In addition, the S transform is closely related to wavelet transform and can be regarded as the deformation of the continuous wavelet transform.
Both S transform and HHT transform can be well applied to the analysis of blasting vibration signals with nonlinearity and nonstationarity, and the main characteristic information of the time-history curve can be extracted.
Based on different time scales, IMF is directly decomposed from the original signal by the EMD method, and it can reflect the inherent characteristics of the original signal itself. Moreover, the IMF component is usually of practical physical meaning. EMD method does not need a preset basic function and has the characteristics of self-adaptation and high efficiency. However, the S transform is decomposed by the basic function, and its decomposition results are closely related to the selected basic function, so it is difficult to select the suitable basic function.
The energy distribution of the nonstationary signals in the time-frequency spectrum could be described by the S transform and Hilbert transform. It should be noticed that as restricted by the Heisenberg uncertainty principle, S transform results in spectrum diffusion and energy leakage, which makes the energy distribution with frequency wider and just a qualitative representation. Hilbert energy spectrum can clearly show the specific distribution of energy at any time and frequency, and most of the energy is concentrated on the limited energy spectrum lines.
S transform and Hilbert transform have the characteristics of adaptive multiresolution. The time resolution of the HHT method is invariable and its accuracy is very high. Meanwhile, its frequency resolution can be adaptively adjusted with the inherent characteristics of the signal. However, the time and frequency resolution of the S transform are limited by the uncertainty principle, and the resolution improvement is at the cost of reducing the resolution of another.
Both the S transform and Hilbert transform have a good effect in nonstationary signal denoising and filtering. Obviously, compared with the S transform, the EMD-based denoising and filtering method is much more flexible and adaptable.
Blasting vibration signal analysis is one of the main techniques to study the blasting vibration effect. Due to the interference of signal detection background and the influence of signal form and processing method, noise exists widely in the signal. Therefore, signal filtering and denoising are an important part of the signal process and detection. It is required that the analyzed signal should be linear, periodic, or stationary when using Fourier transform. Because the time-frequency characteristics of the signal cannot be extracted at the same time, the denoising effect of the nonstationary signal, which may contain many peaks or sudden changes, is not good. In contrast, wavelet analysis can be used for time-frequency analysis at the same time. Moreover, it has the characteristics of time-frequency localization and flexible selection of filter basis. It can effectively distinguish abrupt parts and noises in nonstationary signals. The signal-to-noise ratio could be improved while maintaining a good resolution of abrupt information, but it has not yet gotten rid of the defects based on the Fourier transform.
Furthermore, the EMD method can decompose a time series signal into a set of intrinsic mode functions with different scales, which has the characteristics of multiresolution and self-adaptation. The main purpose of the EMD method is to get IMF. Moreover, the filter and denoise could also be applied. That is, starting from the scale characteristics of the signal itself, the decomposed IMFs are combined to form high-pass, low-pass, and band-pass filters to eliminate noise or carry out the specific analysis. It is better than the wavelet method for filtering and denoising.
6. Conclusions
Through the above analysis and comparison, the following are indicated: (1)It is not scientific to use different segment detonators for the same circle of the auxiliary hole. Under the premise of the same charge quantity, the vibration energy peak value produced by the same segment detonator is less than that of several segments of detonators.(2)Both the S transform and HHT method have the characteristics of adaptive multiresolution. EMD method has the characteristics of self-adaptation and high efficiency; its frequency resolution can be adaptively adjusted with the inherent characteristics of the signal. In contrast, S transform is decomposed by the basic function. Furthermore, the time and frequency resolution of the S transform are limited by the uncertainty principle, and the resolution improvement is at the cost of reducing the resolution of another.
It is believed that the HHT method and S transform have broad application prospects in the study of wave propagation, attenuation law, structural dynamic response characteristics, and blasting vibration damage in rocks.
Data Availability
The data that support the findings of this study are available from the corresponding author upon reasonable request.
Conflicts of Interest
The authors wish to confirm that there are no known conflicts of interest associated with this publication
Acknowledgments
This work was supported by the National Natural Science Foundation of China (Grant Nos. 51522903 and 51774184), Excellent Project Fund in North China University of Technology (Grant No. 216051360020XN199/006), and Scientific Research Fund in North China University of Technology (Grant No. 110051360002).