CN104360381A - 一种地震资料的偏移成像处理方法 - Google Patents

Abstract

本发明涉及地震资料的偏移成像处理方法，野外采集得到炮记录后，对其进行平面波编码得到若干不同射线参数的平面波记录，使用逆时偏移算子对平面波记录进行偏移成像，并在最小二乘反演框架下，通过共轭梯度法使得成像结果向反射率模型收敛，最终得到信噪比更高、振幅属性更为可靠的成像结果，本发明将地震成像看做最小二乘反演问题，能够较好地抑制逆时偏移算法中的低频噪音，减弱稀疏采集造成的震源效应，补偿深层反射同相轴能量，大大改善了常规逆时偏移的成像质量；此外通过平面波编码压缩数据克服了最小二乘偏移计算效率的缺陷，是一种有效的高精度成像方法。

Description

一种地震资料的偏移成像处理方法

技术领域

本发明涉及地震资料处理方法，具体是地震资料的偏移成像处理方法。

背景技术

近年来，基于双程波动方程的逆时偏移(RTM)方法在地震成像领域吸引了越来越多地球物理学家的关注。逆时偏移方法相对于实际油气勘探过程中被广泛应用的Kirchhoff偏移方法和单程波波动方程偏移方法，具有原理简单、成像精确且无倾角限制、适应复杂速度模型等优点。但是RTM方法在成像时会产生成像噪音，在一定程度上湮没了地下构造。虽然国内外一些研究学者应用Laplace滤波或者其改进方法能够一定程度上压制RTM成像产生的低频噪音，但是成像结果的保幅性及地下振幅的均衡性得不到保证。这些弱点使得逆时偏移难以满足高精度成像以及岩性成像的需要。

而基于反演思想的最小二乘偏移(LSM)是保幅成像的最佳方法。但是，将最小二乘偏移看做数据空间的反演问题，其计算量约为常规偏移的2×niter倍(niter为迭代次数)。尤其是当目标区域较大时，LSM的计算成本无法接受。

发明内容

本发明目的在于提供一种地震资料的偏移成像处理方法，为叠前平面波最小二乘逆时偏移(P-LSRTM)方法，将地震成像看做最小二乘反演问题，通过共轭梯度法使成像结果向真实反射率模型收敛，大大改善了常规RTM的成像质量；通过平面波编码压缩偏移数据体，提高了计算效率，实现了高精度地震成像的方法。

本发明包括以下步骤：

(1)采用叠加前单炮数据，输入初始反射率模型、速度场和平面波编码函数；

(2)对炮记录进行平面波编码：

$d_{\mathrm{obs}}(x_g,t;p)=\underset{x_s}{\mathrm{\Σ}}{d}_{\mathrm{obs}}(x_g,t;x_s)*\mathrm{\δ}(t-p{x}_s)$     (式1)

其中，d_obs(x_g，t；x_s)为观测炮记录，δ(t-px_s)为编码函数，编码过程等价于通过tau-p变换合成平面波记录，时移量px_s随震源点位置x_s线性变化，p为射线参数：

$p=\frac{\sin \mathrm{\θ}}{v}$     (式2)

其中，θ是地表入射角，v是地表速度，

(3)通过反偏移算子计算当前反射率模型下的平面波记录，反偏移过程表示为：

d_cal(p，x)＝ω²∫m(x’)W(p，ω)G₀(x’；x_p)G₀(x；x’)dx’    (式3)

其中，m(x’)代表反射系数，W(p，ω)为编码后的平面波震源，G₀(x’；x_p)为从平面波震源到各散射点的格林函数，G₀(x；x’)为从各散射点到接收点的格林函数：

(4)计算当前反射率模型下不同射线参数的平面波记录与真实记录的残差；

(5)通过逆时偏移算子反传残差计算更新梯度(共轭梯度法)：

g^(k+1)＝L^T[Lm^(k)-d_obs]^T

$\mathrm{\β}=\frac{g^{(k+1)}{g}^{(k+1)}}{g^{\left(k\right)}{g}^{\left(k\right)}}$     (式4)

z^(k+1)＝g^(k+1)+βz^(k)

其中，L为反偏移算子的矩阵形式，L^T为逆时偏移算子的矩阵形式，

(6)计算迭代更新的步长：

$\mathrm{\α}=\frac{{\left[z^{(k+1)}\right]}^T{g}^{(k+1)}}{{\left[{\mathrm{Lz}}^{(k+1)}\right]}^T{\mathrm{Lz}}^{(k+1)}}$     (式5)

(7)更新反射率模型m^(k+1)＝m^(k)-αz^(k+1)；

(8)重复步骤3到7，直到当前的反射率模型满足精度；

(9)叠加不同射线参数平面波记录的成像结果，输出当前模型。

本发明有益效果实现了时间域叠前平面波最小二乘逆时偏移方法，能够有效地克服逆时偏移的低频偏移噪声，减弱稀疏采集造成的震源效应，补偿深层反射同相轴能量，是一种保幅的成像方法；将平面波编码应用到最小二乘偏移中，直接对平面波记录偏移，提高了计算效率。本发明具有保幅成像、稳定性高的优点。

附图说明

图1为胜利探区断块速度模型；

图2为通过式1合成的不同射线参数的平面波道集，其中2.(a)～(c)对应的射线参数分别为：p＝-0.1742ms/m、p＝0ms/m、p＝0.1742ms/m；

图3为应用本发明对合成的24个平面波记录偏移后的成像剖面叠加结果，其中3.(a)～(f)分别为迭代5、10、15、20、25、30次的成像结果；

图4为对图3中成像剖面Laplace滤波后的结果；

图5为残差收敛曲线；

图6为几种成像方法结果的对比，其中6.(a)～(c)分别为RTM、LSRTM、P-LSRTM方法成像结果的对比。

具体实施方式

为使本发明的上述和其他目的、特征和优点能更明显易懂，下文特举出较佳实施例，并配合所附图式，作详细说明如下。

(1)对胜利探区断块速度模型进行正演模拟，其中模型参数为：水平方向大小为5.7km，深度方向为2.0km，网格间距10m；计算参数为：震源为主频30Hz的雷克子波，网格间距10m，时间采样间隔0.5ms，采集时长2s，使用时间2阶空间8阶有限差分正演模拟；观测系统的设计：总炮数571炮，571道检波器接收，炮间距和道间距均为10m；

(2)通过式1对由1)得到的炮记录进行平面波编码，射线参数p的取值范围为-0.1742ms/m～0.1742ms/m，呈线性变化，共合成24个平面波记录，图2展示了3个不同射线参数的平面波记录；

(3)通过式3对当前反射率模型进行反偏移(初始模型的值为0，故第一次迭代时模型的平面波记录也为0)；

(4)计算当前反射率模型下的平面波记录与通过编码合成的真实记录的残差；

(5)通过逆时偏移算子反传残差计算更新梯度，根据式4计算共轭梯度法修正后的梯度；

(6)根据式5计算更新步长；

(7)更新反射率模型m^(k+1)＝m^(k)-αz^(k+1)；

(8)重复步骤3到7，直到当前的反射率模型满足精度；

(9)输出每次迭代后的叠加成像结果，如图3所示。

由图3可知：在残差收敛过程中，由成像产生的低频噪音和串扰噪音得到压制；反射同相轴能量增强；深部反射同相轴能量均衡性变好；但对成像结果右侧的强低频噪音不能完全压制，这可能是模型右侧的平层反射系数过大造成的。对于常规平面波逆时偏移，需要90个平面波记录偏移结果的叠加才能完全压制串扰噪音，但是最小二乘偏移本身具有压制噪音的作用，因此本文只叠加了24个平面波记录，串扰噪音已基本压制。

由图4可以知：随着迭代次数增加，可以明显看到红色圆圈中所标示的串扰噪音逐渐收敛；在迭代5次的结果中，红色箭头所指的小断块底部没有成像出来，在迭代30次后，该同相轴已正确成像，且同相轴能量均衡，成像结果精度高。

由图6几种偏移方法的对比可知：RTM成像剖面中含有较强的低频偏移噪音，在一定程度上掩盖了地下的真实构造，尤其是深层构造难以辨认；LSRTM成像质量最高，低频噪音得到压制，深层反射轴能量得到补偿，高陡构造和小断块都准确地成像；P-LSRTM相对于RTM成像质量有了较大改善，主要表现在：成像噪音减弱、剖面均衡性变好、深层反射轴能量得到补偿、小断块形态和分界清晰，其成像质量接近LSRTM，但P-LSRTM的计算时间仅为LSRTM的0.26倍，常规RTM的7.8倍。

Claims (1)

1.一种地震资料的偏移成像处理方法，其特征在于，该地震资料的偏移成像处理方法包括以下步骤：

(1)采用采用叠加前单炮数据输入初始反射率模型、速度场和平面波编码函数；

(2)对炮记录进行平面波编码：

(3)通过反偏移算子计算当前反射率模型下的平面波记录，

(4)计算当前反射率模型下不同射线参数的平面波记录与真实记录的残差；

(5)通过逆时偏移算子反传残差计算更新梯度

(6)计算迭代更新的步长：

(7)更新反射率模型

(8)重复步骤3到7，直到当前的反射率模型满足精度；

(9)叠加不同射线参数平面波记录的成像结果，输出当前模型。