0 引言

地震波速度是影响地震成像和解释结果的重要参数之一^[1]，可以反映地层岩性和构造特征^[2]。速度分析是速度建模的常用手段，速度谱是其最核心的参考依据，通常通过对共中心点叠前道集进行速度扫描得到^[3]，并且需要综合考虑多种因素多次迭代和验证，以获得准确可靠的结果^[4]。

海洋地震数据在速度分析时受崎岖地层、多次波和绕射波等多种因素影响，导致速度谱能量团发散，单一同相轴呈现多能量团的特征，给速度分析增加了难度。优化速度谱能帮助更精准地捕捉地下岩层的速度信息，直接提高速度建模的准确性并降低其多解性，有效提升地震资料的成像品质。

近年来，随着人工智能的迅速发展，不少专家致力于半自动/自动速度拾取的研究工作，没有考虑速度谱本身的精度问题。方法主要有蒙特卡罗法^[5]、非线性函数优化法、共轭梯度法、基于加权相似性^[6-9]等自动寻优速度拾取算法，以及BP神经网络、卷积神经网络、循环神经网络和混合神经网络等更加智能化的速度拾取算法^[10-13]。影响速度谱的精度和分辨率的主要因素包括炮检距分布、覆盖次数、信噪比、切除和地震数据的频带宽度等。刘立彬等^[14]、潘宏勋^[15-16]、Xu等^[17]、谢玉洪等^[18]、王小青等^[19]、周锦明^[20]采用的方法都是从炮检距分布、信噪比、切除等角度改进速度谱品质。

本文提出基于经验模态分解(EMD)的地震数据速度谱优化方法。首先对地震数据进行频移处理，通过Hilbert变换获得地震数据的瞬时振幅，对瞬时振幅进行EMD；然后在分解后的本征模量(IMF)中，选择具有有益表达速度谱信息的本征模态模量，构建新的速度谱数据，最终拓宽共中心点道集的有效频带，加强速度谱超道集的低频成分，从而实现优化速度谱的目的，为后续速度拾取提供高品质数据。

1 频移处理与速度谱优化

本文的目的是通过对地震数据频移处理，突出数据低频能量，进而提升速度谱的精度，为后续偏移成像提供精确可靠的速度基础。图 1是本文速度谱优化的技术方案。

1.1 计算瞬时振幅谱

对于一个实值的地震信号$ x\left(t\right) $，可以构造它的解析信号$ c\left(t\right) $

$ c\left(t\right)=x\left(t\right)+\mathrm{i}H\left(t\right) $

(1)

式中：$ c\left(t\right) $是地震实信号$ x\left(t\right) $的复信号，也称为解析信号；$ \mathrm{i} $=$ \sqrt[]{-1} $是虚数单位；$ H\left(t\right) $是$ x\left(t\right) $的Hilbert变换，即

$ H\left(t\right)=x\left(t\right)\mathrm{*}\frac{1}{\mathrm{\pi }t}=\frac{1}{\mathrm{\pi }}{\int }_{-\mathrm{\infty }}^{+\mathrm{\infty }}\frac{x\left(\tau \right)}{t-\tau }\mathrm{d}\tau $

(2)

式中$ \tau $为积分变量，用于表示在积分中对信号进行卷积的点。Hilbert变换的实质是将原始信号和信号$ 1 / \pi t $做卷积。式(1)还可以写为

$ c\left(t\right)=A\left(t\right){\mathrm{e}}^{\mathrm{i}\phi \left(t\right)} $

(3)

式中$ \phi \left(t\right) $为瞬时相位，可以由下式求取

$ A\left(t\right)={\left[{x}^{2}\left(t\right)+{H}^{2}\left(t\right)\right]}^{\frac{1}{2}} $

(4)

$ \phi \left(t\right)=\mathrm{a}\mathrm{r}\mathrm{g}\left[c\right(t\left)\right] $

(5)

1.2 瞬时振幅EMD

Hilbert变换是地震信号处理中常用的方法。地震信号是实信号，在复信号分析方法中需计算其对应的复信号。因此，可以将待分析的地震信号作为实部，其Hilbert变换作为虚部，即可构成复信号，该复信号可用于计算信号的瞬时属性。基于Hilbert变换获得地震资料的三瞬属性，即瞬时振幅、瞬时频率和瞬时相位，可用于地震数据的频率研究。这些瞬时属性可以揭示地震信号的时变特性，对地震数据的频域分析提供了有价值的信息。EMD是由Huang等^[21]于1998年提出的一种新型自适应信号时频处理方法，特别适用于非线性非平稳信号的分析处理。EMD可将地震信号分解成一系列IMF，每个IMF都代表了信号在不同尺度上的振动模式。EMD有助于更好地理解地震信号的时变特性和频率分布。

EMD可以将非线性非平稳信号分解为几个本征模态分量和一个残差之和，对式(4)求得的瞬时振幅$ A\left(t\right) $进行EMD，即

$ A\left(t\right)=\sum\limits_{k=1}^{K}\mathrm{I}\mathrm{M}{\mathrm{F}}_{k}\left(t\right)+\mathrm{R}\mathrm{e}\mathrm{s}\left(t\right) $

(6)

式中：$ \mathrm{I}\mathrm{M}{\mathrm{F}}_{k}\left(t\right) $为本征模态分量；$ K $为本征模态分量个数。

本征模态分量必须满足以下两个条件：

(1) 整个时间范围内，局部极值点和过零点的数目必须相等，或最多相差一个；

(2) 任意时刻点，局部最大值的包络(上包络线)和局部最小值的包络(下包络线)平均必须为零。

1.3 筛选本征模态分量

根据实际地震数据的频带范围及本征模态分量的频率范围，筛选出有益于提高速度谱品质的本征模态分量$ \mathrm{I}\mathrm{M}{\mathrm{F}}_{1}\left(t\right) $，$ \mathrm{I}\mathrm{M}{\mathrm{F}}_{1}\left(t\right) $，…，$ \mathrm{I}\mathrm{M}{\mathrm{F}}_{L}\left(t\right) $，将它们组合构建新的地震数据$ {x}_{\mathrm{n}\mathrm{e}\mathrm{w}}\left(t\right) $，结果可表示为

$ {x}_{\mathrm{n}\mathrm{e}\mathrm{w}}\left(t\right)=\sum\limits_{k=1}^{L}IM{F}_{k}\left(t\right) \quad L\le K $

(7)

新的残差$ \mathrm{R}\mathrm{e}{\mathrm{s}}_{\mathrm{n}\mathrm{e}\mathrm{w}}\left(t\right) $由以下公式表示

$ \mathrm{R}\mathrm{e}{\mathrm{s}}_{\mathrm{n}\mathrm{e}\mathrm{w}}\left(t\right)=A\left(t\right)-{x}_{\mathrm{n}\mathrm{e}\mathrm{w}}\left(t\right) $

(8)

与原始地震数据$ x\left(t\right) $的振幅谱相比，新的地震数据$ {x}_{\mathrm{n}\mathrm{e}\mathrm{w}}\left(t\right) $的有效频率向低频端移动，移动的尺度取决于$ L $的大小，$ L $越小，频移后频率越低。以上步骤的全过程称为频移处理。

1.4 构建新的速度谱

在CMP道集中，使用频移后的地震数据$ {x}_{\mathrm{n}\mathrm{e}\mathrm{w}}\left(t\right) $计算叠加速度谱的方法如下。

(1) 在待分析速度点附近取出若干CMP道集。

(2) 对这些CMP道集使用一系列常速度值进行动校正和叠加。每个速度值对应一个叠加剖面，这一系列由不同速度得到的叠加剖面称为CVS剖面。CVS剖面的叠加能量随“一系列速度”变化而增大或减小，叠加能量相对速度的变化就构成了速度谱。

(3) 当速度合适时，经过动校正后，反射波同相轴被拉平，此时叠加效果最好，叠加振幅出现极大值，叠加能量最强。通过在CVS叠加能量图像中寻找叠加能量最大的位置，可以解释出理想的叠加速度，即所求动校正速度。

叠加能量计算公式如下

$ \overline{A}=\frac{1}{M}\sum\limits_{m=1}^{M}\sum\limits_{n=1}^{N}x_{\text {new }}^{\prime}\left(t_{m, n}\right) $

(9)

式中：$ M $为CMP道集的横向地震道数；N为计算叠加速度的纵向时窗内的样点数；$ x_{\mathrm{new}}^{\prime}\left(t_{m, n}\right)$为$ {x}_{\mathrm{n}\mathrm{e}\mathrm{w}}\left(t\right) $动校正后的第$ m $道的第$ n $个样点的振幅值；$ \overline{A} $为时窗内的平均振幅即叠加能量。

速度分析是以人工经验干预为主的过程，修改$ L $，重复步骤1.3和本节计算叠加速度谱的步骤，迭代至速度谱寻优区间符合速度分析经验预期，输出优化后的速度谱。依据经验$ L $取值通常为2~4。通过以上步骤，可以利用频移后的地震数据计算叠加速度谱，并通过多次优化迭代的方式改善速度谱的质量。

2 频移技术实验

以一条实际海洋拖揽地震资料的叠加数据为例进行频移处理，该数据记录长度为7 s，采样率为1 ms，道间距为3.125 m，有效频带范围约为4~120 Hz，该叠加数据总体上纵横向分辨率较好。从图 2中可以看出，瞬时振幅(图 2b)主要反应了原始叠加数据(图 2a)的反射振幅强度。由图 3可见，第一次分解出的本征模态分量$ \mathrm{I}\mathrm{M}{\mathrm{F}}_{1}\left(t\right) $(图 3a)具有最高的频率，后续每次EMD的本征模态分量(图 3b~图 3d)频率逐渐降低。第5次分解出的本征模态分量$ \mathrm{I}\mathrm{M}{\mathrm{F}}_{5}\left(t\right) $(图 3e)的频率小于4 Hz，而残差$ \mathrm{R}\mathrm{e}\mathrm{s}\left(t\right) $(图 3f)的振幅很小，需要放大10⁷倍才能被肉眼所见。这里将本征模量$ \mathrm{I}\mathrm{M}{\mathrm{F}}_{1}\left(t\right) $和$ \mathrm{I}\mathrm{M}{\mathrm{F}}_{2}\left(t\right) $相加构建新的地震数据$ {x}_{\mathrm{n}\mathrm{e}\mathrm{w}}\left(t\right) $并计算新的残差$ \mathrm{R}\mathrm{e}{\mathrm{s}}_{\mathrm{n}\mathrm{e}\mathrm{w}}\left(t\right) $(图 4)，即完成频移处理全过程。对比频移处理前(图 2a)、后(图 4a)叠加剖面，频移处理后的叠加数据呈现明显的低频能量增强(图 5); 瞬时振幅的频谱呈现指数形式，且谱低频强；频移处理后，有效频带向低频端移动，与处理后叠加数据显示的低频能量特征相吻合。综上所述，本实验证明频移处理能够有效实现地震数据频带向低频端移动的目的。

3 速度谱优化实例及效果分析 3.1 二维拖缆地震数据

以一条发育崎岖地层的海洋拖缆数据为例，覆盖次数为80，记录长度为12 s，采样率为2 ms，道间距为12.5 m，有效频带相对较窄(约为4~90 Hz)。由图 6可见：处理前，有效反射层较多但缺乏波组特征，呈窄带信号特征；处理后，振幅强弱对比明显，波组关系清晰明确。由图 7可见，频移处理后，有效频带向低频端移动，增加了道集的低频信息。此外，由图 8也可看到，优化前，受崎岖地层的影响，速度谱能量团发散，道集上单一同相轴在速度谱上对应多个能量团，造成速度拾取多解性，进而增加速度误差。优化后，能量团空间上更加聚焦，且与CMP道集主要反射同相轴能一一对应，有效降低了速度分析难度，为精细速度建模提供良好数据基础。该实例证明本文研究方法在崎岖地层发育的拖缆数据中具有显著应用效果。

由图 9可见，速度整体趋势一致，基于优化前的速度谱建立的速度场较平滑(图 9a)，能确定速度的宏观趋势，而经速度谱优化处理后拾取的速度更加精细(图 9b)，能够更清晰地显示局部速度特征。图 10展示了采用图 9的速度场进行偏移的局部结果。整体成像品质略有提升，并在复杂凹陷内幕构造和高倾角基底成像方面有所改善。虽然速度场更加精细，但由于时间偏移对速度信息不够敏感，无法在偏移剖面上呈现显著的效果。

3.2 二维深反射地震数据

低频信息相对于高频信息具有衰减慢、穿透能力强的特点^[22]。受波前扩散和地层吸收衰减作用的影响，地震波在传播至深部时信号变得非常微弱。因此，深部地壳结构和莫霍面的成像主要依赖于低频有效反射能量。深层目标反射通常在双程旅行时间为6~10 s之间，且在基底以下的地层通常具有低信噪比特征，在速度谱上表现为能量团成片的特点^[23]。图 11展示了南海一条深反射地震测线在应用本文速度谱优化技术前(图 11a)、后(图 11b)的速度模型(图 11上)和偏移剖面(图 11下)。由图可见，低频能量突出之后，优化后拾取的叠加速度在目的层(双程旅行时8.5 s附近)与优化前相比提高了1000 m/s以上，更符合壳幔地层的真实速度规律。同时，深部的Moho界面反射也得以凸显，连续性显著提升，这为深部成像特征和南海构造演化研究提供了更可靠的数据基础。

4 小结

本文研究的速度谱优化方法在处理深部低信噪比的壳幔结构成像地震数据速度谱时效果显著，利用频移处理技术将地震数据的有效带宽向低频端拓展，增大了倍频程，提高了速度谱的分辨率。偏移剖面只是从一个角度验证了其效果，沉积层发育较厚的地质条件下无法体现本文速度场的精度优势。另外，地震数据速度谱优化只是频移处理技术应用的一个领域，如果将频移处理应用于地震成果数据，也能有效提高成果数据的解释能力。