An Analysis Method for Time-, Frequency-, and Energy-Domain Characteristics of Downhole Microseismic Signals and Its Application

Based on SOS microseismic (MS) monitoring system of a certain mine, time-, frequency-, and energy-domain characteristics of MS signals are analyzed by comprehensively applying one-dimensional discrete Fourier transform (DFT) and one-dimensional discrete wavelet transform (DWT).,e technology is employed to explore MS signals from incident and transmitted waves on the fault plane. ,e research results show that, within the time domain, the transmitted waves are found to have the rapidly reducing velocity and a few wave crests; within the frequency domain, the fault plane presents a low-pass filtering effect on incident waves; within the energy domain, the fault generally shows the vibration-attenuating effect onMS signals.,e analysis method is used for the No. 3307 working face of a certain mine to realize the prediction on a buried fault in the middle part of the working face, and it is validated to have a favorable effect.


Introduction
Underground coal mines exhibit complex geological environment.
e complex geological structures centered on faults not only influence the continuity of mine production, but also probably trigger some dynamic disasters such as fault-induced rockburst, which greatly threatens the safety production of coal mines [1,2]. e time-, frequency-, and energy-domain curves of a microseismic (MS) signal, as the concomitant signal during the deformation and fracture of coal-rock mass, can indirectly reflect the generation, propagation, and attenuation processes of the signal.
It is widely known that the analysis based on vibration signal is used for detection of mechanical and electrical faults of the motor, and the change of generated vibration signals can help detect the type of fault [3][4][5][6], which can provide a reference for the MS effect analysis in underground coal mining. To date, numerous studies show that MS signals are subjected to local abrupt change or transition when going through a fault. erefore, it is possible to predict a buried fault when excavating the unexposed working face by analyzing the population characteristics of MS signals in an area.
At the beginning of the century, Dou et al. [3][4][5] first introduced the SOS MS monitoring system in China to perform the precursor identification, and monitoring and early warning for rockburst, attaining numerous MS data. Based on signal analysis theory, Cai et al. [6] developed a series of software for postprocessing of MS signals. ese studies provide equipment and software bases for analyzing the time-, frequency-, and energy-domain characteristics of MS signals; however, they focus on exploring the timedependent law of MS signals. e frequency-and energydependent laws of MS signals are maturely applied in the geological exploration field, while their MS signal source generally needs to be actively activated. According to the research in a coalfield in Xinmi, Henan Province, China, Chen et al. [7] found that it is feasible to effectively improve the detection accuracy and efficiency of buried faults after eliminating the influence of high-frequency noises in vibration signals with the aid of the wavelet transform. Liu et al. [8] and Zhou et al. [9] explored the influence law of the fault fracture zones with different depths and dips on propagation of Rayleigh waves. e research reveals the vibration-attenuating effect of faults on Rayleigh surface waves. Jia [10] also employed the technology for time-and frequency-domain characteristic analysis of Rayleigh waves to detect the buried faults.
At present, the signals measured by a MS monitoring instrument in underground coal mines of China are a group of discrete time-dependent signals; by connecting various points, the time-domain curve of MS signals can be drawn. Owing to the time-domain curve corresponding to an amplitude-time function, abundant time-dependent information is contained, while no any frequency-dependent information is included in the curve. e time-domain curve can be reconstructed by utilizing the discreate Fourier transform (DFT) and one-dimensional discrete wavelet transform (DWT). Furthermore, the accurate frequencyamplitude curve of MS signals is attained; however, the curve fails to reflect the change of the frequency of MS signals with time. With the development of the signal analysis technology, the time-frequency representation method of signals, that is, the method for representing signals through the joint function of time and frequency, has been applied to some extent [11][12][13][14]. On this basis, the frequency spectra are further transformed into energy spectra of MS signals, and thus the law of energy attenuation of MS signals in the propagation process can be obtained. erefore, the existing SOS MS monitoring system and monitored MS data during mining of the working faces of mines were applied. On this basis, the time-, frequency-, and energy-domain characteristics of MS signals from incident and transmitted waves within the fault zone were separately analyzed by comprehensively using multiple signal analysis methods. It is expected to determine a quantified threshold to predict whether a fault exists or not, thus providing index reference for predicting the buried faults in working faces.

Transformation between Time-and Frequency-Domain
Curves. Aiming at the characteristic that the measured MS signals appear as the time-amplitude curve, fast Fourier transform (FFT) allows the transformation between timeand frequency-domain information of MS signals. In essence, MS signals are decomposed into many sine waves and high-order harmonics with different frequencies consisting of the basis functions thereof. e basic form of Fourier transform is expressed as follows [15]: where X(t) refers to the function between the wave velocity and time history of MS signals by taking the time t as a variable; X(w) stands for the amplitude-frequency function by taking the angular frequency w of MS signals as a variable. Moreover, the wave files are required to be a continuous function in FFT method, while the measured MS data in underground coal mines are stored as a group of arrays with the fixed frequency. erefore, it is necessary to process the wave files by using the DDT [15]. e DDT delivers the same principle as the FFT method, while it also realizes the data discretization. e DDT method is expressed as where X(n) and w separately represent the wave velocities of various data points and the angular frequency of MS signals.

Transformation between the Time-Domain and Time-Frequency
Curves. However, due to hiding all time information in the transform process, the FFT method cannot show the change law of the frequency with time. Based on the one-dimensional DWT [15], it is possible to conduct the binary decomposition on MS signals; afterwards, the lowfrequency parts are selected to conduct the transform at a half of the original scale. In this way, the high-resolution time information is attained at the high-frequency bands and high-resolution frequency information at the low-frequency bands of signals.
It is supposed that a MS signal X shows the length N, and, therefore, the whole calculation process is completed within log 2 N steps at most. Starting from original signals in the first step, two groups of parameters are generated: one group refers to the approximate signals cA 1 obtained using the lowpass filter, while the other group corresponds to the detailed signals cD 1 attained by applying the high-pass filter. e decomposition result of MS signals X(N) by employing the i-order DWT is expressed as follows: where cA i and cD i separately denote the approximate and detailed coefficients of MS signals at nodes in the ith sublayer; moreover, the bandwidths of the MS signals at the ith node can be separately expressed as follows: e decomposition result of MS signals with the aid of the one-dimensional DWT is displayed in Figure 1: According to the research results by Lin [13] and Zhao et al. [14], the db8 wavelet basis function is selected, which shows favorable adaptability in processing nonstationary vibration signals. e order of the wavelet decomposition equals 5.

Transformation between Frequency-and Energy-Domain
Curves.
e energy of MS signals, a kind of spherical wave, is generally calculated by applying the energy density method [1]. e MS energy W is defined as the known energy flux of seismic waves on the reference spherical surface, which is generated in the vibration process of the wave group s, τ:

2
Shock and Vibration where U and n refer to the energy flux density vector and the one-dimensional normal vector of the microunit dS with the surface area of S, respectively. e energy density of all measuring points is attained by calculating the sum of the densities of P and S waves in various groups. e energy density ε(r) is calculated according to the following equation: where ρ, c, v, and τ denote the media density, propagation velocity of waves, amplitude of vibration velocity, and the duration of vibration of wave groups, with i � P, S and k � 1, 2, ...; According to the above equation, it can be found that the MS signal energy is directly proportional to the square of the amplitude. erefore, the energy distribution of MS signals can be expressed by using the relative energy. In the frequency-domain curve, if the MS signals are partitioned into n frequency bands, the relative energy within the ith frequency band can be expressed as follows: where E i denotes the energy of MS signals at the S i th frequency band, and x i,k (i � 1, 2, . . .n, k � 1, 2, . . .m; m refers to the number of discrete points of MS signals) represents the amplitude at the kth sampling point within the ith frequency band of MS signals. e total energy E of MS signals can be expressed as follows: where E and E i represent the total energy of MS signals and the energy at the S i th frequency band of MS signals. us, the relative energy density in the ith frequency band can be shown as follows: where E 0 and E i denote the total energy of MS signals and the energy at the S i th frequency band of MS signals, respectively.

e Method for Comprehensively Analyzing the Time-, Frequency-, and Energy-Domain
Characteristics. e MS signals recorded by using the SOS MS monitoring instrument are taken as examples. e main steps of the analysis method for time-, frequency-, and energy-domain characteristics of MS signals are summarized as follows: (1) the effective downhole MS signals are monitored by using the SOS MS monitoring instrument, and then data are stored; (2) the MS signals are extracted, and a time window is selected to attain the slices of MS signals with a reasonable length based on the time window; (3) the low-pass filtering is conducted on MS signals to eliminate the noise interference factors such as downhole mechanical equipment; (4) the time-amplitude curve (that is, time-domain curve) of MS signals is drawn and analyzed; (5) the programming is performed by utilizing MATLAB software, and the timedomain curve is subjected to DDT according to Nyquist sampling frequency to draw and then analyze the frequencyamplitude curve (that is, frequency-domain curve) of MS signals; (6) the time-domain curve is subjected to DWT by applying the wavelet analysis toolbox of MATLAB software to draw and analyze the time-frequency curve of MS signals; (7) the time, frequency, and energy data obtained through the aforementioned processing are normalized, and the time-, frequency-, and energy-domain characteristic curves of MS signals are drawn and analyzed by applying the threedimensional drawing function of the MATLAB software. e method for comprehensively analyzing the time-, frequency-, and energy-domain characteristics is presented in Figure 2.
For MS signals from underground coal mines, the seismic waves synchronously experience reflection and transmission due to being influenced by the interfacial dynamic effect at the interface when waves from the seismic source propagate to the structural plane such as faults. e incident and reflected wave surfaces of the same seismic wave, as well as different types of seismic wave surfaces, are superposed, which causes the abrupt change of the spectra of the transmitted waves. It is mainly shown as changes of the time, frequency, and energy spectra of the transmitted waves. Owing to the fact that vibration lasts for a short time, and the seismic source is close to the observation station, it is difficult to attain the hidden law of MS signals through the conventional observation. Based on the method for comprehensively analyzing the time-, frequency-, and energydomain characteristics, it is feasible to partially enlarge different characteristics of MS signals. By analyzing the abrupt change of time and frequency spectra and energy dissipation of MS signals, it is possible to predict the positions of large-scale structural planes within areas. Moreover, MS events generally intensively appeared in the area 10∼60 m in front of the working face. Considering this, the area 70∼120 m in front of the open-off cut was mainly monitored, and the monitoring was constantly performed during the mining of the working face. e condition of the working face is displayed in Figure 3. Some parameters of the MS monitoring system applied in the No. 3307 working face are shown as follows: 16 seismic geophones were used to sample signals at 8000 groups of points at the sampling frequency of 500 Hz, with the time of recording signals for 16 s.

Comparison of the Time-and Frequency-Domain Curves of the Incident and Transmitted Waves.
e MS event occurring at 10 : 57 : 53 a.m. on October 10, 2015, is used as an example. At first, the MS data measured by the 16 groups of geophones are screened to attain the signals measured by geophones 1#, 2#, and 3#. Afterwards, the signals are partially enlarged, and the slices within the time window of 4 s are drawn. e enlarged MS signals are smoothed through the low-pass filtering based on the MATLAB software. e time-domain curves of smoothed signals in three groups of channels are drawn, as shown in Figure 4.
Obviously, it can be seen from Figure 4 that the timedomain curve in the 2# channel shows different characteristics than those in the 1# and 3# channels. As shown in Figure 3, if the F1 fault is taken as the boundary, the 1# and 2# geophones are distributed in front of the F1 fault, which are used to monitor the incident waves; the 3# geophone lies in the rear of the F1 fault, which is applied to monitor the transmitted waves. In addition, considering that the distances of the intensively monitored area to the 2# and 3# geophones are approximated, the data in 2# and 3# channels are separately selected for comparison and analysis. e relationship curves v(t) between the particle vibration velocities in the two channels and time are drawn, as shown in Figure 5. By comparing curves, it can be seen that the peak velocity of incident waves is 4.37 × 10 −4 m/s, with a low reduction rate of wave velocity, and there are about 30 wave crests. Moreover, P wave first motion is clearly found, and vibration lasts for about 1 s; the peak velocity of the transmitted waves equals 6.27 × 10 −4 m/s, showing a large reduction rate of the wave velocity, and about 10 wave crests appear. In addition, the coda waves are developed, and P wave first motion is hard to be identified due to the interference of high-frequency signals. In this case, the duration of the vibration is about 3 s. e digital signals of the incident and transmitted waves are imported into the MATLAB software. By performing the programming based on Equation (2), the measured signals are subjected to DDT, and then the frequency-amplitude distribution curves of incident and transmitted waves are drawn, as shown in Figure 6.
Owing to the fact that the sampling frequency of seismic geophones is 500 Hz, the Nyquist sampling frequency is determined as 250 Hz according to the sampling theorem. Not any amplitude information is contained within the frequency band of 200∼250 Hz. For convenience of analysis, the MS signals in the frequency band of 0∼200 Hz are analyzed, in which 0∼60, 60∼120, and 120∼250 Hz are separately defined as low-, middle-, and high-frequency bands. rough analysis, it is found that (1) the frequencies of incident and transmitted waves are mainly found within the frequency band of 0∼160 Hz, in which the dominant frequency of the incident waves mainly appears in the middleand low-frequency bands of 0∼120 Hz; by contrast, the main frequency of the transmitted waves is transferred from highto low-frequency bands, and their dominant frequency corresponds to a low-frequency band of 0∼20 Hz. It indicates that the presence of the fault blocks the middle-and high-frequency vibration signals and some low-frequency vibration signals (20∼120 Hz) of incident waves, showing a low-pass filtering effect on incident waves; (2) within the dominant frequency bands, the incident waves deliver the uniformly distributed frequencies, and their frequency spectra show significant multipeak characteristics; as for the transmitted waves, the frequency is mainly found to be about 2 Hz, and the frequency spectra exhibit the sharp single-peak characteristic.

Comparison of Time-Frequency Curves of Incident and Transmitted Waves.
e fifth-order DWT is separately conducted on MS signals from incident and transmitted waves to attain the approximate and detailed coefficients of two MS signals, as shown in Figure 7.

Shock and Vibration
According to the above decomposition results, d1 refers to the part with the highest frequency in the detailed coefficient of MS signals, with a small amplitude, low energy, and high discreteness, which is not suitable to be used as the quantified threshold. erefore, it is not intensively analyzed. In the decomposition results, the main energy of MS signals is mainly found in wavelet decomposed signals d2, d3, d4, and d5. By observing the waveform graph, it can be seen that P wave first motions of four waveforms d3, d4, d5, and a5 are clearly observed. e onset amplitudes of the four wavelet decomposed signals of incident and transmitted waves are compared with their corresponding maximum amplitudes, and the change laws of the onset amplitudes and maximum amplitudes are shown in Figure 8. e above analysis results show that (1), among wavelets at different frequency bands, the onset amplitude of incident waves is logarithmically correlated with the frequency band of wavelets. e lower the frequency band, the smaller the onset amplitude of wavelets; (2) the onset amplitudes of transmitted waves at the high-and middle-frequency bands are also logarithmically correlated with the frequency band of wavelets; however, the wavelet presents obviously abrupt change in the amplitude at the low-frequency band a5. e onset amplitude of the transmitted wave is 14.47 times that of the incident wave at the frequency band a5, far larger than those of the wavelets at the other frequency bands; (3) the maximum amplitudes of incident waves are relatively discrete, while those of transmitted waves regularly increase only at the frequency band d5; (4) comparatively speaking, the increment of the onset amplitude of the wavelet at the frequency band a5 is more suitable to be used as the quantified threshold to measure the influence of faults on the energy in MS signals. e DWT result is calculated by performing programming based on MATLAB software, and nephograms of timefrequency characteristics of MS signals from incident and transmitted waves are separately drawn, as shown in Figure 9.
By observing the nephograms, it can be seen that (1), along the timeline, the vibrations of incident and transmitted waves at the dominant frequency band last for about 0.8 and 0.4 s, respectively; the time window of the dominant frequency band of incident waves is narrower than that of the transmitted waves; (2) along the frequency axis, the main vibration frequency band of transmitted waves narrows relative to that of incident waves, showing the transfer of the dominant frequency band from the middle-to low-frequency bands.

Comparison of Energy-Domain Curves of Incident and
Transmitted Waves. According to the directly proportional relationship between the energy and the square of amplitude of MS signals, it is feasible to characterize the energy distribution of MS signals at different frequency bands with the aid of the relative energy density. e frequency-amplitude curve of MS signals in Figure 6 is segmented at first, in which the frequency band of 0 ∼ 160 Hz is divided into eight sections (S 1 ∼S 8 ).
According to Equations (7)-(9), the energy density distributions of incident and transmitted waves at different frequency bands are separately calculated, as shown in Tables 1 and 2. e corresponding energy density distribution curves are drawn, as shown in Figure 10.
According to the above analysis, the energy of incident waves mainly concentrates in the frequency band of 20∼80 Hz, at which the energy accounts for about 67.96% of the total waveform energy. e high-frequency components of transmitted waves significantly reduce, so the energy at the high-frequency band dramatically declines. e waveform energy is mainly found at the low-frequency band of 0∼20 Hz, at which the energy takes up 97.33% of the total waveform energy.
According to results obtained through DDT and DWT and the above analysis results, the time, frequency, and energy data are normalized by utilizing the MATLAB software. e three-dimensional graphs of time-, frequency-, and energy-domain characteristics of incident and transmitted waves are separately shown in Figure 11. By comparing the above graphs, it can be found that (1) the energy of MS signals varies in the leap pattern with time. Globally, the energy of incident waves with high density is    continuously distributed, and the energy is absent only in few frequency bands. e energy of transmitted waves in many vibration frequency bands is lost at the timeline. By discretizing and calculating the curved surface, it can be attained that, at the frequency axis, the energy density of transmitted waves at the low-frequency band of 0∼60 Hz is about 84.7%, while that at the middle-frequency band of 60∼120 Hz is about 6.43%; the energy density at the highfrequency band nearly equals 0. It generally shows the vibration-attenuating effect of faults on MS signals. (2) Locally, the energy of transmitted waves is higher than that of incident waves within some areas, which indicates that faults deliver a vibration-increasing effect on vibration waves in some periods and frequency bands.

Engineering
Verification. According to the above principle, the three-dimensional distribution graphs of time-, frequency-, and energy-domain characteristics of the MS event occurring at 10 : 57 : 53 a.m. on October 10, 2015, in the channels 1#, 2#, and 3# are separately drawn. For convenience of comparison, their time-energy profiles are analyzed, as shown in Figure 12.
As shown in Figure 12, the energy distribution laws in the three channels present an obvious difference: the energy spectra in 1# and 3# channels are distributed in a relatively similar fashion, which greatly differs from that in the 2#    channel. According to the actual condition of the working face, the geophones 1# and 2# are both located in front of the F1 fault, so data collected thereby should show the energy distribution characteristics of incident waves; however, the energy distribution characteristic of data collected by the geophone 1# is similar to that of geophone 3#, showing the energy characteristics of transmitted waves. On this basis, it is inferred that a large-scale structural plane is probably present between the monitored area and the geophone 1#. Also, the amplitude-frequency and time-frequency distribution are extract, as shown in Figures 13 and 14. It can be seen that the amplitude-frequency and time-frequency distribution laws are similar to those of the energy distribution. Afterwards, the buried fault F2 was exposed in the mining process of the working face, which practically verifies the aforementioned analysis results. e position of the measured buried fault F2 is displayed in Figure 15.
As the fault was not explored during the excavation of the No. 3307 working face, the area was not identified at high rockburst risk and subjected to targeted measures when performing region division of the working face according to the rockburst risk. As a result, high-energy mine tremors frequently occurred during the mining of the working face. If the presence of the existing fault had been predicted through the analysis technology for time-, frequency-, and energy-domain characteristics, it is possible to manage the area as a zone with high rockburst risk, thus avoiding the occurrence of potential fault-induced rockburst.

Discussion
e question arises as to whether the aforementioned analysis method for time-, frequency-, and energy-domain characteristics of MS signals could be used to predict the buried faults in underground working faces. Such a method could take full advantage of the SOS MS monitoring system in rock burst coal mines. e advantages of the method are as follows: (1) a single frequency component pulse physical model of vibration wave with Fourier spectrum characteristics is established. e model shows that when the pulse wave with Fourier spectrum passes through the nonlinear joint surface, its waveform changes, and its corresponding spectrum will be modified due to absorption. (2) Discrete Fourier transform (FFT) is used to analyze the incident wave and transmission wave of fault plane. e results show that the nonlinear joint plane shields the high-frequency signal and part of the low-frequency signal of incident wave and acts as a lowpass filter for incident wave. e attenuation range of high-frequency signal is much larger than that of lowfrequency signal. (3) One-dimensional discrete wavelet transform (DWT) is applied to incident wave and transmission wave in fault plane. It is pointed out that the abnormal amplitude increment of a5 wavelet P wave in low-frequency band of transmission wave is more suitable than S-wave to measure the influence of fault on seismic wave propagation. (4) e time-frequency characteristic nephogram and time-energy density distribution diagram of incident wave and transmission wave are calculated and drawn by MATLAB programming, and the "isolation" effect of fault surface on incident wave energy is explained. (5) e application of this technology to the detection of hidden faults in a deep working face can predict a small hidden fault in the middle of the working face. rough the pressure relief work in advance, the potential risk of fault burst can be avoided. However, it is definitely promising to develop an optimization scheme for the detection of buried fault, but depending on the implementation, the resulting method will probably be more complex and time consuming computationally, as a large number of field examples are needed for verification. Despite these remarks, this idea has the potential to deliver interesting results.

Conclusions
(1) e transformation between time-and frequency-domain curves of MS signals is realized by using the DDT, while that between the time-domain and time-frequency curves of MS signals is achieved through DWT; moreover, the relative energy density representation method is employed to realize the transformation between the frequency-and energy-domain curves of MS signals. Eventually, the change laws of the frequency and energy distributions of MS signals with time are attained by utilizing the analysis method for time-, frequency-, and energy-domain characteristics of MS signals. (2) rough research, the differences of time-, frequency-, and energy-domain characteristic curves of incident and transmitted waves of the fault plane are obtained. e differences are mainly shown as follows: within the time domain, the transmitted waves present rapidly decreasing velocity and a few wave crests; within the frequency domain, the structural plane blocks the highfrequency signals and some low-and middle-frequency (20∼120 Hz) signals of incident waves, playing a lowpass filtering effect on incident waves. e vibration frequency of transmitted waves mainly occurs in the low-frequency band; within the energy domain, the fault generally delivers a vibration-attenuating effect on MS signals, while it still shows a vibration-increasing effect on MS signals within several time periods. e onset amplitude of the wavelet at the frequency band a5 is more suitable for being used as the quantified threshold to measure the influence of the fault on propagation of MS signals.
(3) e comprehensive analysis method for time-, frequency-, and energy-domain characteristics of MS signals is applied in a certain working face, and it can predict a buried fault in the middle part of the working face, so it is verified to have a favorable effect. However, considering the complexity of downhole geological structures and nonrepeatability of the case, how to establish a quantitative identification index for a structural plane based on MS signals that is applicable for all mines also needs plenty of data accumulation and remains to be further explored.

Data Availability
Some or all data, models, or code generated or used during the study are available from the corresponding author upon request.

Conflicts of Interest
e authors declare that they have no conflicts of interest.