肖建恩, 白英哲, 李 博, 万城程, 李仲庚

(1.中国石化 石油物探技术研究院，南京 211103;2.中国石化 河南油田分公司 勘探开发研究院，郑州 450000;3.青岛市海产博物馆，青岛 266000)

0 引言

地震偏移成像是地震勘探的关键环节之一，其主要目标是实现反射波的准确归位和恢复反射界面的反射系数、振幅变化以及波场特征。大量研究表明，地下介质中粘滞性广泛存在，采用各向同性成像算法来处理这些探区的数据，会导致成像位置不准确，成像振幅欠估计和分辨率降低等问题。因此，研究基于粘滞性的成像算法，对于获得高精度地震剖面具有重要意义。

叠前成像算法依据理论不同可以分为：①射线类偏移;②波动方程类偏移。射线类偏移基于几何射线理论来计算波场的走时、振幅等信息，进而实现波场延拓成像，灵活性和计算效率比较高。Krichhoff偏移是业界广泛应用的射线类成像算法，主要通过波动方程积分解来进行波场的反向延拓和成像，不仅具有面向目标成像能力，而且对复杂观测系统具有很好的适用性[1-3]。然而其依赖常规射线追踪方法来求取旅行时，存在焦散区、阴影区和多次波至问题。作为积分类偏移算法的改进，高斯束偏移不仅能够处理多次波至的问题，而且对低信噪比数据适用性更好。高斯束偏移算法采用彼此独立的射线束来表征地震波场，不仅继承了射线类偏移的灵活性和高效性，而且包含了波场的动力学信息，易于向各向异性、粘滞性和起伏地表等情况推广[4-6]。

波动方程类偏移通过求解地震波场传播过程中的传播算子，进而采用递归算法来实现偏移成像，不仅能够避免速度剧烈变化引起的焦散问题，而且能够处理多值走时的情况。其依据方程不同主要分为：①单程波偏移;②逆时偏移。单程波偏移通过解耦之后的单程波方程来获取延拓算子，并在此基础上进行延拓，最后采用合适的成像条件进行成像[7-8]。相比射线类偏移，单程波偏移有更高的成像精度，然而其计算效率相对较低，不仅不能进行真振幅成像，而且无法处理陡倾构造成像问题。逆时偏移则是将地表接收到的地震记录进行逆时延拓，进而和震源处正向传播波场进行互相关来求取成像值，不受成像角度的限制，然而其不仅计算量巨大，而且对速度场的精确性要求较高，在实用化方面还有一定的局限性[9-11]。

在Kirchhoff积分的基础上，Popov等[12]提出了采用高斯束加权积分来构建格林函数的逆时偏移算法，能够面向目标成像；黄建平等[13]结合Popov的思想，采用高斯束表征的格林函数来实现地震波场正反向延拓，进而通过互相关成像条件来求取成像值，实现了格林函数高斯束逆时偏移方法。针对地下介质中广泛存在的各向异性，张凯等[14]将弹性参数表征的各向异性射线追踪算法应用于高斯束逆时偏移，实现了角度域各向异性高斯束逆时偏移方法，然而该方法不仅计算效率比较低，同时需要确定地下介质的弹性参数。针对上述问题，肖建恩等[15-16]研究了基于相速度的各向异性声波和转换波高斯束逆时偏移算法，进一步提升了算法的计算效率和适用性。

1 基本原理

各向异性动态聚焦束逆时偏移的实现与高斯束逆时偏移类似，主要包括动态聚焦束算子的构建、格林函数的加权表征、互相关成像以及粘滞性的校正等。

1.1 格林函数的动态聚焦束表征

图1为二维中心射线坐标系，由初始点O(s0,0)到点R(sp,0)的传播矩阵为式(1)。

(1)

由互易原理可知，从点R(sp,0)到初始点O(s0,0)的传播矩阵可以表示为式(2)。

∏(sp,s0)=∏-1(s0,sp)=

(2)

图1 射线中心坐标系Fig.1 Central ray coordinate system

假定高斯射线束在P点聚焦，且宽度为l0，此时可以求得初始点位置处的动力学射线追踪参数：

(3)

引入中间变量：

(4)

此时束参数为式(5)。

(5)

其中：s0为初始点到P点的距离；v0为初始点处速度。通常称由式(5)决定束参数的高斯束为聚焦束。此时束参数的选取是静态的，仅取决于聚焦位置，不随射线路径而变化，不能对整条射线束的能量进行约束。Nowack等[17]进一步修改了复值束参数，构建了动态聚焦束算子，对整条射线束进行约束：

(6)

其中:ωr为参考频率；[p1(s),q1(s)]和[p2(s),q2(s)]分别为动力学射线追踪方程的平面波解和球面波解。

对于动态聚焦束逆时偏移，格林函数可以由动态聚焦束的加权积分来表征：

(7)

其中:θ为出射角;Φ(θ)为权值系数;Uθ(s,n)为动态聚焦束：

(8)

图2 粘声介质正演与偏移示意图Fig.2 Schematic diagram of visco-acoustic media forward and migration(a)正演；(b)偏移

图3 断块模型速度场和Q值场Fig.3 The velocity field and Q value field of fault block model(a)速度场；(b)Q值场

其中:p和q为复值动力学射线追踪参数，具体可以表示为：

(9)

其中:ε为式(6)表征的复值束参数。动态聚焦束本质上是先给出两组实初值来求解动力学射线追踪方程组，然后将计算出的解通过复系数进行线性组合来获取最终的动力学射线追踪参数。

1.2 动态聚焦束逆时偏移

记地表在0≤t≤T时刻内接收到地震波场为U0(x,t)，通过Kirchhoff积分可以求得地下偏移域内任一点x0在t0时刻的反向延拓波场：

(10)

其中:xs为震源点；∂Ω为闭合区间Ω边界；∂/∂nx为外法线方向导数；nx为闭合区间Ω的外法线；G为格林函数。采用动态聚焦束加权积分表征的格林函数(满足Kirchhoff近似：GGB(x,t;x0,t0)|z=0=0)替代式(10)中的格林函数，可以将反向延拓波场进一步简化为式(11)。

(11)

震源处正向传播波场也可以由动态聚焦束表征的格林函数来构建：

(12)

其中:fF(ω)为地震子波的傅里叶变换。

构建完正反向延拓波场，可以依据反射波成像准则，通过两者的互相关来获取成像结果，成像公式为式(13)。

(13)

1.3 粘滞性校正

图2为粘声介质中正演和偏移示意图，exp(-αL1)和exp(-αL2)分别为震源处到反射界面和反射界面到接收点的衰减项，S和R分别为震源处和检波点处的波场，有如下关系：

图4 断块模型试算结果Fig.4 The imaging results of the fault block model(a)声波动态聚焦束逆时偏移；(b)粘声高斯束逆时偏移；(c)本文方法

R(x,t)=S(x,t)exp(-αL1)exp(-αL2)

(14)

其中:α为衰减系数；L1和L2分别为震源点到反射界面和反射界面到检波点的路径。为了得到准确的成像值，需要通过补偿因子exp(αL1)和exp(αL2)对两个波场进行补偿。

地震波在粘滞性介质中传播，传播速度为复值，有如下形式：

(15)

其中:Q为品质因子；v0(x)为声波速度，其实部表征速度频散特性，虚部则体现衰减效应。对于式(15)中一阶项，当1/Q≪1时，射线传播路径保持不变，吸收衰减效应主要由与频率有关的复值旅行时体现：

(16)

其中:t(x)为声波旅行时;t*(x)为沿射线旅行时，具体可以表示为式(17)。

(17)

在粘声动态聚焦束逆时偏移算法中，采用式(16)所示的复值旅行时替换射线追踪的旅行时，并对t*(x)有关项取反号，即可实现粘滞性的校正。

1.4 算法流程

粘声介质动态聚焦束逆时偏移主要分为以下五步：

1)读入速度场、Q值场和炮记录。

2)通过运动学和动力学射线追踪方程，计算中心射线的旅行时、路径和振幅等信息。

3)根据式(6)构建动态聚焦束算子，构建格林函数，进而实现正反向波场的延拓 (式(11)、式(12))。

4)求取正反向延拓波场的互相关，得到单炮成像结果，并根据式(16)进行粘滞性校正。

5)将所有炮成像结果进行叠加，得到最终的成像结果。

2 模型试算

2.1 断块模型

首先采用断块模型对算法有效性进行验证。模型速度场和Q值场如图3所示，横纵向网格分别为1 000和550，间距分别为10 m和6 m。正演地震记录通过有限差分算法获得，共200炮，采用全接收的方式，炮间隔和道间隔分别为50 m和10 m，时间采样点为2 000，间隔为1 ms。

图5 Hess模型速度场和Q值场Fig.5 The velocity field and Q value field of Hess model(a)速度场；(b)Q值场

图6 Hess模型试算结果Fig.6 The imaging results of Hess model(a)声波动态聚焦束逆时偏移；(b)粘声高斯束逆时偏移；(c)本文方法

模型试算结果如图4所示。通过对比可以看出，由于忽略了粘滞性因素，在声波动态聚焦束逆时偏移结果中，深层能量较弱，同相轴聚焦性也比较差(如图4(a)中红色剪头所示)。而在粘声高斯束逆时偏移(图4(b))和动态聚焦束逆时偏移(图4(c))结果中，深层能量得到了比较好地恢复，同相轴也更加聚焦和连续。相比于粘声高斯束逆时偏移，本文方法结果中同相轴的能量分布更为均匀，中深层构造成像质量得到了有效提高，反射界面的刻画也更为清晰，验证了方法的有效性。

2.2 Hess模型

通过Hess模型验证算法的适用性。图5为模型的速度场和Q值场，横向宽度和纵向深度分别为15 km和3.5 km，网格间距均为10 m。合成地震数据共251炮，每炮1501道接收，炮间隔和道间隔分别为60 m和10 m。地震记录时间采样长度为3.0 s，间隔为1 ms。

模型试算结果如图6所示。可以看出在声波聚焦束逆时偏移结果中，由于粘滞性的影响，导致深层能量较弱，高速岩丘下方反射层的同相轴连续性和聚焦性也比较差(如图6(a)红色箭头所示)。在粘声高斯束逆时偏移(图6(b))和本文方法(图6(c))的结果中，深层能量得到了有效的恢复，同相轴的连续性和聚焦性也得到了提升。但相比于粘声高斯束逆时偏移，本文方法结果中同相轴更加聚焦和连续，能量分步也更为均匀，中深层构造的能量聚焦性更为理想，整体剖面更为清晰(如图6(b)和图6(c)中红色箭头所示)。

3 结论

高斯束逆时偏移是一种兼具计算效率和成像精度的深度域成像算法，具有面向目标成像的能力。笔者首先通过构建动态聚焦束算子，并将其应用在高斯束积分中，表征出动态聚焦束，进而计算出动态聚焦束表征的地震波场，解决了高斯射线束宽度随传播距离增加而快速发散的问题；随后引入与频率有关的复值旅行时，分别对炮点和检波点波场的吸收衰减效应进行补偿。模型试算表明,本文算法不仅可以将地震波能量约束在一定范围内，有效提高中深层构造成像质量,而且可以对地层的吸收衰减效应进行有效补偿，对反射界面进行更为清晰的刻画。下一步拟将该方法推广至各向异性介质。