CN104937441A - 一种二维迴折射线积分法叠前深度偏移方法 - Google Patents

Abstract

一种二维迴折射线积分法叠前深度偏移方法，步骤为1)采用地震方法得到成像空间对应的深度域速度场；2)根据速度场划分矩形网格；3)在已经划分的矩形网格上进行射线追踪，而且对射线角度大于90度的射线继续追踪；记录射线在空间位置走过的时间，并将射线没有到达的地方进行插值得到充满空间的旅行时表；4)应用旅行时对叠前炮集数据进行叠前深度偏移处理，得到共成像点道集；5)将共成像点道集数据进行叠加得到该工区深度域的构造成像结果，并将结果显示为地层剖面图像。该方法射线追踪的精度高，为迴折射线旅行时产生创造条件，提高旅行时计算精度。

Description

一种二维迴折射线积分法叠前深度偏移方法

技术领域

本发明涉及反射波地震数据处理过程中的二维叠前深度偏移构造成像技术，具体的说是一种精确的二维迴折射线积分法叠前深度偏移构造成像方法。

背景技术

随着油田勘探程度的提高，复杂地区和复杂油气藏的油气勘探逐步成为勘探的主要目标。根据地表复杂度和地下构造复杂度可以分为以下四类：地表简单地下构造简单的情况用叠后时间偏移就可以得到很好的成像效果；地表复杂地下构造简单的情况用叠前时间偏移可以得到很好的成像效果；地表简单地下构造复杂的情况用叠后深度偏移可以得到很好的成像效果；地表复杂地下构造复杂的情况用叠前深度偏移可以得到很好的成像效果。在国内西北部有许多地区既是复杂地表又是地下构造复杂的情况，对于这种情况必须利用叠前深度偏移加以解决。

叠前深度偏移单程波方法可以分为三大类：积分法、差分法和傅立叶变换法。三类方法各有优缺点，其中积分法的优点是计算效率高、输入输出灵活和陡倾角成像效果好，缺点是整体成像信噪比不如其它两类方法高；差分法的优点是成像精度高，缺点是计算效率低；傅立叶变换法介于以上两种方法之间。需要注意的是叠前深度偏移的计算量是非常大的，大的非常可怕。如果计算机的运算速度足够快，差分法叠前深度偏移是最好的选择，然而目前计算机运算速度还不足以对大块的三维工区进行差分法叠前深度偏移，目前，国内外进行大块工区叠前深度偏移的绝对主力是积分法，至少在未来的几年内是不会改变的。

常规的单程波叠前深度偏移方法能够得到地下复杂构造较好的成像结果，这里常规的单程波叠前深度偏移方法指的是利用下行波场外推成像。但对于逆掩褶皱、逆冲断层和倒转盐体等非常复杂的构造，侧翼成像还存在问题，尤其是大于90度的侧翼部分不能很好的成像。解决这类问题需要利用迴折波成像，在高频近似情况下为迴折射线。所谓的迴折波是没有经过任何地质体的反射而波变换传播方向，由向下传播变为向上传播再到达地表的波(如图1(a)所示)。正是迴折波可以到达盐丘侧翼大于90度的部分，再反射被地面接受，如果能很好的利用迴折波，倒转盐体侧翼大于90度部分的构造成像就不再是问题了。改善上述特定地质体的构造成像质量，需要对常规单程波叠前深度偏移方法进行改进。

将单程波有限差分法和傅立叶变换法改造成带有耦合项单程波方程算法可以对构造倾角大于90度的构造成像(主要针对盐丘侧翼)，当然这种方法的计算量是非常大的，对于大块三维数据进行叠前深度偏移还不现实。如果能利用积分法叠前深度偏移计算效率高的特点，将其改造成具有迴折波成像能力的技术对复杂区大规模叠前深度偏移的推广应用有非常大的意义。目前，在国内还没有出现波积分法对构造倾角大于90度的构造进行成像的技术，所以开展积分法迴折射线叠前深度偏移研究势在必行，本发明在二维积分法叠前深度偏移的基础上，开发出具有迴折波成像能力的积分法叠前深度偏移算法。

发明内容

本发明提供一种二维迴折射线积分法叠前深度偏移方法。

具体具体步骤包括：

1)采用地震的方法取得成像空间对应的深度域速度场；

2)根据速度场划分矩形网格；

3)在已经划分的矩形网格上进行射线追踪，而且对射线角度大于90度的射线继续追踪；记录射线在空间位置走过的时间，并将射线没有到达的地方进行插值得到充满空间的旅行时；

4)应用旅行时对叠前炮集数据进行叠前深度偏移处理，得到共成像点道集；

5)将共成像点道集数据进行叠加即可得到该工区深度域的构造成像结果，并将结果显示为地层剖面图像。

本发明的二维迴折射线积分法叠前深度偏移方法，其具体实施方式为：

利用下面所示的射线追踪公式在已经划分的矩形网格上进行射线追踪：

$x(\mathrm{\τ}+\mathrm{\Δ\τ})=x\left(\mathrm{\τ}\right)+\frac{\mathrm{\Δ\τ}}{6}(K_{x1}+2K_{x2}+2K_{x3}+K_{x4})$

$z(\mathrm{\τ}+\mathrm{\Δ\τ})=z\left(\mathrm{\τ}\right)+\frac{\mathrm{\Δ\τ}}{6}(K_{z1}+2K_{z2}+2K_{z3}+K_{z4})$

$P_x(\mathrm{\τ}+\mathrm{\Δ\τ})=P_x\left(\mathrm{\τ}\right)+\frac{\mathrm{\Δ\τ}}{6}(L_{x1}+2L_{x2}+2L_{x3}+L_{x4})$

$P_z(\mathrm{\τ}+\mathrm{\Δ\τ})=P_z\left(\mathrm{\τ}\right)+\frac{\mathrm{\Δ\τ}}{6}(L_{z1}+2L_{z2}+2L_{z3}+L_{z4})$

K_x1＝f(τ，x(τ)，P_x(τ))K_z1＝f(τ，x(τ)，P_z(τ))

L_x1＝f(τ，x(τ)，P_x(τ))L_z1＝f(τ，x(τ)，P_z(τ))

$K_{x2}=f(\mathrm{\τ}+\frac{\mathrm{\Δ\τ}}{2},x\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{K}_{x1},P_x\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{L}_{x1})$

$K_{z2}=f(\mathrm{\τ}+\frac{\mathrm{\Δ\τ}}{2},z\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{K}_{z1},P_z\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{L}_{z1})$

$L_{x2}=g(\mathrm{\τ}+\frac{\mathrm{\Δ\τ}}{2},x\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{K}_{x1},P_x\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{L}_{x1})$

$L_{z2}=g(\mathrm{\τ}+\frac{\mathrm{\Δ\τ}}{2},z\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{K}_{z1},P_z\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{L}_{z1})$

$K_{x3}=f(\mathrm{\τ}+\frac{\mathrm{\Δ\τ}}{2},x\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{K}_{x2},P_x\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{L}_{x2})$

$K_{z3}=f(\mathrm{\τ}+\frac{\mathrm{\Δ\τ}}{2},z\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{K}_{z2},P_z\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{L}_{z2})$

$L_{x3}=g(\mathrm{\τ}+\frac{\mathrm{\Δ\τ}}{2},x\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{K}_{x2},P_x\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{L}_{x2})$

$L_{z3}=g(\mathrm{\τ}+\frac{\mathrm{\Δ\τ}}{2},z\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{K}_{z2},P_z\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{L}_{z2})$

K_x4＝f(τ+Δτ，x(τ)+ΔτK_x3，P_x(τ)+ΔτL_x3)

K_z4＝f(τ+Δτ，z(τ)+ΔτK_z3，P_z(τ)+ΔτL_z3)

L_x4＝g(τ+Δτ，x(τ)+ΔτK_x3，P_x(τ)+ΔτL_x3)

L_z4＝g(τ+Δτ，z(τ)+ΔτK_z3，P_z(τ)+ΔτL_z3)

其中，其中x、z是位置坐标分量，v是速度矢量，关于x、z的函数。P_x、P_z是慢度矢量在X、Z方向的分量，也是关于X、Z的函数。

本发明能更好的参与旅行时空间插值，采用高阶龙格库塔公式使射线追踪的精度更高，同时在单条射线追踪的过程中对射线方向不做限制，使迴折射线能很好的产生，本发明在较大的射线阴影区不进行旅行时插值的方式代替阴影区填充的方式以提高旅行时计算精度。

附图说明

附图1表示迴折射线形成以及迴折射线反射示意图；

附图2表示V(z)介质纵波速度模型，即V(Z)介质的速度场；

附图3表示V(z)介质中某一点的射线追踪结果图；

附图4表示脉冲响应偏移测试结果，其中，附图4(a)为常规积分法脉冲响应，(b)为迴折波偏移的脉冲响应测试结果；

附图5表示典型的倒转盐体模型纵波速度场；

附图6表示常规积分法叠前深度偏移结果；

附图7表示迴折射线积分法叠前深度偏移结果；

附图8表示常规积分法和迴折射线偏移结果之差；

附图9表示某地区实际资料的积分法叠前深度偏移结果，其中附图9(a)表示常规积分法偏移结果，(b)表示迴折射线偏移结果。

具体实施方式

如图1所示，本发明的二维迴折射线叠前深度偏移方法基于射线方程，通过四阶龙格库塔积分公式进行求解，得到每一条射线的射线路径和旅行时。传统的射线追踪方法是逐个波前面向前追踪，并对每一条射线的传播方向做了不能超过90度的限制。本发明所采用的方法和传统的射线追踪不同之处在于先对全部射线进行追踪然后再进行旅行时插值，而且射线追踪的过程中对每一条射线的传播方向不作明确的限制。本发明这样做的好处在于能实现迴折射线追踪并可以进行迴折射线叠前深度偏移。

实例1：

脉冲响应测试是偏移测试必不可少的环节，如图2所示为V(Z)介质的速度场，图3为某一点的射线追踪图，图4为脉冲响应测试结果，其中(a)为常规积分法脉冲响应，由于不能对迴折波成像，因此成像角度最大能达到90度。(b)为迴折波偏移的脉冲响应测试结果，可以看到对大于90度的构造也能成像。

实例2：

如图5所示为倒转盐体模型纵波速度场，图6为该模型的常规积分法叠前深度偏移的结果，可以看到盐体侧翼大于90度的构造不能很好的成像。图7所示为该模型迴折射线积分法叠前深度偏移结果，相对于常规偏移方法，迴折射线积分法叠前深度偏移方法能更好的对大于90度的构造进行成像，充分体现了该方法的优越性。如图8所示为以上两个偏移结果之差，可以看到迴折射线对盐体侧翼偏移的贡献非常大，对大倾角成像是一个很好的补充。

发明实施例3：

如图9所示为某地区实际资料的积分法深度偏移结果，其中，(a)为该常规积分法叠前深度偏移结果，(b)为迴折射线积分法叠前深度偏移结果。在盐丘侧翼迴折射线积分法叠前深度偏移的结果要好一些，但盐丘围盐构造过于复杂使速度建模非常空难，因而迴折射线的结果不是非常理想。

结合图例说明具体实施方式：

本发明的二维迴折射线积分法叠前深度偏移方法，其具体实施方式为：

1)读取成像空间对应的深度域速度场；

2)根据速度场划分矩形网格(如图2所示)；

3)利用下面所示的射线追踪公式在已经划分的矩形网格上进行射线追踪，而且对射线角度大于90度的射线继续追踪(如图3所示)；记录射线在空间位置走过的时间，并将射线没有到达的地方进行插值得到充满空间的旅行时表。

$x(\mathrm{\τ}+\mathrm{\Δ\τ})=x\left(\mathrm{\τ}\right)+\frac{\mathrm{\Δ\τ}}{6}(K_{x1}+2K_{x2}+2K_{x3}+K_{x4})$

$z(\mathrm{\τ}+\mathrm{\Δ\τ})=z\left(\mathrm{\τ}\right)+\frac{\mathrm{\Δ\τ}}{6}(K_{z1}+2K_{z2}+2K_{z3}+K_{z4})$

$P_x(\mathrm{\τ}+\mathrm{\Δ\τ})=P_x\left(\mathrm{\τ}\right)+\frac{\mathrm{\Δ\τ}}{6}(L_{x1}+2L_{x2}+2L_{x3}+L_{x4})$

$P_z(\mathrm{\τ}+\mathrm{\Δ\τ})=P_z\left(\mathrm{\τ}\right)+\frac{\mathrm{\Δ\τ}}{6}(L_{z1}+2L_{z2}+2L_{z3}+L_{z4})$

K_x1＝f(τ,x(τ),P_x(τ))K_z1＝f(τ,x(τ),P_z(τ))

L_x1＝f(τ,x(τ),P_x(τ))L_z1＝f(τ,x(τ),P_z(τ))

$K_{x2}=f(\mathrm{\τ}+\frac{\mathrm{\Δ\τ}}{2},x\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{K}_{x1},P_x\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{L}_{x1})$

$K_{z2}=f(\mathrm{\τ}+\frac{\mathrm{\Δ\τ}}{2},z\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{K}_{z1},P_z\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{L}_{z1})$

$L_{x2}=g(\mathrm{\τ}+\frac{\mathrm{\Δ\τ}}{2},x\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{K}_{x1},P_x\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{L}_{x1})$

$L_{z2}=g(\mathrm{\τ}+\frac{\mathrm{\Δ\τ}}{2},z\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{K}_{z1},P_z\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{L}_{z1})$

$K_{x3}=f(\mathrm{\τ}+\frac{\mathrm{\Δ\τ}}{2},x\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{K}_{x2},P_x\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{L}_{x2})$

$K_{z3}=f(\mathrm{\τ}+\frac{\mathrm{\Δ\τ}}{2},z\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{K}_{z2},P_z\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{L}_{z2})$

$L_{x3}=g(\mathrm{\τ}+\frac{\mathrm{\Δ\τ}}{2},x\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{K}_{x2},P_x\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{L}_{x2})$

$L_{z3}=g(\mathrm{\τ}+\frac{\mathrm{\Δ\τ}}{2},z\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{K}_{z2},P_z\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{L}_{z2})$

K_x4＝f(τ+Δτ,x(τ)+ΔτK_x3,P_x(τ)+ΔτL_x3)

K_z4＝f(τ+Δτ,z(τ)+ΔτK_z3,P_z(τ)+ΔτL_z3)

L_x4＝g(τ+Δτ,x(τ)+ΔτK_x3,P_x(τ)+ΔτL_x3)

L_z4＝g(τ+Δτ,z(τ)+ΔτK_z3,P_z(τ)+ΔτL_z3)

其中，其中x、z是位置坐标分量，v是速度矢量，关于x、z的函数。

P_x、P_z是慢度矢量在X、Z方向的分量，也是关于X、Z的函数。

4)读取叠前炮集地震数据；

5)应用旅行时对叠前炮集数据进行叠前深度偏移处理，得到共成像点道集；

6)将共成像点道集数据进行叠加即可得到该工区深度域的构造成像结果，并将结果显示为地层剖面图像(如图4所示)。

本发明可以对高陡构造，甚至构造倾角大于90度的倒转盐体进行成像，且具有如下特点：

(1)本发明在旅行时计算过程中采用先追踪后插值的方式，使产生的迴折射线能更好的参与旅行时空间插值，而射线追踪过程采用高阶龙格库塔公式进行，使射线追踪的精度更高，同时在单条射线追踪的过程中对射线方向不做限制，使迴折射线能很好的产生。

(2)本发明在射线的旅行时插值方面采用与众不同的处理方式，在多条交叉射线旅行时插值方面采用加权插值方式代替传统的相邻射线插值方式为迴折射线旅行时产生创造条件，在较大的射线阴影区不进行旅行时插值的方式代替阴影区填充的方式以提高旅行时计算精度。

Claims (4)

1.一种二维迴折射线积分法叠前深度偏移方法，其特征在于：

1)采用地震方法得到成像空间对应的深度域速度场；

2)根据速度场划分矩形网格；

3)在已经划分的矩形网格上进行射线追踪，而且对射线角度大于90度的射线继续追踪；记录射线在空间位置走过的时间，并将射线没有到达的地方进行插值得到充满空间的旅行时表；

4)应用旅行时对叠前炮集数据进行叠前深度偏移处理，得到共成像点道集；

5)将共成像点道集数据进行叠加得到该工区深度域的构造成像结果，并将结果显示为地层剖面图像。

2.根据权利要求1所述的一种二维迴折射线积分法叠前深度偏移方法，其特征在于：步骤3)所述的射线追踪按照以下公式计算：

$x(\mathrm{\τ}+\mathrm{\Δ\τ})=x\left(\mathrm{\τ}\right)+\frac{\mathrm{\Δ\τ}}{6}(K_{x1}+2K_{x2}+2K_{x3}+K_{x4})$

$z(\mathrm{\τ}+\mathrm{\Δ\τ})=z\left(\mathrm{\τ}\right)+\frac{\mathrm{\Δ\τ}}{6}(K_{z1}+2K_{z2}+2K_{z3}+K_{z4})$

$P_x(\mathrm{\τ}+\mathrm{\Δ\τ})=P_x\left(\mathrm{\τ}\right)+\frac{\mathrm{\Δ\τ}}{6}(L_{x1}+2L_{x2}+2L_{x3}+L_{x4})$

$P_z(\mathrm{\τ}+\mathrm{\Δ\τ})=P_z\left(\mathrm{\τ}\right)+\frac{\mathrm{\Δ\τ}}{6}(L_{z1}+2L_{z2}+2L_{z3}+L_{z4})$

K_x1＝f(τ,x(τ),P_x(τ)) K_z1＝f(τ,x(τ),P_z(τ))

L_x1＝f(τ,x(τ),P_x(τ)) L_z1＝f(τ,x(τ),P_z(τ))

$K_{x2}=f(\mathrm{\τ}+\frac{\mathrm{\Δ\τ}}{2},x\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{K}_{x1},P_x\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{L}_{x1})$

$K_{z2}=f(\mathrm{\τ}+\frac{\mathrm{\Δ\τ}}{2},z\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{K}_{z1},P_z\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{L}_{z1})$

$L_{x2}=g(\mathrm{\τ}+\frac{\mathrm{\Δ\τ}}{2},x\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{K}_{x1},P_x\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{L}_{x1})$

$L_{z2}=g(\mathrm{\τ}+\frac{\mathrm{\Δ\τ}}{2},z\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{K}_{z1},P_z\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{L}_{z1})$

$K_{x3}=f(\mathrm{\τ}+\frac{\mathrm{\Δ\τ}}{2},x\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{K}_{x2},P_x\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{L}_{x2})$

$K_{z3}=f(\mathrm{\τ}+\frac{\mathrm{\Δ\τ}}{2},z\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{K}_{z2},P_z\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{L}_{z2})$

$L_{x3}=g(\mathrm{\τ}+\frac{\mathrm{\Δ\τ}}{2},x\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{K}_{x2},P_x\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{L}_{x2})$

$L_{z3}=g(\mathrm{\τ}+\frac{\mathrm{\Δ\τ}}{2},z\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{K}_{z2},P_z\left(\mathrm{\τ}\right)+\frac{1}{2}\mathrm{\Δ\τ}{L}_{z2})$

K_x4＝f(τ+Δτ,x(τ)+ΔτK_x3,P_x(τ)+ΔτL_x3)

K_z4＝f(τ+Δτ,z(τ)+ΔτK_z3,P_z(τ)+ΔτL_z3)

L_x4＝g(τ+Δτ,x(τ)+ΔτK_x3,P_x(τ)+ΔτL_x3)

L_z4＝g(τ+Δτ,z(τ)+ΔτK_z3,P_z(τ)+ΔτL_z3)

其中，其中x、z是位置坐标分量，v是速度矢量，关于x、z的函数。

P_x、P_z是慢度矢量在X、Z方向的分量，也是关于X、Z的函数。

3.根据权利要求1所述的一种二维迴折射线积分法叠前深度偏移方法，其特征在于：步骤3)所述的射线追踪过程中不对射线的方向做任何限制。

4.根据权利要求1所述的一种二维迴折射线积分法叠前深度偏移方法，其特征在于：步骤3)所述的已追踪完的射线沿射线计算相应的旅行时，对射线重叠的区域进行旅行时加权，使上行的迴折射线参与旅行时插值，得到迴折射线的旅行时。