地球物理学报  2013, Vol. 56 Issue (2): 626-634   PDF    
过渡内蕴模态函数对经验模态分解去噪结果的影响研究及改进算法
李月 , 彭蛟龙 , 马海涛 , 林红波     
吉林大学, 长春 130012
摘要: 经验模态分解(Empirical Mode Decomposition, EMD)是一种具有较大应用潜力的去噪算法.目前, 该算法存在的一个较大问题是过渡内蕴模态函数(Intrinsic Mode Function, IMF)中混叠噪声不能有效处理.过渡内蕴模态函数中混叠噪声不易剔除, 限制了该算法的应用.本文针对此问题, 通过研究过渡IMF的特点, 首次提出一种有效去除过渡IMF中混叠噪声的方法.该方法首先对原信号进行一次EMD处理, 得到包含过渡IMF的初步去噪结果, 并将其与合适的余弦信号结合, 改变其包络分布, 然后对其结果再次进行EMD处理, 仿真实验表明该方法在保留有效信号的同时, 可以有效的去除过渡IMF中混叠的噪声, 并将该方法用于实际地震资料随机噪声压制, 处理效果令人满意.
关键词: 地震信号处理      经验模态分解      余弦函数      内蕴模态函数      过渡内蕴模态函数     
Study of the influence of transition IMF on EMD do-noising and the improved algorithm
LI Yue, PENG Jiao-Long, MA Hai-Tao, LIN Hong-Bo     
Jilin University, Changchun 130012, China
Abstract: The Empirical Mode Decomposition is a good method of do-noising which has great application potential. The existing do-noising methods based on EMD (Empirical Mode Decomposition) cannot effectively deal with the transition IMF (Intrinsic Mode Function) from the EMD of the original signal. Because of the noise aliasing in the transition IMF cannot be easily removed, the application of EMD in do-noising is limited. To solve this problem this paper studied the influence of transition IMF on EMD do-noising result, and proposes a method to remove the noise aliasing in the transition IMF. The envelope structure of inchoate do-noising results processe dby single EMD woul dbe changed with an opportune cosine signals, then the signal received is decomposed again by EMD and the noise aliasing in the transition IMF is remored. The do-noising results of analog signal show that it can effectively keep the useful signal and remove the noise aliasing in the transition IMF. Using this method to do-noise the actual seismic signal, it can suppress random noise and improve the SNR (signal to noise ratio) of seismic signal..
Key words: Seismic signal processing      Empirical mode decomposition      Cosine function      Intrinsic mode function      Transition intrinsic mode function     
1 引言

地震勘探中,由于噪声的存在,得到的地震数据受到严重的污染.噪声严重影响着地震记录信噪比,提高地震资料信噪比是油气精细勘探的重要保障.

Huang等[1]提出的EMD算法,它依据信号自身的时间尺度特征将混有噪声的原始信号分解为有限个IMF,经分解后的高频IMF分量通常为包含的噪声,因此,EMD分解方法可以有效地进行滤波[2-3].同Fourier变换和小波变换相比,EMD在客观性和分辨率方面具有明显的优越性,对信号处理后能够提取更多、更接近实际的特性[4-7]. EMD方法自提出以来获得了广泛的应用[8-12],并取得较好的研究成果.在地震信号[13-14]处理方面,Bradley等[10]将EMD和HHT(Hilbert-Huang Transform,希尔伯特-黄变换)用于地震反射信号的分析,评估了EMD和HHT对地震信号在时域和时频域的处理能力,并利用EMD的瞬时特性进行滤波,提高地震资料的信噪比,证明了EMD方法用于地震信号处理的可行性及有效性;Kaslovsky和Meyer[16]对含噪信号进行了EMD研究,分析此时EMD得到的IMF的特点,并得出由噪声主导的IMF和由信号主导的IMF之间存在过渡IMF的结论,并用模拟信号对此结论进行了验证,讨论了过渡IMF对瞬时频率估计的影响.基于EMD的去噪方法如果不能有效处理过渡IMF,则存在去噪效果不理想的问题.国内外学者将EMD和HHT用于地震模型分析[15, 17]、信号分析[5-8, 11, 18-20]以及特征提取[21-23]等方面,但至今对过渡IMF少有研究.

针对这一问题,本文通过分析EMD的特点及含噪信号IMF的特点,研究了过渡IMF对EMD去噪结果的影响,并提出一种有效去除过渡IMF中混叠噪声的方法.将该方法运用于模拟数据和实际地震信号的去噪处理,相比于一次EMD去噪处理后的结果,该方法的结果可以更彻底地去除信号中混杂的噪声,提高信噪比.

2 过渡IMF对EMD去噪结果的影响 2.1 EMD基本原理

经验模态分解的假设条件是所要分解的信号是由许多不同的内蕴模态函数叠加而成的[1].根据信号的局部特征,EMD将信号分解为不同时间尺度的IMF,分解过程是自适应的[24].

EMD是一个筛选的过程,筛选得到的每个IMF需要满足以下两个条件:

1) 整个IMF的极值点和过零点数目相等或者最多相差1;

2) IMF的极大值和极小值对应包络线均值为零.

IMF的第二个条件在实际信号处理中难以完全满足,实际应用中作近似处理.本文要求包络线均值与相应IMF的比例小于一个给定值.

对某一实信号s(t),假设原始残差信号为r0(t)=s(t),通过EMD得到IMF的步骤为:1)找到待分解信号s(t)的所有极值点,然后对所有极大值点和极小值点分别进行三次样条插值拟合得到信号的上包络线和下包络线,并求取上下包络线的均值m(t).

2) 用残差信号r0(t)减去上下包络线的均值m(t),得到h1(t),

(1)

判断h1(t)是否满足IMF的两个条件;若不满足,则将h1(t)作为新的待分解信号重复(1)(2)两步骤;若满足,则认为h1(t)是一个IMF,并令

(2)

3) 去除已得到的IMF,得到残差信号r1(t),

(3)

判断r1(t)是否满足分解结束条件,分解结束条件是剩余信号r1(t)是单调的或者是直流信号.若满足,则结束分解,否则将r1(t)作为新的待分解信号,重复以上1)2)3)步骤,直至筛选出所有的IMF.

分解结束后,分解出的n个IMF包含了原信号从高到低不同频带的成分,残差rn(t)则代表了信号的平均趋势.原信号s(t)可以表达为

(4)

EMD分解得到的一系列IMF代表的是信号的一系列时间尺度的特征模态,是信号内在的模态特征.

2.2 EMD去噪处理过渡IMF存在的问题

实际应用中,有用信号往往处于低频部分[25],而EMD得到的IMF是从高频到低频依次筛选出来的,对于混有随机噪声[26-27]的信号来说,分解得到的高频IMF即认为是噪声.已有方法通过去除高频IMF,用剩余的IMF对原信号进行重构,实现去噪.但是简单通过舍弃某些IMF实现去噪,存在一些不足之处.

Wu和Huang[28]对均匀白噪声序列进行EMD分析,并将结果拓展到一般的白噪声;Flandrin P等[29]对分形噪声进行了EMD研究.从已有噪声EMD的研究结果以及文献[16]可知,简单的去除高频IMF存在噪声主导的IMF和信号主导的IMF之间的过渡IMF不能进行有效处理的问题.

对模拟Ricker子波进行EMD验证已有方法存在的问题. Ricker子波谱峰值频率为35Hz,峰值为1;Ricker子波信号采样时间为ms,数据长度为1024点.

图 1 原信号和含噪信号及其频谱 Fig. 1 Original signal and signal with noise and its spectra

对含噪信号进行EMD,得到各个IMF,并计算各IMF的频谱,结果如图 2所示.EMD的结束条件为残差信号rn(t)是单调的或其能量小于原信号的10%.

图 2 第一次EMD得到的IMF及其频谱 Fig. 2 IMF of first EMD and its spectra

图 2与纯信号频谱对比可以看出EMD得到的第三个IMF中信号和噪声不能有效的分离开.

本文主要解决过渡IMF中混叠噪声的问题.已有基于EMD的去噪方法不能把过渡IMF中有效信号与噪声有效分离,若保留过渡IMF则可以较好的保留有效信号,但噪声也会被保留下来;若去除过渡IMF,则可更好的去除噪声,但同时会损失有效信号.故一次EMD处理得到的包含过渡IMF的结果存在去噪不理想的问题,而且结果中噪声能量较低,不易分离.针对这个问题,本文提出一种有效去除过渡IMF中混叠噪声并保留有效信号的方法.

3 改进的EMD去噪算法

本文方法将包含过渡IMF的第一次EMD处理得到的结果与一个已知的合适的余弦信号相加,并对相加得到的混合信号再次进行EMD处理从而实现有效去除过渡IMF中的混叠噪声,改善地震信号去噪效果,提高信噪比.具体的实现过程如下:

对原始的信号s(t)进行EMD,并将残差信号作为最后一阶IMF,得到一系列IMF,记作c1(t),c2(t),…,cn(t),则有

(5)

根据信号和噪声的特点以及EMD得到的IMF的特征可知,前m个IMFc1(t),c2(t),…,cm(t)(mn)中,噪声对每个IMF的主导作用逐渐降低,而信号对每个IMF的主导作用逐渐增强.因此,首先通过公式(6)计算各IMF与原含噪信号的相关系数确定IMF中信号和噪声的分界点k

(6)

式中ss(t)表示原信号,cici(t)表示第i个IMF,p表示数据长度.

第一次EMD后计算各IMF与原信号的相关系数,确定相关系数的第一个极小值对应的IMF,将此IMF的后一个IMF作为信号和噪声的分界,并认为其为过渡IMF.实验表明,此判断准则是有效的、适用的,能够有效保留原信号.第一次EMD后,将分界点与分界点之后的IMF,即ck(t),ck+1(t),…,cn(t)作为信号部分,用这部分IMF重构信号,则原信号的重构为

(7)

经过第一步去噪后信号的信噪比有所改善,但包含了混叠在过渡IMF中的噪声.第二次EMD处理时,将过渡IMF中频率成分高于有效信号的混叠噪声分离出来并舍弃,即可去除过渡IMF中的混叠噪声.由于EMD是依据信号自身的时间尺度特征将信号分解为有限个IMF,因此在第二次EMD处理之前,如果不对第一次去噪结果进行处理,直接进行EMD处理,则分解得到的IMF就是第一次EMD处理后对重构的各IMF,即为ck(t),ck+1(t),…,cn(t).针对这个问题,将初步去噪结果与一个已知的合适的余弦信号s′(t)=A1cos(2πf1t)相加,得到新的信号

(8)

这样处理后将改变原信号的包络结构,实现信号的有效分解.已知余弦信号参数按以下要求确定.

对于余弦信号的幅值的选取:相加后的信号应该能够突出过渡IMF中混叠的噪声部分,最好满足余弦信号的谱峰值不小于信号的谱峰值.所加余弦信号的频率应该高于有效信号的上限频率,可以根据过渡IMF中混叠噪声的频率成分和过渡IMF的上一阶IMF的频率确定,具体频率范围应包含在这两阶IMF频谱的主要重叠部分.

对于任一实信号c(t)一般可以写成[30]

(9)

其中φ(t)表示信号的相位.若不考虑EMD得到的单调趋势项,即残差信号rn(t),则任一实信号S(t)有

(10)

其中Re表示取实部.

从(10)式可知,EMD实际得到了一个用于信号分解的自适应的广义基,在该分解方法中,基函数是一系列可变幅度和可变频率的正余弦函数,它是由原信号自适应得到的.对s1(t)进行EMD,由于分解过程是根据信号自身的时间尺度特征从高频到低频进行的,因此首先会将加入的相对有效信号为高频的余弦信号s′(t)与过渡IMF中的混叠噪声nIMF(t)分解出来,作为第一阶IMF,即

(11)

通过去除第二次EMD得到的第一阶IMF即可去除过渡IMF中的混叠噪声,改善去噪效果.即

(12)

选择余弦信号改变原信号的包络,是由于余弦信号在时域和频域都较容易辨识,易于去除.若对加入的余弦信号s′(t)选择合适的幅值A1和频率f1,则可以有效去除过渡IMF中的混叠噪声,得到最终去噪后的有用信号s0(t).

4 实验结果及分析

为了验证本文方法对地震信号处理的有效性,首先使用模拟Ricker子波信号进行正演试验.使用的Ricker子波信号如图 1,Ricker子波信号第一次EMD结果如图 2.

第一次EMD之后计算得到的各IMF与原信号的相关系数如表 1.

表 1 各IMF与原信号的相关系数 Table 1 The correlation coefficient of INF and original signal

按照本文分界点判断方法,从表 1可知,IMF3即为过渡IMF.图 2的各IMF波形和频谱也验证了这一结果.因此,使用过渡IMF和信号主导的各IMF对原信号进行重构,得到第一次去噪结果.

将第一次EMD得到包含过渡IMF的结果与一个已知的合适的余弦信号相加,改变信号的包络结构,得到一个新的信号.本文使用的余弦信号幅值和频率分别为0.05和12 Hz.将新的信号进行第二次EMD处理,处理方法和第一次相同,图 3给出第二次EMD得到的IMF的时域图和频谱图(去除余弦信号后的结果).

图 3 第二次EMD得到的IMF及其频谱 Fig. 3 IMF of second EMD and its spectra

图 3可以看出过渡IMF中的噪声部分与信号被分离开,因此可以有效的去除过渡IMF中的混叠噪声.由于此次EMD是为了分解出过渡IMF中混叠的噪声,噪声主导的IMF只有一个,因此各IMF与分解前信号的相关系数可能不满足分界点判断准则.但通过分析信号的特点可以确定第一阶IMF即为第一次EMD分解的过渡IMF中混叠的噪声部分,通过去除第一阶IMF得到去噪结果.

图 4给出了包含过渡IMF的去噪结果、去除过渡IMF的去噪结果以及本文方法去除过渡IMF中混叠噪声之后的结果的对比.

图 4 包含过渡IMF的结果与本文方法的结果 Fig. 4 Result with transition IMF and result of the method of this paper

从图中可以看出,包含过渡IMF的去噪结果可以有效的保留有效信号,但去噪效果不理想;去除过渡IMF的去噪结果虽然能够有效去除噪声,但是信号出现较大衰减;而本文方法在保留有效信号的同时,有更好的去噪效果.原含噪信号的信噪比为0 db,包含过渡IMF的去噪结果信噪比为5.42 db,通过本文方法得到的去噪结果的信噪比为7.36 db,而且本文方法对有效信号具有较好的保持,可见该方法是可行且有效的.

用本文的方法对某地区实际共炮点数据进行处理,此共炮点数据共168道,每道6145个采样点,采样时间间隔为1 ms.图 5给出原始168道地震信号,横向为道序号,纵向为采样点数.图 6是包含过渡IMF的结果,图 7是本文方法处理的结果.

图 5 原地震信号 Fig. 5 Original seismic signal
图 6 包含过渡IMF的处理结果 Fig. 6 Results with transition IMF
图 7 本文方法处理结果 Fig. 7 Results by the method of this paper

图 6中,对于A、B两部分,已有方法处理的结果与原信号相比,由于噪声去除的不彻底,几乎没有改善,但是通过本文方法处理后,更好的去除噪声,轴变得更加清晰了,如图 7中A、B所示.从图 6图 7中的C、D两部分也可以看出本文方法处理后的结果相比于已有方法的结果,由于去除了过渡IMF中的混叠噪声,同相轴变得更加清晰连续.从两幅图中E部分的处理结果的对比明显可以看出本文方法对于随机噪声的压制更彻底,处理效果更好.从图 6图 7的对比,总体上来说本文的方法和已有方法相比,能够在保留有效信号的同时,有效去除过渡IMF中混叠的噪声,实现地震信号的进一步去噪,提高信噪比,具有更优的去噪效果.

5 结论

本文针对已有EMD方法不能有效去除过渡IMF中混叠噪声的问题,研究了过渡IMF对EMD去噪结果的影响,并提出一种改进去噪结果的算法.通过分析EMD以及信号分解所得IMF的特点,将包含过渡IMF的第一次EMD处理得到的结果与一个已知的合适的余弦信号相加得到一个混合信号,再对相加得到的结果进行EMD处理,从而实现在保留有效信号的同时,有效的去除过渡IMF中混叠的噪声,解决了传统EMD方法不能有效去除过渡IMF中混叠噪声问题.理论分析及仿真实验表明本文方法是可行的,并且是有效的.将本文方法应用于实际地震勘探数据处理,去噪效果表明,该方法相比于已有基于EMD的去噪方法有所改进.

参考文献
[1] Huang N E. The empirical mode decomposition and the Hilbert spectrum for nonlinear and nonstationary time series analysis. Proc. R. Soc. Lond A , 1998, 454(4): 903-995.
[2] Boudraa O A, Cexus C J. EMD-Based Signal Filtering. IEEE Transactions on Instrumentation and Measurement , 2007, 56(6): 2196-2202. DOI:10.1109/TIM.2007.907967
[3] 陈凯. 基于经验模式分解的去噪方法. 石油地球物理勘探 , 2009, 44(5): 603–608. Chen K. A new denoising metho dbased on Empirical Mode Decomposition (EMD). Oil Geophysical Prospecting (in Chinese) , 2009, 44(5): 603-608.
[4] Huang N E, Shen Z, Long R S. A new view of nonlinearwater waves:The Hilbert spectrum. Annu Rev. Fluid Mech , 1999, 31: 417-457. DOI:10.1146/annurev.fluid.31.1.417
[5] 杨光亮, 朱元清, 于海英. 基于HHT的地震信号自动去噪算. 大地测量与地球动力学 , 2010, 30(3): 39–42. Yang G L, Zhu Y Q, Yu H Y. Automatic de-noising algorithm of earthquake signal based on HHT decomposition. Journal of Geodesy and Geodynamics (in Chinese) , 2010, 30(3): 39-42.
[6] Kizhner S, Flatley T S, Huang N E, et al. On the Hilbert-Huang Transform Data Processing System Development. IEEE Aerospace Conference Proceedings , 2004, 3: 1961-1979.
[7] Daubechies I, Lu J F, Wu H T. Synchrosqueezed wavelet transforms:An empirical mode decomposition-like tool. Applied and Computational Harmonic Analysis , 2011, 30: 243-261. DOI:10.1016/j.acha.2010.08.002
[8] Mendez M O, Corthout J, Huffel S Van, et al. Automatic screening of obstructive sleep apnea from the ECG based on empirical mode decomposition and wavelet analysis. Physiological Measurement , 2010, 21: 273-289.
[9] Peng Z K, Peter W T, Chu F L. An improved Hilbert-Huang transform and its application in vibration signal analysis. Journal of Sound and Vibration , 2005, 286: 187-205. DOI:10.1016/j.jsv.2004.10.005
[10] Naveed R and Danilo P M. Empirical Mode Decomposition for Trivariate Signals. IEEE Transactions on Signal Processing , 2010, 3(58): 1059-1068.
[11] Saurabh P, Madhuchhanda M. Empirical mode decomposition based ECG enhancement and QRS detection. Computers in Biology and Medicine , 2011, 42: 83-92.
[12] Macelloni L, Battista M B, Knapp C C. Optimal filtering high-resolution seismic reflection data using a weighted-mode empirical mode decomposition operator. Journal of Applied Geophysics , 2011, 75: 603-614. DOI:10.1016/j.jappgeo.2011.09.018
[13] Zhang Z J, Wang G J, Harris J M. Multi-component wavefield simulation in viscous extensively dilatancy anisotropic media. Phys Earth planet Inter , 1999, 114: 25-28. DOI:10.1016/S0031-9201(99)00043-6
[14] Zhang Z J, Lin Y H, Liu E R. Large static correction using a hybrid optimization method in complex terrain:some experience learnt from China. Journal of Seismic Exploration , 2005, 13(4): 337-352.
[15] Battista M B, Knapp C, McGee T, et al. Application of the empirical mode decomposition and Hilbert-Huang transform to seismic reflection data. Geophysics , 2007, 72(2): 29-37.
[16] Kaslovsky N D, Meyer G F. Noise Corruption of empirical mode decomposition and its effect on instantaneous frequency. Adaptive Data Analysis , 2010, 2(3): 373-396. DOI:10.1142/S1793536910000537
[17] 段生全, 贺振华, 黄德济. HHT方法及其在地震信号处理中的应用. 成都理工大学学报(自然科学版) , 2005, 32(4): 396–400.
[18] 徐可君, 秦海勤, 江龙平. 基于EMD和HHT的航空发动机转子-机匣振动信号分析. 振动与冲击 , 2011, 30(7): 237–240. Xu K J, Qin H Q, Jiang L P. Rotor-case vibration signal analysis of an aero engine based on EMD and HHT. Journal of Vibration and Shock (in Chinese) , 2011, 30(7): 237-240.
[19] 谭善文, 秦树人, 汤宝平. Hilbert-Huang变换的滤波特性及其应用. 重庆大学学报 , 2004, 2(27): 237–240. Tan S W, Qin S R, Tang B P. The Filtering Character of Hilbert-Huang Transform and Its Application. Journal of Chongqing University (in Chinese) , 2004, 2(27): 237-240.
[20] 白大为, 底青云, 王光杰. Hilbert-Huang变换与ELF信号处理. 地球物理学进展 , 2009, 24(3): 1032–1038. Bai D W, Di Q Y, Wang G J. Hilbert-Huang transformation and ELF signal processing. Progress in Geophys. (in Chinese) , 2009, 24(3): 1032-1038.
[21] 林玉荣, 王强. 基于一维经验模态分解的图像细节提取方法. 吉林大学学报(工学版) , 2011, 41(6): 1766–1770. Lin Y R, Wang Q. Extracting details from images based on 1-EMD. Journal of Jilin University (Engineering and Technology Edition) (in Chinese) , 2011, 41(6): 1766-1770.
[22] 刘庆敏, 杨午阳, 田连玉, 等. 基于经验模态分解的地震相分析技术. 石油地球物理勘探 , 2010, 45(supp1): 145–149. Liu Q M, Yang W Y, Tian L Y, et al. Seismic facies analysis technique based on the Empirical Mode Decomposition. Oil Geophysical Prospecting (in Chinese) , 2010, 45(supp1): 145-149.
[23] 杨培杰, 印兴耀, 张广智. 希尔伯特-黄变换地震信号时频分析与属性提取. 地球物理学进展 , 2007, 22(5): 1585–1590. Yang P J, Yin X Y, Zhang G Z. Seismic signal time-frequency analysis and attributes extraction based on HHT. Progress in Geophysics. (in Chinese) , 2007, 22(5): 1585-1590.
[24] 王宏禹. 信号处理相关理论综合与统一法. 北京: 国防工业出版社, 2005 : 208 -221. Wang H Y. Synthesis and unification method in correlation theory of signal processing (in Chinese). Beijing: National Defence Industry Press, 2005 : 208 -221.
[25] 陈林, 宋海斌. 地震信号瞬时频率的估算. 地球物理学报 , 2009, 52(1): 206–214. Chen L, Song H B. The estimation of instantaneous frequency of seismic signal. Chinese J. Geophys. (in Chinese) , 2009, 52(1): 206-214.
[26] 李月, 林红波, 杨宝俊, 等. 强随机噪声条件下时窗类型局部线性化对TFPF技术的影响. 地球物理学报 , 2009, 52(7): 1899–1906. Li Y, Lin H B, Yang B J, et al. The influence of limited linearization of time window on TFPT under the strong noise background. Chinese J. Geophys (in Chinese) , 2009, 52(7): 1899-1906.
[27] 鲁来玉, 何正勤, 丁志峰, 等. 华北科学探测台阵背景噪声特征分析. 地球物理学报 , 2009, 52(10): 2566–2572. Lu L Y, He Z Q, Ding Z F, et al. Investigation of ambient noise source in North China array. Chinese J. Geophys. (in Chinese) , 2009, 52(10): 2566-2572.
[28] Wu Z H, Huang N E. A study of the characteristics of white noise using the empirical mode decomposition method. Proc. R. Soc. Lond A , 2004, 454: 1597-1611.
[29] Flandrin P, Rilling G, Goncalves P. Empirical Mode Decomposition as a Filter Bank. IEEE Signal Processing Letters , 2003, 10(20): 1-4.
[30] 孙晖.经验模态分解理论与应用研究[博士论文].浙江:浙江大学, 2005. Sun H. Research on Empirical Mode Decomposition Theoty and Its Application[Ph. D. thesis]. Zhejiang:Zhejiang University, 2005. http://cdmd.cnki.com.cn/article/cdmd-10335-2006073589.htm