2019, Vol. 62
2. 中国科学院大学, 北京 100049;
3. 中国科学院地球科学研究院, 北京 100029;
4. 中国石油集团东方地球物理勘探有限责任公司, 河北涿州 072750
2. University of Chinese Academy of Sciences, Beijing 100049, China;
3. Institutions of Earth Science, Chinese Academy of Sciences, Beijing 100029, China;
4. Research and Develop center of BGP, CNPC, Zhuozhou Hebei 072750, China
多次波问题始终是地震数据处理的重点和难点之一.尽管多次波成像近年来一直是研究的热点(Liu et al., 2011, 2016, 2017; Li et al., 2017a, 2017b)，然而，对于常规的地震勘探，多次波压制仍然是地震数据处理的重要环节.对于海上地震勘探，由于水层和其上、下介质存在较大的物性差异，所以水层多次波被视为最严重的噪声.尽管在理论上，很多方法都可以压制水层多次波(Wiggins, 1988; Wang et al., 2011; 宋家文等，2014; Sun et al., 2014; 张志军和魏天罡，2013; 李列等，2015；徐强等，2015; Li and Li, 2016; 孙维蔷和王华忠，2016)，但在处理实际数据时，仍然会遇到一些问题，比如在浅水情况下，由于记录的数据中缺少近偏移距信息，浅水多次波的压制面临着严峻的挑战.
浅水多次波压制方法可以被分为三类：滤波类方法、基于波动方程的模型域方法以及基于波动方程的数据域方法.滤波类方法根据一次波和多次波的差异性，用滤波的方式直接压制多次波(Lokshtanov, 1995, 1999; Herrmann et al., 2000).预测反褶积就属于这一类方法(Weglein, 1999; Li et al., 2012)，在工业上，它被广泛用来压制浅水多次波.然而，在压制多次波的同时，预测反褶积也会压制跟水层周期接近的一次反射波.
基于波动方程的模型域方法根据已知的水层模型，利用波场传播原理，先预测出水层相关多次波，再将其从数据中去除(Wiggins, 1988; Pica et al., 2005; Spadavecchia et al., 2013, 彭海龙等，2016).水层模型驱动的去多次波方法(Model-based Water-layer Demultiple, 简称MWD)就是典型的模型域方法(Wang et al., 2011).MWD方法先模拟出水底一次反射波的Green函数，再将其与输入数据褶积产生水层相关多次波.水层模型精确时，该方法可以预测较准确的水层多次波.然而，在水很浅的情况下，精确的水层模型通常不容易获得(Yang and Hung, 2012)，多次波模型的精度受到限制.
基于波动方程的数据域方法在工业上被广泛应用于多次波压制.该方法可以被分为两类：第一类是基于波形反演直接估计一次波方法(Wang, 2004; van Groenestijn and Verschuur, 2009; Lin and Herrmann, 2013; Song et al., 2013; Lopez and Verschuur, 2015; 白兰淑等，2017)；第二类是基于预测-相减的多次波衰减方法(Verschuur et al., 1992; Berkhout and Verschuur, 1997; 李鹏等，2007；Liu et al., 2010; 武银婷等，2012；董烈乾等，2013).第一类方法将多次波的压制视为一个直接估计一次波的全波形反演过程.该类方法不但精度高，而且当数据中缺少近偏移距信息时，还可以重建出缺失的信息，适合用于浅水多次波压制.然而，这类方法通常需要较大的计算量和存储量.第二类方法首先预测出多次波，然后通过自适应相减将预测的多次波从数据中去除.在这一类方法中，最具有代表性的是地表相关多次波压制方法(Surface-related multiple elimination，简称SRME) (Verschuur et al., 1992; Berkhout and Verschuur, 1997).
SRME方法不需要任何地下结构信息，可以高效的压制所有类型的地表相关多次波.由于其实用性，SRME方法已得到学术界、工业界的广泛认可.然而，在浅水情况下，SRME却面临问题：(1)由于观测数据中缺少近偏移距信息，导致多次波的预测存在困难(Verschuur, 2006; Hung et al., 2010)；(2)预测的多次波模型中存在一个额外的震源子波，而且高阶多次波的振幅被过度预测，这会导致在首次迭代中，多次波模型产生严重的串扰噪声(Moore and Bisley, 2006; Wang et al., 2011).这不仅会使迭代次数大大增加，还可能会在自适应相减阶段伤害一次反射波.针对SRME方法在浅水情况下存在的问题，Biersteker (2001)提出了利用多道预测算子压制浅水多次波的方法(工业上也将其称为二维预测反褶积).经过后续学者的发展(Hugonnet and Tichatschke, 2005; Hargreaves, 2006; Hung et al., 2010; Yang and Hung, 2012)，该方法逐渐成为工业上压制浅水多次波的重要方法之一.该方法首先从输入数据中估计出含有水底信息的多道预测算子，然后将预测算子与输入数据褶积产生多次波，最后将预测的多次波从输入数据中去除.由于多道预测算子是从数据中估计得到的，所以，预测的多次波与真实的多次波的振幅和相位很接近，在自适应相减阶段，只需要用很短的滤波器就可以达到满意的效果，这有利于保护一次反射波.即使在数据缺失近偏移距信息的情况下，该方法也能估计出可靠的预测算子，这大大缓解了SRME方法在浅水情况下遇到的困难.然而，利用多道预测算子压制多次波方法仍然存在一些问题.由于多道预测算子是从数据中反演得到的，它很容易受噪声污染(Cooper et al., 2015).甚至在有些情况下，很难估计出稳定的多道预测算子(或者受噪声污染很严重)，此时，用受噪声污染的算子来预测多次波，必然会影响多次波模型的精度.
因此，本文提出了一种改进的多道预测算子压制浅水多次波方法.和Biersteker的方法(2001)一样，首先从输入数据中估计出多道预测算子.但由于算子中含有噪声，所以，根据褶积原理改进了多道预测算子，改进的预测算子具有精确的海底一次波走时信息，振幅能够和原始的多道预测算子很好的匹配，而且还不受噪声污染.
由于改进的算子只含有海底的信息，所以，将其与输入数据褶积产生的多次波分为两类：源端多次波和检波点端多次波(Biersteker, 2001; Zhai et al., 2015).然而，这两类多次波含有共同的高阶项：源端-检波点端多次波，它们被预测了两次.所以，在本文中，我们采用了一种两步相减法(Wilkinson and Bale, 2014)，该方法可以避免高阶多次波的过度预测，从而获得较精确的多次波压制结果.
文中首先回顾了利用多道预测算子压制浅水多次波方法；然后，提出了改进的多道预测算子压制浅水多次波方法；最后，用一个合成数据体和一个浅水实际数据体验证提出的方法的可靠性和实用性.
1 方法原理 1.1 利用多道预测算子压制浅水多次波在二维情况下，用矩阵P表示在地表接收到的单频地震波场，其中矩阵的每一行表示一个共炮点道集，每一列表示一个共检波点道集(Berkhout, 1982)，则对于炮检重合的观测系统，SRME方法可以表示为
|
(1) |
式中，P0表示一次波；A为地表算子(可表示为A=RS-1，其中R为地表反射率矩阵，S为震源子波).然而，在浅水情况下，由于原始数据中近偏移距的缺失，以及多次波模型中产生的串扰噪声，SRME方法存在问题(Verschuur, 2006; Hung et al., 2010; Moore and Bisley, 2006; Wang et al., 2011).
为了解决SRME在浅水环境下存在的问题，Biersteker(2001)提出了利用多道预测算子压制浅水多次波的方法.经过Hargreaves(2006)和Hung等(2010)的进一步发展，该方法已成为近年来工业界压制浅水多次波的重要方法.这种方法首先从输入数据中估计出含有海底信息的多道预测算子，再利用估计的算子预测多次波，最后将预测的多次波从输入数据中去除.在本节中，首先介绍如何从输入数据中估计多道预测算子；再介绍多次波的预测和相减过程；最后分析一下该方法存在的问题.
1.1.1 多道预测算子的估计假设地表反射率为-1，方程(1)可以重写为
|
(2) |
其中，F表示一次波Green函数，其余的符号与方程(1)中的表示一致.由于在方程(2)中存在两个未知数(P0和F)，所以需要额外的约束条件.一般而言，假设一次波P0能量最小，则F可以通过多道反褶积原理来获取(Biersteker, 2001).除此之外，为了保证反演的稳定性，F必须是一个在时间和空间方向截断的算子.这样一来，方程(2)不再严格正确，因为其等号右端多次波项不能包含完整的地表相关多次波.但此时，F却包含了来自海底的反射信息，它正好可以用来预测水层相关多次波.
根据Hung等(2010)的研究，截断的F不再是完整的一次波Green函数，所以将其称为多道预测算子.根据方程(2)，对于固定炮点s处的多道预测算子Fs，可以通过求解如下的最小二乘问题来获得：
|
(3) |
式中，r和k分别表示检波点坐标和多次波预测的孔径.P表示输入数据，L2表示2范数，λ表示规则项的权重，*表示时间域褶积.
为了保证反演的稳定性，在求解方程(3)的过程中，Fs必须是一个在时间上有延迟的算子，这个时间延迟一般略小于海底一次反射波的走时(该走时可通过一个近似的海底深度来获取)；而且，算子的长度必须能够包含地下主要反射体信息，但不能过长，否则反演将失去统计性.所以，算子在时间方向的截断是由时间延迟和算子的长度共同确定的，即算子只在该有效长度内有值，在其他样点处值为零.算子在空间方向上的截断由方程(3)中的k来确定，最大的k为待反演的预测算子的总道数，它的位置一般不超过海底一次反射波临界偏移距的位置.对于近偏移距缺失的输入数据，利用方程(3)仍然可以估计出稳定的多道预测算子.这是因为在物理上，Fs代表海底的一次反射波信息，即使在数据中缺失了，但它仍包含在与它相关的多次波路径中，而方程(3)的求解过程就是重建缺失的反射信息的过程.
方程(3)的求解是一个多道反演过程，所以，有必要对该过程的计算量加以说明.求解方程(3)时，先对未知数Fsk逐一求导，并令其导数等于0，再将它们组成一个线性方程组求解即可.方程组的系数矩阵是一个块状Toeplitz矩阵，该方程组可以用递归的莱文森方法快速求解(Levinson, 1947).所以，方程(3)求解主要的计算量在于方程组系数矩阵的构建.构建系数矩阵时，需要将参与反演的地震道Pkr两两之间做有延迟的相关(k和r都是变化的)，系数矩阵的主对角矩阵是Pkr的自相关，其余块矩阵为Pkr的互相关.所以，对于固定炮点s，如果有r个地震道参与反演，待反演的算子总道数为k，则需要做有延迟的相关的次数为rk，相关时的延迟由算子的长度确定.与SRME相比，预测r道多次波，预测孔径为k时，需要做rk次单道褶积，而对于同样的数据，做一次有延迟的互相关耗时小于一次褶积，因此，该反演过程的计算量低于1次SRME多次波预测的计算量.
1.1.2 多次波预测和相减估计的多道预测算子可以直接用来预测水层相关多次波.但值得注意的是，由于估计的算子只含有海底的信息.所以预测的多次波分为两类：源端多次波和检波点端多次波.当地下界面水平时，这两类多次波没有时差；当地下界面倾斜时，他们的走时不一致(如图 1所示).所以，这两类多次波应该分别预测，再将他们从输入数据中减去.然而，这两类多次波中存在共同的高阶项：源端-检波点端多次波，它们被预测了两次.
|
图 1 倾斜海底情况下，源端多次波(实线)和检波点端多次波(虚线)示意图 Fig. 1 A schematic of source-side (solid line) and receiver- side multiples (dash line) when facing an incline seafloor |
因此，为了避免高阶多次波的过度预测，Wilkinson和Bale(2014)提出了两步相减法.该方法首先用预测算子和输入数据做褶积预测出源端多次波，同时将其从输入数据中减去；然后再用去除了源端多次波的数据和预测算子褶积来预测检波点端多次波，再将其从数据中减去，从而获得不含多次波的数据.该方法可以表示为
|
(4) |
式中，A1和A2表示自适应滤波器；P1表示去除了源端多次波的数据，P和P0分别表示输入数据和不含多次波的数据.该方法可以较精确的预测并去除数据中的水层相关多次波.但值得注意的是，对于实际数据，必须对输入数据P进行近偏移距数据重建.
1.1.3 多道预测算子压制方法存在的问题由于多道预测算子是从输入数据中估计出来的，将其用于多次波预测时(见方程(4))，预测的多次波与输入数据有相似的振幅和相位，从而在自适应相减阶段(见方程(4))，只需要用很短的滤波器就可以达到很好的去多次波效果，这有利于保护数据中的一次反射波.
然而，也正是因为多道预测算子是从数据中估计而来这一特点，极大地限制了该方法对数据的适用性.因为通过反演获得的预测算子很容易受噪声污染(Cooper et al., 2015).原因是多方面的，例如，输入数据中的一次波和多次波、多次波和高阶多次波重合；输入数据受直达波、首波污染而被截断；输入数据中的多次波能量不够强.一旦多道预测算子中存在噪声，那么这些噪声就会被带到多次波模型中，从而影响最终的一次反射波的质量.
图 2中展示了6个炮点处的多道预测算子，它们是从本文数值算例部分实际资料数据体中估计得到的.可以看出，图 2(a—c)中所示的算子信噪比较高，而图 2(d—f)中所示的算子噪声较严重.所以，用图 2(d—f)中所示的算子预测多次波，必然会给多次波模型带来噪声.所以，为了获得高精度的多次波模型，有必要对多道预测算子进行改进.
|
图 2 从本文数值算例部分实际资料数据体中估计的6个不用炮点处的多道预测算子 (a—c)中的算子信噪比较高，而(d—f)中的算子信噪比较低 Fig. 2 The multichannel prediction operators at different source locations from a marine field dataset The operators in (a—c) contain less noise, but the operators in (d—f) contain plenty of noise |
为了解决多道预测算子含有噪声的问题，对算子进行了改进.改进的算子不仅不含噪声，而且其振幅和相位和原始的多道预测算子很接近.这保证了多次波模型的振幅和相位与输入数据能很好地匹配，还不会受到噪声影响.
构建改进的多道预测算子主要包括如下步骤：
(1) 构建精确的水层模型.首先，根据方程(3)，可以估计出所有炮点处的多道预测算子.因为这些算子物理上反映的是海底一次反射波信息，所以，把所有的算子叠加，并做常速度偏移(这里，假设水层是均匀的)，即可获得精确的水层模型.此时，就可以在偏移的剖面上拾取出准确的海底.做常速度偏移时，时间偏移或深度偏移均可，因为获得的偏移剖面最终都会被转化到深度域.
通常情况下，这种构建水层模型的方法是很有效的，因为通过方程(3)，即使在输入数据缺失近偏移距信息的情况下，也能获得比较可信的多道预测算子.而且，即使算子含有噪声，在偏移剖面上，海底的形态也能清楚的反映.
(2) 计算精确的水底一次反射波走时.获得了精确的水层模型后，在海底不是特别复杂的情况下，就可以通过下式来计算精确的水底一次反射波走时(Spadavecchia et al., 2013)：
|
(5) |
式中，Kt表示海底一次波时间响应函数；τ表示地表坐标x和y之间的来自海底的一次波走时；τxy′(或τy′y)表示从地表坐标x(或y)到海底坐标y′之间的单程走时.由于水层被假设为均匀介质，所以地震射线从x(或y)到海底y′之间为直线传播，走时很容易求得.因此，在已知海底各点深度情况下，通过方程(5)，可以求出地表任意两点间的海底一次反射波走时.
(3) 统计性估计多道预测算子的振幅信息.假设通过方程(3)估计出的多通道预测算子F可以被当做其走时响应函数和振幅响应函数褶积的结果，而其走时响应函数理论上应该等价于方程(5)中的Kt (因为Kt的走时精确计算的).所以，F算子的振幅响应函数Ka可以近似的通过统计性反演的方式求得：
|
(6) |
式中，各变量右下角的|x表示在特定的炮点坐标x下；y表示检波点坐标.方程(6)的求解很稳定，即便对于噪声较严重的多道预测算子Fs也能估计出稳定的振幅信息Ka|x(t).
(4) 合成改进的多道预测算子.利用褶积模型，改进的多道预测算子K可以表示为海底走时响应函数Kt和振幅响应函数Ka共同作用的结果：
|
(7) |
式中，α表示权系数，在本文中α=1/s0，其中s0表示x和y两点间的海底一次反射波路径的长度.所以，α可以补偿几何扩散，从而反映改进的算子道间振幅的变化.
改进的多道预测算子是根据褶积原理合成的，不受噪声的影响.然而，需要指出的是，改进的算子的合成只能适合于不太复杂的海底，这是因为在计算走时信息过程中(方程(5))，没有考虑到多至走时的情况.另外，改进的算子的振幅信息被假设为道间独立的，而没有考虑到其相位变化.这些效应给多次模型带来的误差可以在自适应相减阶段得到消除.
改进的算子的构建过程不需要太大的计算量.对于走时的计算(方程(5))，每道需要计算多个路径的走时，做比较，取最小；对于振幅信息的估计(方程(6))，每一炮的算子只需要估计出一个波形，所以该反演过程可视为求解一个单道最小二乘问题，计算量很小；对于算子的构建(方程(7))，每道只需要做一次褶积即可.因此，相对于原始的多道预测算子的估计过程(方程(3))，构建改进的算子只增加了极小的计算量.
改进的算子K构建完成后，可以直接用来预测多次波.多次波的预测和相减可以采用1.1.2节中所示的方法(见方程(4)).此时，只需要在方程(4)中，将算子F替换为改进的算子K，即可获得不含水层多次波的数据.同1.1.2节中所述，对于实际数据，在预测多次波之前需要对输入数据进行近偏移距插值.
2 数值算例 2.1 合成数据测试在这个数值算例中，用一个合成的数据体来测试我们所提出的方法的有效性.该数据体是用有限差分方法合成的.图 3是合成数据时所用的速度模型，其中海底深度为60 m.该模型包含5个反射层，从上至下其密度分别为：1.0 g·cm-3、1.8 g·cm-3、2.0 g·cm-3、2.2 g·cm-3、2.5 g·cm-3.合成的数据体一共包括201炮，每炮含有201个检波器.所有的炮点和检波点是固定的，并且位置重合.横向上，这些炮点均匀分布在从0到2000 m的位置，间隔为10 m；纵向上，它们被放置在水下10 m的位置.模拟中所有的震源子波主频均为30 Hz.
|
图 3 用于产生合成数据的速度模型 Fig. 3 The velocity model used for the generation of the synthetic dataset |
图 4a展示了炮点位于x=1000 m的单炮记录，为了符合实际情况，去除了近偏移距±100 m以内的数据.图 4b是从图 4a所示的炮集中估计出的多道预测算子.可以看出，多道预测算子噪声很严重，而且在中远偏移距，信号很弱，这不利于预测大角度入射的水层多次波.图 4c是用我们提出的方法改进的多道预测算子.可以看出，改进的算子几乎不含噪声，而且在各个偏移距都能有较稳定的波形信息.
|
图 4 (a) 一个缺失近偏移距数据的原始炮集；(b)从(a)中估计出的多道预测算子；(c)改进的多道预测算子 Fig. 4 (a) An original shot-gather with missing near offsets; (b) The multichannel prediction operators estimated from (a); (c) The improved multichannel prediction operators |
图 5比较了用多道预测算子和改进的多道预测算子预测的多次波模型.图 5a中展示了图 4a中所示的没有近偏移距缺失的炮集；图 5b和5c分别展示了利用多道预测算子(图 4b)和改进的预测算子(图 4c)预测的源端多次波.可以看出，由于原始的预测算子(图 4b)中存在噪声，所以它预测的多次波模型(图 5b)也存在噪声，见其中白色椭圆所示部分.相比之下，利用改进的预测算子预测的多次波模型几乎不含噪声(见图 5c).
|
图 5 (a) 一个完整的原始炮集；(b)用图 4b中所示算子预测的源端多次波；(c)用图 4c中所示算子预测的源端多次波 Fig. 5 (a) An original shot-gather; (b) The corresponding source-side multiple model predicted with the operators in Fig. 4b; (c) The corresponding source-side multiple model predicted with the operators in Fig. 4c |
图 6和图 7比较了利用多道预测算子方法和改进的多道预测算子方法的去多次波效果.图 6展示了它们在单炮记录上的比较结果；图 7展示了它们在零偏移距剖面上的比较结果.在图 6中，可以看出，利用原始的多道预测算子和改进的预测算子都能够很好压制水层多次波(见图 6b和6c).相比而言，改进的多道预测算子方法的压制效果更好(见图 6c)，因为在浅部，多次波压制更明显，而且噪声更少(见白色椭圆所示).
|
图 6 (a) 图 5a中所示的炮集；(b)用原始的多道预测算子方法压制多次波之后；(c)用改进的多道预测算子方法压制多次波之后 Fig. 6 (a) The original shot gather shown in Fig. 5a; (b) The demultipled shot gather resulted from the multichannel prediction operator method; (c) The demultipled shot gather resulted from the improved multichannel prediction operator method |
|
图 7 (a) 原始的零偏移距剖面；(b)用原始的多道预测算子方法压制多次波之后；(c)利用改进的多道预测算子方法压制多次波之后 Fig. 7 (a) The original zero-offset section; (b) The demultipled zero-offset section resulted from the original multichannel prediction operator method; (c) The demultipled zero-offset section resulted from the improved multichannel prediction operator method |
同样，在零偏移距剖面上，也能发现相似的现象.图 7a展示了原始的零偏移距剖面，图 7b和7c分别展示了利用原始的多道预测算子和改进的多道预测算子方法压制之后的结果.可以看出，两种方法都能很好的压制水层多次波(见图 7(b—c)中白色箭头所示)，即使在浅部，二者的结果中都有一些多次波残留.然而，与原始的多道预测算子方法压制结果(图 7b)相比，改进的多道预测算子方法的结果(图 7c)噪声更少，而且浅部更明显，如图 7(b—c)中白色椭圆所示的区域.
2.2 浅水实际数据测试用一个浅水实际数据体来展示提出的方法的有效性与实用性.这个实际数据体是通过一个海上二维地震测线采集得到的，采集时，炮点到第一个检波点距离是150 m.利用炮检互换原理，可以得到双边炮集记录.在去多次波之前，需要做一系列预处理，如压制随机噪声，对缺失的炮、道数据插值，数据规则化等.预处理之后的数据体含有321炮，每炮含有321个检波器，所有的炮、检位置重合，炮间隔是10 m.
图 8a中展示了一个原始的炮集记录，图 8b和8c分别展示了从图 8a中估计的多道预测算子以及改进的多道预测算子.可以看到，原始的预测算子(图 8b)含有较严重的噪声，而且远偏移距信息不准确；但在改进的多道预测算子中(图 8c)，几乎不含噪声，而且各个偏移距信息都比较可信.
|
图 8 (a) 一个原始的单炮实际数据；(b)从(a)中估计得到的多道预测算子；(c)改进的多道预测算子 Fig. 8 (a) An original field shot-gather; (b) The multichannel operators estimated from (a); (c) The improved multichannel operators from (b) |
图 9a展示了图 8a中的炮集经过近偏移距数据重建之后的结果.数据重建利用Kabir and Verschuur (1995)和Bai等(2017)中的方法均实现的.图 9b和9c分别展示了利用多道预测算子方法和改进的多道预测算子方法压制的效果.可以看出，两种方法都能很好的压制水层相关多次波.然而，在多道预测算子方法的压制结果中(图 9b)存在噪声；而在改进的多道算子方法压制的结果中(图 9c)，不存在明显的噪声，见图 9(b—c)中白色箭头所示.
|
图 9 (a) 图 8a中所示炮集经过近偏移距数据插值后；(b)利用原始的多道预测算子方法压制多次波之后；(c)利用改进的多道预测算子方法压制多次波之后 Fig. 9 (a) The original field shot-gather shown in Figure 8a after near-offset reconstruction; (b) The demultipled shot gather resulted from the multichannel prediction operator method; (c) The demultipled shot gather resulted from the improved multichannel prediction operator method |
图 10a展示了原始的零偏移距剖面，图 10(b—c)分别展示了利用多道预测算子方法和改进的多道预测算子方法压制后的零偏移距剖面.两种方法都表现出了比较满意的压制效果，如图 10a中箭头所示的较强的多次波同相轴，在10(b—c)中都被有效的压制了.然而，相比于多道预测算子方法的压制结果(图 10b)，在改进的多道预测算子方法的结果中(图 10c)含有更少的噪声和更少的多次波残留，见图 10(b—c)中白色椭圆所示的区域.
|
图 10 (a) 原始的零偏移距剖面；(b)利用原始的多道预测算子方法压制多次波之后；(c)利用改进的多道预测算子方法压制多次波之后 Fig. 10 (a) The original zero-offset section; (b) The demultipled section resulted from the multichannel prediction operator method; (c) the demultipled section resulted from the improved multichannel prediction operator method |
在本文中，我们提出了一种利用改进的多道预测算子压制浅水多次波的方法，该方法不需要任何地下结构信息，能够压制所有类型的水层相关多次波.具体而言，该方法首先从输入数据中统计性估计出含有海底信息的多道预测算子，并利用估计的预测算子构建水层模型；然后，通过计算海底走时信息、估计振幅信息、合成新算子三个步骤构建出不含噪声的多道预测算子；最后，利用改进的算子预测水层多次波，并结合自适应相减把这些预测的多次波从数据中去除.
和原始的多道预测算子压制方法相比，改进的方法同样能够压制所有类型的水层多次波，但几乎不受噪声的影响.在本文数值算例部分，合成数据和实际资料的例子都验证了这一结论.所以，改进的多道预测算子方法能够适应更广泛的数据体.
在改进的多道预测算子压制方法中，仍存在两个问题：(1)改进的方法不适用于海底特别复杂的情况，这是因为在计算改进的预测算子的走时的过程中(见方程(5))，没有考虑到多至走时的情况；(2)计算改进的算子振幅信息时，没有考虑道间振幅变化.这些问题有待在以后的工作中作进一步研究.改进的算子精确的走时信息和可靠的振幅信息可以保证预测的多次波模型的精度.由道间振幅不变效应给多次波模型造成的误差，可以在自适应相减阶段得到校正.
另外，改进的方法只能压制水层相关多次波，其他的地表相关多次波，以及层间多次波还需要结合别的方法共同压制，如SRME方法，逆散射级数方法(Weglein et al. 1997)等.
Bai L S, Liu Y K, Lu H Y. 2017. Strategy of multiple elimination by sparsity inversion and its effectiveness analysis. Chinese Journal of Geophysics (in Chinese), 60(12): 4801-4813. DOI:10.6038/cjg20171221 |
Berkhout A J. 1982. Seismic Migration, Imaging of Acoustic Energy by Wave Field Extrapolation. 2nd ed. Amsterdam: Elsevier.
|
Berkhout A J, Verschuur D J. 1997. Estimation of multiple scattering by iterative inversion, part Ⅰ:Theoretical considerations. Geophysics, 62(5): 1586-1595. DOI:10.1190/1.1444261 |
Biersteker J. 2001. MAGIC: Shell's surface multiple attenuation technique.//71st A Ann. Internat Mtg., Soc. Expi. Geophys.. Expanded Abstracts, 1301-1304.
|
Cooper J, Poole G, Wombell R, et al. 2015. Recursive model-based water-layer demultiple.//77th Annual International Conference and Exhibition. Madrid, EAGE, Tu N114 02.
|
Dong L Q, Li Z C, Yang S C, et al. 2013. Non-causal matching filter multiple elimination method based on correlation and iteration. Chinese Journal of Geophysics (in Chinese), 56(10): 3542-3551. DOI:10.6038/cjg20131029 |
Hargreaves N. 2006. Surface multiple attenuation in shallow water and the construction of primaries from multiples.//76th Ann.Internat Mtg., Soc. Expi. Geophys.. Expanded Abstracts, 2689-2693.
|
Herrmann P, Mojesky T, Magesan M, et al. 2000. De-aliased, high-resolution Radon transforms.//70th Ann. Internat Mtg., Soc. Expi. Geophys.. Expanded Abstracts, 1953-1957.
|
Hugonnet P, Tichatschke C. 2005. 2D Deconvolution for OBC data and for internal multiple attenuation-part 1-theory.//67th EAGE Conference & Exhibition. Madrid, EAGE, A026.
|
Hung B, Yang K L, Zhou J, et al. 2010. Shallow water demultiple.//21st Geophysical Conference. Sydney, ASEG, 1-4.
|
Kabir M M N, Verschuur D J. 1995. Restoration of missing offsets by parabolic radon transform. Geophysical Prospecting, 43(3): 347-368. DOI:10.1111/gpr.1995.43.issue-3 |
Levinson N. 1947. The wiener RMS (root mean square) error criterion in filter design and prediction. Journal of Mathematical Physics, 25: 261-278. |
Li L, Xie Y H, Li Z N, et al. 2015. Research progress on offshore multiple attenuation and imaging. Progress in Geophysics (in Chinese), 30(1): 446-453. DOI:10.6038/pg20150166 |
Li P, Liu Y K, Chang X, et al. 2007. Application of the equipoise pseudomultichannel matching filter in multiple elimination using wave-equation method. Chinese Journal of Geophysics (in Chinese), 50(6): 1844-1853. |
Li Z N, Li Z C, Wang P, et al. 2017a. Reverse time migration of multiples based on different-order multiple separation. Geophysics, 82(1): S19-S29. DOI:10.1190/geo2015-0710.1 |
Li Z N, Li Z C, Wang P, et al. 2017b. One-way wave-equation migration of multiples based on stereographic imaging condition. Geophysics, 82(6): S479-S488. DOI:10.1190/geo2016-0478.1 |
Li Z X, Li Z C, Lu W K. 2012. Multichannel predictive deconvolution based on the fast iterative shrinkage-thresholding algorithm. Geophysics, 81(1): V17-V30. |
Li Z X, Li Z C. 2016. Accelerated and robust multichannel predictive deconvolution for shallow water demultiple based on alternating direction method of multipliers.//86th Ann. Internat Mtg., Soc. Expi. Geophys.. Expanded Abstracts, 4529-4533.
|
Lin T T Y, Herrmann F J. 2013. Robust estimation of primaries by sparse inversion via one-norm minimization. Geophysics, 78(3): R133-R150. DOI:10.1190/geo2012-0097.1 |
Liu X J, Liu Y K, Hu H, et al. 2017. Imaging of first-order surface-related multiples by reverse-time migration. Geophysical Journal International, 208(2): 1077-1087. DOI:10.1093/gji/ggw437 |
Liu Y K, Jin D G, Chang X, et al. 2010. Multiple subtraction using statistically estimated inverse wavelets. Geophysics, 75(6): WB247-WB254. DOI:10.1190/1.3494082 |
Liu Y K, Chang X, Jin D G, et al. 2011. Reverse time migration of multiples for subsalt imaging. Geophysics, 76(5): WB209-WB216. DOI:10.1190/geo2010-0312.1 |
Liu Y K, Liu X J, Osen A, et al. 2016. Least-squares reverse time migration using controlled-order multiple reflections. Geophysics, 81(5): S347-S357. DOI:10.1190/geo2015-0479.1 |
Lokshtanov D. 1995. Multiple suppression by single channel and multichannel deconvolution in the tau-p domain.//65th Ann.Internat Mtg., Soc. Expi. Geophys.. Expanded Abstracts, 1482-1485.
|
Lokshtanov D. 1999. Multiple suppression by data-consistent deconvolution. The Leading Edge, 18(1): 115-119. DOI:10.1190/1.1438136 |
Lopez G A, Verschuur D J. 2015. Closed-loop surface-related multiple elimination and its application to simultaneous data reconstruction. Geophysics, 80(6): V189-V199. DOI:10.1190/geo2015-0287.1 |
Moore I, Bisley R. 2006. Multiple attenuation in shallow-water situations.//68th EAGE Conference and Exhibition. Vienna, EAGE, F018.
|
Peng H L, Ren T, He J W, et al. 2016. Application of shallow water de-multiple methods based on water-layer model. Progress inGeophysics (in Chinese), 31(4): 1646-1654. DOI:10.6038/pg20160432 |
Pica A, Poulain G, David B, et al. 2005. 3D surface-related multiple modeling, principles and results.//75th Ann. Internat Mtg., Soc. Expi. Geophys.. Expanded Abstracts, 2080-2083.
|
Song J W, Verschuur E, Chen X H. 2013. Comparing three feedback internal multiple elimination methods. Journal of Applied Geophysics, 95: 66-76. DOI:10.1016/j.jappgeo.2013.05.010 |
Song J W, Verschuur D J, Chen X H. 2014. Research status and progress in multiple elimination. Progress in Geophysics (in Chinese), 29(1): 240-247. DOI:10.6038/pg20140134 |
Spadavecchia E, Lipari E, Bienati V, et al. 2013. Water-bottom multiple attenuation by Kirchhoff extrapolation. Geophysical Prospecting, 61(4): 725-734. DOI:10.1111/gpr.2013.61.issue-4 |
Sun W Q, Hu J T, Wang H Z, et al. 2014. Enhanced model-based water-layer related demultiple.//84st Ann. Internat Mtg., Soc. Expi. Geophys.. Expanded Abstracts, 349-351.
|
Sun W Q, Wang H Z. 2016. Water-layer related multiple suppression based on plane-wave coding. Geophysical Prospecting for Petroleum (in Chinese), 55(4): 516-523. |
Van Groenestijn G J A, Verschuur D J. 2009. Estimating primaries by sparse inversion and application to near-offset data reconstruction. Geophysics, 74(3): A23-A28. DOI:10.1190/1.3111115 |
Verschuur D J, Berkhout A J, Wapenaar C P A. 1992. Adaptive surface-related multiple elimination. Geophysics, 57(9): 1166-1177. DOI:10.1190/1.1443330 |
Verschuur D J. 2006. Seismic Multiple Removal Techniques-Past, Present and Future. Housten: EAGE Publications.
|
Wang P, Hu J Z, Xu S, et al. 2011. Model-based water-layer demultiple.//81st Ann. Internat Mtg., Soc. Expi. Geophys.. Expanded Abstracts, 3551-3555.
|
Wang Y H. 2004. Multiple prediction through inversion:A fully data-driven concept for surface-related multiple attenuation. Geophysics, 69(2): 547-553. DOI:10.1190/1.1707074 |
Weglein A B, Gasparotto F A, Carvalho P M, et al. 1997. An inverse-scattering series method for attenuating multiples in seismic reflection data. Geophysics, 62(6): 1975-1989. DOI:10.1190/1.1444298 |
Weglein B. 1999. Multiple attenuation:An overview of recent advances and the road ahead (1999). The Leading Edge, 18(1): 40-44. DOI:10.1190/1.1438150 |
Wiggins J W. 1988. Attenuation of complex water-bottom multiples by wave-equation-based prediction and subtraction. Geophysics, 53(12): 1527-1539. DOI:10.1190/1.1442434 |
Wilkinson K, Bale R. 2014. Offshore east coast:Model-based water-layer demultiple breathes new life into old data. CSEG Recorder, 39(8): 46-50. |
Wu Y T, Liu Y K, Chang X, et al. 2012. An improved matching method based on L1-norm for multiple elimination and application in deepwater basin of South China Sea. Chinese Journal of Geophysics (in Chinese), 55(3): 1007-1016. DOI:10.6038/j.issn.0001-5733.2012.03.031 |
Xu Q, Wang Z, Shi Z Y, et al. 2015. Shallow water multiple attenuation. Oil Geophysical Prospecting (in Chinese), 50(2): 238-242. |
Yang K L, Hung B. 2012. Shallow water demultiple with seafloor reflection modeling using multichannel prediction operator.//82nd Ann. Internat Mtg., Soc. Expi. Geophys.. Expanded Abstracts, 1-5.
|
Zhai Y, Sen S, Liu Z, et al. 2015. Shallow water multiple elimination(SWME) on broadband data.//77th EAGE Conference and Exhibition. Madrid, EAGE, Th P4 10.
|
Zhang Z J, Wei T G. 2013. Application of joint shallow-water multiple attenuation technique in LD zone of Bohai Sea. China Petroleum Exploration (in Chinese), (1): 59-65. |
白兰淑, 刘伊克, 卢回忆. 2017. 稀疏反演多次波去除策略与效果分析. 地球物理学报, 60(12): 4801-4813. DOI:10.6038/cjg20171221 |
董烈乾, 李振春, 杨少春, 等. 2013. 基于相关迭代的非因果匹配滤波器多次波压制方法. 地球物理学报, 56(10): 3542-3551. DOI:10.6038/cjg20131029 |
李列, 谢玉洪, 李志娜, 等. 2015. 海上多次波压制与成像方法研究进展. 地球物理学进展, 30(1): 446-453. DOI:10.6038/pg20150166 |
李鹏, 刘伊克, 常旭, 等. 2007. 均衡拟多道匹配滤波法在波动方程法压制多次波中的应用. 地球物理学报, 50(6): 1844-1853. DOI:10.3321/j.issn:0001-5733.2007.06.027 |
彭海龙, 任婷, 赫建伟, 等. 2016. 基于水层建模理论的浅水多次波压制方法研究与应用. 地球物理学进展, 31(4): 1646-1654. DOI:10.6038/pg20160432 |
孙维蔷, 王华忠. 2016. 基于平面波编码的水体相关多次波压制方法研究. 石油物探, 55(4): 516-523. DOI:10.3969/j.issn.1000-1441.2016.04.006 |
宋家文, Verschuur D J, 陈小宏. 2014. 多次波压制的研究现状与进展. 地球物理学进展, 29(1): 240-247. DOI:10.6038/pg20140134 |
武银婷, 刘伊克, 常旭, 等. 2012. 改进的1范数匹配滤波法及在南海深水盆地的应用. 地球物理学报, 55(3): 1007-1016. DOI:10.6038/j.issn.0001-5733.2012.03.031 |
徐强, 王征, 史增园, 等. 2015. 浅水多次波衰减. 石油地球物理勘探, 50(2): 238-242. |
张志军, 魏天罡. 2013. 浅水多次波的联合衰减技术在渤海海域LD地区的应用. 中国石油勘探, (1): 59-65. DOI:10.3969/j.issn.1672-7703.2013.01.010 |