地球物理学报  2018, Vol. 61 Issue (3): 1157-1168   PDF    
引用本文
魏亚杰, 韩立国, 单刚义, 等. 2018. 基于脉冲检测的混合震源数据分离. 地球物理学报, 61(3): 1157-1168, doi: 10.6038/cjg2018L0204.  
Wei Y J, Han L G, Shan G Y, et al. 2018. Separation of mixed source seismic data based on the impulse detection. Chinese J. Geophys. (in Chinese), 61(3): 1157-1168, doi: 10.6038/cjg2018L0204.  
基于脉冲检测的混合震源数据分离
魏亚杰, 韩立国 , 单刚义, 朱鹤文     
吉林大学地球探测科学与技术学院, 长春 130026
收稿日期 2017-04-05, 2018-01-26 收修定稿
基金项目: 国家高技术研究发展计划（863）重大项目课题（2014AA06A605，2014AA06A604）资助
第一作者简介: 魏亚杰, 男, 1990年生, 博士研究生, 从事混采数据分离处理. E-mail: weiyj_jlu@foxmail.com
通讯作者: 韩立国, 男, 1961年生, 教授, 主要从事复杂介质波场模拟与成像方法研究. E-mail: hanliguo@jlu.edu.cn
摘要：混合震源采集技术相对于传统地震数据采集具有改善成像质量、提高采集效率的优势.减小混合炮中单炮之间的随机延时范围能够有效的提高采集效率，但这也给之后的混采数据分离带来了影响.混采数据经伪分离后非共炮域数据中的混叠噪声明显更加集中，不利于对混叠噪声进行压制.本文提出基于脉冲检测方法对混采数据进行分离，并且与迭代的多级中值滤波方法作对比，时间延时范围较大时，两种方法都能得到很好的分离结果；时间延时范围较小时，本文方法能更有效的去除混叠噪声，同时也能更好的保留细节信息.实际数据计算结果表明，本文方法一定程度上还能够有效压制其他随机噪声.
关键词: 混合炮分离      脉冲检测      延迟时间范围     
Separation of mixed source seismic data based on the impulse detection
WEI YaJie, HAN LiGuo, SHAN GangYi, ZHU HeWen     
College of Geo-exploration Sciences and Technoloy, Jilin University, Changchun 130026, China
Abstract: Compared with the traditional seismic data acquisition, the mixed source acquisition technology has the advantages of improving the image quality and acquisition efficiency. Reducing the random delay range between single sources in the mixed sources can effectively enhance the efficiency of the acquisition, but it also has a negative impact on the separation of the mixed data. The blending noise in the non-common-source domain data is significantly more concentrated after pseudo-separation, which is difficult to suppress when the random delay time range is small. In this paper, we propose a method to separate the simultaneous source seismic data based on the pulse detection method, and compare it with the iterative multi-level median filtering method. When the time delay range is large, the two methods can both get a good separation result. When the time delay range is small, the method presented in this paper can be more effective, and retain more details of information. The actual data calculation shows that this method can also suppress other random noise to a certain extent.
Key words: Separation of mixed shots    Impulse detection    Range of delay time    
0 引言

混合震源采集技术相对于传统地震数据采集具有改善成像质量、提高采集效率的优势.混合震源采集是不同位置的多个震源以一定的编码方式同时或者延时激发，从而获得相互干涉的混合炮记录，我们把待研究的炮点激发的信号称为有效信号，接收到其他炮激发的信号称为混叠噪声，这些混叠噪声降低了地震记录的信噪比，从而对成像质量有很大的影响.混合震源数据分离就是对混合炮记录中的混叠噪声进行压制，分离出单炮记录的过程.

混合震源采集最早是由Berkhout(2008)提出的，他对混合震源项进行随机线性编码，使得不同震源按照一定的随机延时进行激发；Ikelle(2007)用相位编码将常规的单炮记录压缩为多源记录，并且探讨了混合震源的采集、模拟和处理技术；佘德平(2012)对多源混合炮记录数值模拟进行了研究.

混合震源采集的关键在于压制混叠噪声，分离出高质量的单炮记录.滤波方法是一种压制混叠噪声的有效方法，经过伪分离处理后混叠噪声只有在共炮域是相关的，在非共炮域(共检波点域、共偏移距域等)以脉冲形式随机分布(Hampson et al., 2008).根据这一理论，Huo等(2009)利用矢量中值滤波在共中心点域去除混叠噪声；Mahdad等(2010)进一步设计了迭代算法，在共检波点域压制混叠噪声；韩立国等(2013)提出了基于迭代去噪的混合数据分离技术，实现了高混合度混合震源数据分离；刘强等(2014)将混采数据分离中插值与去噪在同一个框架下同步处理，提高了分离效率.

上述方法在处理混合震源数据时，并没有考虑到随机延迟时间范围，为了提高采集效率，当随机延迟时间范围较小时，混叠噪声在非共炮域并不能随机分布于整个地震记录中，而是集中到了一起，上述方法并不能很好的去除混叠噪声.在随机延迟时间范围较小的情况下，本文使用基于脉冲检测方法对混叠噪声进行压制，同时能够去除其他随机噪声.

1 方法原理 1.1 混合震源数据伪分离

两个或者两个以上单炮混合在共炮域可以表示为

(1)

(2)

其中Ubl表示混合震源记录，U表示待分离的单炮记录，Γbl表示混合震源算子.但是在实际计算中Γbl是一个欠定矩阵，Γbl-1是不存在的，故没有办法直接求出，但是可以通过式(3)构造出伪逆(Mahdad et al., 2012)，公式为

(3)

联合方程(2)和(3)我们可以得到：

(4)

其中〈U〉表示伪分离地震记录.

1.2 三边滤波器(the trilateral filter)

双边滤波(the bilateral filter)是一种非线性滤波方法，能够去除图像上的随机噪声，同时保持图像的边缘细节信息(Tomasi and Manduchi, 1998).三边滤波(the trilateral filter)是在双边滤波原有的‘距离’权重和‘灰度’权重的基础上加入了一个‘脉冲’权重，使得三边滤波能够保持双边滤波特性的同时还能够去除脉冲噪声(Garnett et al., 2005).

假设u=u(x)是被噪声污染的图像矩阵，其中x=(x1, x2)是当前点的位置，令Ω=Ωx(N)是以x点为中心N为半滤波器宽的邻域.令y∈Ω，那么双边滤波修复后的图像点可以表示为

(5)

其中双边滤波权函数为

(6)

(7)

(8)

其中ωS为距离权重，表示邻域Ω内距离所研究点x越远的点对贡献越小；ωR为灰度权重，表示邻域Ω内与所研究点灰度值u(x)相差越大的点对贡献越小.可以看到当σS、σR→∞时，双边滤波无限接近于均值滤波.

三边滤波是基于ROAD(Rank Ordered Absolute Difference)对双边滤波的改进(Garnett et al., 2005)，对于所研究的点x的ROAD值的计算公式为

(9)

其中rk(x)是y∈Ω，|uy-ux|的值按照从小到大排序的第k个值，本文半滤波器宽N的值一般取1，m的值我们一般取4.那么我们可以用式(10)定义脉冲权重ωI，公式为：

(10)

方程(10)可以看到当x为脉冲噪声点时，具有较大的ROAD值，那么脉冲权重ωI(x)相对取较小的值.

三边滤波权函数可以式子表示为

(11)

其中：

(12)

J(x, y)在方程(11)中起到一个开关的作用，当x点和y点都不为脉冲噪声点时，ROAD(x)和ROAD(y)值都比较小，J(x, y)≈0三边滤波的权函数简化为双边滤波权函数形式；当x点和y点至少有一个脉冲噪声点时，ROAD(x)+ROAD(y)值远远大于0，那么J(x, y)≈1，脉冲权重ωI对图像平滑起主要作用.将方程(12)代入方程(5)就得到三边滤波修复后的图像点.实验结果表明σI、σJ取值较小时，对脉冲噪声滤波效果较好，但是对图像细节保持较差.

实际混合震源记录中除了含有混叠噪声，还有其他形式的随机噪声，使用三边滤波对混合震源数据做分离更符合实际，能够在去除数据中混叠噪声和随机噪声的同时还能够极大的保持图像细节信息.

1.3 基于脉冲检测的混合数据分离

基于脉冲检测的混合震源数据分离流程图如图 1所示，大致分为以下几个步骤：

图 1 基于脉冲检测的混合震源数据分离流程图 Fig. 1 Flow chart of mixed source seismic data deblending based on the impulse detection

(a) 对实际混合采集地震数据(一般为两炮混合)进行预处理得到实际混合震源地震数据U；

(b) 首先假设初始分离单炮炮集记录Pi=0(i=0)；

(c) 在频率域通过混合震源编码Γbl将分离单炮记录混合成混合震源记录Pbl；

(d) 计算实际采集地震数据U与分离单炮记录合成的混合震源记录Pbl之间的差值，即剩余混合信号Udif=U-Pbl；

(e) 在频率域通过混合震源编码Γbl对Udif做伪分离，得到伪分离后的单炮地震数据P'；

(f) 将伪分离后的单炮地震数据由共炮域转换到共偏移距域，得到伪分离后的共偏移距域地震记录，此时共偏移距域内有效信号呈连续分布，混叠噪声以随机脉冲形式分布；

(g) 使用三边滤波方法去除共偏移距域地震记录中的混叠噪声；

(h) 将三边滤波后的共偏移距道集地震记录转换回共炮点域，得到剩余混合信号的分离结果Pdif，那么原始混合震源数据U的分离结果Pi+1=Pdif+Pi；

(i) 取i=i+1，逐渐增大σI、σJ取值，重复过程b~h直到最终得到满足精度要求的结果.

2 模拟试算

本文选取了某海域一套实测单炮数据(含有100个震源、100个检波器)对本文方法进行验证，使用前文所述的方法对该数据混合成混合炮数据，混合观测系统参数：每个混合炮中包含两个单震源，第1个震源和第51个震源组成第1个混合震源，第2个震源和第52个震源组成第二个混合震源，以此类推，第50个震源和第100个震源组成第50个混合震源.检波器采样时长2 s，采样间隔0.002 s.

2.1 不同延迟时间范围下的混采数据分离效果

我们分别选用不同的随机延迟时间范围，将迭代多级中值滤波方法与本文方法做对比.

(1) 当第二炮相对第一炮延迟时间为±0.5 s时, 即延迟时间范围为1 s，随机延迟时间如图 2所示.

图 2 随机延时范围较大时，混合震源中两炮的相对延迟时间 Fig. 2 Relative delay time of two shots in mixed sources when the random delay time range is large

图 3和图 4为多级中值滤波和本文方法在随机延时范围较大时分离结果的对比.图 3a为伪分离以后的单炮记录(以第30炮为例)，其中红色箭头指向的是有效信号，蓝色箭头指向的是混叠噪声；图 3b为伪分离数据零偏移距地震剖面；图 3c为多级中值滤波的分离结果；图 3d为本文方法的分离结果；图 4a为多级中值滤波分离后的零偏移距剖面；图 4b为多级中值滤波分离前后零偏移距差剖面，即图 3b与4a的差；图 4c为本文方法分离后的零偏移距剖面；图 4d为本文方法分离前后零偏移距差剖面，即图 3b与4c的差.对比图 3c与3d、图 4a与4c可以看出随机延时范围较大时，两种方法对混合震源数据分离都能得到满意的结果.

图 3 随机延时范围较大时，多级中值滤波与本文方法分离结果对比 Fig. 3 Comparison between separation results of multilevel median filter and the method proposed in this paper when the random delay time range is large
图 4 (a) 多级中值滤波分离后的零偏移距剖面；(b) 图 3b与图 4a的差；(c)本文方法分离后的零偏移距剖面；(d) 图 3b与图 4c的差 Fig. 4 (a) Zero-offset profile based on multilevel median filter; (b) Difference between Fig. 3b and 4a; (c) Zero-offset profile based on the method proposed in this paper; (d) Difference between Fig. 3b and 4c

抽取图 3c和3d中的第30道与第70道分别与原始未混合的数据做对比，如图 5所示(为了显示清晰，本文截取的时窗为400~1000 ms).

图 5 理论数据(黑线)多级中值滤波结果(黄线)和本文方法结果(红线)的波形对比图 (a)第30道数据；(b)第70道数据. Fig. 5 Waveform comparison of theoretical data (black line) based on multilevel median filter (yellow line) and resultant data based on the method proposed in this paper (red line) (a) Data of 30th trace; (b) Data of 70th trace.

(2) 当第二炮相对第一炮延迟时间范围为±0.1 s时, 即随机延迟时间范围为0.2 s，随机延迟时间如图 6所示.

图 6 随机延时范围较小时，混合震源中两炮的相对延迟时间 Fig. 6 Relative delay time of two sources in mixed sources when the random delay time range is small

图 7和图 8为多级中值滤波和本文方法在随机延迟时间范围较小时分离结果的对比.图 7a为伪分离以后的单炮记录(以第30炮为例)，其中红色箭头指向的是有效信号，蓝色箭头指向的是混叠噪声；图 7b为伪分离数据零偏移距地震剖面；图 7c为多级中值滤波的分离结果；图 7d为本文方法的分离结果；图 8a为多级中值滤波分离后的零偏移距剖面；图 8b为多级中值滤波分离前后零偏移距差剖面，即图 7b与8a的差；图 8c为本文方法分离后的零偏移距剖面；图 8d为本文方法分离前后零偏移距差剖面，即图 7b与图 8c的差.对比图 7c与7d，多级中值滤波分离结果图 7c中可以看到明显的残留混叠噪声；对比图 8b与8d，多级中值滤波分离前后零偏移距差剖面图 8b中可以明显的看到更多有效信号被去除.说明在随机延时范围较小时，本文方法能更有效的去除混叠噪声，同时也能更好的保留细节信息.

图 7 随机延时范围较小时，多级中值滤波与本文方法分离结果对比 Fig. 7 Comparison between separation results of multilevel median filter and the method proposed in this paper when the random delay time range is small
图 8 (a) 多级中值滤波分离后的零偏移距剖面；(b) 图 7b与图 8a的差；(c)本文方法分离后的零偏移距剖面；(d) 图 7b与图 8c的差 Fig. 8 (a) Zero-offset profile based on multilevel median filter; (b) Difference between Fig. 7b and 8a; (c) Zero-offset profile based on the method proposed in this paper; (d) Difference between Fig. 7b and 8c

抽取图 7c和7d中的第30道与第70道分别与理论数据做对比，如图 9所示(为了显示清晰，本文截取的时窗为400~1000 ms).

图 9 理论数据(黑线)多级中值滤波结果(黄线)和本文方法结果(红线)的波形对比图 (a)第30道数据; (b)第70道数据. Fig. 9 Waveform comparison of theoretical data (black line) resultant data based on multilevel median filter (yellow line) and resultant data based on the method proposed in this paper (red line) (a)Data of 30th trace; (b) Data of 70th trace.
2.2 含高斯噪声的混采数据分离效果

我们通过在混合数据中加入信噪比为30%的高斯噪声来检验本文分离方法对非脉冲型噪声的压制效果，图 10a为第30炮的伪分离数据，图 10b为第30个单炮的分离结果，可以看到高斯噪声在分离结果中得到了有效的压制.

图 10 (a) 含信噪比为30%高斯噪声的伪分离数据；(b)本文方法的分离结果 Fig. 10 (a) Pseudo-separation data with Gaussian noise (SNR=30%); (b) Separation result based on the method proposed in this paper

抽取图 10b中的第30道与第70道分别与理论数据做对比，如图 11所示(为了显示清晰，本文截取的时窗为400~1000 ms).

图 11 理论数据(黑线)和本文方法结果(红线)的波形对比 (a)第30道数据; (b)第70道数据. Fig. 11 Waveform comparison of theoretical data (black line) and resultant data based on the method proposed in this paper (red line) (a) Data of the 30th trace; (b) Data of 70th trace.
3 海上实际数据处理

海上实际混合数据采集方式为全排列采集，检波线长度为10 km，相邻两个检波器距离为25 m，共400个检波器，炮线长8 km，相邻两炮点水平间距25 m，共320个炮点，两条震源船分布于检波线两侧，距离检波线垂直距离分别为100 m、900 m，两条震源船起始位置分别为第1个炮点和第161个炮点，同时向320炮点方向航行.检波器采样时长为6 s，采样间隔为0.002 s，随机延时范围为0.5 s.

图 12a为混合震源数据，图 12b和12c为使用本文方法分离出的单炮数据，可以看到本文方法不但能够压制混叠噪声，还能够有效压制其他随机噪声.图 13a为实际混合数据伪分离后零偏移距剖面，图 13b为本文方法分离后的零偏移距剖面，进一步验证了本文方法对混叠噪声的压制效果.

图 12 本文方法对实际数据分离结果 (a)原始混合数据(第50炮)；(b)第50个单炮分离结果；(c)第210个单炮分离结果. Fig. 12 Separation results of actual data by using the method proposed in this paper (a) Raw data of 50th shot; (b) Separation result of 50th shot; (c) Separation result of 210th shot.
图 13 零偏移距剖面对比 (a)伪分离数据零偏移距剖面; (b)本文方法分离结果零偏移距剖面. Fig. 13 Comparing of zero-offset profiles (a) Zero-offset profile of pseudo-separation data; (b) Zero-offset profile of separation data by using the method proposed in this paper.
4 结论

本文通过减小混合震源之间的激发时间来进一步提高采集效率.这种情况下，混叠噪声在非共炮域分布更加集中，传统的迭代多级中值滤波方法不能得到很好的处理结果.本文提出了三边滤波的方法去除混采数据中的混叠噪声，同时分离过程中还能一定程度上压制其他的随机噪声.通过迭代多级中值滤波方法与本文方法做对比，从模拟数据分离结果、零偏移距剖面以及单道数据对比都表明本文方法能够最大限度去除混叠噪声的情况下更好的保留有效信号.实际数据的计算结果更表明了该方法的实用性.

对于三维混采数据，由于其混合度高、数据量巨大，应用本文方法还有诸多限制，从二维混采数据向三维混采数据扩展也是接下来主要研究的内容.

参考文献
Berkhout A J. 2008. Changing the mindset in seismic data acquisition. The Leading Edge, 27(2): 924-938.
Doulgeris P, Verschuur D J, Blacquiere G. 2012. Separation of blended data by sparse inversion utilizing surface-related multiples. //74th EAGE Conference and Exhibition. SPE, EAGE.
Garnett R, Huegerich T, Chui C. 2005. A universal noise removal algorithm with an impulse detector. IEEE Transactions on Image Processing, 14(11): 1747-1754. DOI:10.1109/TIP.2005.857261
Hampson G, Stefani J, Herkenhoff F. 2008. Acquisition using simultaneous sources. The Leading Edge, 27(7): 918-923. DOI:10.1190/1.2954034
Han L G, Tan C Q, Lü Q T, et al. 2013. Separation of multi-source blended seismic acquisition data by iterative denoising. Chinese Journal of Geophysics, 56(7): 2402-2412.
Huo S D, Luo Y, Kelamis P. 2009. Simultaneous sources separation via multi-directional vector-median filter. //2009 SEG Annual Meeting. Houston, Texas: SEG.
Ikelle L T. 2007. Coding and decoding: seismic data modeling, acquisition and processing: US 20070274155 A1.
Liu Q, Han L G, Li H J. 2014. Synchronous interpolation and denoising in simultaneous-source data separation. Chinese Journal of Geophysics, 57(5): 1647-1654. DOI:10.6038/cjg20140527
Mahdad A, Doulgeris P, Blacquière G. 2010. Separation of blended data by iterative estimation and subtraction of blending interference noise. Geophysics, 76(3): Q9-Q17.
Mahdad A, Doulgeris P, Blacquière G. 2012. Iterative method for the separation of blended seismic data:discussion on the algorithmic aspects. Geophysical Prospecting, 60(4): 782-801. DOI:10.1111/gpr.2012.60.issue-4
She D P. 2012. A numerical modeling of simultaneous sources. Progress in Geophysics, 27(4): 1533-1540.
Tomasi C, Manduchi R. 1998. Bilateral filtering for gray and color images. //6th International Conference on Computer Vision. Bombay: IEEE, 839-846.
韩立国, 谭尘青, 吕庆田, 等. 2013. 基于迭代去噪的多源地震混合采集数据分离. 地球物理学报, 56(7): 2402–2412.
刘强, 韩立国, 李洪建. 2014. 混采数据分离中插值与去噪的同步处理. 地球物理学报, 57(5): 1647–1654.
佘德平. 2012. 多震源地震正演数值模拟技术. 地球物理学进展, 27(4): 1533–1540. DOI:10.6038/j.issn.1004-2903.2012.04.029