# Energy Characteristics of Seismic Waves on Cardox Blasting Tube

Authors: Xuejiao Cui; Bo Ke; Songtao Yu; Ping Li; Mingsheng Zhao
Journal: Geofluids
 (2021)

---

## Abstract

In order to study the energy characteristics of seismic waves on the liquid CO2 blasting system, the blasting seismic wave signal of liquid CO2 blasting was obtained by on-site microseismic monitoring tests. The adaptive optimal kernel time-frequency analysis method was used to study the basic time-frequency properties of the seismic wave signal. Combining wavelet packet transform decomposition and reconstruction and adaptive optimal kernel time-frequency analysis method, the liquid CO2 energy distribution of the seismic wave signal was further analyzed. And the energy regression model of seismic wave source of liquid CO2 blasting system was discussed. The results show that the vibration velocity is at a low level, and the main frequency range is between 30 and 70 Hz, and the duration is about 20-30 ms. The energy is mainly distributed in 0-125 Hz, which is composed of two main regions. The power function model can be used to describe the attenuation law of the seismic wave energy. The energy conversion coefficient and characteristic coefficient of the source of liquid CO2 blasting system were defined and analyzed. Combined with the empirical formula of the Sadovsky vibration velocity, the energy regression model of the seismic wave source of liquid CO2 blasting system was obtained.


---

## Body

## 1. Introduction

The cyclic nature of the drill-and-blast method, increasing excavation demands, environmental concerns, and trends towards safety, necessitates the development and exploration of the potential of new and improved concepts of rock excavation. One of the promising and nonexplosive concepts that warrants in-depth evaluation as a tool for rock fragmentation is the “Penetrating Cone Fracture” (PCF) method [1]. Liquid CO2 blasting system is one of the nonexplosive blasting technologies, which is based on liquid carbon dioxide being converted to high-pressure carbon dioxide gas or fluid with the ignition. The gas spreads through fissures and microcracks in the rock and breaks it in tension, rather than compression as with explosives, and the damage at lower tensile stress levels is more efficient in the utilization of energy, less vibration, pollution-free, reduce the damage of the surrounding rock mass, and destruction of the environment. Figure 1 shows a schematic diagram of the components which make up the cartridge. The chemical energizer is activated by a small electrical charge which causes the blasting [2].Figure 1
Schematic diagram of a liquid CO2 blasting system [2].In the past, researches on liquid CO2 blasting technology mainly concentrated on the blasting equipment and its applications. The gas pressures and the velocities of blast waves that travel through the sandy shale were determined by field tests [3]. Other important applications of this technology are in the fields of environment-friendly blasting, such as urban underground construction [4], neighbor rock breaking of the forest [5], highway construction for the cold region [6], and controlled blasting [7, 8]. Recently, with the wide application of this technology, many scholars have done a lot of research on coal seam permeability improvement [9–17], energy calculation [2, 18], pressure characteristics [19, 20], and mechanism of fracture [13, 21–25] and its application.Liquid CO2 blasting fracturing technology is a new green blasting excavation method, which can effectively reduce the vibration effect compared with the traditional explosive blasting. Blasting vibration signal is the carrier and physical manifestation of blasting seismic waves, which determines the crack geometrical form and propagation mode of the rock mass. However, the time-frequency characteristics and energy distribution of seismic wave signals of liquid CO2 blasting are not clear from above literature review.The blasting seismic waves of liquid CO2 blasting are nonstationary signals, and time-frequency analysis/distribution is an important mean of analyzing such nonstationary signals. The time-frequency analysis technique can be divided into two categories: linear time-frequency analysis and nonlinear time-frequency analysis. The linear time-frequency analysis includes short-time Fourier transform (STFT), Gaber expansion, wavelet transform (WT), and S transform. The nonlinear time-frequency analysis includes bilinear time-frequency analysis and adaptive optimal kernel (AOK) [26, 27] time-frequency analysis [28]. Moreover, based on the simulation signal of blasting vibrations and seismic waves, the STFT, WT, S transformation, Wigner-Ville distribution, smooth pseudo-Wigner distribution, cone-shaped kernel time-frequency distribution, and the distribution of the AOK time-frequency analysis were used for the signal processing in the MATLAB software. The results of the comparative analysis show that the kernel function of AOK is adaptive with the change of time, which can effectively suppress the cross term and realize the time-frequency localization with the best precision. Published researches that mainly focused on the signals of blasting vibrations and seismic waves processed by the AOK time-frequency analysis found in the literature were those of Zhao et al. [28], Wang et al. [29], Sejdić et al. [30], and Sun et al. [31].In the analysis method of signal decomposition and reconstruction, compared with wavelet analysis, wavelet packet transform analysis provides a more refinement analysis method for signals [32–34]. The frequency band can be divided into different levels; according to the characteristics of the analyzed signal, the corresponding frequency band is adaptively selected to match with the signal spectrum, which can improve the time-frequency resolution. The signal was decomposed and reconstructed in multiscale and multiresolution by the wavelet packet, which can better express the energy distribution characteristics of the frequency band. However, the signal reconstruction of time-frequency analysis still uses the linear summation calculation method to analyze the different frequency band energy, instead of nonlinear calculation on the double integral of time and frequency to analyze the signal energy, which resulting in a relatively low accuracy of calculation results.In this work, the research is conducted by microseismic tests of blasting seismic signals of liquid CO2 blasting as the foundation after reviewing the signal processing and analysis methods. Firstly, the basic time-frequency characteristics of the seismic signals were analyzed through the adaptive optimal kernel time-frequency analysis method. Secondly, combining wavelet packet transform decomposition and reconstruction and adaptive optimal kernel time-frequency analysis method, the energy distribution of the seismic wave signals from liquid CO2 blasting was further analyzed. Finally, an energy regression model of the seismic wave source of the liquid CO2 blasting system was obtained. These findings can provide a theoretical basis for the propagation law of shock waves of supercritical CO2 jet in the rock mass and for the design of liquid CO2 blasting.


## 2. Experimental Methodology

### 2.1. Test Scheme and Design

Since the seismic waves from liquid CO2 blasting (principle and equipment can be referred to in reference [1]) are weak, the microseismic monitoring system was adopted to monitor the blasting seismic waves. The tests were divided into two parts: one part is for the analysis and research of time-frequency characteristics of the blasting source signals, and the other part is for the analysis of seismic wave attenuation law, as shown in Figure 2, the layout plan of the microseismic test station.(1)
Analysis of Time-Frequency Characteristics. From Figure 2, the monitoring points are arranged as a circle with a total of 8 stations. There are five blasting points; one blasting point is arranged in the center of the circle, and the other four blasting points are arranged in a square at a distance of 0.5 meters from the center of the circle. The distance between the 8 monitoring points and the center is 12 meters(2)
Study on the Attenuation Law of Seismic Waves. Monitoring points are arranged in a “straight line” pattern, with a total of 4 monitoring points. The first monitoring point is 12 meters away from the blasting source, and the distance between each monitoring point is 12 metersFigure 2
Microseismic monitoring schematic diagram of the test scheme.


### 2.2. Implementation of the Test Plan

#### 2.2.1. Test Equipment and Installation

As shown in Figure3, the monitoring point equipment used in the experiment consists of three parts: sensor, collector, and battery. The sensor is installed by deep drilling, with a hole depth of 1 meter and a diameter of 90 mm. After the installation of the whole monitoring point, the monitoring software can be used to check its background noise. Generally, the background noise should be no more than E-5 V. The sensor and the hole wall adopt the yellow mud coupling.Figure 3
Equipment installation at site monitoring points.


#### 2.2.2. Program Implementation

A total of 5 blasting tests were carried out in this test, and the phase change blasting parameters of liquid CO2 are shown in Table 1, which were obtained from the laboratory at Central South University (see the reference [2]). Test 1, test 2, and test 3 were carried out to discuss the time-frequency characteristics, and test 4 and test 5 were carried out to analyze the attenuation law of seismic waves. The recording sampling rate of the test data was 4000 Hz. A total of 8 microseismic instruments were used to record seismic wave data, and the sensitivity of the sensor was set at 200 V/m/s.Table 1
Phase change blasting parameters of liquid CO2.
Serial numberCO2 qualityStomatal directionBurst disc pressureBlasting rod10.452 kgEast and west150 MPa1 bar20.526 kgSouth and north150 MPa1 bar30.454 kgEast and west150 MPa2 bar40.498 kgEast and west150 MPa1 bar50.448 kgEast and west150 MPa2 bar


## 2.1. Test Scheme and Design

Since the seismic waves from liquid CO2 blasting (principle and equipment can be referred to in reference [1]) are weak, the microseismic monitoring system was adopted to monitor the blasting seismic waves. The tests were divided into two parts: one part is for the analysis and research of time-frequency characteristics of the blasting source signals, and the other part is for the analysis of seismic wave attenuation law, as shown in Figure 2, the layout plan of the microseismic test station.(1)
Analysis of Time-Frequency Characteristics. From Figure 2, the monitoring points are arranged as a circle with a total of 8 stations. There are five blasting points; one blasting point is arranged in the center of the circle, and the other four blasting points are arranged in a square at a distance of 0.5 meters from the center of the circle. The distance between the 8 monitoring points and the center is 12 meters(2)
Study on the Attenuation Law of Seismic Waves. Monitoring points are arranged in a “straight line” pattern, with a total of 4 monitoring points. The first monitoring point is 12 meters away from the blasting source, and the distance between each monitoring point is 12 metersFigure 2
Microseismic monitoring schematic diagram of the test scheme.


## 2.2. Implementation of the Test Plan

### 2.2.1. Test Equipment and Installation

As shown in Figure3, the monitoring point equipment used in the experiment consists of three parts: sensor, collector, and battery. The sensor is installed by deep drilling, with a hole depth of 1 meter and a diameter of 90 mm. After the installation of the whole monitoring point, the monitoring software can be used to check its background noise. Generally, the background noise should be no more than E-5 V. The sensor and the hole wall adopt the yellow mud coupling.Figure 3
Equipment installation at site monitoring points.


### 2.2.2. Program Implementation

A total of 5 blasting tests were carried out in this test, and the phase change blasting parameters of liquid CO2 are shown in Table 1, which were obtained from the laboratory at Central South University (see the reference [2]). Test 1, test 2, and test 3 were carried out to discuss the time-frequency characteristics, and test 4 and test 5 were carried out to analyze the attenuation law of seismic waves. The recording sampling rate of the test data was 4000 Hz. A total of 8 microseismic instruments were used to record seismic wave data, and the sensitivity of the sensor was set at 200 V/m/s.Table 1
Phase change blasting parameters of liquid CO2.
Serial numberCO2 qualityStomatal directionBurst disc pressureBlasting rod10.452 kgEast and west150 MPa1 bar20.526 kgSouth and north150 MPa1 bar30.454 kgEast and west150 MPa2 bar40.498 kgEast and west150 MPa1 bar50.448 kgEast and west150 MPa2 bar


## 2.2.1. Test Equipment and Installation

As shown in Figure3, the monitoring point equipment used in the experiment consists of three parts: sensor, collector, and battery. The sensor is installed by deep drilling, with a hole depth of 1 meter and a diameter of 90 mm. After the installation of the whole monitoring point, the monitoring software can be used to check its background noise. Generally, the background noise should be no more than E-5 V. The sensor and the hole wall adopt the yellow mud coupling.Figure 3
Equipment installation at site monitoring points.


## 2.2.2. Program Implementation

A total of 5 blasting tests were carried out in this test, and the phase change blasting parameters of liquid CO2 are shown in Table 1, which were obtained from the laboratory at Central South University (see the reference [2]). Test 1, test 2, and test 3 were carried out to discuss the time-frequency characteristics, and test 4 and test 5 were carried out to analyze the attenuation law of seismic waves. The recording sampling rate of the test data was 4000 Hz. A total of 8 microseismic instruments were used to record seismic wave data, and the sensitivity of the sensor was set at 200 V/m/s.Table 1
Phase change blasting parameters of liquid CO2.
Serial numberCO2 qualityStomatal directionBurst disc pressureBlasting rod10.452 kgEast and west150 MPa1 bar20.526 kgSouth and north150 MPa1 bar30.454 kgEast and west150 MPa2 bar40.498 kgEast and west150 MPa1 bar50.448 kgEast and west150 MPa2 bar


## 3. Decomposition and Reconstruction of Signal and Time-Frequency Analysis Techniques

### 3.1. Wavelet Packet Decomposition and Reconstruction

The time-frequency local property of a seismic wave is the most fundamental and key property of a nonstationary signal. The traditional Fourier transform signal analysis cannot well describe this property. Therefore, in order to analyze and process a nonstationary signal, researchers have improved or created new signal analysis theories based on the Fourier transform. Wavelet analysis and wavelet packet analysis are widely used in seismic signal analysis. The understanding of wavelet analysis and wavelet packets analysis is illustrated by three layers structure diagram of the wavelet analysis tree and wavelet packet analysis tree shown in Figures4(a) and 4(b).Figure 4
Structure diagram of wavelet analysis and wavelet packet analysis decomposition.
(a)
Three-layer wavelet analysis tree structure diagram(b)
Three-layer wavelet packet analysis tree structure diagramWavelet transform has the characteristics of multiresolution analysis and can represent the local features of signals in both time and frequency domains. In Figure4, A represents the low frequency parts, and D represents the high frequency parts, and they are followed by the number of layers (i.e., scale number) of decomposition. As can be seen from Figure 4(a), the low-frequency part is continuously decomposed by wavelet analysis, while the high-frequency part is not considered. The final signal consists of A3+D3+D2+D1. Figure 4(b) shows the structure diagram of the wavelet packet analysis tree. The decomposition relation of signal [34] is S=AAA3+DAA3+ADA3+DDA3+AAD3+DAD3+ADD3+DDD3.From the previous analysis, it can be seen that the wavelet packet analysis can provide a more detailed analysis for the signals, by dividing the frequency band into multiple levels. Decompose the high-frequency part without subdivision of the multiresolution analysis, and adaptively select the corresponding frequency band according to the characteristics of the analyzed signal, so that it can match with the signal spectrum, thus improving the time-frequency resolution. In this paper, wavelet packet analysis is selected for signal decomposition and reconstruction.


### 3.2. Time-Frequency Analysis of Nonstationary Signals

Adaptive optimal kernel time-frequency analysis (AOK) is a nonlinear time-frequency distribution analysis method proposed by Jones and Baraniuk [26], which uses short-time fuzzy function and time-varying adaptive kernel function to distinguish the details of multicomponent signals in the time-frequency distribution. AOK is an optimal time-frequency method in time-frequency matching. Its kernel function changes in an adaptive manner with the change of time, which is characterized by the optimal accuracy of both cross term suppression and time-frequency localization. This ability of good localization in both time and frequency domains is suitable for the analysis of time-frequency characteristics of seismic signals. The AOK is used to analyze signal in the time-frequency domain.The AOK distribution of signals at asst can be expressed as follows.(1)PAOKt,f=12π∫−∞+∞∫−∞+∞At;θ,τΦoptt;θ,τe−jθτ−jτωdθdτ.In the equation,At;θ,τ is the short fuzzy function of the signal; Φoptt;θ,τ is the corresponding optimal kernel function.(2)At;θ,τ=∫−∞+∞hu−τ2wu−t−τ2hu+τ2fu−t+τ2ejθudu.In the formula,wu is the symmetric window function, and t is the central position of wu. When wu=0, only the signal in the range of t−T,t+T can calculate its kernel function. For any detail part of the signal, the short-time ambiguity function can be accurately described. With the definition of the short-time fuzzy function, it is easy to calculate the corresponding Φoptt;θ,τ. The short-time fuzzy function varies with time, so the optimal kernel also varies with time.Φoptt;θ,τ can be obtained by solving the following optimization problems:(3)maxΦ∫02π∫0+∞At;r,φΦt;r,φ2rdrdφ.The constraint conditions are(4)Φt;r,φ=exp−r22σ2φ12π∫02π∫0+∞Φt;r,φ2rdrdφ=12π∫02πσ2φdφ≤α,α≥0.In which,σφ is the extension of the radial Gaussian function in the direction of the radial angle φ, which is called the expansion function. The φ is the angle between the radial and the horizontal φ=arctanτ/θ, r=θ2+τ2. θ is the polar coordinate angle of ambiguity function, τ is the time interval, and f is frequency, which is the energy volume of radial Gaussian kernel function. If α is too small, the kernel function will filter out some self-components. If α is too large, the kernel function can not effectively remove the influence of cross components. The proper selection of α is taken according to the actual signal, and the range of values is generally 1≤α≤5 Gaussian window, and radial Gaussian function was adopted in time-frequency analysis. The energy volume variable 128×128 is 2 [35, 36]; the output resolution is 512.Nonlinear time-frequency representation (TFR) analysis method of a nonstationary signalst has the following properties [37]:(1)
TFR is a real value and positive, indicating the change of energy(2)
TFR gives the signal energy of the double integration of time and frequency, namely(5)E=∫−∞+∞∫−∞+∞TFRt,fdfdt.According to equation (5), time and frequency distribution (TFR) is the two-dimensional spatial distribution of signal energy in time and frequency, which has a clear physical meaning of joint distribution of signal energy in time and frequency domain.


### 3.3. Detailed Analysis of Energy Distribution Characteristics of Blasting Seismic Waves

Although AOK time-frequency analysis can obtain the main frequency information and the frequency range of energy concentration, it cannot show the detailed information of energy distribution in different frequency bands. Wavelet packet transform can be achieved by signal decomposition, and reconstruction can show good energy distribution of different frequency bands after detail, but the energy calculation is the signal amplitude linear summation, and the time-frequency distribution (TFR) is the signal energy in time and frequency; compared to 2 d space distribution, it does not have clear physical meaning when the signal energy in frequency domain on the joint distribution. In this paper, seismic waves are decomposed and reconstructed in different frequency bands by wavelet packet transformation, and then, the time-frequency characteristics and energy of reconstructed signals in each frequency band are analyzed in detail by using AOK time-frequency analysis method. The method combining wavelet packet transform and AOK time-frequency analysis can make up for each other’s shortcomings and accurately analyze the magnitude and distribution rule of different frequency band energy of liquid CO2 blasting seismic wave signals.Wavelet packet decomposition and reconstruction analysis of signals can be realized directly on the MATLAB platform. AOK time-frequency analysis technology adopts the AOK time-frequency analysis toolbox (TFTB, Time-Frequency Toolbox); a time-frequency analysis toolbox developed by professor (CNRS, The National Center for Scientific Research) François Auger [38] of France National Center for Scientific Research was used to calculate the band energy of reconstructed signals. The above analysis is all independent individual analysis, which needs to be further compiled on MATLAB software platform to complete the wavelet packet decomposition and reconstruction of seismic wave signal, AOK time-frequency analysis, and energy calculation. The flow chart of the MATLAB program for signal decomposition and reconstruction, AOK time-frequency analysis, and energy calculation is shown in Figure 5. Firstly, the signal is decomposed and reconstructed by the wavelet packets, and then, the energy calculation formula (5) is added to the AOK time-frequency analysis program to calculate the energy of each frequency band. Finally, the energy normalization, calculation, and analysis of the energy distribution in different frequency bands are carried out.Figure 5
Flow chart of the MATLAB program for signal reconstruction and AOK time-frequency analysis.


## 3.1. Wavelet Packet Decomposition and Reconstruction

The time-frequency local property of a seismic wave is the most fundamental and key property of a nonstationary signal. The traditional Fourier transform signal analysis cannot well describe this property. Therefore, in order to analyze and process a nonstationary signal, researchers have improved or created new signal analysis theories based on the Fourier transform. Wavelet analysis and wavelet packet analysis are widely used in seismic signal analysis. The understanding of wavelet analysis and wavelet packets analysis is illustrated by three layers structure diagram of the wavelet analysis tree and wavelet packet analysis tree shown in Figures4(a) and 4(b).Figure 4
Structure diagram of wavelet analysis and wavelet packet analysis decomposition.
(a)
Three-layer wavelet analysis tree structure diagram(b)
Three-layer wavelet packet analysis tree structure diagramWavelet transform has the characteristics of multiresolution analysis and can represent the local features of signals in both time and frequency domains. In Figure4, A represents the low frequency parts, and D represents the high frequency parts, and they are followed by the number of layers (i.e., scale number) of decomposition. As can be seen from Figure 4(a), the low-frequency part is continuously decomposed by wavelet analysis, while the high-frequency part is not considered. The final signal consists of A3+D3+D2+D1. Figure 4(b) shows the structure diagram of the wavelet packet analysis tree. The decomposition relation of signal [34] is S=AAA3+DAA3+ADA3+DDA3+AAD3+DAD3+ADD3+DDD3.From the previous analysis, it can be seen that the wavelet packet analysis can provide a more detailed analysis for the signals, by dividing the frequency band into multiple levels. Decompose the high-frequency part without subdivision of the multiresolution analysis, and adaptively select the corresponding frequency band according to the characteristics of the analyzed signal, so that it can match with the signal spectrum, thus improving the time-frequency resolution. In this paper, wavelet packet analysis is selected for signal decomposition and reconstruction.


## 3.2. Time-Frequency Analysis of Nonstationary Signals

Adaptive optimal kernel time-frequency analysis (AOK) is a nonlinear time-frequency distribution analysis method proposed by Jones and Baraniuk [26], which uses short-time fuzzy function and time-varying adaptive kernel function to distinguish the details of multicomponent signals in the time-frequency distribution. AOK is an optimal time-frequency method in time-frequency matching. Its kernel function changes in an adaptive manner with the change of time, which is characterized by the optimal accuracy of both cross term suppression and time-frequency localization. This ability of good localization in both time and frequency domains is suitable for the analysis of time-frequency characteristics of seismic signals. The AOK is used to analyze signal in the time-frequency domain.The AOK distribution of signals at asst can be expressed as follows.(1)PAOKt,f=12π∫−∞+∞∫−∞+∞At;θ,τΦoptt;θ,τe−jθτ−jτωdθdτ.In the equation,At;θ,τ is the short fuzzy function of the signal; Φoptt;θ,τ is the corresponding optimal kernel function.(2)At;θ,τ=∫−∞+∞hu−τ2wu−t−τ2hu+τ2fu−t+τ2ejθudu.In the formula,wu is the symmetric window function, and t is the central position of wu. When wu=0, only the signal in the range of t−T,t+T can calculate its kernel function. For any detail part of the signal, the short-time ambiguity function can be accurately described. With the definition of the short-time fuzzy function, it is easy to calculate the corresponding Φoptt;θ,τ. The short-time fuzzy function varies with time, so the optimal kernel also varies with time.Φoptt;θ,τ can be obtained by solving the following optimization problems:(3)maxΦ∫02π∫0+∞At;r,φΦt;r,φ2rdrdφ.The constraint conditions are(4)Φt;r,φ=exp−r22σ2φ12π∫02π∫0+∞Φt;r,φ2rdrdφ=12π∫02πσ2φdφ≤α,α≥0.In which,σφ is the extension of the radial Gaussian function in the direction of the radial angle φ, which is called the expansion function. The φ is the angle between the radial and the horizontal φ=arctanτ/θ, r=θ2+τ2. θ is the polar coordinate angle of ambiguity function, τ is the time interval, and f is frequency, which is the energy volume of radial Gaussian kernel function. If α is too small, the kernel function will filter out some self-components. If α is too large, the kernel function can not effectively remove the influence of cross components. The proper selection of α is taken according to the actual signal, and the range of values is generally 1≤α≤5 Gaussian window, and radial Gaussian function was adopted in time-frequency analysis. The energy volume variable 128×128 is 2 [35, 36]; the output resolution is 512.Nonlinear time-frequency representation (TFR) analysis method of a nonstationary signalst has the following properties [37]:(1)
TFR is a real value and positive, indicating the change of energy(2)
TFR gives the signal energy of the double integration of time and frequency, namely(5)E=∫−∞+∞∫−∞+∞TFRt,fdfdt.According to equation (5), time and frequency distribution (TFR) is the two-dimensional spatial distribution of signal energy in time and frequency, which has a clear physical meaning of joint distribution of signal energy in time and frequency domain.


## 3.3. Detailed Analysis of Energy Distribution Characteristics of Blasting Seismic Waves

Although AOK time-frequency analysis can obtain the main frequency information and the frequency range of energy concentration, it cannot show the detailed information of energy distribution in different frequency bands. Wavelet packet transform can be achieved by signal decomposition, and reconstruction can show good energy distribution of different frequency bands after detail, but the energy calculation is the signal amplitude linear summation, and the time-frequency distribution (TFR) is the signal energy in time and frequency; compared to 2 d space distribution, it does not have clear physical meaning when the signal energy in frequency domain on the joint distribution. In this paper, seismic waves are decomposed and reconstructed in different frequency bands by wavelet packet transformation, and then, the time-frequency characteristics and energy of reconstructed signals in each frequency band are analyzed in detail by using AOK time-frequency analysis method. The method combining wavelet packet transform and AOK time-frequency analysis can make up for each other’s shortcomings and accurately analyze the magnitude and distribution rule of different frequency band energy of liquid CO2 blasting seismic wave signals.Wavelet packet decomposition and reconstruction analysis of signals can be realized directly on the MATLAB platform. AOK time-frequency analysis technology adopts the AOK time-frequency analysis toolbox (TFTB, Time-Frequency Toolbox); a time-frequency analysis toolbox developed by professor (CNRS, The National Center for Scientific Research) François Auger [38] of France National Center for Scientific Research was used to calculate the band energy of reconstructed signals. The above analysis is all independent individual analysis, which needs to be further compiled on MATLAB software platform to complete the wavelet packet decomposition and reconstruction of seismic wave signal, AOK time-frequency analysis, and energy calculation. The flow chart of the MATLAB program for signal decomposition and reconstruction, AOK time-frequency analysis, and energy calculation is shown in Figure 5. Firstly, the signal is decomposed and reconstructed by the wavelet packets, and then, the energy calculation formula (5) is added to the AOK time-frequency analysis program to calculate the energy of each frequency band. Finally, the energy normalization, calculation, and analysis of the energy distribution in different frequency bands are carried out.Figure 5
Flow chart of the MATLAB program for signal reconstruction and AOK time-frequency analysis.


## 4. Results and Discussions

### 4.1. Basic Characteristics and Energy Distribution Analysis of Blasting Seismic Wave Signal

The microseismic monitors have been arranged as a circle in Figure2 of the test plan. Three blasting tests have been carried out. 24 sets of signals from blasting earthquakes have been received by the monitoring equipment of 8 stations, each of which has component data in three directions (D, B, and Z directions). In this paper, the time-frequency analysis of liquid CO2 blasting seismic wave signal is carried out by taking the typical test data obtained from the No. 1 test as an example.As shown in Figure6, the seismic wave signal of the No. 1 test is shown. It can be seen from the diagram that the direct wave energy of liquid CO2 blasting signal is strong, and the attenuation speed is fast. No obvious S wave is found in the blasting process, and the attenuation law of seismic signal waveform accords with the characteristics of the explosive blasting signal. The amplitude of blasting vibration of seismic wave signals in three directions is 0.01-0.04 cm/s, and the vibration velocity is at a lower level. In order to analyze and study seismic wave signals from time domain and frequency domain, this paper selects three components of M1 seismic wave signal in the No. 1 test as typical test data and uses AOK time-frequency analysis technology to analyze seismic wave signal based on signal analysis MATLAB processing platform. AOK time-frequency analysis method adopts program and toolbox compiled by Jones et al. [26] of (Rice University) of Rice University. The time-frequency distribution (TFR) contours and three-dimensional diagrams of the three-component signals of the M1 station of the No. 1 test station are shown in Figure 7 by using the toolbox.Figure 6
Seismic wave of blasting test 1.Figure 7
Time-frequency distribution (TFR) contours and three-dimensional maps of the three component signals for station M1 (in turn,D, B, and Z components, respectively).As can be seen from Figure7, the AOK time-frequency analysis method can localize the time domain and frequency domain of the three-component signals and can obtain the main frequency information and duration of the seismic wave signal. The main frequency of the three components of the M1 station in test 1 is 60 Hz, 54 Hz, and 47 Hz, and the frequency band with large energy is between 30 and 70 Hz. The duration is about 20-30 ms. The frequency range of the D component signal is between 0 and 250 Hz, the frequency range of the B component signal is between 0 and 125 Hz, the frequency range of the Z component signal is between 0 and 200 Hz, and the duration of three-component signals is about 0.1 s.The sampling frequency of the microseismic monitoring system is 4000 Hz. According to the sampling theorem, the sampling frequency of Nyquist is 2000 Hz. When the wavelet packets are used to process the seismic wave signal, the selection of wavelet basis function directly determines the accuracy of signal processing and analysis results. Daubechies wavelet (DB wavelet) basis function series can better reflect the unstable change process of the seismic wave signal in time and frequency distribution. db8 wavelet is often used to transform blasting vibration signal by wavelet packet transform. In the study of blasting seismic wave signals analysis, most of them adopted 8, 16, 32, 64 as the minimum decomposition frequency band. Combined with the Nyquist sampling frequency, the wavelet packet decomposition and reconstruction signal are decomposed by using the db8 wavelet series as the basis function. A total of 256 subbands is obtained, the minimum frequency band is 0.0000 Hz~7.8125 Hz, the highest frequency band is 1992.1875 Hz~2000 Hz, and the intermediate frequency band is increased with 7.8125 Hz as the equal difference sequence.The AOK time-frequency domain analysis of 256 seismic waves with different frequency bands after wavelet packet reconstruction is carried out. The margtfr function in the time-frequency analysis toolbox is used to calculate the time-frequency distribution TFR, timet, and frequency f after AOK time-frequency analysis. The return value is the energy of the signal, and the function expression is margt,margf,E=margtfrtfr,t,f, E as the energy of the signal.Based on the above analysis flow and analysis method, combined with the time-frequency analysis toolbox, the three components recorded in the M1 test station were decomposed and reconstructed by the wavelet packets, and then, AOK time-frequency analysis and energy calculation were carried out, and 256 band energy distribution percentages are obtained as shown in Table2 and Figure 8.Table 2
The percentage of energy distribution for three-component signals in different frequency bands.
Serial numberFrequency band (Hz)DBZFrequency bandDistribute (Hz)DBZDistribute (%)Distribute (%)Distribute (%)Distribute (%)Distribute (%)Distribute (%)10-7.81250.420.4801070.3125-78.1250.30.450.0427.8125-15.6250.560.830.041178.125-85.93751.481.20.18315.625-23.43755.862.594.231285.9375-93.750.650.340.11423.4375-31.251.411.140.741393.75-101.56259.033.180.38531.25-39.062511.3220.51.8214101.5625-109.3757.265.581.54639.0625-46.8756.0610.0513.9415109.375-117.18751.390.520.31746.875-54.68754.8530.6957.716117.1875-1254.225.781.19854.6875-62.532.7614.5717.517-32125-25011.041.50.24962.5-70.31250.180.520.0133-256250-20001.220.080.04Figure 8
The percentage of energy distribution for three-component signals of station M1 in different frequency bands.As can be seen from Table2 and Figure 8, the frequency band energy distribution of liquid CO2 blasting seismic wave signal is as follows:(1)
It can have arranged from the figure that the energy of the signal is mainly concentrated in the 1-16 frequency band (0-125 Hz). The energy distribution of the three signal components in the first 16 frequency bands is, respectively, 87.74%, 98.42%, and 99.72%, and the maximum value occurs in the frequency band where the main frequency is located. The energy distribution of 125 Hz-250 Hz is 11.04%, 1.5%, and 0.24%, respectively. After the frequency of 250 Hz, the energy distribution of the three components is 1.22%, 0.08%, and 0.04%, respectively. It indicates that the energy of the seismic wave signal is very concentrated, and the energy of the high-frequency part decays very rapidly, although the high-frequency part of the energy decays rapidly but still occupies a certain proportion of the distribution(2)
In the 16 frequency bands, the energy distribution is divided into two regions, the 1-8 frequency band region and the 11-16 frequency band region. The energy of the three-component signals distributed in the 1-8 frequency is 63.24%, 80.85%, and 95.97%. The energy of the three-component signals distributed in the 11-16 frequency is 24.02%, 16.59%, and 3.7%, respectively.


### 4.2. Energy Calculation of Measuring Point

According to the linear arrangement of the test scheme, a total of 4 monitoring points and two blasting tests were set up, and a total of 8 groups of test data were set up. Each group of test data included the vibration signal components inD, B, and Z directions. Because of the poor signal-to-noise ratio (SNR) of the signal recorded by station 2, a total of 18 data amounts were obtained in this experiment. According to the waveform information of the blasting signal recorded by the monitoring network, the PPV value of each component of each station is calculated. Table 3 is the PPV table of three directional components of the seismic wave in two explosion tests.Table 3
PPV of three direction components of seismic waves in two explosion tests.
Serial numberRecording stationD-PPV (mm/s)B-PPV (mm/s)Z-PPV (mm/s)4a10.0640.1110.127a20.03950.0670.039a40.0090.00550.00011595a10.040.10.107a20.020.04610.033a40.005230.00350.00012Figure9 shows the trend diagram that the PPV of the three directional components of seismic waves decreases with the increase of distance. It can be seen from the diagram that the whole vibration velocity is at a lower level. With the increase of distance, the vibration velocity decays rapidly, in which the attenuation speed of vibration velocity in D and B directions is first fast and then slow; 24 meters is the inflection point of velocity attenuation rate, and the attenuation rate of vibration velocity in Z direction shows a trend of rapid decrease. When the vibration velocity decreases to the micron level at 48 meters, the vibration velocity in three directions tends to be similar.Figure 9
Relationship between PPV and distance of three directional components of seismic waves in two explosion tests.As can be seen from Table4, the AOK time-frequency distribution method can localize the time domain and frequency domain of the three components of the signal, which can obtain the main frequency information and duration of the seismic wave signal. The main frequency of the three-component signals of the 4 test a1, a2, and a4 stations is mainly distributed between 30 and 53 Hz, the duration is about 20-50 ms, all the seismic wave component signals are concentrated in the frequency domain, and the frequency range is between 0 and 250 Hz. The duration of the signal is about 0.1 s. With the increase of distance, the main frequency decreases, and the energy decreases gradually.Table 4
Main frequency of three directional components of seismic wave in explosion test 4.
Serial numberRecording stationD (Hz)B (Hz)Z (Hz)4a1534147a2413950a4303935In the process of seismic wave propagation, it is very difficult to accurately calculate the seismic wave energy density at each geophone position. In engineering, the maximum peak vibration velocity and the square (discrete signal) of the amplitude of each sampling point in the seismic wave duration are usually used to represent the energy of the measuring point. Table5 shows the seismic wave energy parameters of the liquid CO2 blasting system, and Figure 10 shows the correlation curve between the peak vibration velocity of the measuring point and the double integral energy of TFR to time and frequency. It can be seen from the diagram that the correlation curve of the two parameters has a high linear correlation. The correlation coefficient is R2=0.90. Therefore, it can be considered that the peak velocity of the measuring point can reflect the energy of the measuring point for liquid CO2 and explosive blasting seismic wave, and it is not necessary to calculate the energy obtained by the double integration of time and frequency of TFR in the whole time period of the event. In the process of analyzing the law of seismic wave energy attenuation in liquid CO2 blasting system, on the one hand, it can automatically pick up and calculate the peak velocity and reduce the calculation workload; on the other hand, it can improve the efficiency of evaluating the damage and attenuation degree of seismic wave energy.Table 5
Seismic wave energy parameters of the liquid CO2 blasting system.
Measure pointDistance (m)PPV (cm/s)Energy (cm2/s2)DBZDBZ4-a1120.00640.01110.01270.01560.0292810.0228424-a2240.003950.00670.00390.0040702860.0111590.0060344-a4480.00090.000550.000011590.0004817860.0002054.99E-08Figure 10
Regression curve between peak vibration velocity and energy.


### 4.3. Energy Attenuation Law

The regression models that can be used for seismic wave energy attenuation can be divided into two categories [39]: (1) exponential function form, see formula (6), and (2) power function form, see formula (7).(6)E=E0e−αr,(7)E=E0r−α.In the formula,r is the distance between the measuring point and the source, E0 is the initial energy of the source, E is the energy at r, and α is the attenuation coefficient. Here, E0 is only the initial energy of the source obtained by the regression curve or the initial energy of the virtual source.Use the data in Table5 to compare the exponential function with the fitting curve of the power function according to the above formula. As shown in Figure 11, the parameter coefficients and correlation coefficients of the exponential function model and the power function are shown in Table 6. It can be seen from the table that both the exponential function model and the power function model are suitable to describe the energy attenuation law of seismic waves produced by the liquid CO2 blasting system. From the point of view of correlation coefficient alone, the exponential function model is more suitable for seismic wave attenuation law of liquid CO2 blasting system, but the initial energy E0 of the source is small or has little difference with the energy at the measuring point of 12 meters and is in the same order of magnitude. Although the correlation coefficient of the power function model is smaller than that of the exponential function model, it is closer to the real value from the initial energy E0 of the source. To sum up, the power function is more suitable to describe this typical attenuation law of explosive seismic wave: with the increase of the distance from the source, the early attenuation is rapid, and the late attenuation is slow. The attenuation of explosion seismic wave energy in lithe quid CO2 blasting system is similar to that of explosive explosion seismic wave energy.Figure 11
Regression model of energy with distance.
(a)
Exponential function model(b)
Power function modelTable 6
Parameter of exponential function model and power function model from the regression curves in Figure11.
FormulaE=E0e−αrFormulaE=E0r−αR20.9990.9920.999R20.9960.9490.986DE00.05929DE02.3994α0.1112α2.0253BE00.08088BE01.79948α0.08434α1.65352ZE00.08777ZE04.00521α0.11212α2.078Although the data are few here, while the relationship between energy and distance in this paper conforms to the general law of energy characteristics of blasting seismic waves. Therefore, in a sense, we just proved this relationship, if data points are more in the future work, the other analytical method such as the analysis of variance (ANOVA) can be conducted to deep investigate the influence of distance on energy.


### 4.4. Source Energy Regression Model

#### 4.4.1. The Empirical Formula of Vibration Velocity of Sadovsky

The velocity attenuation law of particle vibration in blasting engineering is commonly expressed by Sadovsky’s empirical formula [39], that is,(8)V=KQ3rα=V0r−α.In the formula, the particle vibration velocity is the particle vibration velocity, the cm/s;V is the site coefficient; K is the charge, kg; r and α are the same as above; V0 is the initial vibration velocity of the source; and the same V0 of cm/s is the vibration velocity of the virtual source here. For the data in Table 5, combined with the formula (8), the average energy of liquid CO2 blasting tube in reference [2] is 0.030 kg TNT. The regression analysis shows that K=4.924 and α=1.0526.


#### 4.4.2. Energy Conversion Factor

It is not convenient for Sadovsky’s formula to be directly used to calculate the energy of the microseismic source. The site coefficientK is related to the source medium and blasting parameters. The seismic wave energy corresponding to explosive quantity Q also needs to be converted. In solid media, only a very small part of explosive explosion energy is converted into the seismic wave. Here is a concept of energy conversion coefficient η that needs to be explained and defined. The ratio of seismic wave energy to total energy produced by the explosive explosion is defined as the seismic wave energy conversion coefficient of explosive blasting. According to this theory, this paper also represents the ratio of seismic wave energy Ec to TNT equivalent total energy produced by liquid CO2 blasting according to this theory that liquid CO2 blasting system has a similar seismic wave energy conversion coefficient ET[39].(9)ηc=EcET.Through theoretical analysis and a large number of experimental data, the energy conversion coefficient of explosive was discussed and analyzed in detail and put forward the formula for calculating the energy conversion coefficient of explosive earthquake based on the statistical parametersK and α based on explosive seismic effect (the unit of Q is m of kg, rand the unit of V is cm/s) and the formula for the energy conversion coefficient of explosive as follow [39, 40]:(10)ηc=K∙10−23/α∙10−3.The seismic wave energy conversion coefficientC=1.88E−07 of the liquid CO2 blasting system is obtained by using A=4.924 and B=1.0526 generation (15).


#### 4.4.3. Source Energy Characteristic Coefficient

The total energy generated by the explosion of TNT with a mass of Q is(11)ET=Q∗QV.QV is the explosion heat of TNT, 4150.2 kJ/kg.The seismic wave energy is obtained by replacing formula (11) with (9).(12)Ec=Q∗QV∗ηc.Then, theQ, QV, and ηc parameters are substituted into (12) to obtain the elastic wave energy of liquid CO2 blasting source Ec=0.0232 J. Replace formula (12) with (8) to get(13)V=KQVηcα/3Ec3rα=V0Ecα/3Ec3rα.Assume(14)Kc=KQVηcα/3=V0Ecα/3.TheKc is defined as the energy characteristic coefficient of the source, and formula (10) is obtained by replacing formula (14) with the energy characteristic coefficient of the source.(15)Kc=102+αQVα/3.From formulas (14) and (15), it can be seen that the value of Kc is mainly related to the dynamic characteristics of the source medium but independent of the specific blasting mode. That is to say, Kc reflects the proportional relationship between the source elastic wave energy Ec and the initial peak vibration velocity V0 of the source. Under the same experimental conditions, the stability ofKcis better than that of site coefficient. The energy characteristic coefficient Kc of the source of liquid CO2 blasting system is 5.376 by replacing K, QV, ηc, and (14).


#### 4.4.4. Source Energy Regression Model

The relationship between the peak vibration velocity of particles and the energy of the seismic wave is obtained by replacing formula (14) with (13).(16)V=KcEc3rα.After determining the energy characteristic coefficientKc, formula (16) is used as the regression model of seismic wave energy attenuation in the liquid CO2 blasting system. Taking the peak vibration velocity V and attenuation coefficient α as parameters, the source energy of liquid CO2 blasting system is directly obtained by the least square method, and the seismic wave source energy Ec=0.0234 J is obtained by regression of Table 5 experimental data. The deviation of elastic wave energy between liquid CO2 blasting source and formula (12) is 8.55%.


## 4.1. Basic Characteristics and Energy Distribution Analysis of Blasting Seismic Wave Signal

The microseismic monitors have been arranged as a circle in Figure2 of the test plan. Three blasting tests have been carried out. 24 sets of signals from blasting earthquakes have been received by the monitoring equipment of 8 stations, each of which has component data in three directions (D, B, and Z directions). In this paper, the time-frequency analysis of liquid CO2 blasting seismic wave signal is carried out by taking the typical test data obtained from the No. 1 test as an example.As shown in Figure6, the seismic wave signal of the No. 1 test is shown. It can be seen from the diagram that the direct wave energy of liquid CO2 blasting signal is strong, and the attenuation speed is fast. No obvious S wave is found in the blasting process, and the attenuation law of seismic signal waveform accords with the characteristics of the explosive blasting signal. The amplitude of blasting vibration of seismic wave signals in three directions is 0.01-0.04 cm/s, and the vibration velocity is at a lower level. In order to analyze and study seismic wave signals from time domain and frequency domain, this paper selects three components of M1 seismic wave signal in the No. 1 test as typical test data and uses AOK time-frequency analysis technology to analyze seismic wave signal based on signal analysis MATLAB processing platform. AOK time-frequency analysis method adopts program and toolbox compiled by Jones et al. [26] of (Rice University) of Rice University. The time-frequency distribution (TFR) contours and three-dimensional diagrams of the three-component signals of the M1 station of the No. 1 test station are shown in Figure 7 by using the toolbox.Figure 6
Seismic wave of blasting test 1.Figure 7
Time-frequency distribution (TFR) contours and three-dimensional maps of the three component signals for station M1 (in turn,D, B, and Z components, respectively).As can be seen from Figure7, the AOK time-frequency analysis method can localize the time domain and frequency domain of the three-component signals and can obtain the main frequency information and duration of the seismic wave signal. The main frequency of the three components of the M1 station in test 1 is 60 Hz, 54 Hz, and 47 Hz, and the frequency band with large energy is between 30 and 70 Hz. The duration is about 20-30 ms. The frequency range of the D component signal is between 0 and 250 Hz, the frequency range of the B component signal is between 0 and 125 Hz, the frequency range of the Z component signal is between 0 and 200 Hz, and the duration of three-component signals is about 0.1 s.The sampling frequency of the microseismic monitoring system is 4000 Hz. According to the sampling theorem, the sampling frequency of Nyquist is 2000 Hz. When the wavelet packets are used to process the seismic wave signal, the selection of wavelet basis function directly determines the accuracy of signal processing and analysis results. Daubechies wavelet (DB wavelet) basis function series can better reflect the unstable change process of the seismic wave signal in time and frequency distribution. db8 wavelet is often used to transform blasting vibration signal by wavelet packet transform. In the study of blasting seismic wave signals analysis, most of them adopted 8, 16, 32, 64 as the minimum decomposition frequency band. Combined with the Nyquist sampling frequency, the wavelet packet decomposition and reconstruction signal are decomposed by using the db8 wavelet series as the basis function. A total of 256 subbands is obtained, the minimum frequency band is 0.0000 Hz~7.8125 Hz, the highest frequency band is 1992.1875 Hz~2000 Hz, and the intermediate frequency band is increased with 7.8125 Hz as the equal difference sequence.The AOK time-frequency domain analysis of 256 seismic waves with different frequency bands after wavelet packet reconstruction is carried out. The margtfr function in the time-frequency analysis toolbox is used to calculate the time-frequency distribution TFR, timet, and frequency f after AOK time-frequency analysis. The return value is the energy of the signal, and the function expression is margt,margf,E=margtfrtfr,t,f, E as the energy of the signal.Based on the above analysis flow and analysis method, combined with the time-frequency analysis toolbox, the three components recorded in the M1 test station were decomposed and reconstructed by the wavelet packets, and then, AOK time-frequency analysis and energy calculation were carried out, and 256 band energy distribution percentages are obtained as shown in Table2 and Figure 8.Table 2
The percentage of energy distribution for three-component signals in different frequency bands.
Serial numberFrequency band (Hz)DBZFrequency bandDistribute (Hz)DBZDistribute (%)Distribute (%)Distribute (%)Distribute (%)Distribute (%)Distribute (%)10-7.81250.420.4801070.3125-78.1250.30.450.0427.8125-15.6250.560.830.041178.125-85.93751.481.20.18315.625-23.43755.862.594.231285.9375-93.750.650.340.11423.4375-31.251.411.140.741393.75-101.56259.033.180.38531.25-39.062511.3220.51.8214101.5625-109.3757.265.581.54639.0625-46.8756.0610.0513.9415109.375-117.18751.390.520.31746.875-54.68754.8530.6957.716117.1875-1254.225.781.19854.6875-62.532.7614.5717.517-32125-25011.041.50.24962.5-70.31250.180.520.0133-256250-20001.220.080.04Figure 8
The percentage of energy distribution for three-component signals of station M1 in different frequency bands.As can be seen from Table2 and Figure 8, the frequency band energy distribution of liquid CO2 blasting seismic wave signal is as follows:(1)
It can have arranged from the figure that the energy of the signal is mainly concentrated in the 1-16 frequency band (0-125 Hz). The energy distribution of the three signal components in the first 16 frequency bands is, respectively, 87.74%, 98.42%, and 99.72%, and the maximum value occurs in the frequency band where the main frequency is located. The energy distribution of 125 Hz-250 Hz is 11.04%, 1.5%, and 0.24%, respectively. After the frequency of 250 Hz, the energy distribution of the three components is 1.22%, 0.08%, and 0.04%, respectively. It indicates that the energy of the seismic wave signal is very concentrated, and the energy of the high-frequency part decays very rapidly, although the high-frequency part of the energy decays rapidly but still occupies a certain proportion of the distribution(2)
In the 16 frequency bands, the energy distribution is divided into two regions, the 1-8 frequency band region and the 11-16 frequency band region. The energy of the three-component signals distributed in the 1-8 frequency is 63.24%, 80.85%, and 95.97%. The energy of the three-component signals distributed in the 11-16 frequency is 24.02%, 16.59%, and 3.7%, respectively.


## 4.2. Energy Calculation of Measuring Point

According to the linear arrangement of the test scheme, a total of 4 monitoring points and two blasting tests were set up, and a total of 8 groups of test data were set up. Each group of test data included the vibration signal components inD, B, and Z directions. Because of the poor signal-to-noise ratio (SNR) of the signal recorded by station 2, a total of 18 data amounts were obtained in this experiment. According to the waveform information of the blasting signal recorded by the monitoring network, the PPV value of each component of each station is calculated. Table 3 is the PPV table of three directional components of the seismic wave in two explosion tests.Table 3
PPV of three direction components of seismic waves in two explosion tests.
Serial numberRecording stationD-PPV (mm/s)B-PPV (mm/s)Z-PPV (mm/s)4a10.0640.1110.127a20.03950.0670.039a40.0090.00550.00011595a10.040.10.107a20.020.04610.033a40.005230.00350.00012Figure9 shows the trend diagram that the PPV of the three directional components of seismic waves decreases with the increase of distance. It can be seen from the diagram that the whole vibration velocity is at a lower level. With the increase of distance, the vibration velocity decays rapidly, in which the attenuation speed of vibration velocity in D and B directions is first fast and then slow; 24 meters is the inflection point of velocity attenuation rate, and the attenuation rate of vibration velocity in Z direction shows a trend of rapid decrease. When the vibration velocity decreases to the micron level at 48 meters, the vibration velocity in three directions tends to be similar.Figure 9
Relationship between PPV and distance of three directional components of seismic waves in two explosion tests.As can be seen from Table4, the AOK time-frequency distribution method can localize the time domain and frequency domain of the three components of the signal, which can obtain the main frequency information and duration of the seismic wave signal. The main frequency of the three-component signals of the 4 test a1, a2, and a4 stations is mainly distributed between 30 and 53 Hz, the duration is about 20-50 ms, all the seismic wave component signals are concentrated in the frequency domain, and the frequency range is between 0 and 250 Hz. The duration of the signal is about 0.1 s. With the increase of distance, the main frequency decreases, and the energy decreases gradually.Table 4
Main frequency of three directional components of seismic wave in explosion test 4.
Serial numberRecording stationD (Hz)B (Hz)Z (Hz)4a1534147a2413950a4303935In the process of seismic wave propagation, it is very difficult to accurately calculate the seismic wave energy density at each geophone position. In engineering, the maximum peak vibration velocity and the square (discrete signal) of the amplitude of each sampling point in the seismic wave duration are usually used to represent the energy of the measuring point. Table5 shows the seismic wave energy parameters of the liquid CO2 blasting system, and Figure 10 shows the correlation curve between the peak vibration velocity of the measuring point and the double integral energy of TFR to time and frequency. It can be seen from the diagram that the correlation curve of the two parameters has a high linear correlation. The correlation coefficient is R2=0.90. Therefore, it can be considered that the peak velocity of the measuring point can reflect the energy of the measuring point for liquid CO2 and explosive blasting seismic wave, and it is not necessary to calculate the energy obtained by the double integration of time and frequency of TFR in the whole time period of the event. In the process of analyzing the law of seismic wave energy attenuation in liquid CO2 blasting system, on the one hand, it can automatically pick up and calculate the peak velocity and reduce the calculation workload; on the other hand, it can improve the efficiency of evaluating the damage and attenuation degree of seismic wave energy.Table 5
Seismic wave energy parameters of the liquid CO2 blasting system.
Measure pointDistance (m)PPV (cm/s)Energy (cm2/s2)DBZDBZ4-a1120.00640.01110.01270.01560.0292810.0228424-a2240.003950.00670.00390.0040702860.0111590.0060344-a4480.00090.000550.000011590.0004817860.0002054.99E-08Figure 10
Regression curve between peak vibration velocity and energy.


## 4.3. Energy Attenuation Law

The regression models that can be used for seismic wave energy attenuation can be divided into two categories [39]: (1) exponential function form, see formula (6), and (2) power function form, see formula (7).(6)E=E0e−αr,(7)E=E0r−α.In the formula,r is the distance between the measuring point and the source, E0 is the initial energy of the source, E is the energy at r, and α is the attenuation coefficient. Here, E0 is only the initial energy of the source obtained by the regression curve or the initial energy of the virtual source.Use the data in Table5 to compare the exponential function with the fitting curve of the power function according to the above formula. As shown in Figure 11, the parameter coefficients and correlation coefficients of the exponential function model and the power function are shown in Table 6. It can be seen from the table that both the exponential function model and the power function model are suitable to describe the energy attenuation law of seismic waves produced by the liquid CO2 blasting system. From the point of view of correlation coefficient alone, the exponential function model is more suitable for seismic wave attenuation law of liquid CO2 blasting system, but the initial energy E0 of the source is small or has little difference with the energy at the measuring point of 12 meters and is in the same order of magnitude. Although the correlation coefficient of the power function model is smaller than that of the exponential function model, it is closer to the real value from the initial energy E0 of the source. To sum up, the power function is more suitable to describe this typical attenuation law of explosive seismic wave: with the increase of the distance from the source, the early attenuation is rapid, and the late attenuation is slow. The attenuation of explosion seismic wave energy in lithe quid CO2 blasting system is similar to that of explosive explosion seismic wave energy.Figure 11
Regression model of energy with distance.
(a)
Exponential function model(b)
Power function modelTable 6
Parameter of exponential function model and power function model from the regression curves in Figure11.
FormulaE=E0e−αrFormulaE=E0r−αR20.9990.9920.999R20.9960.9490.986DE00.05929DE02.3994α0.1112α2.0253BE00.08088BE01.79948α0.08434α1.65352ZE00.08777ZE04.00521α0.11212α2.078Although the data are few here, while the relationship between energy and distance in this paper conforms to the general law of energy characteristics of blasting seismic waves. Therefore, in a sense, we just proved this relationship, if data points are more in the future work, the other analytical method such as the analysis of variance (ANOVA) can be conducted to deep investigate the influence of distance on energy.


## 4.4. Source Energy Regression Model

### 4.4.1. The Empirical Formula of Vibration Velocity of Sadovsky

The velocity attenuation law of particle vibration in blasting engineering is commonly expressed by Sadovsky’s empirical formula [39], that is,(8)V=KQ3rα=V0r−α.In the formula, the particle vibration velocity is the particle vibration velocity, the cm/s;V is the site coefficient; K is the charge, kg; r and α are the same as above; V0 is the initial vibration velocity of the source; and the same V0 of cm/s is the vibration velocity of the virtual source here. For the data in Table 5, combined with the formula (8), the average energy of liquid CO2 blasting tube in reference [2] is 0.030 kg TNT. The regression analysis shows that K=4.924 and α=1.0526.


### 4.4.2. Energy Conversion Factor

It is not convenient for Sadovsky’s formula to be directly used to calculate the energy of the microseismic source. The site coefficientK is related to the source medium and blasting parameters. The seismic wave energy corresponding to explosive quantity Q also needs to be converted. In solid media, only a very small part of explosive explosion energy is converted into the seismic wave. Here is a concept of energy conversion coefficient η that needs to be explained and defined. The ratio of seismic wave energy to total energy produced by the explosive explosion is defined as the seismic wave energy conversion coefficient of explosive blasting. According to this theory, this paper also represents the ratio of seismic wave energy Ec to TNT equivalent total energy produced by liquid CO2 blasting according to this theory that liquid CO2 blasting system has a similar seismic wave energy conversion coefficient ET[39].(9)ηc=EcET.Through theoretical analysis and a large number of experimental data, the energy conversion coefficient of explosive was discussed and analyzed in detail and put forward the formula for calculating the energy conversion coefficient of explosive earthquake based on the statistical parametersK and α based on explosive seismic effect (the unit of Q is m of kg, rand the unit of V is cm/s) and the formula for the energy conversion coefficient of explosive as follow [39, 40]:(10)ηc=K∙10−23/α∙10−3.The seismic wave energy conversion coefficientC=1.88E−07 of the liquid CO2 blasting system is obtained by using A=4.924 and B=1.0526 generation (15).


### 4.4.3. Source Energy Characteristic Coefficient

The total energy generated by the explosion of TNT with a mass of Q is(11)ET=Q∗QV.QV is the explosion heat of TNT, 4150.2 kJ/kg.The seismic wave energy is obtained by replacing formula (11) with (9).(12)Ec=Q∗QV∗ηc.Then, theQ, QV, and ηc parameters are substituted into (12) to obtain the elastic wave energy of liquid CO2 blasting source Ec=0.0232 J. Replace formula (12) with (8) to get(13)V=KQVηcα/3Ec3rα=V0Ecα/3Ec3rα.Assume(14)Kc=KQVηcα/3=V0Ecα/3.TheKc is defined as the energy characteristic coefficient of the source, and formula (10) is obtained by replacing formula (14) with the energy characteristic coefficient of the source.(15)Kc=102+αQVα/3.From formulas (14) and (15), it can be seen that the value of Kc is mainly related to the dynamic characteristics of the source medium but independent of the specific blasting mode. That is to say, Kc reflects the proportional relationship between the source elastic wave energy Ec and the initial peak vibration velocity V0 of the source. Under the same experimental conditions, the stability ofKcis better than that of site coefficient. The energy characteristic coefficient Kc of the source of liquid CO2 blasting system is 5.376 by replacing K, QV, ηc, and (14).


### 4.4.4. Source Energy Regression Model

The relationship between the peak vibration velocity of particles and the energy of the seismic wave is obtained by replacing formula (14) with (13).(16)V=KcEc3rα.After determining the energy characteristic coefficientKc, formula (16) is used as the regression model of seismic wave energy attenuation in the liquid CO2 blasting system. Taking the peak vibration velocity V and attenuation coefficient α as parameters, the source energy of liquid CO2 blasting system is directly obtained by the least square method, and the seismic wave source energy Ec=0.0234 J is obtained by regression of Table 5 experimental data. The deviation of elastic wave energy between liquid CO2 blasting source and formula (12) is 8.55%.


## 4.4.1. The Empirical Formula of Vibration Velocity of Sadovsky

The velocity attenuation law of particle vibration in blasting engineering is commonly expressed by Sadovsky’s empirical formula [39], that is,(8)V=KQ3rα=V0r−α.In the formula, the particle vibration velocity is the particle vibration velocity, the cm/s;V is the site coefficient; K is the charge, kg; r and α are the same as above; V0 is the initial vibration velocity of the source; and the same V0 of cm/s is the vibration velocity of the virtual source here. For the data in Table 5, combined with the formula (8), the average energy of liquid CO2 blasting tube in reference [2] is 0.030 kg TNT. The regression analysis shows that K=4.924 and α=1.0526.


## 4.4.2. Energy Conversion Factor

It is not convenient for Sadovsky’s formula to be directly used to calculate the energy of the microseismic source. The site coefficientK is related to the source medium and blasting parameters. The seismic wave energy corresponding to explosive quantity Q also needs to be converted. In solid media, only a very small part of explosive explosion energy is converted into the seismic wave. Here is a concept of energy conversion coefficient η that needs to be explained and defined. The ratio of seismic wave energy to total energy produced by the explosive explosion is defined as the seismic wave energy conversion coefficient of explosive blasting. According to this theory, this paper also represents the ratio of seismic wave energy Ec to TNT equivalent total energy produced by liquid CO2 blasting according to this theory that liquid CO2 blasting system has a similar seismic wave energy conversion coefficient ET[39].(9)ηc=EcET.Through theoretical analysis and a large number of experimental data, the energy conversion coefficient of explosive was discussed and analyzed in detail and put forward the formula for calculating the energy conversion coefficient of explosive earthquake based on the statistical parametersK and α based on explosive seismic effect (the unit of Q is m of kg, rand the unit of V is cm/s) and the formula for the energy conversion coefficient of explosive as follow [39, 40]:(10)ηc=K∙10−23/α∙10−3.The seismic wave energy conversion coefficientC=1.88E−07 of the liquid CO2 blasting system is obtained by using A=4.924 and B=1.0526 generation (15).


## 4.4.3. Source Energy Characteristic Coefficient

The total energy generated by the explosion of TNT with a mass of Q is(11)ET=Q∗QV.QV is the explosion heat of TNT, 4150.2 kJ/kg.The seismic wave energy is obtained by replacing formula (11) with (9).(12)Ec=Q∗QV∗ηc.Then, theQ, QV, and ηc parameters are substituted into (12) to obtain the elastic wave energy of liquid CO2 blasting source Ec=0.0232 J. Replace formula (12) with (8) to get(13)V=KQVηcα/3Ec3rα=V0Ecα/3Ec3rα.Assume(14)Kc=KQVηcα/3=V0Ecα/3.TheKc is defined as the energy characteristic coefficient of the source, and formula (10) is obtained by replacing formula (14) with the energy characteristic coefficient of the source.(15)Kc=102+αQVα/3.From formulas (14) and (15), it can be seen that the value of Kc is mainly related to the dynamic characteristics of the source medium but independent of the specific blasting mode. That is to say, Kc reflects the proportional relationship between the source elastic wave energy Ec and the initial peak vibration velocity V0 of the source. Under the same experimental conditions, the stability ofKcis better than that of site coefficient. The energy characteristic coefficient Kc of the source of liquid CO2 blasting system is 5.376 by replacing K, QV, ηc, and (14).


## 4.4.4. Source Energy Regression Model

The relationship between the peak vibration velocity of particles and the energy of the seismic wave is obtained by replacing formula (14) with (13).(16)V=KcEc3rα.After determining the energy characteristic coefficientKc, formula (16) is used as the regression model of seismic wave energy attenuation in the liquid CO2 blasting system. Taking the peak vibration velocity V and attenuation coefficient α as parameters, the source energy of liquid CO2 blasting system is directly obtained by the least square method, and the seismic wave source energy Ec=0.0234 J is obtained by regression of Table 5 experimental data. The deviation of elastic wave energy between liquid CO2 blasting source and formula (12) is 8.55%.


## 5. Conclusions

In order to study the time-frequency and energy characteristics of seismic waves in the liquid CO2 blasting system, the liquid CO2 blasting seismic wave signal is obtained by field microseismic monitoring test. According to the nonstationary characteristics of the blasting seismic wave signal, the basic time-frequency characteristics of the seismic wave signal are studied by adaptive optimal kernel time-frequency analysis method. On this basis, combined with wavelet packet transform decomposition and reconstruction technology and adaptive optimal kernel time-frequency analysis method, the energy distribution of the seismic wave signal is analyzed in detail. The main conclusions are as follows:(1)
The direct wave energy of the blasting signal is strong, and the attenuation speed is fast. There is no obvious S wave in the blasting process. The amplitude of blasting vibration of seismic wave signals in three directions is 0.01-0.04 m/s at about 12 m, and the vibration velocity is at a lower level. The main frequency range is between 30 and 70 Hz, and the duration is about 20-30 ms. The energy is mainly distributed in 0-125 Hz, the main frequency appears in the frequency band with the maximum energy, and the energy in the high-frequency part decays rapidly, and there are two main regions in the energy distribution, which indicates that two different peak pressures will be produced in the process of liquid CO2 blasting, which is consistent with the experimental data of the pressure response of the free explosion field in reference [2](2)
The attenuation of explosion seismic wave energy in the liquid CO2 blasting system is similar to that of explosive explosion seismic wave energy. The power function model can be used to describe the attenuation law of seismic wave energy: with the increase of distance from the source, the early attenuation is rapid, and the later attenuation is slow(3)
The energy conversion coefficient and characteristic coefficient of the source of liquid CO2 blasting system are defined and analyzed. Combined with the empirical formula of Sadovsky vibration velocity, the energy regression model of seismic wave source of liquid CO2 blasting system is given



---

*Source: 9919764-2021-10-13.xml*