CN106897704A - 一种微震信号降噪方法 - Google Patents

Abstract

本发明公开了一种微震信号降噪方法，属于信息处理技术领域，本发明首先对微震信号进行NLMS自适应滤波预处理；其次对预处理后的微震信号进行EEMD分解；再次对分解得到的IMF分量进行阈值去噪；最后重构出降噪后的微震信号。本发明在EEMD分解之前消除脉冲干扰，有效地抑制了由此产生的模态混叠，具有抗噪性能好、计算精度高和算法实时性较强等优点，具有很好的技术价值和应用前景。

Description

一种微震信号降噪方法

技术领域

本发明属于信息处理技术领域，具体涉及一种微震信号降噪方法。

背景技术

微震技术已广泛应用于煤岩动力灾害监测预警、油田压裂、边坡稳定性分析等领域。由于监测环境的复杂性及仪器自身局限性，实际监测到的微震信号在采集、传输过程中受到了大量的噪声干扰，严重影响了微震信号的分析与处理，也给后续的震源定位、能量计算、灾变预警等分析研究工作带来极大困难。另一方面由于微震信号具有类似噪声的频带特性，传统的频谱分析方法并不适用，所以研究微震信号的去噪方法具有重要的理论和实际意义。

目前，针对微震信号所具有的随机性、非平稳性等特点，学者提出了很多信号降噪方法，其中有较大影响的两类方法是：小波阈值降噪法、EMD降噪法，这两类方法的降噪性能已经被证实。其中，小波阈值去噪法具有良好的时频分析能力，计算量小、适用于消除非平稳信号中的噪声，但小波基和分解层数的选择仍没有确定的方法，导致其推广应用困难；EMD去噪法具有很好的自适应性，可以将非线性信号分解为一系列按频率由高到低排列的固有模态函数(Intrinsic Mode Function，IMF)，较好地反映了信号在时间局部的频率特征，它作为一种新的自适应信号时频处理方法，已在地震信号去噪、故障诊断、振动信号分析等多个领域得到广泛应用。但EMD方法在分解的过程中容易产生模态混叠问题，大大影响了EMD方法对微震信号的降噪效果。

发明内容

针对现有技术中存在的上述技术问题，本发明提出了一种微震信号降噪方法，设计合理，克服了现有技术的不足，具有良好的效果。

为了实现上述目的，本发明采用如下技术方案：

一种微震信号降噪方法，采用NLMS自适应滤波器，包括如下步骤：

步骤1：设微震信号的时序序列为x(t)，对微震信号的时序序列x(t)进行NLMS自适应滤波预处理，得到信号x'(t)，其中，t＝1，2，…，n，n为微震信号的采样点个数；

步骤2：给信号x'(t)叠加一组高斯白噪声信号m(t)，得到一个总体信号X(t)，即X(t)＝x′(t)+m(t)；

步骤3：根据公式(1)对总体信号X(t)进行EEMD分解，得到m阶IMF分量；

其中，c_j表示第j个IMF分量，r_n为分解残余项；在EEMD分解过程中，进行多次迭代直到符合某种停止准则而结束，依次得到有限个IMF及残余分量。引入评估函数为:

其中，e_max(t)与e_min(t)分别为极大值包络曲线和极小值包络曲线。设定3个门限值θ₁(默认值为0.05)、θ₂(默认值为0.5)与α(默认值为0.95)，规定当δ(t)小于θ₁的比例达到α，且不存在大于θ₂的值时，筛选迭代终止。

步骤4：根据公式(3)对步骤3中的m阶IMF分量进行阈值处理；

其中，表示经阈值处理降噪后的IMF；T_k为通用阈值，其表达式如公式(4)所示；

其中，C为阈值系数；为第i个IMF分量的估计能量；

步骤5：重构出降噪后的微震信号y'(t)，即

优选地，NLMS自适应滤波器的输出信号x'(t)为输入信号x(t)与滤波器系数的卷积，即为：

x′(t)＝w^T(t)*x(t) (5)；

其中，w(t)为滤波器系数；

所述的滤波器系数w(t)的权值使用输入信号x(t)和误差信号e(t)来更新，计算公式为：

其中，μ为步长因子，用于控制自适应滤波的收敛速度；γ是为了避免x^T(t)x(t)过小而设定的参数；w(t+1)为更新后的滤波器系数；x^T()为信号x(t)的转置。

优选地，所述γ取值范围为[0，1]。

优选地，所述阈值系数C设为0.7。

本发明首先通过NLMS自适应滤波对信号进行预处理，滤除信号中的脉冲干扰；其次，对预处理后的信号进行EEMD分解，EEMD分解将含噪信号的各谐波分量依次分解成噪声主导的IMF集合、噪声与有效信息共存的IMF集合以及定义为低频有效信息的余项，最后对IMF分量进行阈值估计和阈值去噪，从而重构出降噪后的微震信号。

本发明的有益效果为：

本发明在EEMD分解之前消除脉冲干扰，有效地抑制了由此产生的模态混叠，具有抗噪性能好、计算精度高和算法实时性较强等优点，具有很好的技术价值和应用前景。

附图说明

图1为本发明一种微震信号降噪方法的流程框图。

图2为微震信号波形及相应的频谱示意图。

图3为降噪后的微震信号波形及相应的频谱示意图。

具体实施方式

下面结合附图以及具体实施方式对本发明作进一步详细说明：

本发明首先对微震信号x(t)进行NLMS自适应滤波预处理，滤除信号中的脉冲干扰；其次，对预处理后的信号进行EEMD分解，EEMD分解将含噪信号的各谐波分量依次分解成噪声主导的IMF集合、噪声与有效信息共存的IMF集合以及定义为低频有效信息的余项，因此需要在m个IMF分量进行阈值估计和阈值去噪，从而重构出降噪后的微震信号。

如图1所示，一种微震信号降噪方法，具体包括如下步骤：

步骤1：NLMS自适应滤波，

所述微震信号时序序列设为x(t)，t＝1，2，…，n，其中n为所述微震信号的采样点个数，对所述微震信号x(t)进行NLMS自适应滤波预处理，得到信号x'(t)；

x′(t)＝w^T(t)*x(t) (1)

公式(1)中w(t)为滤波器权系数。

所述的滤波器系数w(t)的权值使用输入信号x(t)和误差信号e(t)来更新；

公式(2)中，μ为步长因子，用于控制自适应滤波的收敛速度；γ是为了避免x^T(t)x(t)过小而设定的参数，γ取值范围为[0，1]。

步骤2：通过给信号x'(t)叠加一组高斯白噪声信号m(t)获得一个总体信号X(t)；

X(t)＝x′(t)+m(t) (3)

步骤3：对信号X(t)进行EEMD分解，得到m阶IMF分量；

其中，c_j表示第j个IMF分量，r_n为分解残余项；在EEMD分解过程中，进行多次迭代直到符合某种停止准则而结束，依次得到m个IMF及残余分量。引入评估函数:

其中，e_max(t)与e_min(t)分别为极大值包络曲线和极小值包络曲线。设定3个门限值θ₁(默认值为0.05)、θ₂(默认值为0.5)与α(默认值为0.95)，规定当δ(t)小于θ₁的比例达到α，且不存在大于θ₂的值时，筛选迭代终止。

步骤4：对m个IMF分量进行阈值处理。

基于EMD阈值去噪的通用阈值；

公式(7)中，C为阈值系数；为第i个IMF分量的估计能量，本实施例阈值系数为0.7。

步骤5：重构出降噪后的微震信号，即y'(t)。

微震信号波形及相应的频谱示意如图2所示，以时刻t为横轴，振幅Amplitude为纵轴，微震信号时序序列表示为x(t)，t＝1，2，…，500，微震信号采样点数据见表1。

表1微震信号采样点数据(可以存储于Excel中)

序号	时刻t/ms	振幅Amplitude/V
			1	1	1.11E-06
2	2	2.76E-06
			3	3	4.46E-06
4	4	5.29E-06
			…	…	…
499	499	1.34E-06
			500	500	1.58E-06

按照步骤1提供的算法编程对微震信号时序序列x(t)进行NLMS自适应滤波，得到信号x'(t)，t＝1，2，…，500。

按照步骤2提供的算法给信号x'(t)叠加一组高斯白噪声信号m(t)获得一个总体信号X(t)，t＝1，2，…，500。

按照步骤3对X(t)进行EEMD分解，得到m阶IMF分量。

按照步骤4运用所设阈值系数0.7计算出通用的阈值，然后将通用的阈值和m个IMF分量进行比较，当IMF分量大于通用的阈值时，则保留。

降噪后的微震信号波形及相应的频谱示意如图3所示，根据步骤4得到的m个IMF分量，按照步骤5重构出降噪后的微震信号y'(t)。

为解决脉冲干扰带来的模态混叠，本发明引入NLMS-EEMD方法对微震信号进行降噪处理，即首先通过NLMS自适应滤波对信号进行预处理，滤除信号中的脉冲干扰；其次，对预处理后的信号进行EEMD分解，EEMD分解将含噪信号的各谐波分量依次分解成噪声主导的IMF集合、噪声与有效信息共存的IMF集合以及定义为低频有效信息的余项，因此需要在m个IMF分量进行阈值估计和阈值去噪，从而重构出降噪后的微震信号。

当然，上述说明并非是对本发明的限制，本发明也并不仅限于上述举例，本技术领域的技术人员在本发明的实质范围内所做出的变化、改型、添加或替换，也应属于本发明的保护范围。

Claims (6)

1.一种微震信号降噪方法，采用NLMS自适应滤波器，其特征在于，包括如下步骤：

步骤1：设微震信号的时序序列为x(t)，对微震信号的时序序列x(t)进行NLMS自适应滤波预处理，得到信号x'(t)，其中，t＝1，2，…，n，n为微震信号的采样点个数；

步骤2：给信号x'(t)叠加一组高斯白噪声信号m(t)，得到一个总体信号X(t)，即X(t)＝x′(t)+m(t)；

步骤3：根据公式(1)对总体信号X(t)进行EEMD分解，得到m阶IMF分量；

$X\left(t\right)=\underset{j=1}{\overset{n}{\Σ}}{c}_j+r_n---\left(1\right);$

其中，c_j表示第j个IMF分量；r_n为分解残余项；

步骤4：根据公式(3)对步骤3中的m阶IMF分量进行阈值处理；

$\stackrel{\&OverBar;}{{\mathrm{imf}}_k}\left(t\right)=\left\{\begin{array}{cc}{\mathrm{imf}}_k\left(t\right),& |{\mathrm{imf}}_k\left(t\right)\&GreaterEqual;T_k|\\ 0,& |{\mathrm{imf}}_k\left(t\right)<T_k|\end{array}\right.,(k=2,...7)---\left(3\right);$

其中，表示经阈值处理降噪后的IMF；T_k为通用阈值，其表达式如公式(4)所示；

$T_k=C{\widehat{\mathrm{\&sigma;}}}_i\sqrt{(2\&times;\ln (n\left)\right)}---\left(4\right);$

其中，C为阈值系数；为第i个IMF分量的估计能量；

步骤5：重构出降噪后的微震信号y'(t)，即

2.根据权利要求1所述的微震信号降噪方法，其特征在于，NLMS自适应滤波器的输出信号x'(t)为输入信号与滤波器系数的卷积，即为：

x′(t)＝w^T(t)*x(t) (5)；

其中，x(t)为输入信号；w(t)为滤波器系数；

根据公式(6)对滤波器系数w(t)的权值进行更新；

$w(t+1)=w\left(t\right)+\frac{\mathrm{\&mu;}}{\mathrm{\&gamma;}+x^T\left(t\right)x\left(t\right)}x\left(t\right)e\left(t\right)---\left(6\right);$

其中，μ为步长因子；γ是为避免x^T(t)x(t)过小而设定的参数；w(t+1)为更新后的滤波器系数；x(t)为输入信号；x^T(t)为信号x(t)的转置；e(t)为误差信号。

3.根据权利要求2所述的微震信号降噪方法，其特征在于，所述γ取值范围为[0，1]。

4.根据权利要求1所述的微震信号降噪方法，其特征在于，所述阈值系数C设为0.7。

5.根据权利要求1所述的微震信号降噪方法，其特征在于，在步骤3中，具体包括如下步骤：

步骤3.1：引入评估函数，表达式如(2)所示；

$\mathrm{\&delta;}\left(t\right)=\frac{|e_{max}\left(t\right)+e_{min}\left(t\right)}{|e_{max}\left(t\right)-e_{min}\left(t\right)|}---\left(2\right);$

其中，e_max(t)与e_min(t)分别为极大值包络曲线和极小值包络曲线；

步骤3.2：设定3个门限值θ₁、θ₂与α；

步骤3.3：设置停止准则，规定在δ(t)中小于θ₁的个数在δ(t)中的比例达到α，且不存在大于θ₂的值时，筛选迭代终止；

步骤3.4：得到有限个IMF及残余分量。

6.根据权利要求5所述的微震信号降噪方法，其特征在于，门限值θ₁、θ₂、α的默认值分别为0.05、0.5、0.95。