利用最小平方基准面延拓压制自由表面多次波II：海洋拖缆观测数据

朱峰, 程玖兵, 王腾飞, 耿建华

同济大学海洋地质国家重点实验室，上海 200092

0 引言

主动源拖缆采集是迄今在海洋地质调查、油气与天然气水合物勘探中占主导地位的地震数据观测方式，其中大多数是用水听器记录标量声压场信号.为了给鬼波或表面多次波压制、数据规则化等提供有利的数据条件，同时记录声压和质点速度垂直分量信号的观测方式(俗称拖缆双检)在海洋地震勘探中越来越受到重视(Haggerty, 1956; White, 1965; 杜向东, 2018; 徐进师等, 2018; 耿建华等, 2018).

海面和海底等强阻抗界面致使拖缆地震数据中广泛发育多次波.根据下行反射产生的位置不同，多次波分为表面多次波和层间多次波.在大多数海洋地震数据中，表面多次波能量比层间多次波要强得多，对一次反射波构成极大的干扰，因此需要在数据处理过程中加以有效压制.如图1a，表面多次波包括在海面和地层界面之间多次震荡形成的各阶多次波，以及由于震源和检波器沉于海水中形成的鬼波(又称水层虚反射).源端鬼波常与震源信号混在一起，使震源子波发生改变；检波器端鬼波则常伴生在上行反射波之后，构成形如“双眼皮”的复合波形(图1b).鬼波会导致地震数据出现陷波效应，是提高成像分辨率需要克服的问题.对于以下行波形式记录的检波器端鬼波，一般采用与水陆双检、上-下缆以及斜缆或变深度缆等特殊采集方式对应的方法进行压制(Carlson et al., 2007; Özdemir et al., 2008; Soubaras, 2012; 张威等, 2017; 李洪建等, 2018).

图1 海洋拖缆数据中的自由表面多次波

拖缆地震数据表面多次波压制方法有很多，但都有特定的应用条件或局限性.例如，基于多次波和一次波在时-空和其他变换域分布差异的滤波类方法(Ryu, 1982; Foster and Mosher, 1992; 胡天跃等, 2000; 金德刚等，2008)在二者差异不大时(如在近偏移距)遇到的挑战还值得研究(李列等, 2015).在多次波预测相减方法中，模型驱动的波场外推法要求已知震源子波、海水深度和海底反射系数等信息，当海底起伏剧烈时，应用效果不好(张广利等, 2016)；数据驱动类的表面多次波压制方法，如反馈迭代法(SRME)和逆散射级数法，要求震源和检波器规则同构分布在同一深度上，因此在实际应用中需要精细的数据插值或规则化处理(Wang et al., 2009; Shao et al., 2021).在浅水环境下，随着炮检距的增加浅层迅速达到临界反射，同时近偏移距数据缺失和浅层噪声干扰等进一步增大了多次波的预测和压制难度(周小鹏等，2019; 徐强等，2015).此外，当多次波和一次波交叉或重叠时，基于最小能量假设的自适应匹配相减过程通常会损伤一次波信号(van Groenestijn and Verschuur, 2009a；李鹏等，2007).为了避免匹配相减，SRME被升级为估计一次波的稀疏反演方法(EPSI)，它基于稀疏约束同时反演震源子波和一次波信号(van Groenestijn and Verschuur, 2009b; Lin and Herrmann, 2013).实践表明，上述方法在压制源端表面多次波(含源端鬼波)时的针对性和有效性仍需要提升(Amundsen, 2001; Tu and Herrmann, 2015).

为了压制源端表面多次波，本文基于地震干涉原理(Schuster, 2009)和最小平方基准面延拓理论(Zhu and Cheng, 2022; 朱峰和程玖兵, 2022)提出一种估计一次反射波的约束优化反演方法.它利用处于拖缆沉降深度的基准面上的上、下行波数据，借助地震干涉合成在基准面上虚拟观测的一次反射波信号，通过迭代反演消除上覆介质和自由表面效应的影响.本文首先，基于波场互易定理和地震干涉原理，建立最小平方基准面延拓(LSR)正问题的控制方程，提出求解相应反问题的稀疏约束优化方法.然后，分别针对拖缆双检单检两种观测方式，论述基于LSR压制表面多次波的技术流程及其假设条件.最后，通过合成数据和南黄海实际地震资料，验证本文方法的应用效果与实用化潜力.

1 最小平方基准面延拓理论

最小平方基准面延拓(LSR)以波场互易定理为理论基石，通过线性波形反演重构基准面上的虚拟观测信号.波场互易定理指出，在同一封闭无源区域内两个相互独立的波场状态下，通量归一化的入射波场和出射波场相互褶积或相关沿封闭区域的边界积分是守恒的(Wapenaar and Grimbergen,1996).它对封闭区域以外的震源分布和介质特性无特殊要求.对于勘探地震，当选定封闭区域边界为地下两个无限延伸的水平界面时，沿垂直边界的积分可以忽略，此时的入射波和出射波就是穿过水平界面的下行波与上行波(Wapenaar et al., 2004).利用波场互易定理可以建立两种状态下波场之间的相互转化关系.如图2，本文一种状态为实际观测波场状态(即含有自由表面效应)，而另一种状态是为不受自由表面影响的理想波场状态.从理想状态向实际状态的波场转化过程为正过程(或正问题)，反之则为反过程(或反问题).在LSR理论框架下压制表面多次波，首先要构建正问题的控制方程.在本文姊妹篇“利用最小平方基准面延拓压制表面多次波I: 海底观测地震数据”中(朱峰和程玖兵, 2022)，由褶积型单向波互易定理推导了LSR的正演方程，并且给出了相应的物理诠释.下文针对海洋拖缆数据表面多次波压制的LSR方法将直接引用该方程.

1.1 褶积型地震干涉正演

图2 海洋拖缆地震数据基准面延拓示意图

(1)

1.2 相关型地震干涉基准面延拓

参照Berkhout(1982)，将公式(1)写成矩阵形式：

D-=D+R∪.

(2)

等号两侧同乘下行波信号的共轭D+*，得到

D+*D-=D+*D+R∪，

(3)

式中D+*D-是互相关函数的集合，D+*D+是海森算子(Wapenaar et al., 2011).从数学上看，可通过如下多维反褶积计算虚拟观测信号，即(Wapenaar et al., 2000; Amundsen, 2001)：

(4)

上式中分子对应的互相关运算相当于褶积型地震干涉的伴随过程，它对应基于相关型地震干涉的虚源法(Schuster and Zhou, 2006; Bakulin and Calvert, 2006).如图3，相关型干涉一定程度上消除检波器B实际记录的上行波信号中下行经过检波器A的源端多次波路径的影响.由于虚源法采用的多维互相关不是正过程(即多维褶积)的逆运算，获得的虚拟数据仅是理想观测信号的模糊滤波形式，含有一些假信号干扰.海森算子反映了观测系统、震源子波以及上覆介质中的传播与自由表面效应的影响.就单个炮检对而言，空间-时间域的海森算子对应点扩散函数(PSF).按(4)式由海森算子或点扩散函数参与的多维反褶积运算，正好达到去模糊滤波作用，从而更精确地重构虚拟观测信号.

图3 褶积型地震干涉示意图

对于实际规模的地震数据处理而言，海森算子的显式计算、存储以及多维反褶积运算都对计算机资源条件提出过高要求，而且这样操作不便于对模型参数(这里指虚拟观测信号)施加恰当的约束或正则化处理.

1.3 虚拟观测信号约束优化反演

本文采用最小二乘反演方法求解上述反问题.以干涉法基准面延拓估计的虚拟观测信号为初始模型，按方程(2)描述的多维褶积合成上行波数据，在L2范数意义下使其与实际观测的上行波数据的残差(δD-=D+R∪-D-)最小化.同时，由于理想的虚拟观测信号具有脉冲响应特征，故对其施加基于L1范数的正则化处理，以提升反演的鲁棒性.于是针对该稀疏反演问题构建如下目标泛函：

(5)

其中ε为权重因子.

采用共轭梯度算法迭代求解该反问题，即有：

(6)

其中k表示迭代次数，α表示更新步长，zk是由泛函梯度gk修正得到的共轭梯度方向.参照本文姊妹篇(Part I)附录B(朱峰和程玖兵, 2022)，通过伴随状态方法推导得到泛函梯度的计算式：

(7)

在上述最小平方基准面延拓(LSR)的迭代求解过程中，正则化处理通过Split Bregman方法进行施加(Osher et al., 2005; Yang et al., 2016)，具体流程见Part1算法1.

2 用于基准面延拓的上、下行波数据

为了实现数据驱动的基准面延拓，无论采用虚源法还是LSR方法，都需要输入在基准面上提取的上行波和下行波数据(图4).下面分别就海洋拖缆双检与单检观测方式，讨论如何准备这些输入数据.

图4 海洋拖缆地震数据波路径示意图

2.1 拖缆双检观测

在拖缆双检观测条件下，声压场p和垂直速度vz在上行波数据D-中极性相同，在下行波数据D+中极性相反，通常可按如下公式分离出通量归一化的上、下行波信号(Amundsen, 1993; Wapenaar, 1998)：

(8)

2.2 拖缆单检观测

在单分量观测条件下，仅由水听器记录到标量的声压场信号，缺少质点振动矢量场信号，很难实施准确的上、下行波分离处理.通常，震源和检波器均离海面较近(气枪震源深度一般不超过15 m，拖缆沉放深度一般不超过30 m).如图4，假设拖缆沉放深度非常接近海面(即Δh≈0)，上行波与其经海面一次反射的下行波到达水听器的走时差较小，主要表现为相位或极性差异.实际观测的上行波信号由一次反射波M0与多次波Mi叠加而成，即

(9)

其中Mi(i≥1)代表i阶表面多次波.设多次波经海面最后一次反射对应的反射系数为R，则可把观测数据中的下行波信号(图4中虚线波路径)表示为

(10)

将(9)和(10)式代入(2)式，得到

D-=(S+RD-)R∪.

(11)

该式与SRME或EPSI采用的正演方程是一致的(van Groenestijn and Verschuur 2009a).因此，当震源和检波器处于相同水深时，本文LSR退化为EPSI方法.有所不同的是，后者需要同时估计R∪和震源函数S.由于公式(2)左右两侧的上、下行波数据均含有震源子波信息，LSR不依赖于震源函数，故对拖缆数据震源一致性要求较低，反演过程中仅迭代估计虚拟观测信号.

3 数值实验

这里通过三组实验验证本文方法的有效性.首先，采用合成的拖缆双检地震数据测试基于LSR的多次波压制算法.然后，通过合成的拖缆单检地震数据分析上、下行波信号近似处理对多次波压制效果的影响.最后，将该方法用于处理南黄海拖缆双检地震数据，检验算法模块的实用性.在前两组实验中，为了便于控制多次波的阶数，采用全波场模拟方法(Berkhout, 2014)合成叠前地震数据，其中震源采用主频为20Hz的雷克子波.

3.1 拖缆双检合成数据

针对图5所示的速度模型，拖缆双检观测系统如下：在海面以下10 m深度从横向4.0 km到12.3 km以25 m的间隔均匀激发333炮；拖缆沉放在海面以下20 m深度，缆长为4.0 km，检波器间隔为12.5 m.将基准面设置在20 m深度上，假设基准面上虚拟震源和检波器在横向上分布与原始观测系统保持一致.对于合成的声压-垂直速度信号，基于PZ叠加/相减分离出上、下行波数据，从中均可观察到明显的表面多次波(图6a与6b).除直达波之外，上行波与对应的下行波极性正好相反.图6c为相关型干涉法基准面延拓重构出来的虚拟炮记录(作为LSR的初始模型).可见，一次反射波信号得到了加强，但存在较强的串扰同相轴，垂向分辨率也有所下降.经过50次迭代的LSR重构的虚拟炮记录如图6d所示，几乎看不到表面多次波和源端鬼波的残余，而且一次反射波的垂向分辨率也明显提高.

图5 速度模型与拖缆观测示意图

基于偏移距为1.0 km的共偏移距道集直观地展示多次波压制效果.如图7a所示，受自由表面影响，原始声压场信号中每个上行波均立即伴生一个下行波，二者的复合波形呈“双眼皮”现象.经上、下行波分离之后的，上行波数据不再受下行波(主要是检波器端鬼波)的干扰，但保留较强的源端多次波能量，其中的源端鬼波使等效源子波畸变，导致整个上行波剖面的垂向分辨率降低(图7b).经相关型干涉重构的虚拟数据存在大量串扰噪声干扰(图7c)，而经50次迭代的LSR，这些串扰同相轴几乎完全消除，共偏移距剖面上表面多次波(包括源端鬼波)得到了较彻底的压制，垂向分辨率也得到了提升.

3.2 拖缆单检合成数据

对于拖缆单分量观测，在水听器记录的声压场数据中，特定上行波与其对应的下行波之间除了存在极性翻转，还存在因缆深引起的随偏移距变化的时移量.由于缺少垂直速度波场记录，上、下行波很难准确分离.因此，下面尝试两种近似分离方法(方法A与方法B).在方法A中，先将声压场数据去除直达波和检波器端鬼波，视为上行波信号(图8a).将其在极性翻转后与直达波叠加可近似估计下行波信号(图8b).作为数据驱动的波形匹配反演方法，这样的上、下行波近似处理给LSR带来了明显的影响，导致虚拟观测信号仍然残留一些断断续续的多次波能量，而且深层一次反射波信号在小偏移距未能有效重构(图8d).这主要是因为下行波与其对应的上行波时移量在近偏移距较大，导致上述近似影响地震干涉与数据拟合效果.

图8 方法A对应震源位于x=8.15 km处的拖缆单检炮道集与基准面延拓结果

在方法B中，为了避免缆深引入的走时误差，将含检端鬼波的原始声压场数据直接视为下行波信号，进一步去除直达波后视为上行波信号.此时，上行波与下行波分别给对方构成干扰，尤其在近偏移距出现明显的“双眼皮”现象.随着偏移距增加，上行波与对应下行波时差逐渐减小并最终重合在一起.对于近似的下行波数据，复合同相轴的下部以有效的下行波信号为主，上部基本都是上行波串扰.对于近似的上行波数据，现象恰好相反.由于上述串扰噪声的影响，相关干涉法基准面延拓重构的炮记录多次波能量仍然很强，串扰噪声也很多，垂向分辨率较低(图9c).经LSR重构的虚拟炮记录基本消除了表面多次波，近偏移距反射波重构效果略有提升，但所有的一次反射波仍然呈现复合波形特征(图9d).

图9 方法B对应震源位于x=8.15 km处的拖缆单检炮道集与基准面延拓结果

对于拖缆单检观测数据，上述两种对上、下行波的近似处理会一定程度上降低LSR压制多次波的作用.图10中偏移距为1.0 km的共偏移距剖面对比发现，基于近似方法A的LSR处理一定程度较合理地重构出了一次反射波，但因缆深导致的上、下行波之间的时差影响了表面多次波的压制效果(图10a和10c)；基于近似方法B的LSR压制多次波效果更好，但受检波器端鬼波影响，重构的一次反射波波形存在失真，表现出“双眼皮”特征(图10b和10d).与观测的声压信号(图10e)相比，两种方法均明显压制了源端鬼波和震源子波的影响.在实际拖缆单检观测条件下，当拖缆沉放较浅时，建议先压制检波器端鬼波，然后由LSR压制源端多次波.当拖缆沉放较深时，可先由LSR压制源端多次波，再对重构数据施加去检端鬼波处理.压制检波器端鬼波的方法可参考(王芳芳等，2013；Wang et al., 2013).

图10 近偏移距剖面对比

3.3 南黄海水陆双检实际地震数据

针对南黄海地区某二维拖缆双检观测地震数据，开展基于LSR的多次波压制处理.该测线海底深度介于54.6～63.6 m之间(属于浅水环境)，气枪震源深度为10 m，拖缆沉放深度为25 m.震源间隔为25 m，检波器间隔为12.5 m，缆长为8.3 km，最小偏移距为150 m.经过远偏移距(折射波)切除、振幅校正(包括几何扩散补偿、三维至二维传播效应校正)、基于炮检互易补充负偏移距地震记录以及近偏移距插值等一系列前处理之后，采用PZ叠加/相减分别获得上、下行波数据(图11a与11b).由于自由表面和硬海底的存在，这些上、下行波数据均含有很强的表面多次波.由于新生界与中古生界之间存在强阻抗界面，其能量屏蔽作用导致中古生界内幕反射信号很弱，几乎完全被多次波所掩盖(耿建华等，2018).经仔细核对，发现近偏移距相邻阶次多次波的延时约为80 ms，与该区平均水深存在直接对应关系.经相关干涉法基准面延拓处理，多次波并未得到有效压制，同相轴振幅关系也不太合理(图11c).然而，经本文LSR方法处理之后，多次波能量得到显著削弱，一次反射波同相轴连续性增强，振幅特征趋于合理(图11d).

图11 南黄海拖缆双检地震数据上-下行波分离与多次波压制

从多炮道集分选的最小偏移距剖面及其自相关谱可以更直观地评估整条测线多次波压制效果.如图12a，在双检PZ叠加分离出的上行波信号中，80 ms附近海底界面的反射和0.9 s附近强阻抗界面的反射之下均伴生一些较强的水层震荡相关的表面多次波(图12a右侧箭头).相关干涉法基准面延拓仅压制了少部分多次波，反射同相轴振幅一致性较低(图12b).经过多轮迭代之后，LSR压制了绝大多数多次波能量，相对合理地保持了一次反射波的振幅与相位特征(图12c).图13展示了近偏移距上行波数据和LSR重构数据的自相关谱.前者在零延迟两侧存在一系列对称分布的旁瓣，间隔时间约为80 ms，反映了水层中震荡相关表面多次波的强烈影响；后者能量集中分布在零延迟附近，旁瓣能量非常弱，说明表面多次波得到了有效压制.

图12 近偏移距剖面对比

图13 近偏移距剖面自相关谱对比

4 讨论

图14 震源位于x=8.15 km处的预测波场

LSR和EPSI均属于在反演框架下压制自由表面多次波的方法，其区别和联系在于：前者的正演方程由褶积型互易定理推导得到，描述了有自由表面效应的实际观测状态和无自由表面效应的理想观测状态之间的波场转化关系.后者与SRME方法一样，都是借助WRW波场传播模型(Berkhout, 1982)描述多次波和一次波之间的数学物理关系.LSR方法基准面选取更加灵活，同时适用于拖缆和海底观测数据表面多次波压制问题，对子波一致性无明确要求；而EPSI方法要求震源和检波器处于相同水深，仅适用于拖缆观测数据，且对子波一致性有要求.此外，EPSI要求震源和检波器位置尽量耦合排布，而LSR则无此要求.对于震源子波已知，且已去除检波器端鬼波的拖缆单检数据，LSR方法自动退化为EPSI方法.LSR和EPSI在每轮迭代中的正演和求取泛函梯度的核心运算为多维褶积和多维互相关，二者计算量相当，均高于SRME方法.

基于LSR压制表面多次波要求输入在基准面上分离或估计的上、下行波场信号，这依赖于地震数据的采集方式.在本文姊妹篇中，海底四分量检波器同时记录了声压场信号和质点矢量振动三个速度分量信号，在完成仪器响应和耦合矫正以及水-陆检信号标定之后，当转换横波(S波)较弱时，联合声压和垂直速度分量的PZ叠加/相减能分离出上、下行纵波(P波)信号；当转换横波能量较强时，结合海底弹性参数估计的弹性波场分离方法能同时分离上、下行P波和S波(刘学义等，2021).以这些分离后的上、下行波数据作为本文LSR的输入，可以重构不受表面多次波影响的PP波和PS波反射数据.在本文拖缆观测条件下，对于双检数据，完成仪器响应矫正和水-陆检数据标定之后，同样是按PZ叠加/相减分离上、下行波.而对于拖缆单检数据，由于仅记录声压场信号，观测数据的上、下行波分离存在困难，在LSR中对上、下行波数据的近似处理会降低重构信号的质量或多次波压制效果.对于较大的拖缆沉放深度，上、下行波之间的时移量随偏移距变小而增大，可由静校正或基准面延拓予以消除(Tu and Herrmann，2015)，从而降低对单检拖缆数据LSR处理的影响.拖缆观测数据缺失近偏移距信号会给LSR处理效果带来影响，需提前借助炮检互易性和插值技术加以弥补(如van Groenestijn and Verschuur, 2009b; Wang et al.,2009).为了改善深层弱反射信号重构质量，有必要结合曲波变换等对数据的稀疏表征，在LSR中采用更有效的约束反演算法(如Lin and Herrmann，2013).

本质上，LSR是一种数据驱动的波形反演方法，其目的是消除基准面上覆介质中的波场传播效应，包括自由表面效应以及上覆地层中的多次散射(如层间多次波)效应.因此，如果有条件实现水平井观测，则可以把基准面设定在水平井深度，利用本文方法可以重构在基准面虚拟观测的地震信号，从而克服水平井之上复杂近地表和地质体给深部目标探测造成的不良影响.这是今后值得研究的课题.

5 结论

本文回顾了基于波场互易定理和地震干涉原理的最小平方基准面延拓方法，籍此提出了压制海洋拖缆地震数据表面多次波的新方法.对于拖缆双检观测方式，基于声压场和垂直速度信号可以在拖缆深度(即基准面)上准确分离上、下行波数据，因此LSR能够有效地压制表面多次波，同时消除震源子波和源端鬼波效应，高分辨率重构基准面之下的反射或背向散射响应.对于拖缆单检观测方式，只能根据声压场信号近似获得上、下行波数据，会一定程度影响LSR重构数据的质量.如果检波器端鬼波提前得到压制，则LSR处理效果随缆深变浅会得到改善.如果检波器端鬼波没有提前消除，则对上、下行波数据形成干扰，导致LSR重构数据残余鬼波脚印，最终获得的反射或散射波信号出现“双眼皮”现象.南黄海二维浅水拖缆双检地震数据处理结果表明本文LSR方法具有实用性.今后将进一步完善相应的算法流程，配套更精细的近偏移距插值等前处理技术，使之具备大规模应用前景.