WO2019019565A1 - 一种基于能量分布特征的矿山微震信号辨识方法 - Google Patents

一种基于能量分布特征的矿山微震信号辨识方法 Download PDF

Info

Publication number
WO2019019565A1
WO2019019565A1 PCT/CN2018/072533 CN2018072533W WO2019019565A1 WO 2019019565 A1 WO2019019565 A1 WO 2019019565A1 CN 2018072533 W CN2018072533 W CN 2018072533W WO 2019019565 A1 WO2019019565 A1 WO 2019019565A1
Authority
WO
WIPO (PCT)
Prior art keywords
signal
microseismic signal
microseismic
energy distribution
identified
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Ceased
Application number
PCT/CN2018/072533
Other languages
English (en)
French (fr)
Inventor
张杏莉
卢新明
贾瑞生
彭延军
赵卫东
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Shandong University of Science and Technology
Original Assignee
Shandong University of Science and Technology
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Shandong University of Science and Technology filed Critical Shandong University of Science and Technology
Publication of WO2019019565A1 publication Critical patent/WO2019019565A1/zh
Anticipated expiration legal-status Critical
Ceased legal-status Critical Current

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
    • G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
    • G01V1/30—Analysis
    • G01V1/307—Analysis for determining seismic attributes, e.g. amplitude, instantaneous phase or frequency, reflection strength or polarity
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V2210/00—Details of seismic processing or analysis
    • G01V2210/60—Analysis
    • G01V2210/61—Analysis by combining or comparing a seismic data set with other data
    • G01V2210/616—Data from specific type of measurement
    • G01V2210/6161—Seismic or acoustic, e.g. land or sea measurements

Definitions

  • the invention belongs to the field of signal analysis and identification, and particularly relates to a mine microseismic signal identification method based on energy distribution characteristics.
  • Microseismic monitoring is an advanced and effective monitoring method for coal and rock dynamic disasters developed in recent years. It can monitor the microseismic activities of coal and rock mass in real time, continuously and online, and form microseismic monitoring data. Due to the complex environment of the mine, there are a lot of interference signals such as background noise and blasting vibration, which makes the microseismic monitoring system unable to accurately identify and record the effective microseismic events. Later, it is necessary to manually identify the effective microseismic events by technicians, which seriously affects the identification of the microseismic monitoring system. effectiveness. Because coal mine blasting operations often occur, and the micro-seismic and blasting vibration waveforms of coal and rock mass are very similar, the manual identification method often causes mishandling and is difficult to identify.
  • common time-frequency analysis methods for waveform identification of mine microseismic signals include Fourier transform, wavelet transform, wavelet packet transform, frequency slice wavelet transform and EMD.
  • Traditional Fourier transform is mainly used to analyze periodic stationary signals, including spikes and The random and non-stationary microseismic signal analysis of the mutation is not effective; the wavelet analysis can simultaneously perform time-frequency analysis, but the appropriate wavelet base needs to be selected to achieve better decomposition effect; EMD can process random non-stationary signals well.
  • the EMD method has boundary effects and modal aliasing, which leads to instability and non-uniqueness of EMD. These methods have a certain degree of disadvantages in signal analysis, which increases the difficulty of signal identification and high false positive rate.
  • the present invention proposes a mine microseismic signal identification method based on energy distribution characteristics, and uses time division mode decomposition (VMD) to perform time-frequency analysis on signals.
  • VMD is a new signal decomposition method. Compared with other modal decomposition techniques, it has a solid theoretical foundation, eliminates the modal aliasing problem, overcomes the shortcomings of the prior art, and has a good frequency domain adaptive decomposition. effect.
  • Step 2 Perform VMD decomposition on the identified microseismic signal x(t) to obtain K variable-variant modal components ⁇ u 1 ,...u k ,..., u K ⁇ arranged in descending order of frequency:
  • the differential seismic signal x(t) is decomposed into K variable modal components by VMD.
  • the constraint is to minimize the sum of the estimated bandwidths of the modalities, and the sum of the modalities is equal to the microseismic signal x(t) to be identified.
  • the constrained variational model is described as equations (1) and (2):
  • x(t) represents the microseismic signal to be identified
  • ⁇ (t) is a Dirac function
  • * represents a convolution
  • j 2 -1; in equation (2), To sum all the variational modes;
  • ⁇ is a quadratic penalty factor and ⁇ (t) is a Lagrangian multiplication operator
  • Step 2.1 defining the value of the number K of the variational modal component and the value of the penalty factor ⁇ ;
  • Step 2.2 Initialize
  • Step 2.4 Execute the first loop of the inner layer, and update u k according to formula (4);
  • Step 2.6 Execute the second loop of the inner layer, and update ⁇ k according to equation (5);
  • Step 2.8 Perform an outer loop to update ⁇ according to equation (6);
  • ⁇ is an update step size parameter of the Lagrangian multiplication operator ⁇ (t);
  • Step 2.9 Repeat steps 2.3 to 2.8 until the iterative stop condition is satisfied, as shown in equation (7), ending the entire loop, and outputting the result to obtain K variable modal components;
  • Step 3 Calculate the energy distribution vector P of the microseismic signal x(t) to be identified
  • the energy percentage value of the modal component u k can be obtained.
  • Step 4 Calculate the center of gravity coefficient cx of the energy distribution of the X-axis of the energy distribution of the microseismic signal x(t) to be identified;
  • Step 5 Identify the microseismic signal x(t) to be identified according to the identification threshold T. If cx>T is the microseismic signal of the mine rock mass rupture, cx ⁇ T is the blasting vibration signal;
  • Step 6 adaptively update the value of the identification threshold T
  • the identification threshold T is updated according to the system of equations (10):
  • W 1 is a set of cx values of the microseismic signal of the coal rock mass in the training concentrated
  • W 2 is a set of cx values of the vibration signal of the training concentrated blasting.
  • the present invention utilizes the characteristic difference of the energy distribution of the two microseismic signals, firstly reads the microseismic signal to be identified and performs the VMD decomposition, and obtains K according to the frequency from the high.
  • To the low-ordered variational modal components calculate the band energy of each modal component, extract the energy percentage value of each modal component from the original signal to form the energy distribution vector P; calculate the energy distribution based on the energy distribution vector P
  • the X-axis center-of-gravity coefficient cx identifies the mine microseismic signal according to the identification threshold T.
  • the detected microseismic signal is the mine rock mass rupture microseismic signal. If cx ⁇ T, the microseismic signal is detected as the blasting vibration signal.
  • the method can effectively identify the microseismic signals and blasting vibration signals of coal rock mass rupture.
  • the present invention adopts the above technical solutions, and has the following advantages compared with the prior art:
  • the invention automatically divides the micro-seismic signal of the mine, and according to the significant difference of the energy distribution of the micro-seismic signal and the blasting vibration signal of the coal-rock mass in different frequency bands, the X-axis center-of-gravity coefficient of the energy distribution of the microseismic signal is calculated.
  • the effective identification of the microseismic signals of the two types of mines is realized.
  • the method has the characteristics of simple algorithm, adaptability and real-time performance, and has good technical value and application prospect.
  • 1 is a flow chart of a mine microseismic signal identification method based on energy distribution characteristics.
  • FIG. 2 is a schematic diagram of a microseismic signal x(t) to be identified and a time-frequency diagram thereof.
  • FIG. 3 is a schematic diagram of six variational modal components obtained by decomposing the microseismic signal x(t) by VMD and its time-frequency diagram.
  • Figure 4 is a histogram of the energy distribution of the microseismic signal x(t) to be identified.
  • Figure 5 is a graph showing the energy vector, center of gravity coefficient and identification results of the 15 sets of coal rock mass fracture microseismic test signals.
  • Figure 6 is a graph showing the energy vector, center of gravity coefficient and identification results of 15 groups of blasting vibration test signals.
  • Figure 7 shows the classification and recognition results of the microseismic signals in the test group.
  • a mine microseismic signal identification method based on energy distribution characteristics the flow of which is shown in Figure 1, which specifically includes the following steps:
  • Step 2 Perform VMD decomposition on the identified microseismic signal x(t) to obtain K variable-variant modal components ⁇ u 1 ,...u k ,..., u K ⁇ arranged in descending order of frequency:
  • the differential seismic signal x(t) is decomposed into K variable modal components by VMD.
  • the constraint is to minimize the sum of the estimated bandwidths of the modalities, and the sum of the modalities is equal to the microseismic signal x(t) to be identified.
  • the constrained variational model is described as equations (1) and (2):
  • x(t) represents the microseismic signal to be identified
  • ⁇ (t) is a Dirac function
  • * represents a convolution
  • j 2 -1; in equation (2), To sum all the variational modes;
  • ⁇ is a quadratic penalty factor and ⁇ (t) is a Lagrangian multiplication operator
  • Step 2.1 defining the value of the number K of the variational modal component and the value of the penalty factor ⁇ ;
  • Step 2.2 Initialize
  • Step 2.4 Execute the first loop of the inner layer, and update u k according to formula (4);
  • Step 2.6 Execute the second loop of the inner layer, and update ⁇ k according to equation (5);
  • Step 2.8 Perform an outer loop to update ⁇ according to equation (6);
  • ⁇ is an update step size parameter of the Lagrangian multiplication operator ⁇ (t);
  • Step 2.9 Repeat steps 2.3 to 2.8 until the iterative stop condition is satisfied, as shown in equation (7), ending the entire loop, and outputting the result to obtain K variable modal components;
  • Step 3 Calculate the energy distribution vector P of the microseismic signal x(t) to be identified
  • the energy percentage value of the modal component u k can be obtained.
  • Step 4 Calculate the center of gravity coefficient cx of the energy distribution of the X-axis of the energy distribution of the microseismic signal x(t) to be identified;
  • Step 5 Identify the microseismic signal x(t) to be identified according to the identification threshold T. If cx>T is the microseismic signal of the mine rock mass rupture, cx ⁇ T is the blasting vibration signal;
  • Step 6 adaptively update the value of the identification threshold T
  • the identification threshold T is updated according to the system of equations (10):
  • W 1 is a set of cx values of the microseismic signal of the coal rock mass in the training concentrated
  • W 2 is a set of cx values of the vibration signal of the training concentrated blasting.
  • FIG. 4 is a histogram of the energy distribution of the microseismic signal, and the black solid in the figure The circle is the position of the center of gravity of the energy distribution plane of the microseismic signal.
  • test group 15 sets of coal rock mass rupture microseismic signals and 15 sets of blasting vibration microseismic signals are given respectively.
  • the energy vector, center of gravity coefficient and identification result of the coal rock rupture microseismic test signal are shown in Fig. 5; blasting vibration microseismic test
  • the energy vector, center of gravity coefficient and its identification result of the signal are shown in Fig. 6.
  • the test group has a total of 30 sets of microseismic signals, of which 29 groups are correctly identified, 1 group is identified incorrectly, and the correct rate of identification is 96.67%.
  • the classification and recognition results of the microseismic signals of the test group are shown in Fig. 7.
  • the microseismic signal is a non-stationary random signal, and its frequency distribution is relatively scattered.
  • the energy distribution of different types of microseismic signals is significantly different in different frequency bands. Therefore, according to this feature, the energy distribution vector of the microseismic signal can be extracted, and the center of gravity coefficient of the energy distribution can be calculated and identified. By comparing the threshold values, the classification and identification of the microseismic signals to be measured can be realized.

Landscapes

  • Engineering & Computer Science (AREA)
  • Remote Sensing (AREA)
  • Physics & Mathematics (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Acoustics & Sound (AREA)
  • Environmental & Geological Engineering (AREA)
  • Geology (AREA)
  • General Life Sciences & Earth Sciences (AREA)
  • General Physics & Mathematics (AREA)
  • Geophysics (AREA)
  • Geophysics And Detection Of Objects (AREA)
  • Investigating Or Analyzing Materials By The Use Of Ultrasonic Waves (AREA)
  • Measurement Of Mechanical Vibrations Or Ultrasonic Waves (AREA)
  • Other Investigation Or Analysis Of Materials By Electrical Means (AREA)

Abstract

一种基于能量分布特征的矿山微震信号辨识方法，属于信号分析及识别领域，包括如下步骤：读取待辨识微震信号x(t)；对x(t)进行VMD分解，得到K个按照频率从高到低顺序排列的变分模态分量；计算出各模态分量的频带能量，提取各模态分量占原信号的能量百分比值构成能量分布向量P；以能量分布向量P为基础计算出能量分布X轴重心系数cx；根据辨识阈值T识别矿山微震信号，若cx>T为矿山煤岩体破裂微震信号，若cx≤T为爆破震动信号；最后对辨识阈值T的值进行自适应更新。该方法能有效区分煤岩体破裂微震信号和爆破震动信号，具有自适应性强、准确性高等特点。

Description

一种基于能量分布特征的矿山微震信号辨识方法 技术领域
本发明属于信号分析及识别领域，具体涉及一种基于能量分布特征的矿山微震信号辨识方法。
背景技术
微震监测是近年来发展起来的先进且行之有效的煤岩动力灾害监测手段，它能够对煤岩体微震活动实时、连续、在线监测，形成微震监测数据。由于矿山环境复杂，存在现场背景噪声、爆破震动等大量干扰信号，使得微震监测系统无法准确识别并记录有效微震事件，后期需要依靠技术人员人工识别出有效微震事件，严重影响了微震监测系统的识别效率。由于煤矿爆破作业经常发生，而煤岩体微震和爆破震动波形又极为相似，采用人工识别方式，经常出现误处理，识别难度大。
目前，针对矿山微震信号波形识别的常用时频分析法包括傅立叶变换、小波变换、小波包变换、频率切片小波变换和EMD等，传统傅立叶变换主要用于分析周期性平稳信号，对包含有尖峰和突变的随机性、非平稳性微震信号分析效果欠佳；小波分析能同时进行时频分析，但需要选择合适的小波基才能达到较好的分解效果；EMD能较好地处理随机非平稳信号，但EMD方法存在边界效应及模态混叠现象，导致EMD具有不稳定性和不唯一性。这些方法用于信号分析时均存在一定程度的弊端，为信号辨识增加了难度，误判率高。
发明内容
针对现有技术中存在的上述问题，本发明提出了一种基于能量分布特征的矿山微震信号辨识方法，采用变分模态分解(VMD)对信号进行时频分析。VMD是一种新的信号分解方法，相比于其它模态分解技术，它具有坚实的理论基础，消除了模态混叠问题，克服了现有技术的不足，具有良好的频域自适应分解效果。
为了实现上述目的，本发明采用如下技术方案：
一种基于能量分布特征的矿山微震信号辨识方法，包括如下步骤：
步骤1：读取待辨识微震信号x(t)，其中，t＝1，2，…，N，N为微震信号的采样点个数；
步骤2：对待辨识微震信号x(t)进行VMD分解，得到K个按照频率从高到低顺序排列的变分模态分量{u 1,…u k,…,u K}：
对待辨识微震信号x(t)采用VMD分解为K个变分模态分量，约束条件为使各个模态的估计带宽之和最小，且各模态之和等于待辨识微震信号x(t)，约束变分模型描述为式(1)和式(2)：
Figure PCTCN2018072533-appb-000001
s.t.∑ ku k＝x(t)  (2)；
其中，x(t)表示待辨识的微震信号，{u k}:＝{u 1,…,u K}代表分解得到的K个有限带宽的变分模态分量，{ω k}:＝{ω 1,…,ω K}表示各分量的频率中心，δ(t)为狄拉克(Dirac)函数，*表示卷积，j 2＝-1；式(2)中，
Figure PCTCN2018072533-appb-000002
为对所有的变分模态求和；
为求解式(1)和式(2)的最优解，引入扩展的Lagrange将约束变分问题变为非约束变分问题，其表达式为式(3)：
Figure PCTCN2018072533-appb-000003
其中，α为二次惩罚因子，λ(t)为拉格朗日乘法算子；
求解该变分问题的具体步骤如下：
步骤2.1：定义变分模态分量个数K值与惩罚因子α的值；
步骤2.2：初始化
Figure PCTCN2018072533-appb-000004
步骤2.3：令n＝n+1，执行整个循环；
步骤2.4：执行内层第一个循环，根据式(4)更新u k；
Figure PCTCN2018072533-appb-000005
其中，
Figure PCTCN2018072533-appb-000006
为待辨识微震信号x(t)的傅立叶变换，
Figure PCTCN2018072533-appb-000007
步骤2.5：令k＝k+1，重复步骤2.4，直到k＝K，结束内层第一个循环；
步骤2.6：执行内层第二个循环，根据式(5)更新ω k；
Figure PCTCN2018072533-appb-000008
步骤2.7：令k＝k+1，重复步骤2.6，直到k＝K，结束内层第二个循环；
步骤2.8：执行外层循环，根据式(6)更新λ；
Figure PCTCN2018072533-appb-000009
其中，τ为拉格朗日乘法算子λ(t)的更新步长参数；
步骤2.9：重复步骤2.3至步骤2.8，直到满足迭代停止条件如式(7)所示，结束整个循环，输出结果，得到K个变分模态分量；
Figure PCTCN2018072533-appb-000010
其中，ε为求解精度；
步骤3：计算待辨识微震信号x(t)的能量分布向量P；
根据公式(8)计算各模态分量u k对应的能量E k；
Figure PCTCN2018072533-appb-000011
其中，x ik(t)(i＝1,2,…N；k＝1,2,…,K；N为采样点个数，K为变分模态个数)表示模态分量u k时序序列的离散点幅值；
根据每个模态分量u k的能量以及待辨识微震信号x(t)的总能量，可以得到模态分量u k的能量百分比值
Figure PCTCN2018072533-appb-000012
从而得到该微震信号的能量分布向量P，即P＝[P(1),…,P(k),…,P(K)]；
步骤4：计算待辨识微震信号x(t)的能量分布X轴的重心系数cx；
根据公式(9)计算能量分布X轴重心系数cx：
Figure PCTCN2018072533-appb-000013
步骤5：根据辨识阈值T识别待辨识微震信号x(t)，若cx>T为矿山煤岩体破裂微震信号，cx≤T为爆破震动信号；
步骤6：自适应更新辨识阈值T的值；
根据方程组(10)更新辨识阈值T：
Figure PCTCN2018072533-appb-000014
其中，W 1为训练集中煤岩体破裂微震信号的cx值集合，W 2为训练集中爆破震动信号的cx值集合。
本发明原理如下：
为实现煤岩体破裂微震信号和爆破震动信号的有效分类辨识，本发明利用两种微震信号能量分布差异显著的特点，首先读取待辨识微震信号并进行VMD分解，得到K个按照频率从高到低顺序排列的变分模态分量；计算出各模态分量的频带能量，提取各模态分量占原信号的能量百分比值构成能量分布向量P；以能量分布向量P为基础计算出能量分布X轴重心系数cx；根据辨识阈值T识别矿山微震信号，若cx>T时，检测微震信号为矿山煤岩体破裂微震信号，若cx≤T时，检测微震信号为爆破震动信号。该方法可以实现对煤岩体破裂微震信号和爆破震动信号的有效辨识。
本发明采用以上技术方案，与现有技术现比，具有以下优点：
本发明依据VMD良好频谱分解特征对矿山微震信号进行自适合剖分，依据煤岩体破裂 微震信号和爆破震动信号在不同频段上能量分布的显著差异，通过计算微震信号能量分布X轴重心系数，实现对两类矿山微震信号的有效辨识，该方法具有算法简单、自适应性和实时性强的特点，具有很好的技术价值和应用前景。
附图说明
图1为本发明一种基于能量分布特征的矿山微震信号辨识方法的流程图。
图2为待辨识微震信号x(t)的示意图及其时频图。
图3为待辨识微震信号x(t)经VMD分解后得到的6个变分模态分量示意图及其时频图。
图4为待辨识微震信号x(t)的能量分布直方图。
图5为15组煤岩体破裂微震测试信号的能量向量、重心系数及其辨识结果图。
图6为15组爆破震动测试信号的能量向量、重心系数及其辨识结果图。
图7为测试组微震信号分类识别结果。
具体实施方式
下面结合附图以及具体实施方式对本发明作进一步详细说明：
一种基于能量分布特征的矿山微震信号辨识方法，其流程如图1所示，具体包括如下步骤：
步骤1：读取待辨识微震信号x(t)，其中，t＝1，2，…，N，N为微震信号的采样点个数；
步骤2：对待辨识微震信号x(t)进行VMD分解，得到K个按照频率从高到低顺序排列的变分模态分量{u 1,…u k,…,u K}：
对待辨识微震信号x(t)采用VMD分解为K个变分模态分量，约束条件为使各个模态的估计带宽之和最小，且各模态之和等于待辨识微震信号x(t)，约束变分模型描述为式(1)和式(2)：
Figure PCTCN2018072533-appb-000015
s.t.Σ ku k＝x(t)    (2)；
其中，x(t)表示待辨识的微震信号，{u k}:＝{u 1,…,u K}代表分解得到的K个有限带宽的变分模态分量，{ω k}:＝{ω 1,…,ω K}表示各分量的频率中心，δ(t)为狄拉克(Dirac)函数，*表示卷积，j 2＝-1；式(2)中，
Figure PCTCN2018072533-appb-000016
为对所有的变分模态求和；
为求解式(1)和式(2)的最优解，引入扩展的Lagrange将约束变分问题变为非约束变分问题，其表达式为式(3)：
Figure PCTCN2018072533-appb-000017
其中，α为二次惩罚因子，λ(t)为拉格朗日乘法算子；
求解该变分问题的具体步骤如下：
步骤2.1：定义变分模态分量个数K值与惩罚因子α的值；
步骤2.2：初始化
Figure PCTCN2018072533-appb-000018
步骤2.3：令n＝n+1，执行整个循环；
步骤2.4：执行内层第一个循环，根据式(4)更新u k；
Figure PCTCN2018072533-appb-000019
其中，
Figure PCTCN2018072533-appb-000020
为微震信号x的时序序列x(t)的傅立叶变换，
Figure PCTCN2018072533-appb-000021
步骤2.5：令k＝k+1，重复步骤2.4，直到k＝K，结束内层第一个循环；
步骤2.6：执行内层第二个循环，根据式(5)更新ω k；
Figure PCTCN2018072533-appb-000022
步骤2.7：令k＝k+1，重复步骤2.6，直到k＝K，结束内层第二个循环；
步骤2.8：执行外层循环，根据式(6)更新λ；
Figure PCTCN2018072533-appb-000023
其中，τ为拉格朗日乘法算子λ(t)的更新步长参数；
步骤2.9：重复步骤2.3至步骤2.8，直到满足迭代停止条件如式(7)所示，结束整个循环，输出结果，得到K个变分模态分量；
Figure PCTCN2018072533-appb-000024
其中，ε为求解精度；
步骤3：计算待辨识微震信号x(t)的能量分布向量P；
根据公式(8)计算各模态分量u k对应的能量E k；
Figure PCTCN2018072533-appb-000025
其中，x ik(t)(i＝1,2,…N；k＝1,2,…,K；N为采样点个数，K为变分模态个数)表示模态分量u k时序序列的离散点幅值；
根据每个模态分量u k的能量以及待辨识微震信号x(t)的总能量，可以得到模态分量u k的能量百分比值
Figure PCTCN2018072533-appb-000026
从而得到该微震信号的能量分布向量P，即P＝[P(1),…,P(k),…,P(K)]；
步骤4：计算待辨识微震信号x(t)的能量分布X轴的重心系数cx；
根据公式(9)计算能量分布X轴重心系数cx：
Figure PCTCN2018072533-appb-000027
步骤5：根据辨识阈值T识别待辨识微震信号x(t)，若cx>T为矿山煤岩体破裂微震信号，cx≤T为爆破震动信号；
步骤6：自适应更新辨识阈值T的值；
根据方程组(10)更新辨识阈值T：
Figure PCTCN2018072533-appb-000028
其中，W 1为训练集中煤岩体破裂微震信号的cx值集合，W 2为训练集中爆破震动信号的cx值集合。
如图2所示，步骤1获取以时间(s)为横轴，振幅为纵轴，采样频率fs＝1000Hz的待辨识微震信号x(t)，t＝1*1/fs,2*1/fs,…,5000*1/fs，微震信号采样点数据见表1。
表1 监测信号采样点数据(可以存储于Excel表中)
序号 采样点(N) 振幅
1*1/fs 1 4.34E-08
2*1/fs 2 1.69E-07
3*1/fs 3 1.41E-07
4*1/fs 4 -2.43E-07
5*1/fs 5 -6.50E-07
… …
4999*1/fs 4999 -6.10E-06
5000*1/fs 5000 -6.25E-06
按照步骤2的VMD算法对待辨识微震信号x(t)进行变分模态分解，取K＝6，二次惩罚因子α＝2000，分解后的6个变分模态分量及其时频谱如图3所示。
按照步骤3的方法，计算出各模态分量的频带能量，提取各模态分量占原信号的能量百分比值，得到该微震信号的能量分布向量P，即P＝[0.02，0.09，3.01，1.45，13.41，82.02]。
按照步骤4的方法，应用公式(9)计算待辨识微震信号x(t)的能量分布X轴重心系数cx，得到cx＝0.957，图4为该微震信号的能量分布直方图，图中黑色实心圆为该微震信号的能量分布平面重心位置。
按照步骤5的方法，根据辨识阈值T＝0.56以及待检测微震信号x(t)的cx＝0.957对该待测微震信号进行辨识，因为cx>T，所以，该待检测微震信号x(t)辨识为煤岩体破裂微震信号。
测试组中分别给出了15组煤岩体破裂微震信号和15组爆破震动微震信号，煤岩体破裂微震测试信号的能量向量、重心系数及其辨识结果如图5所示；爆破震动微震测试信号的能 量向量、重心系数及其辨识结果如图6所示。根据检测结果，该测试组共30组微震信号，其中29组辨识正确，1组辨识错误，辨识正确率为96.67％。测试组微震信号分类识别结果如图7所示。
按照步骤6的方法，将测试组数据加入已有训练组中，对辨识阈值自适应更新为T＝0.61，以继续提高辨识正确率。
微震信号是非平稳随机信号，其频率分布较为分散，不同种类的微震信号在不同频带能量分布差异显著，因此可根据这一特点，提取微震信号的能量分布向量，通过计算能量分布重心系数并与辨识阈值进行对比，即可实现对待测微震信号的分类辨识。
当然，上述说明并非是对本发明的限制，本发明也并不仅限于上述举例，本技术领域的技术人员在本发明的实质范围内所做出的变化、改型、添加或替换，也应属于本发明的保护范围。

Claims (1)

  1. 一种基于能量分布特征的矿山微震信号辨识方法，其特征在于：包括如下步骤：
    步骤1：读取待辨识微震信号x(t)，其中，t＝1，2，…，N，N为微震信号的采样点个数；
    步骤2：对待辨识微震信号x(t)进行VMD分解，得到K个按照频率从高到低顺序排列的变分模态分量{u 1,…u k,…,u K}：
    对待辨识微震信号x(t)采用VMD分解为K个变分模态分量，约束条件为使各个模态的估计带宽之和最小，且各模态之和等于待辨识微震信号x(t)，约束变分模型描述为式(1)和式(2)：
    Figure PCTCN2018072533-appb-100001
    s.t.∑ ku k＝x(t) (2)；
    其中，x(t)表示待辨识的微震信号，{u k}:＝{u 1,…,u K}代表分解得到的K个有限带宽的变分模态分量，{ω k}:＝{ω 1,…,ω K}表示各分量的频率中心，δ(t)为狄拉克(Dirac)函数，*表示卷积，j 2＝-1；式(2)中，
    Figure PCTCN2018072533-appb-100002
    为对所有的变分模态求和；
    为求解式(1)和式(2)的最优解，引入扩展的Lagrange将约束变分问题变为非约束变分问题，其表达式为式(3)：
    Figure PCTCN2018072533-appb-100003
    其中，α为二次惩罚因子，λ(t)为拉格朗日乘法算子；
    求解该变分问题的具体步骤如下：
    步骤2.1：定义变分模态分量个数K值与惩罚因子α的值；
    步骤2.2：初始化
    Figure PCTCN2018072533-appb-100004
    n＝0；
    步骤2.3：令n＝n+1，执行整个循环；
    步骤2.4：执行内层第一个循环，根据式(4)更新u k；
    Figure PCTCN2018072533-appb-100005
    其中，
    Figure PCTCN2018072533-appb-100006
    为待辨识微震信号x(t)的傅立叶变换，
    Figure PCTCN2018072533-appb-100007
    步骤2.5：令k＝k+1，重复步骤2.4，直到k＝K，结束内层第一个循环；
    步骤2.6：执行内层第二个循环，根据式(5)更新ω k；
    Figure PCTCN2018072533-appb-100008
    步骤2.7：令k＝k+1，重复步骤2.6，直到k＝K，结束内层第二个循环；
    步骤2.8：执行外层循环，根据式(6)更新λ；
    Figure PCTCN2018072533-appb-100009
    其中，τ为拉格朗日乘法算子λ(t)的更新步长参数；
    步骤2.9：重复步骤2.3至步骤2.8，直到满足迭代停止条件如式(7)所示，结束整个循环，输出结果，得到K个变分模态分量；
    Figure PCTCN2018072533-appb-100010
    其中，ε为求解精度；
    步骤3：计算待辨识微震信号x(t)的能量分布向量P；
    根据公式(8)计算各模态分量u k对应的能量E k；
    Figure PCTCN2018072533-appb-100011
    其中，x ik(t)(i＝1,2,…N；k＝1,2,…,K；N为采样点个数，K为变分模态个数)表示模态分量u k时序序列的离散点幅值；
    根据每个模态分量u k的能量以及待辨识微震信号x(t)的总能量，可以得到模态分量u k的能量百分比值
    Figure PCTCN2018072533-appb-100012
    从而得到该微震信号的能量分布向量P，即P＝[P(1),…,P(k),…,P(K)]；
    步骤4：计算待辨识微震信号x(t)的能量分布X轴的重心系数cx；
    根据公式(9)计算能量分布X轴重心系数cx：
    Figure PCTCN2018072533-appb-100013
    步骤5：根据辨识阈值T识别待辨识微震信号x(t)，若cx>T为矿山煤岩体破裂微震信号，cx≤T为爆破震动信号；
    步骤6：自适应更新辨识阈值T的值；
    根据方程组(10)更新辨识阈值T：
    Figure PCTCN2018072533-appb-100014
    其中，W 1为训练集中煤岩体破裂微震信号的cx值集合，W 2为训练集中爆破震动信号的cx值集合。
PCT/CN2018/072533 2017-07-26 2018-01-13 一种基于能量分布特征的矿山微震信号辨识方法 Ceased WO2019019565A1 (zh)

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
CN201710615340.9 2017-07-26
CN201710615340.9A CN107505652B (zh) 2017-07-26 2017-07-26 一种基于能量分布特征的矿山微震信号辨识方法

Publications (1)

Publication Number Publication Date
WO2019019565A1 true WO2019019565A1 (zh) 2019-01-31

Family

ID=60689480

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/CN2018/072533 Ceased WO2019019565A1 (zh) 2017-07-26 2018-01-13 一种基于能量分布特征的矿山微震信号辨识方法

Country Status (2)

Country Link
CN (1) CN107505652B (zh)
WO (1) WO2019019565A1 (zh)

Cited By (10)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111307277A (zh) * 2020-03-20 2020-06-19 北京工业大学 基于变分模态分解和预测性能的单模态子信号选择方法
CN111413588A (zh) * 2020-03-31 2020-07-14 陕西省地方电力（集团）有限公司咸阳供电分公司 一种配电网单相接地故障选线方法
CN111751134A (zh) * 2020-06-22 2020-10-09 西安科技大学 一种基于vmd与rls的采煤机振动信号降噪方法
CN111915865A (zh) * 2020-07-29 2020-11-10 东北大学 一种基于采动震源参数的煤矿复合地质灾害预警方法
CN113985481A (zh) * 2021-10-26 2022-01-28 长江大学 一种基于再约束的变分模态降噪方法及装置
CN116009084A (zh) * 2022-12-15 2023-04-25 深圳市中金岭南有色金属股份有限公司凡口铅锌矿 观测耦合模拟的金属矿深井开采诱发矿震机制分析方法
CN117092692A (zh) * 2023-07-20 2023-11-21 东北大学 适用于钻爆法隧道开挖过程的即时型岩爆时间预警方法
CN120908867A (zh) * 2025-10-10 2025-11-07 长沙飞翼智联科技有限公司 一种矿井下微震智能监测与可视化方法
CN121325274A (zh) * 2025-12-16 2026-01-13 西安奥华电子仪器股份有限公司 用于中子发生器的数据采集与时序同步处理方法及系统
CN121559607A (zh) * 2026-01-22 2026-02-24 中国海洋大学 一种随钻地震钻柱参考信号特征提取方法

Families Citing this family (10)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107505652B (zh) * 2017-07-26 2018-11-13 山东科技大学 一种基于能量分布特征的矿山微震信号辨识方法
CN107941513B (zh) * 2017-12-25 2019-08-06 北京建筑大学 一种列车走行部轴承非平稳运维的时频阶比跟踪方法
CN108846307B (zh) * 2018-04-12 2021-12-28 中南大学 一种基于波形图像的微震与爆破事件识别方法
CN109633566B (zh) * 2019-01-25 2023-08-15 西安电子科技大学 基于vmd算法的电子侦察信号预处理方法
CN109682561B (zh) * 2019-02-19 2020-06-16 大连理工大学 一种自动检测高速铁路桥梁自由振动响应以识别模态的方法
WO2020220416A1 (zh) * 2019-04-28 2020-11-05 山东科技大学 一种基于深度学习的微震信号分类辨识方法
CN110161125A (zh) * 2019-06-17 2019-08-23 哈尔滨工业大学 基于加速度与声发射感知技术相结合的航空发动机智能监测方法
CN110333530B (zh) * 2019-06-25 2020-11-10 广东石油化工学院 一种微震事件检测方法和系统
CN111044619A (zh) * 2019-12-26 2020-04-21 辽宁工程技术大学 一种基于vmd的煤岩体破裂声发射信号处理方法
CN120448783B (zh) * 2025-04-25 2026-03-27 中国矿业大学(北京) 一种矿用微震信号智能识别方法及系统

Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2013169937A1 (en) * 2012-05-08 2013-11-14 Octave Reservoir Technologies, Inc. Microseismic event localization using both direct-path and head-wave arrivals
CN105740840A (zh) * 2016-02-29 2016-07-06 中南大学 一种岩体破裂信号与爆破振动信号的非线性识别方法
CN105956526A (zh) * 2016-04-22 2016-09-21 山东科技大学 基于多尺度排列熵的低信噪比微震事件辨识方法
CN106483563A (zh) * 2015-08-25 2017-03-08 中国石油天然气股份有限公司 基于互补集合经验模态分解的地震能量补偿方法
CN106556865A (zh) * 2016-11-25 2017-04-05 成都理工大学 一种串联型地震信号优化时频变换方法
CN106814396A (zh) * 2017-03-13 2017-06-09 山东科技大学 一种基于vmd的矿山微震信号的降噪滤波方法
CN107505652A (zh) * 2017-07-26 2017-12-22 山东科技大学 一种基于能量分布特征的矿山微震信号辨识方法

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104266894B (zh) * 2014-09-05 2016-12-07 中国矿业大学 一种基于相关性分析的矿山微震信号初至波时刻提取方法

Patent Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2013169937A1 (en) * 2012-05-08 2013-11-14 Octave Reservoir Technologies, Inc. Microseismic event localization using both direct-path and head-wave arrivals
CN106483563A (zh) * 2015-08-25 2017-03-08 中国石油天然气股份有限公司 基于互补集合经验模态分解的地震能量补偿方法
CN105740840A (zh) * 2016-02-29 2016-07-06 中南大学 一种岩体破裂信号与爆破振动信号的非线性识别方法
CN105956526A (zh) * 2016-04-22 2016-09-21 山东科技大学 基于多尺度排列熵的低信噪比微震事件辨识方法
CN106556865A (zh) * 2016-11-25 2017-04-05 成都理工大学 一种串联型地震信号优化时频变换方法
CN106814396A (zh) * 2017-03-13 2017-06-09 山东科技大学 一种基于vmd的矿山微震信号的降噪滤波方法
CN107505652A (zh) * 2017-07-26 2017-12-22 山东科技大学 一种基于能量分布特征的矿山微震信号辨识方法

Cited By (13)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111307277A (zh) * 2020-03-20 2020-06-19 北京工业大学 基于变分模态分解和预测性能的单模态子信号选择方法
CN111307277B (zh) * 2020-03-20 2021-10-01 北京工业大学 基于变分模态分解和预测性能的单模态子信号选择方法
CN111413588A (zh) * 2020-03-31 2020-07-14 陕西省地方电力（集团）有限公司咸阳供电分公司 一种配电网单相接地故障选线方法
CN111751134A (zh) * 2020-06-22 2020-10-09 西安科技大学 一种基于vmd与rls的采煤机振动信号降噪方法
CN111751134B (zh) * 2020-06-22 2021-12-14 西安科技大学 一种基于vmd与rls的采煤机振动信号降噪方法
CN111915865A (zh) * 2020-07-29 2020-11-10 东北大学 一种基于采动震源参数的煤矿复合地质灾害预警方法
CN113985481A (zh) * 2021-10-26 2022-01-28 长江大学 一种基于再约束的变分模态降噪方法及装置
CN113985481B (zh) * 2021-10-26 2023-07-18 长江大学 一种基于再约束的变分模态降噪方法及装置
CN116009084A (zh) * 2022-12-15 2023-04-25 深圳市中金岭南有色金属股份有限公司凡口铅锌矿 观测耦合模拟的金属矿深井开采诱发矿震机制分析方法
CN117092692A (zh) * 2023-07-20 2023-11-21 东北大学 适用于钻爆法隧道开挖过程的即时型岩爆时间预警方法
CN120908867A (zh) * 2025-10-10 2025-11-07 长沙飞翼智联科技有限公司 一种矿井下微震智能监测与可视化方法
CN121325274A (zh) * 2025-12-16 2026-01-13 西安奥华电子仪器股份有限公司 用于中子发生器的数据采集与时序同步处理方法及系统
CN121559607A (zh) * 2026-01-22 2026-02-24 中国海洋大学 一种随钻地震钻柱参考信号特征提取方法

Also Published As

Publication number Publication date
CN107505652A (zh) 2017-12-22
CN107505652B (zh) 2018-11-13

Similar Documents

Publication Publication Date Title
WO2019019565A1 (zh) 一种基于能量分布特征的矿山微震信号辨识方法
CN106814396B (zh) 一种基于vmd的矿山微震信号的降噪滤波方法
CN114779330B (zh) 一种基于微震监测的采掘工作面主裂隙方位分析预测方法
CN107991706B (zh) 基于小波包多重阈值和改进经验模态分解的煤层水力压裂微震信号联合降噪方法
Zhang et al. Identification of blasting vibration and coal-rock fracturing microseismic signals
CN110196448B (zh) 一种滑坡次声信号识别方法
CN110737023B (zh) 一种矿用微震监测信号的处理方法
Chen et al. Joint sound denoising with EEMD and improved wavelet threshold for real-time drilling lithology identification
WO2020220416A1 (zh) 一种基于深度学习的微震信号分类辨识方法
CN113887398A (zh) 一种基于变分模态分解和奇异谱分析的gpr信号去噪方法
CN109799531B (zh) 一种基于地震分频相干属性的裂缝储层预测方法
Zhong et al. A study on the stationarity and Gaussianity of the background noise in land-seismic prospecting
CN107255563A (zh) 实现齿轮箱混合故障信号盲源分离方法
CN106886044A (zh) 一种基于剪切波与Akaike信息准则的微地震初至拾取方法
CN111007569A (zh) 一种集成经验模态分解的样本熵阈值微地震信号降噪方法
Tian et al. A novel identification method of microseismic events based on empirical mode decomposition and artificial neural network features
CN117471528A (zh) 一种基于多峰森林优化算法的hvsr反演方法
Zhang et al. Microseismic source location based on improved artificial bee colony algorithm: performance analysis and case study
CN106771598A (zh) 一种自适应谱峭度信号处理方法
CN106548031A (zh) 一种结构模态参数识别方法
CN116304638B (zh) 一种基于改进lcd的管道泄漏孔径识别方法
CN103913771A (zh) 地震数据处理方法、装置和系统
Wu et al. Statistical significance test of intrinsic mode functions
CN120974303A (zh) 一种基于深度学习与自适应时窗的微震p波初至拾取方法
CN112526611A (zh) 表层地震波品质因子的提取方法及装置

Legal Events

Date Code Title Description
121 Ep: the epo has been informed by wipo that ep was designated in this application

Ref document number: 18837916

Country of ref document: EP

Kind code of ref document: A1

NENP Non-entry into the national phase

Ref country code: DE

122 Ep: pct application non-entry in european phase

Ref document number: 18837916

Country of ref document: EP

Kind code of ref document: A1