多方向矢量中值滤波在多分量地震数据中的应用

寻　超,汪　超,王　赟

(1.中国科学院地球化学研究所,贵州贵阳550081;2.中国科学院大学,北京100049;3.中国地质大学(北京)地球物理与信息技术学院,北京100083)



多方向矢量中值滤波在多分量地震数据中的应用

寻超1,2,汪超1,王赟3

(1.中国科学院地球化学研究所,贵州贵阳550081;2.中国科学院大学,北京100049;3.中国地质大学(北京)地球物理与信息技术学院,北京100083)

多方向矢量中值滤波在压制噪声的同时能够有效地保护多分量地震信号的矢量特性。将多方向矢量中值滤波扩展应用于多分量地震数据去噪处理,把多分量地震数据当作一个矢量波场进行处理,每个时间采样点均为一个矢量。基于地震数据沿同相轴变化缓慢的特点,在矢量空间中,以当前滤波点为中心沿不同试验方向截取相邻道矢量信号片段,根据这些矢量信号片段之间的相关性自适应地选取最佳矢量中值滤波窗。扩展后的方法在压制噪声的同时能够有效保护多分量地震数据各分量之间的相对振幅信息。理论模型和实际数据测试结果表明,多方向矢量中值滤波应用于多分量地震数据去噪处理是有效的。

矢量中值滤波;多方向;多分量地震数据;去噪

压制相干或非相干噪声一直是地震数据处理的关键。在P波处理中,有许多压制噪声的技术,如F-X反卷积[1]、Radon变换[2-3]、保边界滤波[4]、中值滤波[5-6]、曲波变换[7-8]等。这些技术虽然也可用于多分量数据处理,但都是将多分量数据的各个分量当作标量场来处理,容易破坏多分量地震数据的矢量特性。中值滤波是一种常用的非线性滤波,其基本功能是消除非平稳信号中的峰值噪声,对数据起平滑作用。在中值滤波的基础上,ASTOLA等[9]提出了矢量中值滤波(VMF)方法。该方法将彩色图像的每个像素点当作一个三维矢量,利用矢量中值滤波法进行去噪处理,取得了很好的去噪效果;LIU等[10]将矢量中值滤波引入到地震资料处理中。HUO等[11]将矢量中值滤波扩展到多方向矢量中值滤波(MD-VMF),并成功对多炮连续采集的多源地震数据进行了波场分离;LIU等[10]和HUO等[11]将VMF用于处理地震信号时,仍是对多分量地震数据中的各个分量单独处理,而不是将其当作一个矢量波场进行处理。目前,针对多分量地震数据联合处理的技术较少,NAGHIZADEH[12]利用在F-X域建立向量自回归(VAR)模型将传统的F-X去噪方法应用于三分量地震数据。VAR模型能够有效识别多分量地震数据中各分量间的相关性,但在实际数据测试时会损失有效信号。RODRIGUEZ等[13]基于组稀疏约束时频变换提出了一种三分量去噪方法。基于稀疏约束的去噪方法能利用矢量场的振幅信息,但不能利用矢量的方向信息,且计算成本高。WANG等[14]提出多通道顺序统计滤波方法。

我们将多方向矢量中值滤波扩展到多分量地震数据应用中,将多分量地震数据看作一个矢量波场,利用多方向矢量中值滤波同时处理多分量地震数据的各个分量。该方法在压制噪声的同时能够有效保护多分量地震信号的矢量特性。既克服了矢量中值滤波的单方向特性,降低对数据过分拉平的要求,又能够在有效去除噪声的同时很好地保留多分量地震数据所包含的丰富的运动学与动力学信息。

1　方法原理

1.1中值滤波(SMF)

中值滤波方法的实质就是中位数法,可用排序进行定义。取定中值滤波的窗长度为N(N为奇数),在信号序列中以当前滤波点为中心顺次截取N个数据{xi|i=1,2,…,N},将其按照数值大小递增(或递减)的顺序进行排列,中间位置的值,即为中值,并将其作为当前滤波点的输出值。也可根据数据序列中一个数到其它所有数的距离来定义,距离和最小的数即为中值。当采用最小距离定义时,一个数到其它所有数的距离和可用(1)式表示:

(1)

式中:j=1,2,…,N;l表示范数的阶数(l=1或2)。中值xm可由(2)式定义:

(2)

当l=1时,为L1范数,用排序法与最小距离法定义的中值完全相同。然而,当采用非L1范数时,基于两种定义的中值有所不同。例如,设数据序列为{1,2,3,4,5000},当最小距离采用L1范数时,两种方法的中值均为3;而当最小距离采用L2范数时,中值为4[15]。可根据实际情况选取不同的范数,故在实际应用中用最小距离法定义中值较基于排序法的定义更加灵活[9]。

1.2矢量中值滤波(VMF)

矢量中值滤波是对一个矢量集求取中值矢量。

将方程(1)扩展为矢量形式:

(3)

方程(3)与方程(1)相同,只是将标量xi变为矢量Xi,Xi=[xi1,xi2,…,xiM]T,N为矢量集中矢量的个数,M为矢量的维数,即多分量地震数据的分量个数(M=2,3,4或9)。中值矢量Xm定义为:

(4)

因此,中值滤波(SMF)可以看作矢量中值滤波(VMF)的特例。

将矢量中值滤波应用于多分量地震数据时,所选取的滤波矢量集为沿同相轴用长度为W的窗截取的W个矢量点,每个矢量的维数为地震数据的分量个数。以三分量地震数据为例(即M=3),截取的矢量点如(5)式所示(这里假设同相轴水平):

A=

(5)

式中:j表示时间采样点;i表示空间地震道;w为半窗长度,有w=(W-1)/2,W为窗长度,一般为奇数,且W≥3。矩阵A中每一列元素即为一个矢量。

zji,xji,yji分别为三分量地震数据在Z,X,Y分量上相同时间、相同道号上的采样值。常规单分量地震数据去噪采用矢量中值滤波[10,16]时,是将单个分量上沿时间方向截取的n个点看成一个n维向量。而本文将矢量中值滤波扩展应用于多分量地震数据时,是将多分量地震数据中各个分量同一时间同一地震道处的地震数据组成一个矢量,从而能够很好地保持弹性波场在每个物理点上的振动信息。

1.3多方向矢量中值滤波(MD-VMF)

由于矢量中值滤波将地震记录中相邻地震道在相同时间点上的矢量作为待滤波的矢量集,过分依赖于数据的拉平作用,且实际地震数据经校正后,很难实现同相轴真正意义上的拉平。因此,本文将MD-VMF应用到多分量地震数据中。

多方向矢量中值滤波在上述矢量中值滤波的基础上,考虑了实际地震记录中同相轴形态的多样性,即除水平同相轴外,还有倾斜同相轴和弯曲同相轴等,克服了矢量中值滤波所要求的水平同相轴单方向特性。即多方向矢量中值滤波可以根据局部同相轴的形状自适应地调整滤波窗口,从而保证矢量点均是沿局部同相轴选取。

为了能够确定局部同相轴的形状,可沿不同倾斜方向从相邻地震道截取矢量信号片段,根据这些矢量信号片段之间的相似度来确定最优的局部同相轴轨迹。如图1所示,图中的几条倾斜线即为用滤波窗截取相邻道矢量信号时时窗在每道的中心位置。搜索局部同相轴时所需的参数为:滤波窗长度W(一般为奇数且W≥3),最大倾角pmax,最小倾角pmin,倾角增量pi,时窗长度N(一般为奇数且N≥3)。滤波窗长度W是指每个矢量集所涉及的地震道数,时窗长度N是指每个矢量信号片段所包含的时间采样点数。最大倾角pmax与最小倾角pmin共同确定同相轴倾角搜索范围,倾角单位为相邻两道波至时差;倾角增量pi与倾角范围共同决定矢量集的个数(即搜索方向的个数)[16-17]。

图1　滤波时窗划分示意

以当前滤波点为中心,沿p方向截取的矢量信号片段,以三分量地震数据为例(即M=3),如(6)式所示:

(6)

式中:Sji(p)=[zji(p),xji(p),yji(p)]为一个矢量点,p=pmin,…,pmax,ns=(N-1)/2,为半时窗长度。矩阵A(p)中每一列为一个滤波矢量片段。

矢量信号片段之间的相似性由距离函数(7)给定:

(7)

(8)

由(7)式可知,矢量信号片段之间的距离是方向p的函数。距离最小时对应的方向pbest即为同相轴的局部方位,也即是最佳的滤波方向。沿着该最佳滤波方向pbest从相邻的W道地震信号上选取出W个矢量采样点(每道只选取一个矢量点),然后根据(3)式和(4)式对所选取的W个矢量点进行矢量中值滤波,得到的中值矢量即为当前滤波点滤波后的值。

用上述方法逐道逐点进行处理,即可完成对整个地震剖面的滤波。给定滤波窗长度W和时窗长度N后,对每个点进行多方向矢量中值滤波处理的具体处理步骤如下:

1) 由pmin,pmax确定倾角范围,倾角增量pi与倾角范围共同决定矢量信号片段集的个数(即搜索方向的个数);

2) 沿着方向pk截取矢量信号片段集,根据公式(7),求取该方向上的矢量信号片段之间的距离;

3) 进行p方向循环,即方向由原先的pk变为pk+1,遍历完所有的方向后,由(8)式求取出最佳滤波方向pbest;

4) 沿着最佳滤波方向pbest从相邻地震道选取

W个矢量点,进行矢量中值滤波,其输出即为当前滤波点的最终取值。

2　合成数据测试

为比较多分量多方向矢量中值滤波与单分量多方向矢量中值滤波及低通滤波,本文利用合成地震数据进行测试。图2a和图2b分别为添加了均值为0,标准差为0.5的高斯噪声的Z分量和X分量数据。多分量多方向矢量中值滤波与单分量多方向矢量中值滤波的滤波参数取值为:滤波窗长度W=7,时窗长度N=7,倾角范围为-4～4(采样点/道),倾角增量为0.05(采样点/道)。低通滤波的截频为75Hz。图2c和图2d分别为经低通滤波处理的Z分量和X分量数据,由图可见,经低通滤波处理后,地震剖面上仍存在低频干扰,其中弱同相轴处较为严重。图3a和图3b分别为经多分量MD-VMF处理的Z分量和X分量,图3c和图3d分别为经单分量MD-VMF处理的Z分量和X分量。由图3a到图3d可见,两种方法都有效压制了随机噪声,但其中多分量MD-VMF的去噪效果较单分量MD-VMF好些。如图3a到图3d中黑色框所示,在地震同相轴相交处和弱同相轴处,多分量MD-VMF较单分量MD-VMF好些;且不同波型在Z方向和X方向上的投影保持得较好。图3e和图3f分别为多分量MD-VMF压制的Z分量和X分量的噪声,图3g和图3h分别为单分量MD-VMF压制的Z分量和X分量的噪声。为进一步展示多分量MD-VMF的应用效果,对部分多分量地震数据片段进行矢端图分析,如图4所示。图4a 和图4b中蓝色虚线表示合成地震数据,红色实线表示含噪地震数据;图4c和图4d为单分量MD-VMF滤波后的矢端图;图4e和图4f为多分量MD-VMF滤波后的矢端图。由图4可知,与单分量MD-VMF相比较,多分量MD-VMF能够更好地保持矢量信号的矢量特性。

实验表明:当滤波窗长度W较小时,能够较好地保护信号的细节信息,但不能有效地滤除随机噪声;而当滤波窗长度W较大时,能更好地抑制噪声,但却会严重损失有效信息。因此,需根据噪声性质调节滤波窗口的长度,从而达到滤除随机噪声的同时能够更好地保护有效信号以及减少信号泄漏。

图2　合成地震数据Z分量和X分量低通滤波处理效果分析a 添加高斯噪声的Z分量; b 添加高斯噪声的X分量; c 低通滤波处理后的Z分量; d 低通滤波处理后的X分量

图3　合成地震数据Z分量和X分量的多分量MD-VMF和单分量MD-VMF处理效果分析及其压制的噪声a 经多分量MD-VMF处理后的Z分量; b 经多分量MD-VMF处理后的X分量; c 经单分量MD-VMF处理后的Z分量; d 经单分量MD-VMF处理后的X分量; e 多分量MD-VMF压制的Z分量噪声; f 多分量MD-VMF压制的X分量噪声; g 单分量MD-VMF压制的Z分量噪声; h 单分量MD-VMF压制的X分量噪声

图4　合成地震数据Z分量和X分量的矢端图示a,b 含噪地震数据; c,d 经单分量MD-VMF滤波后的数据; e,f 经多分量MD-VMF滤波后的数据

3　应用实例

为了验证多方向矢量中值滤波方法在多分量地震数据处理中的有效性,将其应用于叠前多分量地震数据随机噪声压制处理。图5a和图5b分别为实际地震数据的Z分量与X分量,由图可见,原始资料信噪比较低,存在随机干扰,有效波的连续性受到干扰。采用本文方法对其进行滤波,滤波参数的取值为:滤波窗长度W=7,时窗长度N=7,倾角范围为-5～5(采样点/道),倾角增量为0.05(采样点/道)。图5c和图5d分别为多分量MD-VMF处理后的Z分量和X分量;图5e和图5f 分别为多分量MD-VMF压制的Z分量和X分量的噪声。由图5c和图5d可见,经多分量MD-VMF处理后的资料信噪比较高,随机干扰得到了有效压制,有效波的连续性得到了提高。

为了更好地比较多分量MD-VMF与单分量MD-VMF的滤波效果,图6a和图6b分别给出了实际地震数据X分量的多分量MD-VMF与单分量MD-VMF部分滤波结果。由于信噪比较低,将X分量单独处理时,很难有效区分噪声与有效信号。因此,经单分量MD-VMF处理的X分量地震剖面上短同相轴的方向不准确,经多分量MD-VMF处理后的X分量地震剖面上此种现象有所缓解。

图5　实际地震数据Z分量和X分量经多分量MD-VMF滤波后的结果及其压制的噪声a 实际地震数据的Z分量; b 实际地震数据的X分量; c 经多分量MD-VMF滤波后的Z分量; d 经多分量MD-VMF滤波后的X分量; e 多分量MD-VMF压制的Z分量噪声; f 多分量MD-VMF压制的X分量噪声

图6　实际地震数据X分量的部分滤波结果a 多分量MD-VMF滤波结果; b 单分量MD-VMF滤波结果

4　结论

中值滤波在实际地震资料处理中的应用已较为成熟,但常规的中值滤波技术只应用于单分量地震资料处理,并未涉及多分量地震资料。本文将多方向矢量中值滤波应用于多分量地震资料处理,既克服了常规中值滤波单方向、单分量的缺点,又能够很好地保持多分量地震数据所包含的丰富的运动学与动力学信息。实际数据处理结果表明,该方法在有效压制噪声的同时,对同相轴起到保护作用,并且不破坏多分量地震数据的矢量波场,从而能够有效地保护矢量信号的矢量特性。该方法特别适合于快、慢横波旋转分离之前的地震数据,由于介质各向异性原因使得单独的R分量和T分量上某些同相轴呈现出时弱时强的变化,而从两分量合成的矢量空间看这些同相轴则是均衡的。

[1]MARFURT K J.Robust estimates of 3D reflector dip and azimuth[J].Geophysics,2006,71 (4):29-40

[2]SACCHI M D,PORSANI M.Fast high resolution parabolic Randon tranform[J].Expanded Abstracts of 69thAnnual Internat SEG Mtg,1999:1477-1480

[3]SACCHI M D,VERSCHUUR D J,ZWARTJES P M.Data reconstruction by generalized deconvolution[J].Expanded Abstracts of 74thAnnual Internat SEG Mtg,2004:1989-1992

[4]LUO Y,DOSSARY S M,ALFARAJ M N.Edge-preserving smoothing and applications[J].The Leading Edge,2002,21(2):136-141

[5]SCHUSSTER G T.Fermat's interometric principle for target-oriented traveltime tomography[J].Geophysics,2005,70(4):47-50

[6]BIONDI B,KOSTOV C.High-resolution velocity spectra using eigenstructure methods[J].Expanded Abstracts of 68thAnnual Internat SEG Mtg,1998:890-893

[7]彭才,常智,朱仕军.基于曲波变换的地震数据去噪方法[J].石油物探,2008,47 (5):461-464

PENG C,CHANG Z,ZHU S J.Noise elimination method based on curvelet transform[J].Geophysical Prospecting for Petroleum,2008,47(5):461-464

[8]张之涵,孙成禹,姚永强,等.三维曲波变换在地震资料去噪处理中的应用研究[J].石油物探,2014,53 (4):421-430

ZHANG Z H,SUN C Y,YAO Y Q,et al.Research on the application of 3D curvelet transform to seismic data denoising[J].Geophysical Prospecting for Petroleum,2014,53(4):421-430

[9]ASTOLA J,HAAVISTO P H,NEUVO Y.Vector median filters[J].Proceedings of the IEEE,1990,78:678-689

[10]LIU Y K,LUO Y,WANG Y,et al.Vector median filter and its applications in geophysics[J].Expanded Abstracts of 79thAnnual Internat SEG Mtg,2009:3342-3346

[11]HUO S D,LUO Y,KELAMIS P,et al.Simultaneous sources separation via multi-directional vector-median filter[J].Expanded Abstracts of 79thAnnual Internat SEG Mtg,2009:31-35

[12]NAGHIZADEH M A S.Multicomponent f-x seismic random noise attenuation via vector autoregressive operators[J].Geophysics,2012,77 (2):91-99

[13]RODRIGUEZ I V,BONAR D,SACCHI M.Microseismic data denoising using a 3C group sparsity constrained time-frequency transform[J].Geophysics,2012,77 (2):21-29

[14]WANG C,WANG Y,XUN C.Multicomponent seismic noise attenuation with multivariate order statistic filters[J].Journal of Applied Geophysics,2016,133:70-81

[15]LIU Y.Noise reduction by vector median filtering[J].Geophysics,2013,78 (3):79-86

[16]HUO S D,LUO Y,KELAMIS P G.Simultaneous sources separation via multidirectional vector-median filtering[J].Geophysics,2012,77 (4):123-131

[17]曹务祥,孙哲,刘丽丽.方向中值滤波去噪方法[C].SPG/SEG深圳2011国际地球物理会议论文集.深圳:中国石油学会物探专业委员会&美国勘探地球物理家学会,2011:996-1000

CAO W X,SUN Z,LIU L L.Direction Median filtering denoising method,[C].SPG/SEG Shenzhen 2011 International Geophysical Conference Technical Program Expanded Abstracts.Shenzhen:SPG & SEG,2011:996-1000

(编辑：陈杰)

The application of multi-directional vector median filtering in multi-component seismic data

XUN Chao1,2,WANG Chao1,WANG Yun3

(1.InstituteofGeochemistry,ChineseAcademyofSciences,Guiyang550081,China;2.UniversityofChineseAcademyofSciences,Beijing100049,China;3.SchoolofGeophysicsandInformationTechnology,ChinaUniversityofGeosciences,Beijing100083,China)

The multi-directional vector median filtering (MD-VMF) can attenuate random noise of multi-component seismic data,and effectively protect the vector characteristics of multi-component seismic wavefields simultaneously.We deployed a vector-median filtering (VMF),along lines with different dips to relax the flattening requirement,to suppress noise contained in seismic data represented by multi-dimensional,multi-component vector fields.Although the VMF can be applied to each component of a vector field individually,the MD-VMF is applied to all components simultaneously.The MD-VMF treats multi-component seismic data as a vector wavefield rather than a set of scalar fields and every time samples is a vector.According to the characteristics of seismic signals,we implement the MD-VMF along local events.Then,the optimal local events are selected according to the cross-correlation between the vector signals,which are windowed from neighboring seismic traces with a sliding time window along each trial direction.The MD-VMF can attenuate noise of multi-component seismic data while preserving the relative amplitude information effectively than the single-channel filtering.We demonstrated the effectiveness of the MD-VMF on simulated data and real multi-component seismic data.

vector median filtering,multi-direction,multi-component seismic data,noise attenuation

2016-03-12;改回日期：2016-05-05。

寻超(1990—),女,硕士在读,主要从事多分量地震数据去噪方法研究。

汪超(1982—),男,博士,助理研究员,主要从事地球物理矢量场数字信号处理研究。

国家自然科学基金项目(41425017,41504107)和贵州省科技基金项目(黔科合J字[2014]2165号)联合资助。

P631

A

1000-1441(2016)05-0703-08

10.3969/j.issn.1000-1441.2016.05.009

This research is financially supported by the National Natural Science Foundation of China (Grant Nos.41425017,41504107) and the Guizhou Province Science and Technology Fund (Grant No.J[2014]2165).