Electronic detonator misfire identification method based on comparison of multiple characteristic information of blasting vibration

Abstract

The present invention discloses a method for identifying misfired electronic detonators based on comparison of multiple characteristic information of blasting vibration. The method performs EMMD signal decomposition based on the original blast source information, establishes an AR model to extract multiple characteristic information, and uploads the model and multiple characteristic information to the data center to form a basic database of blasting vibration. At the same time, the original blast source information of multiple blasting sites is collected and its characteristics are extracted to establish a large database of blasting vibration signals under different delay intervals, typical hole network parameters, and charge structures. Then, the actual blasting vibration signal is collected and the actual multiple characteristics are extracted, and they are compared with the multiple characteristics in the large database of blasting vibration signals. When there are characteristic differences, it is determined that there are misfired detonators in this blasting project. The present invention improves the blasting construction efficiency, reduces the great danger in the process of misfired detonator treatment, and improves the safety of using digital electronic detonators at the same time.

Description

Technical Field

[0001] The present invention relates to the technical field of misfire identification, and particularly to an intelligent misfire identification method for digital electronic detonators based on the comparison of multiple characteristic information of blasting vibration and the technology of acoustic-vibration fusion. Background Art

[0002] Digital electronic detonators utilize integrated circuit technology and have the function of precise delay initiation, realizing the refinement of blasting technology. With the increase in labor costs and the decrease in the cost of microelectronic technology, and the mass production of high-precision electronic detonators, digital electronic detonators have gradually replaced detonator tubes in fields such as controlled blasting in complex environments, large-scale mine exploitation, and rare mineral development, and have the characteristics of online detectability of the blasting network, precise delay, high reliability, and good safety.

[0003] During the use of digital electronic detonators, safety accidents caused by misfires occur, reducing the blasting construction efficiency, and the process of dealing with misfires is extremely dangerous.

[0004] For example, the patent application with the application number CN201310296680.1 discloses a misfire identification method in blasting engineering. A complete blasting vibration signal is obtained through a blasting vibration tester. The blasting vibration signal is analyzed by the time-energy density method based on wavelet transform to obtain the blasting vibration time-energy density curves at different decomposition scales. By comparing to determine the blasting vibration time-energy density curve at the optimal decomposition scale among them, a window function is used to smooth the blasting vibration time-energy density curve, and the window function bandwidth is selected by the window closing method. Then, the smoothed time-energy density curve that can more clearly reflect the actual initiation time of each section of detonators is selected. By comparing the peak points of the curve with the blasting design, it can be determined whether a misfire occurs and the location where the misfire occurs. This solution analyzes the misfire location based on the on-site blasting vibration signal, and its signal is greatly interfered by the environment, cannot accurately analyze whether a misfire occurs, and the analysis time is long, resulting in low blasting construction efficiency. Summary of the Invention

[0005] The purpose of the present invention is to overcome the deficiencies of the prior art and provide an electronic detonator misfire identification method based on the comparison of multiple characteristic information of blasting vibration. By designing a wavelet function to establish a prediction model and extract multiple characteristics of blasting information, and constructing a large database of blasting vibration based on this. During blasting construction, only by extracting the multiple characteristics of the actual on-site blasting vibration signal for characteristic comparison can it be determined whether a misfire occurs in the current blasting project. This method improves the safety of using digital electronic detonators, reduces the misfire analysis time, and improves the blasting construction efficiency.

[0006] The purpose of the present invention is achieved through the following technical solutions:

[0007] An electronic detonator misfire identification method based on the comparison of multiple characteristic information of blasting vibration, including:

[0008] Step 1: Export the original blast source information of the blasting vibration signal of digital electronic detonators from the digital electronic detonator handheld device, perform EMMD signal decomposition based on the original blast source information, establish an AR model to extract multiple characteristic information, and upload the AR model and multiple characteristic information to the blasting vibration measurement data center to form a blasting vibration basic database;

[0009] Step 2: Based on the blasting vibration basic database, collect the original blast source information of multiple blasting sites and perform feature extraction, so as to establish a large database of blasting vibration signals under different delay intervals, typical hole network parameters and charge structures;

[0010] Step 3: Design a misfire monitoring vibration test plan for digital electronic detonators, collect the actual blasting vibration signal, and extract the multiple characteristics of the actual blasting vibration;

[0011] Step 4: Compare the multiple characteristics of the actual blasting vibration with the multiple characteristics in the large database of blasting vibration signals. When there are characteristic differences, it is determined that there is a misfire in this blasting project.

[0012] Specifically, the first step specifically includes the following sub-steps:

[0013] S1. Automatically export the number of blast holes, charge amount per hole, delay time and initiation time of each blast hole detonator from the digital electronic detonator handheld device, and use the above information as the original blast source information of the blasting vibration signal of digital electronic detonators;

[0014] S2. Perform EMMD signal decomposition on the original blast source information to obtain several IMF components and establish an AR model, and use this model to extract the waveform multi-peak characteristic, frequency spectrum characteristic, amplitude characteristic and lifting wavelet packet relative energy spectrum distribution characteristic from the original blast source information;

[0015] S3. Upload the AR model and the waveform multi-peak characteristic, frequency spectrum characteristic, amplitude characteristic and lifting wavelet packet energy spectrum characteristic to the blasting vibration measurement data center as the blasting vibration basic database.

[0016] Specifically, the S2 specifically includes: for the original explosion source information, using EMMD as a signal processing tool to decompose the blasting vibration signal in the original explosion source information to obtain N IMF components with different characteristics; constructing an AR model, and using the waveform, spectrum, amplitude, and lifting wavelet packet relative energy spectrum distribution of the blasting vibration signal as the model parameters of the AR model respectively; using the AR model to extract the initial eigenvector matrix composed of initial eigenvectors from the N IMF components, performing singular value decomposition on the initial eigenvector matrix to obtain the singular values of the initial eigenvector matrix, normalizing each IMF component, and obtaining the singular value entropy of the initial eigenvector matrix as the eigenvalue corresponding to the model parameter, so as to extract the waveform multi-peak feature, spectrum feature, amplitude feature, and lifting wavelet packet relative energy spectrum distribution feature from the original explosion source information.

[0017] Specifically, the feature difference in step four is as follows: compared with the multiple features in the blasting vibration signal large database, the number of waveform multi-peak features in the actual multiple features becomes smaller, and there are obvious differences in the spectrum and the distribution law of the lifting wavelet energy spectrum and the amplitude data.

[0018] The beneficial effects of the present invention:

[0019] The present invention extracts multiple features from the collected original explosion source information to form a basic database by constructing a prediction model for the blasting vibration waveform time history curve, and collects multiple original explosion source information under different blasting conditions for feature extraction, thereby forming a large database of blasting vibration signals. Then, the multiple features of the blasting vibration signals collected on the actual site are extracted and compared with the multiple features in the large database. When there are obvious feature differences, it is determined that a misfire occurs in the current blasting process. This method improves the blasting construction efficiency, reduces the great danger in the misfire treatment process, and improves the safety of using digital electronic detonators at the same time. Brief Description of the Drawings

[0020] Figure 1 is the method flow chart of the present invention;

[0021] Figure 2 is the technical roadmap of the present invention;

[0022] Figure 3 is the schematic diagram of lifting wavelet decomposition and reconstruction;

[0023] Figure 4 is the curve graph of the SGW(10, 10) scaling function and wavelet function;

[0024] Figure 5 is the typical blasting vibration time history curve graph;

[0025] Figure 6 is the power spectrum waveform graph. Detailed implementation mode

[0026] In order to have a clearer understanding of the technical features, objectives, and beneficial effects of the present invention, the technical solutions of the present invention are described in detail below. Obviously, the described embodiments are part of the embodiments of the present invention, rather than all embodiments, and should not be construed as limiting the scope of implementation of the present invention. Based on the embodiments of the present invention, all other embodiments obtained by those of ordinary skill in the art without creative efforts belong to the protection scope of the present invention.

[0027] Embodiment 1:

[0028] In this embodiment, as Figure 1 shown, the electronic detonator misfire recognition method based on the comparison of multi-source feature information of blasting vibration includes:

[0029] Step 1: Export the original blast source information of the digital electronic detonator blasting vibration signal from the digital electronic detonator handheld device, perform EMMD signal decomposition based on the original blast source information, establish an AR model to extract multi-source feature information, and upload the AR model and multi-source feature information to the blasting vibration data center to form a blasting vibration basic database;

[0030] Step 2: Based on the blasting vibration basic database, collect the original blast source information of multiple blasting sites and perform feature extraction, so as to establish a large database of blasting vibration signals under different delay intervals, typical hole network parameters, and charge structures;

[0031] Step 3: Design a misfire monitoring vibration test plan for digital electronic detonators, collect actual blasting vibration signals, and extract multi-source features of the actual blasting vibration;

[0032] Step 4: Compare the multi-source features of the actual blasting vibration with the multi-source features in the large database of blasting vibration signals. When feature differences occur, it is determined that misfires occur in this blasting project.

[0033] In this embodiment, Step 1 specifically includes the following sub-steps:

[0034] S1. Automatically export the number of blast holes, charge amount per hole, detonator delay time and initiation time of each blast hole from the digital electronic detonator handheld device, and use the above information as the original blast source information of the digital electronic detonator blasting vibration signal;

[0035] S2. Perform EMMD signal decomposition on the original blast source information, obtain several IMF components and establish an AR model, and use this model to extract waveform multi-peak characteristics, frequency spectrum characteristics, amplitude characteristics, and lifting wavelet packet relative energy spectrum distribution characteristics from the original blast source information;

[0036] S3. Upload the AR model, along with the waveform multi-peakness feature, spectral feature, amplitude feature, and lifted wavelet packet energy spectrum feature, to the blasting vibration data center as the basic database for blasting vibration.

[0037] In this embodiment, S2 specifically includes: for the original blast source information, using EMMD as the signal processing tool, decompose the blasting vibration signal in the original blast source information to obtain N IMF components with different characteristics; construct an AR model, and use the waveform, spectrum, amplitude, and relative energy spectrum distribution of the lifted wavelet packet of the blasting vibration signal as the model parameters of the AR model respectively; use the AR model to extract the initial feature vectors from the N IMF components to form an initial feature vector matrix, perform singular value decomposition on the initial feature vector matrix to obtain the singular values of the initial feature vector matrix, perform normalization processing on each IMF component, and obtain the singular value entropy of the initial feature vector matrix as the eigenvalue corresponding to the model parameter, so as to extract the waveform multi-peakness feature, spectral feature, amplitude feature, and relative energy spectrum distribution feature of the lifted wavelet packet from the original blast source information.

[0038] In this embodiment, a feature extraction method based on Extremum Field Meanmodedecomposition (EMMD) and AR (autoregressive) singular value entropy is used to extract the features of the blasting vibration signal. The Extremum Field Meanmode Decomposition (EMMD) method is an improvement of the Empirical Mode Decomposition (EMD), which can well process non-stationary signals and reduce the interference of noise. This feature extraction method uses EMMD as the signal processing tool to decompose the signal into several Intrinsic Mode Functions (IMFs), and then uses the characteristic that the autoregressive parameters of the AR model are sensitive to the state change law, selects the AR model parameter vector of the IMF component of the blasting vibration signal as the initial feature vector of the fault, and then obtains the singular value entropy H of the AR model parameter vectors of several IMF components; finally, compare the singular value entropy H to judge whether there is a feature difference. This method can be effectively applied to the feature extraction of nonlinear and non-stationary fault signals.

[0039] In this embodiment, the process of using the above EMMD algorithm to decompose the blasting vibration signal is as follows:

[0040] 1. For the original blasting vibration signal x(t), first initialize: let r0(t) = x(t), i = 1.

[0041] 2. By calculating the local mean at each extreme point of the signal and fitting the local mean curve, obtain the i-th IMF component. The specific process is as follows:

[0042] (1) Initialize: h j-1 (t) = ri-1(t); let j = 1;

[0043] (2) Find h j-1 (t) local extreme points;

[0044] (3) Perform definite integral mean value processing on the local extrema of h j-1 (t) to obtain the local mean m corresponding to each point j-1 (t n ), n = 1, 2, ..., N j-1 ;

[0045] (4) Perform cubic spline fitting on the local mean m j-1 (t n ), n = 1, 2, ..., N j-1 to obtain the mean curve m j-1 (t);

[0046] (5) Let h j (t) = h j-1 (t) - m j-1 (t);

[0047] (6) If the termination screening criterion 0.2 ≤ Δ SD ≤ 0.3, then c i (t) = h j (t), otherwise j = j + 1, return to step (2);

[0048] 3. r i (t) = r i-1 (t) - c i (t).

[0049] 4. If r i (t) has no less than 2 extreme points, then i = i + 1, go back to step 2; otherwise the decomposition ends and r i (t) is the residual component.

[0050] For each IMF component c obtained by the EMMD method i , adopt the FPE criterion to establish the following autoregressive model AR(n):

[0051]

[0052] In the formula, and n are the model parameters and model order of the autoregressive parameter model AR(n) respectively; η(t) is a white noise sequence with a mean of zero and a variance of η 2 . According to the properties of the autoregressive parameter model AR(n), and η 2 can be selected as the feature vector to characterize c iCharacteristic indicators. The initial eigenvector matrix C composed of the initial eigenvectors extracted from N IMF components is expressed as

[0053] C = [C1 C2…C N T Perform singular value decomposition on C to obtain the singular values σ1, σ2,..., σ of the initial eigenvector matrix C N . Normalize each component: p i = σ i / P, where σ1 + σ2 +... + σ N , so, there is p1 + p2 + … p N = 1. Thus, using the following formula

[0054]

[0055] The singular value entropy H of the initial eigenvector matrix C can be calculated and used as the feature of the corresponding parameter for the next step of analysis.

[0056] In this embodiment, the feature difference in step four is: compared with the multi - features in the blasting vibration signal large database, the number of waveform multi - peak features in the actual multi - features becomes smaller, and there are obvious differences in the spectral and lifting wavelet energy spectrum distribution rules and amplitude data.

[0057] Embodiment 2:

[0058] In this embodiment, compared with Embodiment 1, this embodiment uses a curve prediction model to extract multi - feature information, and the method includes the following steps:

[0059] Step 1: Export the original blast source information of the digital electronic detonator blasting vibration signal from the digital electronic detonator handheld device, and establish a curve prediction model based on the original blast source information to extract multi - feature information, and upload the curve prediction model and multi - feature information to the blasting vibration data center to form a blasting vibration basic database;

[0060] Step 2: Based on the blasting vibration basic database, collect the original blast source information of multiple blasting sites and perform feature extraction, so as to establish a large database of blasting vibration signals under different delay intervals, typical hole - net parameters and charge structures;

[0061] Step 3: Design a vibration test plan for monitoring misfires of digital electronic detonators, collect the actual blasting vibration time - history signal, and extract the multi - features of the actual blasting vibration;

[0062] Step 4: Compare the multi - features of the actual blasting vibration with the curve prediction model and multi - features in the blasting vibration signal large database, and determine that there is a misfire in this blasting project when there is a feature difference. ​

[0063] In this embodiment, as shown in Figure 2, the technical implementation process of the method is as follows:

[0064] (1) Automatically export information such as the number of blast holes, charge amount per hole, detonator delay time of each blast hole, initiation time, etc. from the digital electronic detonator handheld device, and the above information serves as the original blast source information of the digital detonator blasting vibration signal;

[0065] (2) Based on the original blast source information, establish a prediction model for the blasting vibration waveform time history curve, and extract waveform multi-peak characteristics, frequency spectrum characteristics, amplitude characteristics, and the relative energy spectrum distribution characteristics of the lifting wavelet packet;

[0066] (3) Upload the above prediction model for the blasting vibration waveform time history curve and the waveform multi-peak characteristics, frequency spectrum characteristics, amplitude characteristics, and the energy spectrum characteristics of the lifting wavelet packet to the blasting vibration measurement data center as the basic database for blasting vibration; when the original blast source parameters are changed, regenerate the corresponding prediction model for the blasting vibration time history curve and the waveform multi-peak characteristics, frequency spectrum characteristics, amplitude characteristics, and the energy spectrum characteristics of the lifting wavelet packet.

[0067] (4) Design a vibration test plan for monitoring misfires of digital electronic detonators;

[0068] (5) After issuing the initiation command, collect the actual blasting vibration time history signal;

[0069] (6) Extract the multi-peak characteristics, frequency spectrum characteristics, amplitude characteristics, and the relative energy spectrum distribution characteristics of the lifting wavelet packet of the actual blasting vibration time history signal;

[0070] (7) Establish a large database of blasting vibration signals under conditions of different delay intervals, typical hole pattern parameters, charge structures, and other blast source parameters;

[0071] (8) Compare the extracted characteristic values with the characteristic values in the prediction model for the blasting vibration waveform time history curve and the large database of blasting vibration signals;

[0072] (9) Compare the multi-source characteristic information. When there is a decrease in the number of multi-peak characteristics, obvious differences in the frequency spectrum and the distribution law of the lifting wavelet energy spectrum, and amplitude, it is determined that a misfire has occurred in this blasting project.

[0073] Among them, the prediction model for the blasting vibration waveform time history curve is mainly implemented based on the construction method of the lifting wavelet basis function. The construction method of the lifting wavelet basis function is realized by using the lifting algorithm. The core idea of the lifting algorithm is to decompose the existing wavelet filter into basic construction modules. Decompose the wavelet transform process into: decomposition (Split), prediction (Predict), and update (Update). Its decomposition and reconstruction principle is as Figure 3 shown.

[0074] Assume a sampling sequence S = {s k , k ∈ Z, s k ∈ R}, and assume that the length of the signal is 2 k . After one wavelet decomposition, the obtained approximation sequence and detail sequence are represented by s k-1 and d k-1 respectively, and their lengths are both 2 k-1 ; similarly, if the approximation sequence and detail sequence of the signal are s k-1 and d k-1 respectively, the signal sequence obtained after the inverse wavelet transform is represented by s k , and its length is 2 k . The signal decomposition process is as follows:

[0075] (1) Decomposition. The sampling sequence {s(k), k ∈ Z} is partitioned into an odd sample sequence S o = {s o (k), k ∈ Z} and an even sample sequence S e = {s e (k), k ∈ Z}, and this decomposition is called lazy wavelet transform. Among them:

[0076] s o (k) = s(2k + 1) k ∈ Z (1)

[0077] s e (k) = s(2k) k ∈ Z (2)

[0078] (2) Prediction, also called dual lifting. Based on the correlation of the original data, the predicted value P[s e (k)] of the even sequence S e is used to predict (or interpolate) the odd sequence s o , that is, the result of applying the filter P(·) to the even signal is used as the predicted value of the odd signal, and the difference between the actual value and the predicted value of the odd signal is used to obtain the residual signal. In practice, although s o (k) cannot be accurately predicted from s e (k), P[s e (k)] is very close to s e (k). Therefore, the difference between P[s e (k)] and s e (k) can be used to replace the original s e (k). The s e (k) generated in this way contains less information than the original s e (k). Define the prediction deviation as the detail signal:

[0079] d(k) = s o (k) - P[s e(k)] k ∈ Z (3)

[0080] (3) Update, also called original lifting. Update is also referred to as lifting. In order to make some global characteristics of the original signal set continue to be maintained in its subsets and eliminate the frequency aliasing effect generated during the dissection process, an update must be performed. The basic idea of the update is to find a better subset that keeps the mean and high-order vanishing moments of the original signal unchanged. Let U(.) be the updater, and based on the detail signal sequence D = {d(k), k ∈ Z}, update s e (k), and the result is defined as the approximation signal c(k)

[0081] c(k) = s e (k) + U(D) k ∈ Z (4)

[0082] It can be seen from the above basic principle that the lifting algorithm repeats through the iteration of the sampled data to achieve the multi-resolution decomposition of the multi-level transformation.

[0083] From the perspective of the frequency domain analysis, predicting P[s e (k)] means smoothing, which can be regarded as low-frequency; the detail signal d(k) means the error between the signal in the local area and its own low-frequency component, reflecting the high-frequency component in the signal s(k), also known as the wavelet coefficient. The smaller the wavelet coefficient, the more accurate the prediction smoothing and the better the prediction effect; the approximation signal c(k) reflects the low-frequency component in the signal. Continuing to decompose c(k), the approximate signal and detail signal of the next layer can be obtained, and so on.

[0084] The reconstruction process of the lifting wavelet transform is the inverse process of the decomposition process. As Figure 3 can be seen, the signal reconstruction consists of three steps: restoring the update, restoring the prediction, and merging the odd and even sampling sequences. By performing the same prediction operator and update operator operations on the wavelet decomposition result of the lifting format and changing the operation signs, we can get:[[]]

[0085] (1) Anti-update. s e (k) = c(k) - U(D) (5)

[0086] (2) Anti-prediction. s o (k) = d(k) + P[s e (k)] (6)

[0087] (3) Merge. s(k) = merge(s o (k), s e (k)) (7)

[0088] When processing signals with the lifting algorithm and the Mallat algorithm, only the filter coefficients are required to achieve signal decomposition and reconstruction, without the need to know the specific expressions of the scaling function and wavelet function. When analyzing the influence of the characteristics of the wavelet function on the signal analysis results, the specific forms of the wavelet function and the scaling function are studied again.

[0089] In this embodiment, based on the above principle process, the construction process of the model is as follows:

[0090] 1. Construct wavelets in lifting format

[0091] Each sampling sequence s j,k (k) corresponds to a scaling function Correspondingly, given an arbitrary δ sequence such that s j,k = δ j,k , when interpolating and subdividing it until infinity, a scaling function can be generated. When sampling at equal intervals, and the 0,0 generated by interpolating s 0,0 satisfy the two-scale equation of the first-generation wavelet

[0092]

[0093] Similarly, each d j,k corresponds to a wavelet function ψ j,k (x). Given an arbitrary sequence of d j,k = δ j,k , after performing an inverse wavelet lifting transformation on it and then interpolating and subdividing until infinity, a wavelet function ψ j,k (x) can be generated; when sampling at equal intervals, ψ j,k (x) and the wavelet mother function ψ 0,0 generated by d 0,0 = δ 0,0 ψ j,k (x) = ψ(x) satisfy the two-scale equation of the first-generation wavelet

[0094] ψ j (x) = ψ(2 N -k) (9)

[0095] 2. Design of P and U coefficients

[0096] It can be seen from the basic principle of constructing the second-generation wavelet from the boosting algorithm that the selection of the prediction coefficient P and the update coefficient U is the key to the construction of the second-generation wavelet, which is related to the effect of signal analysis. The key to the selection of the prediction coefficient P is to make the prediction accurate. The more accurate the prediction, the better the smoothing effect, and the more effectively the high-frequency signal can be decomposed. The key to the selection of the update operator U is to make the mean, high-order moments, etc. of the updated low-frequency component remain unchanged from the original signal. When the requirement for operation real-time performance is not high, the value of P can be directly obtained through the Neville interpolation algorithm, and then U can be determined by constructing an auxiliary sequence. When the requirement for data processing real-time performance is relatively high, the coefficients of P and U can be obtained in advance, and then the lifting wavelet decomposition and reconstruction can be completed according to Equations (1) to (7).

[0097] 3. Construction of Interpolation Subdivision Wavelet

[0098] Given N + 1 distinct points x0, x1, …, x N at which the function values are y0, y1, …, y N , and y i = f(x i ), i = 0, …, N, then there exists a unique polynomial L n (x) of degree not greater than N such that L n (x i ) = f(x i ), then

[0099]

[0100] where

[0101]

[0102] Each time of subdivision, take N (N = 2D, D ∈ Z + ) known samples y j,k-D+1 , … y j,k , … y j,k+D , and assume that these samples are sampled at equal time intervals, and their corresponding sampling times are x k + 1, x k + 2, …, x k + N, x k is an arbitrary starting time. The new sampling values generated by subdivision are at the intermediate positions of these known samples. The interpolation point (or prediction point) is: x = x k + (N + 1) / 2. In this way, the prediction coefficient can be determined by Equation (12), that is

[0103]

[0104] The operation of obtaining the coefficients in the above formula is actually independent of x k . For the sake of simplicity, xk = 0。

[0105] When N = 10, the prediction coefficients are obtained according to Equation (12).

[0106] Let the signal s be a δ sequence, that is

[0107] s = [0 0 0 1 0 0 0] (13)

[0108] See Figure 3 In the wavelet reconstruction part, the reconstruction relationship is simplified to

[0109] s j+1,2k = s j,k , s j+1,2k+1 = p(x j+1,2k+1 ) (14)

[0110] Interpolate and subdivide s according to the above formula, and use the values of s at four points k - 1, k, k + 1, k + 2 to predict the value of s j,k value, and the interpolation boundary is processed by padding zeros. Use the above interpolation and subdivision method to compile a calculation program through Matlab. The scaling function waveform generated after 30 iterations is shown in j+1,2k+1 See Figure 3 .

[0111] To construct a wavelet function, it is necessary to first design an update operator U. The role of the update operator is to ensure that the mean value and low - order vanishing moments of the signal remain unchanged before and after wavelet transform. When (N is the number of predictors, is the number of updaters), the predictor coefficients can be directly divided by 2 as the update coefficients.

[0112] Assume the structure of the update operator is

[0113]

[0114] where U and D are the operation coefficients of the update operator and the detail signal sequence respectively. Similar to the construction of the scaling function, let the detail signal be a δ sequence and the approximate signal s be a zero sequence. Refer to Figure 2 the wavelet reconstruction part, perform an inverse lifting - scheme transform once, and the even sequence of the signal is

[0115] s e = [0 - u6 - u5 - u4 - u3 - u2 - u1] (16)

[0116] This sequence undergoes a recovery prediction operation through the prediction operator P(·) to obtain the odd sequence of the signal as

[0117] s o = D + P(s​e ) (17)

[0118] From equations (16) to (17), s can be obtained o , for s o Interpolating and subdividing according to equation (17), the corresponding wavelet function can be obtained after 30 iterations, as Figure 4 shown.

[0119] The scaling function and wavelet function based on the lifting scheme are symmetric, compactly supported, and have an impulse shape. When N and are small, the support intervals of the scaling function and wavelet function are small; conversely, the support intervals are large and have better continuity. Generally, wavelet functions with small support intervals are suitable for processing non-stationary signals, and wavelet coefficients can effectively describe the transient components of signals, while wavelets with large and well-connected support intervals are suitable for describing stationary signals. Blasting vibration signals have typical short-time non-stationary characteristics, and it is proposed to use SGW(10, 10) for blasting vibration signal analysis.

[0120] 4 Construction of the lifted db wavelet

[0121] Daubechies wavelets have good compact support, smoothness, and approximate symmetry, and have been successfully applied to the analysis of non-stationary signals including blasting vibration signals. The wavelet coefficients have different sequences (dbN) according to positive integers N. In this invention, the lifting scheme (Daubechies lifting schemes) liftingdb1 provided in Matlab is used, and the second-generation wavelet SGW(10, 10) based on interpolating subdivision is used for denoising of blasting vibration signals.

[0122] 2. Method for extracting the energy spectrum of the lifted wavelet packet

[0123] Let the blasting vibration signal {s[i]} (i = 1, 2,... 2 N , N ∈ Z + ) After being decomposed by the lifted wavelet packet, the coefficient of the (j, n) node is denoted as

[0124] Define the normalized energy:[[]]

[0125]

[0126] In the above formula, E is the total energy of the signal, that is, E is used to normalize the energy in each frequency band, and the corresponding eigenvector:[[]]

[0127] V = (E1, E2,..., E j ) (19)

[0128] V is called the normalized eigenvector.

[0129] Figure 5 It is the blasting vibration time history curve collected by using the single-hole and single-section initiation method in the field experiment. From the power spectrum Figure 6 it can be intuitively obtained that the main vibration frequencies of this signal are in three frequency bands of (38 - 80) Hz, (80 - 122) Hz, and (150 - 180) Hz. The power spectral density contained in the frequency part above 200 Hz is small, and the resolution in the spectrogram is low. Therefore, the power spectral density distribution in the range of 0 - 256 Hz was analyzed during misfire identification.

[0130] SGW(10, 10) is used for analyzing the blasting vibration signal, and the formed blasting time history curve is closer to the waveform shape of the actual blasting project. The analysis effect is better than that of SGW(6, 6) and SGW(8, 8)

[0131] The above shows and describes the basic principles, main features and advantages of the present invention. Those skilled in the art should understand that the present invention is not limited by the above embodiments. What is described in the above embodiments and the specification only illustrates the principles of the present invention. Without departing from the spirit and scope of the present invention, the present invention will have various changes and improvements, and these changes and improvements fall within the scope of the present invention claimed. The scope of the present invention claimed is defined by the appended claims and their equivalents.

Claims

1. An electronic detonator misfire recognition method based on comparison of multiple characteristic information of blasting vibration, characterized in that Including: Step 1: Export the original source information of the blasting vibration signal of digital electronic detonators from the digital electronic detonator handheld device, perform EMMD signal decomposition based on the original source information, establish an AR model to extract multi-source feature information, and upload the AR model and multi-source feature information to the blasting vibration measurement data center to form a basic database of blasting vibration; the AR model is specifically an autoregressive model; Step 2: Based on the basic database of blasting vibration, collect the original source information of multiple blasting sites and perform feature extraction, so as to establish a large database of blasting vibration signals under different delay intervals, typical hole network parameters and charge structures; Step 3: Design a vibration test plan for monitoring misfires of digital electronic detonators, collect actual blasting vibration signals, and extract multi-source features of the actual blasting vibration; Step 4: Compare the multi-source features of the actual blasting vibration with the multi-source features in the large database of blasting vibration signals. When there are feature differences, it is determined that there are misfires in this blasting project: when the number of multi-peak features in the actual multi-source features decreases, and there are obvious differences in the distribution laws of the spectrum and the lifting wavelet energy spectrum and the amplitude, it is determined that there are misfires in this blasting project.

2. The electronic detonator misfire identification method based on blasting vibration multi - feature information comparison according to claim 1, wherein, The specific steps of Step 1 include the following sub-steps: S1. Automatically export the number of blast holes, the charge amount per hole, the delay time and the initiation time of each blast hole detonator from the digital electronic detonator handheld device, and use the above information as the original source information of the blasting vibration signal of the digital electronic detonator; S2. Perform EMMD signal decomposition on the original source information to obtain several IMF components and establish an AR model, and use this model to extract the waveform multi-peak feature, spectrum feature, amplitude feature and relative energy spectrum distribution feature of the lifting wavelet packet from the original source information; S3. Upload the AR model and the waveform multi-peak feature, spectrum feature, amplitude feature and energy spectrum feature of the lifting wavelet packet to the blasting vibration measurement data center as the basic database of blasting vibration.

3. The electronic detonator misfire identification method based on blasting vibration multi - feature information comparison according to claim 2, characterized in that The specific content of S2 includes: For the original source information, using EMMD as the signal processing tool, decompose the blasting vibration signal in the original source information to obtain N IMF components with different features; construct an AR model, and use the waveform, spectrum, amplitude and relative energy spectrum distribution of the lifting wavelet packet of the blasting vibration signal as the model parameters of the AR model respectively; use the AR model to extract the initial feature vectors from the N IMF components to form an initial feature vector matrix, perform singular value decomposition on the initial feature vector matrix to obtain the singular values of the initial feature vector matrix, perform normalization processing on each IMF component, and obtain the singular value entropy of the initial feature vector matrix as the eigenvalue of the corresponding model parameter, so as to extract the waveform multi-peak feature, spectrum feature, amplitude feature and relative energy spectrum distribution feature of the lifting wavelet packet from the original source information.

4. The electronic detonator misfire recognition method based on blasting vibration multi-characteristic information comparison according to claim 1, characterized in that The feature difference in Step 4 is as follows: Compared with the multi-source features in the large database of blasting vibration signals, in the actual multi-source features, the number of waveform multi-peak features decreases, and there are obvious differences in the distribution laws of the spectrum and the lifting wavelet energy spectrum and the amplitude data.

Images