Abstract
The limitation associated with Fourier transform and wavelet analysis that they often fail to produce satisfactory resolution simultaneously in time and frequency when dealing with nonlinear and nonstationary signals is frequently encountered. Therefore, this paper aims at using the HHT (Hilbert–Huang transform) method, which is built on the basis of the EMD- (empirical mode decomposition-) based wavelet threshold denoising technique and the Hilbert transform, to analyze the blasting vibration signals in a south China lead-zinc mine. The analysis is conducted in terms of three-dimensional Hilbert spectrum, marginal spectrum, and instantaneous energy spectrum. The results indicate that the frequencies of the blasting vibration signals lie mainly within 0∼200 Hz, which consists of more than 90% of the total signal energy. At the onset of the blasting, the vibration frequency tends to be low, with the frequency that is less than 50 Hz being dominant. By using instantaneous energy spectrum, which can reveal the condition of energy release for detonator explosion, the initiation moments of detonators with 7 time-lag levels are accurately identified. This accurate identification demonstrates the superiority of the HHT method in coping with nonlinear and nonstationary signals. Additionally, the HHT method that is characterized by adaptivity, completeness, strong reconfigurability, and high accuracy provides an opportunity for reflecting signals’ change features with regard to time domain, frequency domain, and energy irrespective of the limitation of the Heisenberg uncertainty principle.
1. Introduction
The blasting technique has been widely applied to engineering constructions of railway, mining, tunneling, and so on [1, 2]. During blasting, the explosion energy will apply work to surrounding rock mass, and then, part of them propagates away in the form of wave [3–5]. Considering that abundant information about blasting vibration is included in the blasting vibration wave, the extraction and analysis of the blasting vibration signals provide effective tools that can be used in studying the blasting efficiency.
Since Fourier published the theory of heat conduction analysis in 1822, the Fourier transform has been extensively used as a tool of analysis in the signal processing field [6, 7]. This tool transforms signals from the time domain to the frequency domain and describes the variation of signals in the frequency domain by using the overall frequency components that are included in the signals. Because of the inability to account for the instantaneous variation in some signal frequency, the tool possesses limitations in dealing with nonlinear and nonstationary signals. Under the restriction of the Heisenberg uncertainty principle [8, 9], this tool also fails to generate satisfactory resolution simultaneously in time and frequency [10].
On the contrary, since the introduction of the concept of wavelet in the 1980s by Morlet, a French geophysicist, the wavelet analysis theory has begun being gradually established [11]. At present, wavelet transform is serving as a tool that is extensively used in the analysis of nonstationary signals. In the effort to analyze vibration signals by using the wavelet technique, Newland [12] extended the engineering applications of this technique. Based on the usage of Mexican hat wavelets, Zhou and Adeli [13] developed a method for the analysis of time-frequency signals in earthquake records. Suárez and Montejo [14] proposed a wavelet-based procedure to produce an accelerogram with a response spectrum that was compatible with the target spectrum. This procedure contributes to the application of wavelet technique in the seismic wave analysis. Zhong et al. [15] analyzed the wavelet packet energy spectra for blasting vibration signals through utilizing the technique of the wavelet packet analysis. They also investigated the characteristics of attenuation of blasting vibration wave for various explosion parameters.
However, wavelet transform is still not free from the Heisenberg uncertainty principle. Under certain scales, wavelet transform cannot achieve high accuracy at both time and frequency. In addition, the wavelet basis is difficult to choose [16].
The HHT method, which was proposed by Huang et al. [17] in 1998, is a new alternative for dealing with nonlinear and nonstationary signals. This method can reveal accurate time-frequency information for signals. Compared to traditional methods of signal processing and failure testing [18–20], the HHT method is more accurate [21–23]. Because of this, it seems preferable to introduce the HHT method to the analysis of blasting vibration signals. Consequently, the objective of this paper is to use the HHT method to analyze the blasting vibration signals in a south China lead-zinc mine in terms of three-dimensional Hilbert spectrum, marginal spectrum, and instantaneous energy spectrum.
2. Basic Principle
Taking into account the inherent characteristics of signals, the HHT method decomposes signals into a series of intrinsic mode functions (IMFs) using the empirical mode decomposition (EMD) method. And then, Hilbert transform is applied to these IMF components, leading to the derivation of energy distribution spectrogram on the time-frequency plane. The derived energy distribution spectrogram conveys accurately various information with respect to time, frequency, and energy. It has been well known that the HHT method is mainly composed of EMD and Hilbert transform.
2.1. The EMD Method
The EMD method hypothesizes that a signal is a compound signal consisting of different IMFs. Each of these IMF components must satisfy the following: (1) over the entire time series, the number of extrema and zero crossing point must be equal or differ at most by one; and (2) at any point, the mean value of the envelope defined by the local maxima and the envelope defined by the local minima is zero. In this way, any signal can be decomposed into the addition of a finite number of IMF components.
The detailed implementation procedure of EMD can be summarized as follows:(1)Firstly, find out all the maximum and the minimum points for the signal . These points are then fitted by using cubic spline functions to generate two envelopes, which are defined, respectively, by the maxima and the minima. Calculate the mean value of the two envelopes, . The subtraction of from the original series results in deriving a new series without low frequency, :(2)Repeat Step 1 k times, until satisfying the requirement defined for IMF: the derived mean value approaches to zero. Thus, the first IMF component , which represents the highest frequency of signal , is obtained.(3)Separate from to derive a signal without high-frequency component, :(4)Treat as the original data, and repeat the above steps. Then, the second IMF component is derived. Continue repeating times, leading to the derivation of IMF components with a total number of . At the end of the repetition, the following relationship is obtained:where is the residual function, representing the average trend of the original signal. Note that the frequency bands are varied for different IMF components. And the sequence of is ranked in the descending order.
2.2. Hilbert Transform
For the time series , its Hilbert transform iswhere is Cauchy’s principal value. When the relationship of and is the complex conjugate, then the following analytic signal is derived:where and . Meanwhile, it is requisite to define instantaneous frequency as . Therefore, the Hilbert transform provides a specific function which can be used to calculate the instantaneous frequency and the amplitude.
Applying the Hilbert transform to each of the IMF components gives rise to the derivation of the Hilbert spectrum:where represents taking the real component. Equation (6) indicates the distribution of amplitude on the frequency-time plane. The time-domain integral of the Hilbert spectrum in (6) generates the marginal spectrum :
Furthermore, the instantaneous energy spectrum is derived through
3. Engineering Background and Data Sources
Field testing of blasting vibration signals was implemented in the S5#S stope, middle section of SH-455m, a south China lead-zinc deep mine. The all-band wave velocity and the vibration data that were collected during the 4th field blasting have been chosen to analyze using the HHT method.
The plane layout in field is shown in Figure 1, in which point B represents blasting source, and points 1, 2, and 3 are three vibration monitoring points (MP). Points 1 and 2 locate, respectively, on two sides of the chamber. Point 3 locates on one side of the roadway. These three MPs are collinear which brings convenience to study the blasting effect at different distances from the blasting source.
The Blastmate III manufactured in the USA, a vibration monitoring instrument as shown in Figure 2, has been used to monitor the blasting vibration effect and to collect data. The detector was selected as the Triaxial Geophone, which is capable of collecting the radial, the normal, and the tangential wave velocities in rock mass during blasting. The component of the instrument includes microphones, geophone, recorder, microcomputer, and microprinter, which have various functions such as vibration signal acquisition, signal analysis and processing, and report printing. The main parameters of the instrument are shown in Table 1. (The sampling rate of this experiment is 4026.)
The blasting was triggered by a high-precision millisecond detonator. The time lag levels of the detonator and the corresponding emulsion explosive weight for this blasting are presented in Table 2.
4. Signal Processing and Analysis
The preliminary analysis of the collected signal data has been conducted by applying the Fourier transform. The comparisons among the three monitoring points in terms of peak vibration velocities (PVV) and main vibration frequencies (MVF) along different directions are presented in Table 3. It is suggested from Table 3 that the peak vibration velocities and the main vibration frequencies in the vertical direction are almost consistently greater than that in radial and tangential directions, so the vertical vibration signals have been selected to analyze based on the HHT method. For MP 2, its original vibration waveform in the vertical direction is shown in Figure 3.
Since signal collection is often disturbed by numerous factors associated with the complex mining environment, noise signals inevitably exist in the original signals. However, the noise signals can pollute the original vibration signals, reducing greatly the accuracy of the signal analysis. Therefore, denoising before the signal analysis is very significant.
4.1. EMD-Based Wavelet Threshold Denoising
The EMD-based wavelet threshold denoising technique has been widely used in the signal denoising area, thanks to the advantages it possesses [24]. This technique is characterized by the multiresolution analysis, good time-frequency localization, and flexible threshold selection. It can also detect the transient state of the normal signals. Considering that the collected vibration signals are composed of low-frequency blasting vibration signals and high-frequency noise signals, and that the EMD method is capable of decomposing the blasting vibration wave into a series of IMF components that are ranked in descending order of frequency, in this paper, denoising has been conducted for the high-frequency IMF components using this technique.
The treatment of the vertical vibration signals at MP 2 by using the EMD method produces 11 IMF components. The waveforms of these IMF components are presented in Figure 4. It can be seen that the included frequencies are varied for different IMF components. From IMF1 to IMF10, the frequency and the amplitude decrease, while the time period increases. For IMF11, because no complete period is observed over the entire time-domain and the amplitude is relatively small, this component is the residual component representing the average trend of the original signals. It can also be indicated that the EMD method can analyze directly the signals taking considerably into account the signal characteristics, while without introducing any limitations. This method decomposes adaptively the signals into IMF components of finite number. And loss or emitting of IMF components is not omission. Consequently, the advantages of the EMD method, that is, adaptivity, completeness, and strong reconfigurability, are demonstrated.
Among all the IMF components, the frequencies of IMF1, IMF2, and IMF3 are distinctly greater, indicating that high-frequency noise signals are included in them. Meanwhile, the amplitudes of these three IMFs are relatively greater than those of other IMF components. This means that they are the dominant components that take up most of the total energy in the original signals. If wavelet-forced denoising is applied to IMF1, IMF2, and IMF3, then the valuable information included in the original signals is possibly removed which may lead to distortion of signals. On the contrary, the valuable information can be effectively extracted if the EMD-based wavelet threshold denoising technique is used. The detailed procedure of denoising using this technique is presented below.
Firstly, obtain the removed noise thresholds for IMF1, IMF2, and IMF3 through the wavelet function. Secondly, using the db4 wavelet basis function, decomposition of 5 levels is applied to the three IMFs to isolate the high-frequency coefficients. Thirdly, quantitative processing of the high-frequency coefficients is performed using the removed-noise thresholds. Fourthly, reconstruct the signal to complete denoising of the IMF components. Finally, the IMF1, IMF2, and IMF3 components after denoising are combined with the remaining 8 IMF components to derive the signals without noise. The signals without noise and the removed noise signals are shown in Figure 5. It can be seen that the vibration curve of the signals after denoising is smoother, and the trend of the vibration curve is clearer compared to that of the original signals.
(a)
(b)
4.2. Analysis of Three-Dimensional Hilbert Spectrum
Reapply the EMD method to the signals after denoising to obtain new IMF components. These components are then subjected to the Hilbert transform, leading to the derivation of three-dimensional Hilbert spectrum, marginal spectrum, and instantaneous energy spectrum for the signals after denoising.
Figure 6 shows the derived three-dimensional Hilbert spectrum, which can reflect visually the instantaneous characteristics of the vibration signals and reveal clearly the distribution of signal energy on the frequency-time plane [25]. In this figure, each of the bars of different colors represents the normalized instantaneous energy at some specific frequency and time. It is recognized that the bars are mainly distributed within the range defined by the time of 0.2–0.4 s and frequency of 0–200 Hz, while the energy corresponding to frequency greater than 200 Hz can be neglected. This means that most of the frequencies of the blasting vibration signals are less than 200 Hz. The published studies of blasting vibration signals are mainly performed through either Fourier transform or wavelet transform, with mere consideration of the influence of frequency or amplitude while without analyzing the combining effect of frequency, energy, and vibration duration [26]. However, this defect can be compensated by the introduction of three-dimensional Hilbert spectrum.
4.3. Analysis of Marginal Spectrum
The marginal spectrum is the time-domain integral of the Hilbert spectrum and thus represents the addition of the amplitude for each of the frequencies over the time domain [27]. Therefore, marginal spectrum can reflect the condition of energy concentration for the frequencies. The marginal spectrum of the signals is shown in Figure 7. From this figure, it can be seen that the energy is mainly located in the low-frequency band where the frequency is less than 200 Hz. In particular, within the band of 0∼10 Hz, the accumulated amplitude is relatively greater with a maximum of 1043 mm/s when compared to other frequency bands. The concentration of most of the energy in this frequency band indicates that the frequency tends to decay to the low-frequency band (less than 10 Hz) during the blasting process. After reaching 10 Hz, the amplitudes for the frequencies drop dramatically and then fluctuate around 100 mm/s. When the frequency becomes greater than 200 Hz, the amplitude begins to decay to zero.
For the convenience of quantitative description of the energy included in different frequency bands, 6 frequency bands, that is, 0∼50, 50∼100, 100∼200, 200∼300, 300∼400, and >400 Hz, are selected. According to (8), the integrals of frequency over these frequency bands are performed to calculate the ratios of energy taken up by these frequency bands. To ensure the accuracy of the calculated results, the data for the three MPs are analyzed parallelly. The ratios of energy taken up by these frequency bands are listed in Table 4. It is indicated from Table 4 that the energy taken up by the frequency band of 0∼50 Hz is the maximum, being more than 50% of the total energy. This also means that the frequency band of 0∼50 Hz is the dominant frequency band. The ratio of energy taken up by the <200 Hz frequency band reaches 93.1%, which further indicates that the frequency of the blasting vibration wave is almost totally within the range of 0∼200 Hz.
4.4. Analysis of Instantaneous Energy Spectrum
The instantaneous energy spectrum can reveal the accumulation of blasting vibration energy over the time domain and its variational characteristics. Figure 8 shows the waveform and the instantaneous energy spectrum of the signals. It can be seen that the distribution of instantaneous energy agrees well with the vibration curve. And the peaks of the instantaneous energy are reached when the mutations of vibration waveform take place. This further validates the ability of the HHT method in recognizing signals and its better resolution. For the instantaneous energy spectrum in Figure 8, it is indicated that the energy begins to appear at 0.22 s and approximately disappears at 0.45 s, showing that the blasting vibration at this MP keeps being observable during this duration. Also, 7 peaks of instantaneous energy are clearly observed during this duration. These peaks, however, represent the release of energy brought by emulsion explosive explosion. Therefore, the peaks in the instantaneous energy spectrum can be used to determine the initiation instant of time of the detonator and to check if blasting initiates at the accurate instant of time for each of the time-lag levels.
The instant of time corresponding to these 7 peaks are 0.240, 0.266, 0.295, 0.326, 0.353, 0.426, and 0.450 s, which represent, respectively, the initiation instant of time of time-lag level 1, 2, 3, 4, 5, 6, 7, and 8. Note that the third peak is the maximum, indicating that the explosive quantity for the time-lag level 3 is relatively large. Similarly, the magnitude of the sixth and the seventh peaks shows that their explosive quantities are relatively small. The comparison between the theoretical and the measured time lags for these time lag levels is presented in Table 5, which shows that the measured time lags are all within the theoretical ranges of time lag. If an assumption is made that the time-lag level 1 detonator explodes punctually at an instant of time 0, then the real initiation instant of time for these detonators of different time-lag levels can be determined according to the measured time lags. It can be indicated that the time-lag level 8 detonator explodes with 15 ms earlier compared to the theoretical initiation instant of time. Despite this, the overall blasting effect is comparatively satisfactory.
5. Conclusions
Based on the blasting (4 times) vibration data collected from the S5#S stope, the middle section of SH-455m, a south China lead-zinc mine, the HHT method has been used to analyze the blasting vibration waves in terms of time-frequency and energy characteristics. The following conclusions can be drawn:(1)The marginal spectrum shows that the accumulative amplitude is relatively great in the frequency domain of less than 10 Hz. In this domain, the maximum amplitude is 1043 mm/s. The fact that most of the signal energy is concentrated in this frequency domain indicates a direction to which the frequency tends to develop when explosion initiates.(2)Most of the blasting vibration frequencies are concentrated within the range of 0∼200 Hz. This range takes up over 90% of the total energy. Additionally, in this range, half of the total energy is taken up by the 0∼50 Hz frequency band, demonstrating that the 0∼50 Hz frequency band is the dominant frequency band.(3)The calculated initiation instant of time of the detonators by recognition of the instantaneous energy spectrum agrees well with the measured results, demonstrating the feasibility of using this method in recognizing delay blasting. Also, good resolution and excellent recognition capability are validated for the HHT method in signal processing.
Conflicts of Interest
The authors declare that there are no conflicts of interest regarding the publication of this paper.
Acknowledgments
This research work was funded by Open Research Fund Program of Hunan Province Key Laboratory of Safe Mining Techniques of Coal Mines (Hunan University of Science and Technology, Grant no. 201505), National Natural Science Foundation of China (Grant nos. 51004007 and 51774107), the Open Program of State Key Laboratory of Explosion Science and Technology (Beijing Institute of Technology, KFJJ17-12M), the Fundamental Research Funds for the Hefei Key Project Construction Administration (2013CGAZ0771), and the Fundamental Research Funds of the Housing and Construction Department of Anhui Province (2013YF-27).