0 引言

在最小二乘偏移出现之前，最高级的偏移算法是逆时偏移技术，虽然逆时偏移能够进行高精度的成像，但是无法反演出地层的弹性参数^[1-6]，而正确估计地下介质的弹性参数是寻找岩性油气藏的基础。为解决这一难题，Tarantola^[7]于20世纪80年代在反演框架下提出了最小二乘偏移。该方法通过将正演数据和实际数据的误差达到最小来估计地下介质弹性参数。

最小二乘偏移相当于对常规偏移结果校正，校正的桥梁就是Hessian矩阵。常规偏移结果是Hessian矩阵作用后的最小二乘偏移结果，是模糊了的、振幅畸变了的最小二乘偏移结果，因此Hessian矩阵也被称为模糊化算子。将Hessian矩阵的逆作用于常规偏移结果就可获得最小二乘偏移结果。最小二乘偏移可以对常规偏移结果进行照明补偿，提高偏移成像的分辨率和保真度，成像结果更接近地层的真实反射系数。

最小二乘偏移方法分为两大类，即数据域方法和成像域方法。前者需要多次迭代，使正演与实际炮集数据的误差达到最小^[8-11]。该类方法运算量巨大、运行时间长，但收敛比较稳定、结果更加精确。后者需要显式计算Hessian矩阵^[12-13]，计算难度大，但计算量小、运行时间相对短。两类方法的用时都是常规偏移的数十倍^[14-17]。本文采用的是数据域的最小二乘偏移方法。

在数据域最小二乘偏移中，子波和处理频带是两个非常重要的参数。子波的选取对偏移效果的影响至关重要^[18-22]；频带的选取同样影响偏移效果，同时对偏移耗时有重大影响。

首先通过Sigsbee模型数据试验，论证了当子波频带和数据的频带一致时，最小二乘偏移结果的频带最宽，分辨率最高，而且使用从数据提取子波的偏移结果比使用Ricker理论子波的偏移结果有更宽的频带和更高的分辨率。通过对比陆地与海上数据的最小二乘偏移结果，探讨了前者的最小二乘偏移效果不够理想的原因是子波变化剧烈。最后，通过一块实际数据的试验性处理获得了偏移运算时间与处理频带、数据量之间的经验关系，可为实际生产中最小二乘偏移所需时间的估算提供参考。

1 最小二乘偏移的基本原理

地震散射波理论认为地震记录是地震子波和格林函数的卷积（图 1）

$ d(\boldsymbol{r}, t;\boldsymbol{s})=\int m\left(\boldsymbol{x}\right)W\left(t\right)*G(\boldsymbol{r}, t;\boldsymbol{x})\mathrm{*}G(\boldsymbol{x}, t;\boldsymbol{s})\mathrm{d}\boldsymbol{x} $

(1)

式中：$ d(\boldsymbol{r}, t;\boldsymbol{s}) $表示点s激发、点r接受的记录；m（x）为反射系数；W（t）为震源子波；$ G(\boldsymbol{x}, t;\boldsymbol{s}) $为从震源s到地下散射点x的格林函数；$ G(\boldsymbol{r}, t;\boldsymbol{x}) $为从散射点x到接收点r的格林函数。将式（1）写成矩阵形式

$ \boldsymbol{d}=\boldsymbol{L}\boldsymbol{m} $

(2)

式中L表示正演算子。

在常规偏移中，由于L的逆很难求取，所以用它的转置近似逆，有

$ \boldsymbol{m}={\boldsymbol{L}}^{\mathrm{T}}\boldsymbol{d} $

(3)

用正演算子的转置代替正演算子的逆会使偏移结果的误差增大，但是直接求取正演算子的逆并不现实。为了解决这一问题，可以构建一个误差函数，然后通过迭代的方式，使正演数据d^cal和采集数据d^obs的误差达到最小，获得更准确的偏移结果。这就是最小二乘偏移的核心思想。误差函数的定义为

$ F\left(\boldsymbol{m}\right)={‖{\boldsymbol{d}}^{\mathrm{o}\mathrm{b}\mathrm{s}}-\boldsymbol{L}\boldsymbol{m}‖}^{2} $

(4)

上式的最小二乘解可表示为

$ \boldsymbol{m}=({\boldsymbol{L}}^{\mathrm{T}}{\boldsymbol{L})}^{-1}\boldsymbol{L}{\boldsymbol{d}}^{\mathrm{o}\mathrm{b}\mathrm{s}}={\boldsymbol{H}}^{-1}{\boldsymbol{L}}^{\mathrm{T}}{\boldsymbol{d}}^{\mathrm{o}\mathrm{b}\mathrm{s}} $

(5)

式中：H=L^TL定义为Hessian矩阵，也称模糊化算子；H^-1称为照明算子或去模糊化算子；L^Td^obs为常规偏移的结果。常规偏移结果相当于对最小二乘偏移结果做了Hessian矩阵滤波，是模糊了的最小二乘偏移成像结果，也使振幅发生了畸变。最小二乘偏移结构相当于用Hessian矩阵的逆算子作用于常规偏移结果，是一个去模糊化过程，理论上可以补偿照明不足，提高成像的分辨率和保幅性。

数据域的最小二乘偏移的主要步骤（图 2）如下：

（1）对第n步迭代结果m_n做正演，计算正演数据和实际观测数据的误差

$ {\boldsymbol{d}}^{\mathrm{e}\mathrm{r}\mathrm{r}}=\boldsymbol{L}{\boldsymbol{m}}_{n}-{\boldsymbol{d}}^{\mathrm{o}\mathrm{b}\mathrm{s}} $

(6)

（2）对误差进行偏移，得到成像结果的更新梯度

$ \mathrm{\Delta }\boldsymbol{m}={\boldsymbol{L}}^{\mathrm{T}}{\boldsymbol{d}}^{\mathrm{e}\mathrm{r}\mathrm{r}}{\boldsymbol{m}}_{n} $

(7)

（3）利用更新梯度和正演算子求取更新步长

$ {k}_{n}=\frac{{\left(\mathrm{\Delta }\boldsymbol{m}\right)}^{\mathrm{T}}\mathrm{\Delta }\boldsymbol{m}}{{\left(\boldsymbol{L}\mathrm{\Delta }\boldsymbol{m}\right)}^{\mathrm{T}}\left(\boldsymbol{L}\mathrm{\Delta }\boldsymbol{m}\right)} $

(8)

（4）利用更新梯度和更新步长修正偏移结果

$ {\boldsymbol{m}}_{n+1}={\boldsymbol{m}}_{n}+{k}_{n}\mathrm{\Delta }\boldsymbol{m} $

(9)

重复以上四步，迭代更新偏移结果，直到误差达到设定的门槛值或迭代次数达到最大。

2 子波对最小二乘偏移效果的影响

在最小二乘的迭代过程中，偏移和正演运算都有子波参与。子波对最小二乘偏移的质量影响很大。

2.1 Sigsbee模型数据的子波影响分析

使用Sigsbee模型数据研究不同子波对最小二乘偏移结果的影响。由于Sigsbee模型的炮集数据的主频是20 Hz，因此选取了10、20、30 Hz三个主频的Ricker子波作为输入参数进行最小二乘偏移，分别低于、一致于、高于炮集数据的主频。三个子波的最小二乘偏移迭代次数均设为10次，以避免迭代次数对结果的影响，其结果如图 3所示。与单程波偏移结果相比，三个最小二乘偏移结果都提高了分辨率，能量更均衡，盐下照明也得到了补偿。其中图 3c的最小二乘偏移结果的能量最聚焦、畸变最小、成像质量最好。而图 3c选取的子波主频和炮集数据的主频一致。对不同主频子波偏移结果做频谱分析，结果如图 4所示，差异更加明显。由频谱对比可以看出：

（1）如果选取的子波主频低于炮集数据主频，迭代后低频损失，高频抬升，成像质量较差（图 4a）；

（2）如果选取的子波主频高于炮集数据主频，迭代后高频损失，低频抬升，成像质量也较差（图 4b）；

（3）如果选取的子波主频等于炮集数据主频，迭代后高、低频均抬升，频带得到拓宽，成像质量最好（图 4a和图 4b）。

2.2 实际数据的子波影响分析

实验选用国外M工区的地震资料。首先从实际数据中提取一个子波，其主频为35 Hz；然后使用该子波和主频为35 Hz的Ricker子波进行最小二乘偏移，迭代次数均为10次，结果如图 5所示，其频谱如图 6所示。主频为35Hz的Ricker子波的最小二乘偏移结果同相轴聚焦较差、剖面画弧现象较为严重（图 5a）；数据提取子波的最小二乘偏移结果成像更聚焦、画弧等偏移假象也得到较好压制、剖面的分辨率更高（图 5b）。从数据提取子波（红线）比主频为35 Hz的Ricker子波（紫线）的偏移结果频带更宽。综合以上两方面的分析对比，本文推荐使用从实际数据提取的子波做最小二乘偏移，而不是理论Ricker子波。

2.3 陆地与海上资料的最小二乘偏移结果对比

通过Sigsbee模型数据试验和实际生产数据试验可知，最小二乘偏移相对于常规偏移结果能够拓宽频带，提高成像分辨率和振幅的均衡性。

对于数据域的最小二乘偏移算法，根据最小二乘偏移理论，子波会参与最小二乘偏移的正演和反演过程，子波的选取对最小二乘偏移的效果至关重要。对于海上资料，子波比较稳定，可以提取到一个和实际数据非常符合的子波，最小二乘偏移效果较好。对于陆地资料，由于低降速带的影响，子波在空间上变化剧烈，想提取一个与整体数据相符的子波非常困难或不可能，只能在局部子波比较稳定处达到相对较好效果。

图 7为海上地震资料应用最小二乘偏移前、后的剖面对比，图 8为其频谱对比。可以发现，最小二乘偏移剖面相对于常规单程波偏移剖面频带更宽、分辨率更高，提高了振幅一致性。

图 9为陆地地震资料应用最小二乘偏移前、后的剖面对比，图 10为其频谱对比。虽然应用最小二乘偏移后拓宽了频带、提高了分辨率，如果只看红框内的局部，横向振幅较均衡。但由于提取子波只适用于局部，并不适应整个工区，横向整体均衡性并没有得到明显改善。所以，如果希望最小二乘偏移在陆地资料应用中取得理想的效果，就需在地震子波的一致性处理方面多下功夫，否则只能一定程度上拓宽频带、提高分辨率，但不会改善整体成像的振幅均衡性。

3 处理频带和数据量对偏移运算时间的影响

在实际项目的运作中发现，最小二乘偏移的运算时间与处理频带的宽度以及数据量都呈正相关关系。做生产项目时，通常会先选取一小块数据进行试验，以验证节点的可靠性和数据准备的正确性。如果知道最小二乘偏移的运算时间与处理频带以及数据量之间的关系，就能够通过试验数据大概估算后续整体的运算时间。由于最小二乘偏移是一个非常耗时的过程，运算时间的正确估算对项目运作的时间安排和计划有重要的实际意义，可以帮助项目组更好地把控项目的整体进度。

3.1 处理频带与最小二乘偏移用时的关系

如果处理频带选择得比较窄，运算时间会很短，有可能导致有效信号丢失，不能得到理想的偏移结果；如果处理频带过宽，虽然可以保证不会丢失有效信号，但会大大提高时间成本。另外，在实际资料的处理中，有时候需要快速得到一个试验结果验证并行节点是否运行正常、偏移数据的准备是否正确。为缩短程序的运算时间，通常会选择一个将主频包括在内、窄于实际资料频带的宽度，以便得到一个成像效果尚可且花费较少运算时间的结果。如果能知道频带宽度和程序运算时间的定量关系，就可以根据试验大概估算整体生产数据的处理时间。

试验数据是从国外三维陆地工区抽取的一小块数据，有效频带为2~50 Hz。试验数据规模为3212.9 MB，偏移频带为0-1-50-55 Hz，使用60个节点，每个节点24核，单节点内存为48 GB。最小二乘偏移的一次迭代各步骤用时统计如表 1所示。

使用同样的数据和同样的处理节点，将偏移频带改为0-1-85-95 Hz，再进行最小二乘偏移，其一次迭代的用时统计如表 2所示。

首先计算两次最小二乘偏移试验的运算时间比为116.39/64.15=1.814；再计算两次试验频带宽度比为95/55=1.7272。可见两次试验的用时比与频带宽度比基本相等，由此可知，最小二乘偏移的运算时间与处理频宽大致成正比。

3.2 数据量与最小二乘偏移用时的关系

使用上述国外三维陆地数据抽取规模为6160.31 MB的一块数据做试验。选择的处理频带为0-1-50-55 Hz，其他参数和试验条件保持不变，一次迭代用时统计如表 3所示。

首先计算两次最小二乘偏移的运算时间比为173.17/64.15=2.699；然后计算数据量比为6160.31/3212.91=1.917；两个比值关系为2.699∶1.917≈1.41≈$ \sqrt[]{2} $∶1。可见最小二乘偏移的运算时间与数据量的关系并不是简单的正比关系。如果数据量增大N倍，运算时间大概会增加$ \sqrt[]{2}N $倍，即数据量对偏移时间的影响更为显著。

4 结论

最小二乘偏移作为近几年发展的地震数据处理新技术，虽然计算成本巨大，但是其优势也是显而易见的，不仅能够拓宽频带，提高成像分辨率，还能提高成像的振幅均衡性和补偿盐下成像照明。

（1）结合模型数据和实际生产数据讨论了子波对最小二乘偏移结果的影响，即子波的频带和数据的频带一致时偏移结果最佳，在实际数据处理中最好从数据中提取子波。

（2）对比海上与陆地数据的成像结果发现，陆地资料子波的横向不稳定性降低了最小二乘偏移结果的横向整体均衡性。要想使最小二乘偏移在陆地资料中得到较好的应用效果，需多对子波一致性处理下功夫。

（3）通过试验总结了最小二乘偏移的运算时间与频带宽度以及数据量的经验关系，其中数据量对最小二乘偏移运算时间的影响更大，为实际生产中估算整个数据体的最小二乘偏移所需时间提供了参考。