基于连续小波变换的面波衰减方法研究

岳　龙,刘怀山,尹燕欣,刘　凯

(中国海洋大学海底科学与探测技术教育部重点实验室,山东青岛266100)



基于连续小波变换的面波衰减方法研究

岳龙,刘怀山,尹燕欣,刘凯

(中国海洋大学海底科学与探测技术教育部重点实验室,山东青岛266100)

摘要:面波压制是陆上地震资料处理的主要任务之一。根据面波和有效波在小波域能量和分布区域的不同,提出利用连续小波变换自动压制面波的方法。首先利用大尺度的小波系数自动识别面波区,然后计算每道面波的频带范围,最后对面波区小波系数进行滤波处理并估计出有效信号,根据计算的面波频带范围对估计出的有效信号进行高通滤波得到去除面波后的地震记录。实际地震资料的处理结果表明,与常用的高通滤波方法相比,该方法自动化程度高,压制面波的同时可以减少对有效信号的影响。

关键词：连续小波变换;面波压制;面波区域;尺度范围;高通滤波

面波是陆上地震勘探中普遍存在的规则干扰波[1-2]。面波在单炮记录上一般呈扇形分布,其特点是频率低、视速度低、能量强,具有一定的频带宽度。面波干扰会降低地震资料的信噪比,特别是对深部弱反射波影响较大。因此,如果不能很好地压制面波,会严重影响地震资料处理结果的质量。

针对面波低频、强能量、低视速度等特点,前人已经提出了对应的压制方法。如基于面波和有效波频率差异的高通滤波法,该方法的缺点是对没有受面波影响区域的有效信号也进行了滤波,损害了部分有效信号[3-4];基于面波视速度差异的F-K滤波[5]的面波压制方法,其缺点是切除区的有效信号被同时切除,并且产生了“炕席”状的相干噪声;余波等[6]提出利用径向道滤波方法压制面波等线性干扰;谢金娥等[7]提出了自适应波束压制面波的方法,它是利用有效波与面波在视速度方面的差异,不断优化滤波函数,在面波区实现压制面波的目的;孔庆丰[8]提出利用Hilbert-Huang变换来进行面波压制;胡江涛等[9]提出利用拉东谱约束的射线束实现面波压制;张恒磊等[10]提出在曲波域实现面波压制。由于小波变换具有多尺度特性,也被广泛应用于面波压制。其中,利用一维小波变换压制面波[11]的方法是将面波区小波系数充零,再重构实现面波压制,但是其没有实现面波区域的自动识取,而且没有考虑面波频宽随偏移距变化的问题;张华等[12]通过二维小波变换将地震记录变换到时间、频率、空间和波数四维域中,利用面波低频和高波数的特性,将低频、高波数的分量充零,实现面波压制;陈文超等[13]利用连续小波变换,根据相邻道有效信号在小波时频域分布的相似性,提出自适应面波压制的方法;包乾宗等[14]联合利用连续小波变换和脊波变换实现面波衰减,即首先根据面波强振幅、低频的特征,利用连续小波变换压制面波,然后在脊波域压制剩余的面波;王德营等[15]利用改进的S变换压制面波取得了很好的效果;李庆忠[4]提出的“内切除滤波”是根据面波的强能量特征,自动拾取面波区域,仅在面波区进行高通滤波,自动化程度高。

本文基于前人的研究成果,改进了基于连续小波变换压制面波的算法。首先利用连续小波变换得到地震记录的小波系数;然后利用几个大尺度的小波系数,自动计算每道面波出现的时间;同时考虑到面波频带宽度随偏移距的变化,计算每道面波的频带范围并利用不同的滤波系数对面波区、面波频带范围内的小波系数进行滤波;将压制掉的面波区的小波系数进行重构得到面波估计并从原始数据中减去估计的面波信号得到有效信号估计;根据计算的面波频带范围对有效信号估计进行高通滤波,实现面波压制。最后给出本方法在实际资料中的应用效果,并与高通滤波进行了对比。本文方法的特点是:自动识别面波区,只对面波频带宽度内的小波系数进行处理,减小了对有效信号的影响。

1方法实现

1.1连续小波变换

连续小波变换具备时域和频域双重良好的局部性和随尺度变化的自动调焦功能[16],而自动调焦可以通过母小波的伸缩和平移来实现。首先给出母小波函数ψ(t)的定义,把满足(1)式的小波称为母小波函数:

(1)

式中:Ψ(ω)是ψ(t)的傅里叶变换。

对ψ(t)进行平移和伸缩,得到一簇小波函数ψa,b(t)。

(2)

式中:a为尺度参数;b为平移参数。

对于给定的信号f(t)∈L2(R),其连续小波变换的定义为:

(3)

在计算过程中,规定连续小波变换尺度变量的表达式为:

(4)

式中:scal为尺度的个数;s=scal,scal-1,…,1;fw为母小波的中心频率[17],其表达式为:

(5)

其中,ξ(f)是母小波的傅里叶变换,f1和f2分别为母小波频带的上、下限。

每个尺度下的小波函数都有一个中心频率,其表达式为:

(6)

式中:fa为尺度a对应的频率;dt为待分析信号的采样间隔。在(4)式到(6)式的约束下,尺度变量对应的频率fa的采样是等间隔的。某个尺度下小波变换的结果即是待分析信号与具有某个确定中心频率的小波函数的相似性,反映出待分析信号在该频率下的能量分布特征。

利用(3)式的小波系数,可以重构原信号:

(7)

为了便于表示,文中提到的尺度都是(4)式中尺度变量a的序号。

1.2面波区的确定

一般情况下,有效波与面波的频带范围有差异[18],而且在小波域的有限个大尺度范围内面波能量要明显强于有效波[19]。因此可以利用有限个大尺度的小波系数来确定面波区,具体步骤如下:

1) 选择中远偏移距的含面波的共偏移距道集或者共炮点道集,进行频谱分析。由于面波频带一般为4～18Hz[20],计算18Hz以内的振幅谱最大值对应的频率即为面波主频F0,依据(4)式到(6)式计算面波主频F0对应的尺度S0以及18Hz对应的尺度S18。

2) 对单道地震记录进行Hilbert变换,得到其解析信号;然后对该解析信号进行连续小波变换得到尺度S18到尺度Smax的小波系数,其中Smax是尺度变量a的最大值的序号。

3) 根据步骤1)计算的面波主频对应的尺度S0,依次提取每道记录的尺度S0到Smax的小波系数的模,并计算得到平均值Mf(b),如(8)式所示;然后以Mf(b)的平均值的1/2为门槛,刚刚大于该门槛值且在Mf(b)的峰值之前距离Mf(b)的峰值最近的点就是面波存在的起始时间。

(8)

4) 重复步骤2)到步骤3)完成所有地震道的计算,得到面波区的范围曲线。为了减小小波反变换产生的边界效应,对每道的面波时间点进行修正,在计算的面波时间点附近向前搜索第1个过零点,并将面波时间点修正为该过零点位置;以中间放炮为例,震源两侧曲线变化最大的两个点所对应的地震道就是面波存在的边界道,曲线时间点以下是面波区域。

5) 当工区地表地质条件变化比较大或者施工参数变化时,需要每隔一定的炮数,更新面波主频F0,然后再重复步骤2)到步骤3),继续面波区的拾取。

当同一工区的面波视速度变化不大且相邻炮号在空间上的分布也是相邻的情况下,为了提高计算效率,可以在邻域内每10炮计算一次面波区。如果为了提高面波压制的准确性而需要计算每炮记录的面波区,可以在上一炮的基础上计算下一炮面波区,即在上一炮面波区的边界道两边各扩充10道来计算下一炮的面波区,减少面波区的计算道数来提高效率。

1.3面波尺度范围确定

地震波能量因地层吸收衰减,随着传播距离的增大逐渐减小。而高频成分衰减快,低频成分衰减慢,因此面波主频变化小[21],而频宽随偏移距的变化较大。鉴于此,需要计算不同偏移距面波的尺度范围(频带范围)。假定面波频带范围为f1～f2,根据频率和尺度的反比关系,其对应的尺度范围为S2～S1,其中,f1对应S1,f2对应S2。为了更好地压制面波,在计算过程中指定S1为最大的尺度Smax。面波尺度范围计算步骤如下:

1) 通过1.2节计算出每道面波出现的时间TS,将TS之后的小波系数Wf(a,b)的模沿时间方向叠加在一起,得到面波区小波系数随尺度的变化Pf(a):

(9)

2) 为了计算面波频带高频端对应的尺度,以Pf(a)最大值的1/4为门槛值,刚刚大于该值且在Pf(a)峰值之前的点所对应的尺度为S2,其计算示意图如图1所示。根据实际资料的不同,需要对该门槛值进行试验,选取合适的门槛值。

图1　面波尺度范围计算示意图解

3) 相较于远道,近道的反射波和面波能量差别较小甚至强于面波能量,导致在时频域有效波和面波在尺度方向的界限不明显。因此拾取的面波尺度S2会比较小,超过面波的频带范围。由于只计算尺度S18(频率18Hz对应的尺度)到Smax的小波系数,所以当步骤2)计算的门槛值曲线与Pf(a)在Pf(a)峰值之前没有交点时,将S2定为S18,避免出现S2过小的情况。

4) 依次完成面波区所有道的计算,得到面波的尺度范围。

1.4面波压制

在计算出面波区和不同偏移距面波频宽之后,需要在时频域进行面波压制。李庆忠[4]提出的“内切滤波法”在压制面波的过程中,没有将低频面波完全充零,而是去除低频面波的98%。本文也采用类似的方法,没有对面波区的小波系数完全充零,而对面波区小波系数采用不同的滤波系数进行处理。滤波系数的设计原则是:面波主频对应的滤波系数为0.01,面波高截频对应的滤波系数为0.03,然后通过二次多项式拟合得到面波其它尺度的滤波系数。

为了得到去除面波后的记录,陈文超等[13]提出通过原始记录减去估计的有效反射波得到估计的面波成分,再对估计的面波成分进行低通滤波并从原始记录中减去低通滤波之后的面波成分实现面波压制。本文借鉴了这种思路,具体步骤如下:

1) 根据面波主频和滤波系数设计原则,确定面波各尺度下小波系数的滤波系数;

2) 按照不同的滤波系数,对炮集记录面波区、面波频带范围内的小波系数进行滤波处理,即将小波系数与不同的滤波系数相乘;

3) 从原始时频谱中减去步骤2)滤波后的时频谱,得到面波时频谱估计;

4) 逆变换步骤3)得到的时频谱,获得面波信号的估计;

5) 从原始记录中减去面波信号的估计得到有效信号估计,并对估计的有效信号成分进行高通滤波,高通滤波的参数根据计算得到的面波频带宽度确定,得到去除面波后的记录;

6) 在面波压制区和非面波区的间断点,为了应对因为滤波产生的相邻数据点突变,在面波压制后,采用在面波时间点附近7点平滑的处理方法,处理数据的突变。

因为工区和施工参数的差异,实际资料面波主

频可能不同。因此需要根据实际资料的情况,实时更新面波主频和滤波函数。假定面波主频为4Hz,频带1～18Hz的情况下,利用多项式拟合得到各个尺度的滤波函数为:

ffilt=0.0001×1.02x2-0.0514x+6.49

(10)

式中:ffilt为滤波函数;x为面波尺度范围。

2实际资料处理

2.1关键步骤与参数

压制面波过程中,首先需要确定小波变换的两个关键参数:母小波和尺度变量。本次研究采用高斯函数的二阶导数作为母小波,尺度变量的表达式如(4)式所示,scal设定为256。

选取某油田陆上二维地震资料进行面波压制试验。该资料采集参数为:中间放炮,记录道数为180道,最小偏移距为150m,最大偏移距为4600m,道间距为50m,采样间隔为2ms,记录长度为3s,共有215炮。

压制面波前,首先对面波在时频域的分布和能量特征进行初步分析。由图2a可知,大尺度范围内的时频谱不存在异常能量,而图2b在2000ms之后,几个大尺度范围内存在比较强的面波能量。

图2　某油田陆上二维地震资料单道记录及其时频谱a 不含面波; b 含面波

图3　单炮地震记录面波区

按照1.2节的步骤确定可靠的面波区如图3 所示。对曲线变化最大的两个点所对应的道之外的区域不作处理,仅在三角形区域内进行面波压制。确定面波区之后,计算得到每道面波的尺度范围(图4)。其中面波主频(红线)对应的尺度为252,蓝色曲线是计算的面波频带的高频端对应的尺度,绿色直线是尺度238。第90道为中间道,从图4可以看出,随着偏移距的增大,面波的频带宽度逐渐减小。

对地震记录分别采用固定的面波频带宽度和随偏移距变化的频带宽度(图4)进行面波压制,以第75道和第110道为例,两者压制的面波区如图5所示。图5中红线是采用随偏移距变化的频宽的面波压制区边界,绿线是固定频宽的面波压制区的边界。随着偏移距变大,面波频宽固定的处理方式由于处理的频带宽度要大于计算的面波频带宽度,会压制掉一部分有效信号,如图4中A,B部分即为被压制掉的面波频带宽度之外的部分尺度的小波系数。

图4　面波尺度范围变化情况

经过上述步骤,得到了面波区和面波尺度范围。再根据设计的面波区小波系数的滤波函数完成面波压制。

图5　地震记录第75道(a)和第110道(b)压制面波区

2.2应用效果分析

为了检验该方法的有效性,选择单道记录和单炮记录来检验面波压制效果。图6是单道记录面波压制效果图。图6a是原始地震记录及其时频谱,在2000ms以后,尺度250附近时频谱能量比较强,该区域即是面波时频谱分布区域;图6b是采用本文的思路,压制面波后的地震记录及其时频谱,可见2000ms以后的强能量被压制掉,时间域信号也显示面波被压制;图6c是压制掉的面波及其时频谱。图7a是人工识别面波区域并在该区域内利用高通滤波压制面波的结果;图7b是采用本

文方法自动识别面波区并压制面波的结果。由图7a 和图7b可见,高通滤波的滤波效果与选择的高截频有很大关系,大的高截频在压制面波的同时会损失有效信号。图7c和图7d分别为图7a和图7b 中红色线框区域的放大显示。图7d面波区的反射波同相轴能量要比图7c反射波同相轴能量强,说明本文方法在一定程度上保护了有效信号。

对高通滤波和本文方法压制面波后的数据进行常规叠加处理,处理过程中没有进行其它去噪处理,得到的叠后时间剖面如图8所示。由图8可见,面波压制后,信噪比得到了提高。高通滤波高截频固定,选取不当会损害有效信号或者造成面波压制不彻底,同时高通滤波人工划定面波区工作效率较低。图9为采用本文方法压制面波前、后的单炮记录频谱分析。由图9可知,本文方法很好地压制了低频段面波,而其它频段的有效波基本没有损伤。

图6　单道记录面波压制效果a 面波压制前记录及其时频谱; b 面波压制后记录及其时频谱; c 压制的面波及其时频谱

图7　面波压制效果对比a 高截频18Hz高通滤波去噪结果; b 本文方法的去噪结果; c,d分别为a,b中红色线框区域的放大显示

图8　未压制面波(a)和经高通滤波(b)以及本文方法(c)处理后的叠加剖面

图9　采用本文方法压制面波前、后频谱对比(红线代表面波压制前,蓝线代表面波压制后)

实际资料的应用结果表明,利用连续小波变换确定的面波区不受面波视速度空间变化的影响,拾取的面波区域准确。同时考虑到面波能量的衰减,计算出不同偏移距处面波的频带范围,在面波区和面波频带范围内压制面波,减小了对有效信号的影响。

3结束语

连续小波变换压制面波的改进算法,关键在于确定面波区和面波的尺度范围。本文利用面波与有效波在时频域能量分布的差异,提取几个大尺度的小波系数,利用小波系数均值的一半作为阈值,确定面波出现的时间点;然后再计算面波的尺度范围,在相应尺度范围内压制面波。实际资料应用结果表明:该方法可以自动识别面波区域,确定面波的尺度范围,压制面波的同时尽量减小对有效信号的影响。

参考文献

[1]夏洪瑞,朱海波,邹少峰,等.山前带地震资料噪声消除中存在的问题与对策[J].石油物探,2012,51(6):562-569

XIA H R,ZHU H B,ZOU S F,et al.Actuality and countermeasure of denoising for seismic data in foothill area[J].Geophysical Prospecting for Petroleum,2012,51(6):562-569

[2]吕景贵,刘振彪,管叶君,等.压制叠前相干噪音的速度变换域滤波方法[J].石油物探,2001,40(4):94-99

LV J G,LIU Z B,GUAN Y J,et al.Prestack coherent noise suppression via filtering in different velocity domains[J].Geophysical Prospecting for Petroleum,2001,40(4):94-99

[3]刘志刚,谢言光,陈峰.地震数据处理中噪声衰减方法的探讨[J].石油地球物理勘探,2009,44(增刊1):67-71

LIU Z G,XIE Y G,CHEN F.Discussion on noise attenuation methods in seismic data processing[J].Oil Geophysical Prospecting,2009,44 (S1):67-71

[4]李庆忠.走向精确勘探的道路[M].北京:石油工业出版社,1993:107-112

LI Q Z.The way to obtain a better resolution in seismic prospecting[M].Beijing:Petroleum Industry Press,1993:107-112

[5]张立,龚立,肖敏.CXZB地区浅层地震资料处理实例[J].石油物探,2002,41(4):479-483

ZHANG L,GONG L,XIAO M.Shallow seismic data processing:a case study in CXZB area[J].Geophysical Prospecting for Petroleum,2002,41(4):479-483

[6]余波,黄中玉,谈大龙,等.径向道滤波法去线性干扰[J].石油物探,2005,44(2):109-112

YU B,HUANG Z Y,TAN D L,et al.RT filter attenuating the linear noise[J].Geophysical Prospecting for Petroleum,2005,44(2):109-112

[7]谢金娥,郭全仕,刘财,等.高密点地震资料面波压制的自适应波束方法[J].石油物探,2009,48(2):110-114

XIE J E,GUO Q S,LIU C,et al.Self-adaptive beam method for suppressing ground roll in high-density seismic data[J].Geophysical Prospecting for Petroleum,2009,48(2):110-114

[8]孔庆丰.基于Hilbert-Huang变换的面波压制方法研究[J].石油物探,2012,51(5):446-450

KONG Q F.Surface wave suppressing method based on Hilbert-Huang transform[J].Geophysical Prospecting for Petroleum,2012,51(5):446-450

[9]胡江涛,王华忠,王雄文.拉东谱约束的射线束形成方法及应用[J].石油物探,2013,52(4):335-338

HU J T,WANG H Z,WANG X W.Beam forming by Radon spectrum constraint and its application[J].Geophysical Prospecting for Petroleum,2013,52(4):335-338

[10]张恒磊,刘天佑,李红巧.Curvelet域面波衰减方法研究[J].中南大学学报(自然科学版),2011,42(8):2372-2378

ZHANG H L,LIU T Y,LI H Q.Attenuation of surface wave in Curvelet domain[J].Journal of Central South University(Science and Technology),2011,42(8):2372-2378

[11]ANDREW D J,DOYLE W R.Ground-roll suppression using the wavelet transform[J].Geophysics,1997,62(6):1896-1903

[12]张华,潘冬明,张兴岩.二维小波变换在去除面波干扰中的应用[J].石油物探,2007,46(2):147-150

ZHANG H,PAN D M,ZHANG X Y.Application of 2-D wavelet transformation in eliminating surface wave interference[J].Geophysical Prospecting for Petroleum,2007,46(2):147-150

[13]陈文超,高静怀,包乾宗.基于连续小波变换的自适应面波压制方法[J].地球物理学报,2009,52(11):2854-2861

CHEN W C,GAO J H,BAO Q Z.Adaptive attenuation of ground roll via continuous wavelet transform[J].Chinese Journal of Geophysics,2009,52(11):2854-2861

[14]包乾宗,李庆春,陈文超,等.联合连续小波变换和脊波变换的面波衰减方法及应用[J].石油物探,2011,50(4):367-372,397

BAO Q Z,LI Q C,CHEN W C,et al.Continuous wavelet transform and ridgelet transform joint ground-roll attenuation method and its application[J].Geophysical Prospecting for Petroleum,2011,50(4):367-372,397

[15]王德营,李振春,王姣.S域时频滤波分析与改进[J].石油物探,2015,54(4):396-403

WANG D Y,LI Z C,WANG J.The analysis and improvement on time-frequency filtering in S-transform domain[J].Geophysical Prospecting for Petroleum,2015,54(4):396-403

[16]岳龙,刘怀山,刘凯,等.基于时频分析的初至拾取方法研究[J].石油物探,2015,54(5):508-520

YUE L,LIU H S,LIU K,et al.First-break picking based on time frequency analysis[J].Geophysical Prospecting for Petroleum,2015,54(5):508-520

[17]樊计昌,刘明军,海燕,等.计算尺度函数和小波函数

中心频率的GUI及其应用[J].科技导报,2007,25(24):36-39

FAN J C,LIU M J,HAI Y,et al.GUI for computing center-frequency of the scale and wavelet functions and its application[J].Science Technology Review,2007,25(24):36-39

[18]熊永学,黄青宇,晁元萍,等.FOCUS系统在低信噪比资料处理中的去噪方法探讨[J].石油物探,2002,41(增刊1):290-294

XIONG Y X,HUANG Q Y,CHAO Y P,et al.Denoising method discussion of low signal to noise ratio seismic data processing using FOCUS system[J].Geophysical Prospecting for Petroleum,2002,41(S1):290-294

[19]罗国安,杜世通.小波变换及信号重建在压制面波中的应用[J].石油地球物理勘探,1996,31(3):337-349

LUO G A,DU S T.Application of wavelet transform and signal reconstruction in surface wave elimination[J].Oil Geophysical Prospecting,1996,31(3):337-349

[20]夏洪瑞,吕秋玲,陈胜红,等.小波域爬山拟合法消除面波[J].石油地球物理勘探,2013,48(3):345-350

XIA H R,LV Q L,CHEN S H,et al.A fitting method via mountain climbing in wavelet domain and ground roll elimination[J].Oil Geophysical Prospecting,2013,48(3):345-350

[21]李文杰,魏修成,宁俊瑞,等.基于频率衰减特性的面波压制方法[J].石油物探,2008,47(3):225-227

LI W J,WEI X C,NING J R,et al.A method for suppressing surface wave by considering the charcateristic of frequency attenuation[J].Geophysical Prospecting for Petroleum,2008,47(3):225-227

(编辑：陈杰)

Attenuation of ground roll based on continuous wavelet transform

YUE Long,LIU Huaishan,YIN Yanxin,LIU Kai

(KeyLabofSubmarineGeosciencesandProspectingTechniques,MinistryofEducation,OceanUniversityofChina,Qingdao266100,China)

Abstract:Ground roll attenuation is one of the main tasks of land seismic data processing.According to the energy and distribution characteristics of ground roll and effective waves in wavelet domain,a method is presented to attenuate ground roll using continuous wavelet transform.Firstly,wavelet coefficients of large scale were used to recognize the ground roll region automatically.Then,the ground roll frequency band of each trace was calculated.Finally,the ground roll wavelet coefficients of ground roll region were filtered and the effective signal was approximately obtained.Moreover,the effective signal was filtered using the precalculated ground roll frequency band.Real data examples show that our method has high-degree automation and can attenuate the ground roll effectively with less effective signal loss compared with conventional high-pass filtering methods.

Keywords:continuous wavelet transform,ground roll attenuation,ground roll region,scale range,high-pass filtering

文章编号：1000-1441(2016)02-0214-09

DOI:10.3969/j.issn.1000-1441.2016.02.007

中图分类号：P631

文献标识码：A

基金项目：国家自然科学基金项目(41176077,41230318)和国家高技术研究发展计划(863计划)项目(2013AA092501)共同资助。

作者简介：岳龙(1988—),男,博士在读,主要从事地震勘探资料处理方法研究工作。通讯作者：尹燕欣(1984—),女,讲师,从事地震资料解释工作。

收稿日期：2015-04-13;改回日期：2015-10-09。

This research is financially supported by the National Natural Science Foundation Project (Grant Nos.41176077 and 41230318) and the National High-tech R&D Program (863 Program) (Grant No.2013AA092501).