莺歌海盆地位于印支地块与华南地块之间，天然气资源丰富。盆地包括莺东斜坡、莺西斜坡、中央坳陷和临高凸起四大构造单元，其中，中央坳陷发育了大量串珠状排列的、由泥和热流体幕式上侵活动而形成的“泥—流体底辟构造”，由于底辟运动，该区域断裂发育，天然气沿着断裂从深部层系垂直运移至浅层层段中聚集成藏，形成浅层气^[1-3]。莺歌海盆地地震资料受浅层气和底辟构造的影响，其下覆地层在地震剖面上表现为模糊区，即地震反射连续性变差或中断，内部反射较杂乱，甚至为空白反射，局部见同相轴上拱或下拉现象^[4-6]。图 1为该区典型地震剖面，地震资料最明显的特点是多次波非常发育，且种类多，能量强，目标层微弱的有效信号完全被由于鸣震、自由表面和由浅层气产生的强多次波掩盖，严重影响有效地层的识别，因此多次波问题是该区勘探必须解决的问题。

目前多次波衰减主要有三大类方法^[7-12]，一是基于多次波周期性的预测类方法，如预测反褶积；一是基于视速度差异的多次波压制方法，如Radon变换；还有一种是基于波动方程理论的预测与匹配相减方法，如自由表面多次波衰减方法(SRME)。相对而言，基于波动方程理论的预测与匹配相减方法更适应于复杂构造区域的多次波剔除处理，是目前应用最为广泛的多次波剔除技术。

基于波动方程理论的预测与匹配相减方法中多次波模型与原始数据的匹配相减是决定多次波压制效果的重要环节之一。自适应匹配相减的目的就是改善所预测多次波模型的振幅、相位、频率和旅行时信息，借助于匹配滤波器，使预测多次波与原始地震数据的多次波模型更好的匹配^[8]。最常用的匹配滤波方法是最小二乘自适应匹配，该方法是基于L2范数准则，此方法可以有效、快捷地对多次波模型进行匹配相减，但是当一次波和多次波的相关性高时，滤波算子也会对一次波进行匹配，造成多次波残留或损伤一次波能量^[13-14]，因此常规基于L2范数准则的自适应匹配衰减难以适应莺歌海盆地低信噪比地区的多次波剔除。

Sergi^[15]提出了一种结合稀疏变化及相关的匹配滤波器的复小波域的联合多次模型自适应减法，其在连续小波框架下使用非固定Wiener匹配滤波器进行多道处理。在该子波框架下，灵活冗余和多道应用可以更好的管理模型错位误差的时间变化，因此联合多次波模型可以在保护一次波不受损伤的前提下，最大限度地压制多次波，是低信噪比地区的理想多次波剔除方法。本文将该方法引入到莺歌海盆地的多次波剔除处理中，有效压制了该区域地震数据的多次波，显著提高了该弱信号区域的信噪比，也为该类地区的多次波高精度剔除处理提供指导和借鉴。

1 复小波域的联合多次模型自适应减法

小波分析方法具有时频变换的局部性和自适应性，是信号处理的有效算法之一。小波变换通过伸缩和平移运算对信号进行多尺度分解，来逼近原始信号，根据逼近程度能够有效地从信号中获取各种时频信息，它在时域和频域同时具有良好的局部化性质，克服了短时傅里叶变换分析时窗固定的缺陷，具有多分辨率能力，且计算效率高^[16-20]。

复小波变换是将信号沿两正交空间同时做实小波变换，同时得出2个空间的分量信息，然后将这2个空间的信息联合即可得出复小波变换的幅值和相位信息。复小波变换时频分析能力较强，幅值、相位信息对信号突变更为敏感，能检测出信号突变的时间信息^[21-27]。考虑到复Morlet小波采用时频窗面积最小的高斯窗函数，且时频域局部化性能都很好，并且对称性较好，故通常选用复Morlet小波作为母小波进行振动信号小波分析。

该方法具体为：在复Morlet小波框架下模拟在特定尺度下计算复数导数，然后在小波域中的每个尺度上使用自适应一维滤波器。该方法分为两步：首先，在变换中选择适当的冗余，以确保对非相干绕动的稳定性，它可以在具有不同冗余度和噪声水平的合成数据上使用；其次，将单模型复小波自适应多次波减去法扩展到联合多次波模型方法。该方法利用小波框架的辨别力简化了一个长匹配滤波器，重新设计了全局和局部复合一维滤波器，最大限度地减少多次波同向轴与模型匹配之间的误差^[28-31]。

1.1 复小波变换分解

经典地震道描述模型如下：

$ d[n] = p[n] + m[n] + w[n]。$

(1)

其中：d[n], p[n], m[n]和w[n]分别代表的是记录数据、一次波、多次波和背景噪音，n为离散时间。

用连续小波框架的离散近似来执行每个数据d[n]和多次波模型x_k[n]的时间尺度分解。复Morlet小波在转换域中可表征出振幅和相位延迟信息，且其可在小波尺度上模拟复地震道，因此本文选择复小波进行小波分析。其公式为：

$ \psi \left( t \right) = {\pi ^{ - 1/4}}{e^{ - i{\omega _0}t}}{e^{ - {t^2}/2}}。$

(2)

其中：ω₀是调制高斯中心频率，t是连续时间变量。

相关离散函数被定义为母小波的采样：

$ \psi _{r,j}^\nu[n] = \frac{1}{{\sqrt {{2^{j + \frac{\nu}{V}}}} }}\psi \left( {\frac{{nT - r{2^j}{b_0}}}{{{2^{j + \frac{\nu}{V}}}}}} \right)。$

(3)

其中：T是时间采样，V为每倍频程的音数。指数r, j∈Ζ, ν∈[0, …, V－1]，分别表示离散时间，倍频程音数。最后，b₀代表零值时的采样周期，接近$\frac{{2v}}{{{b_0}}}$的整体冗余控制着计算效率和稳定性之间的平衡。

在时间尺度上，d[n]的表达式可写为：

$ \begin{array}{l} d = d_{r,j}^\nu = \left\langle {d[n],\psi _{r,j}^\nu[n]} \right\rangle = \\ \sum\limits_n {d[n]} \overline {\psi _{r,j}^\nu[n]}。\end{array} $

(4)

假设从d的一些时间尺度处理得到结果 $\widehat d $，在此进行模型的匹配和相减，然后把应用双框架的滤波后结果合成转换回时间域。

$ \hat d[n] = \sum\limits_r {\sum\limits_{j,\nu } {\hat d_{r,j}^\nu \widetilde \psi _{r,j}^\nu [n]} } 。$

(5)

$\widetilde \psi _{r, j}^\nu$[n]作为双框架组件的总和，$\widehat d_{r, j}^\nu$适应多次波分解加权。实际上，双框架较好地近似于一个常数因子的分析框架，V≥3和b₀≤1.5。

1.2 单模型一维滤波估计

当模型和多次波序列之间的延迟差异小于所有尺度的半个周期时，单个多次波模型x₁可以进行最小二乘误差(LSE)方法时间校正。最优一维滤波器是指在一个给定的小波尺度下，道的局部或全局被定义为复合体标量a₁，乘以时间尺度分解的多次波模型x₁，使滤波数据与滤波模型正交，如公式(6)表示：

$ {a_{opt}} = \mathop {\arg \;\min }\limits_{{a_1}} {\left\| {d - {a_1}{x_1}} \right\|^2}。$

(6)

其中复数标量a_opt用于补偿局部延迟和振幅不匹配^[7]。

1.3 联合多次波模型一维滤波估计

当有几种不同的多次波模型x_k可用时，不同的延迟和振幅可能会影响每个可用模型和实际的多次波序列之间的耦合。公式(6)相应修改为：

$ {a_{opt}} = \mathop {\arg \;\min }\limits_{\{ {a_k}\} \left( {k \in K} \right)} {\left\| {d - {a_k}{x_k}} \right\|^2}。$

(7)

LSE的意义是使最佳滤波信号$\widehat d = d-\sum {_k{a_k}{x_k}} $与k个滤波模型a_kx_k正交，使用公式是：

$ \left\langle {d - \sum\limits_k {{a_k}{x_k},{a_m}{x_m}} } \right\rangle = 0,\forall {a_m} \ne 0。$

(8)

应用第一个参数线性和第二个参数共轭线性的性质，可得：

$ \left\langle {d,{x_m}} \right\rangle = \sum\limits_k {{a_k}\left\langle {{x_k},{x_m}} \right\rangle } 。$

(9)

即复合信号的向量维纳方程。在实践中，由于一些多次波模型是局部相似的，互相关矩阵常常接近于奇异值。因此需保证互相关矩阵特征值对应的特征向量高于规定的阈值。

2 模型数据测试

为了了解上述方法的适用性，首先建立了如图 2所示的含有浅层气特殊异常体的层状介质模型，海水层速度为1 520 m/s、密度为1.5 g/cm³、Q值为6 000；第一层地层速度层速度为2 500 m/s、密度为2.2 g/cm³、Q值为100；第二层地层速度层速度为2 300 m/s、密度为2.1 g/cm³、Q值为270；第三层地层速度层速度为3 000 m/s、密度为2.4 g/cm³、Q值为200；第四层地层速度层速度为2 400 m/s、密度为2.15 g/cm³、Q值为120；在浅层地层中设置了1个浅层气异常体，异常体的层速度为1 580 m/s、密度为0.8 g/cm³、Q值为10；根据地质模型进行了粘弹性的波动方程正演模拟，然后在模型数据上进行测试。

图 3(a)和图 4(a)为模型正演后得到的炮集和叠加剖面，图 3(b)和图 4(b)为使用常规方法衰减多次波后的炮集和叠加剖面，图 3(c)和图 4(c)为使用本文介绍的方法衰减多次波后的炮集和叠加剖面。在(a)图炮集和叠加剖面中，可以明显看到海底、各地层、浅层气的一阶、二阶多次波，这些都是有关自由表面的多次波，这种多次波的特征是形态和产生源一致，能量强，多次波之间时差为海底的倍数；同时，尤其在浅层气下覆能看到能量较弱的层间多次波，这些多次波的存在严重干扰着有效反射。(b)图中各种多次波有明显的残留，尤其对夹杂在浅层气下覆地层中的层间多次波衰减效果较弱，(c)图中多次波衰减更加干净，残留不明显，尤其对夹杂在浅层气下覆地层中的层间多次波衰减效果更加理想。

3 莺歌海盆地地震数据多次波剔除分析

本文以莺歌海盆地L气田资料为例，深入分析本文介绍方法对于莺歌海盆地地震数据多次波剔除的有效性。L气田某测线实际炮集和CMP道集如图 5(a), 6(a)所示，应用常规方法剔除多次波后的炮集和CMP道集如图 5(b)和6(b)所示，而应用本文介绍方法剔除多次波后的炮集和CMP道集见图 5(c)和6(c)。

由图 6可知，该测线实际CMP中发育有强能量短周期多次波及层间多次波。应用常规方法衰减多次波后的CMP中，尚能见到明显的多次波残余，而应用本文介绍方法剔除多次波后的CMP中，残余多次波已基本剔除干净，充分显示出本文介绍方法在该区域多次波剔除中的有效性。

图 7(a)为多次波衰减前的叠加剖面，图 7(b)使用常规方法衰减多次波后的叠加剖面，图 7(c)为使用本文方法衰减多次波后的叠加剖面。从叠加剖面的效果分析来看，常规多次波自适应相减法处理后，大部分多次波得到有效衰减，但仍然可以清晰地看到残留了较多多次波，尤其在浅层气下覆，依然影响有效地层的识别。应用本文介绍方法处理后，鸣震和由浅层气产生的强短周期多次波衰减效果更加彻底，模糊区的信噪比提高，微弱的有效信号更加突出，连续性明显，为之后精确成像提供可靠基础数据。

4 结语

莺歌海盆地由于地质原因，产生的多次波种类复杂且能量较强，常规多次波自适应减去法难以取得理想的多次波剔除效果，多次波残留明显。本文将复小波域的联合多次模型自适应相减方法引入莺歌海盆地的多次波剔除处理中，处理结果显示鸣震和由浅层气产生的强短周期多次波衰减效果更加彻底，模糊区的信噪比提高，微弱的有效信号更加突出，连续性明显，为之后精确成像提供可靠基础数据。