波动方程偏移方法

射线偏移是一种近似的几何偏移，虽然地震波的运动学特点得以恢复，但波的动力学特点(如振幅、波形、相位等)却受到畸变，因此，射线偏移已逐渐被高精度的波动方程偏移所代替。波动方程偏移是以波动理论为基础的偏移处理方法，团液轿其基本思路是，当地表产生d性波向下传播(称为下行波)，遇到反射界面时将产生反射，这时可将反射界面看作新的波源，又有新的波以波动理论向上传播(称为上行波)，在地表接收到的地震记录就可看作反射界面产生的波场效应。偏移就是将地表接收到的波场按波动方程的传播规律反向向下传播，通常称为波场反向延拓，当波场反向延拓到反射界面时成像(成像剖面为偏移剖面)，从而找到了真实反射界面，达到了偏移处理的目的。可见波动方程偏移主要由波场延拓和成像两部分组成。波场延拓可用多种不同的方法实现，随之形成了多种不同的波动方程偏移方法。成像也有成像的原理，叠前和叠后偏移各有不同的成像条件。

3.4.3.1 波动方程偏移的成像原理

波动方程成像原理分叠后偏移成像原理和叠前偏移成像原理。

3.4.3.1.1 爆炸反射界面成像原理

该原理属叠后偏移成像原理。叠加剖面相当自激自收剖面，若将剖面中时间除2，或将传播速度减一半，就可将自激自收剖面看作在反射界埋明面上同时激发的地震波沿界面法线传播到地表所接收的记录，即可将界面看作爆炸源，称为爆炸反射界面。若用波动方程将地表接收的波场(叠加剖面)作反时间方向传播(向下延拓)，当波场延拓到时间t为零(t=0)时，该波场的所在位置就是反射界面位置。因此，t=0成为叠后波动方程偏移的成像条件。从延拓的结果(地下各点的波场)中取出地下各点处零时刻的波场值组成的剖面就为成像剖面，该剖面为叠后波动方程偏移结果。

3.4.3.1.2 波场延拓的时间一致性成像原理

图3-22 时间一致性成像原理示意图

时间一致性成像原理适用于叠前偏移。此成像原理可描述塌肆为：在地下某一深度存在一反射界面R(如图3-22(a))，在地面S点激发的下行波D到达界面R时产生反射上行波U，到达G点被接收，下行波D到达界R面的时间(或空间位置)与上行波U产生的时间(或空间位置)是一致的，即称为时间(或空间位置)一致性。设波从S点到R的传播时间为ts，从R至G的传播时间为tg，从S到G的总时间为tsg=ts+tg。在叠前偏移中，若模拟一震源函数D自S点正向(向下)延拓，而将G点接收到的上行波U反向延拓，当D和U延拓深度为Z1时，D的正向传播时间和U的反向传播时间分别为ts1和tg1，因Z1＜ZR(ZR为反射点深度)，tsg-tg1＞ts1，说明上行波和下行波所在的时间(或空间位置)不一致(如图3-22(b))，当D和U延拓深度为zz=ZR时，下行波正向传播时间为ts1=ts，上行波反向传播时间为tg2=tg，即有tsg-tg2=ts2，或tsg-tg=ts，这时上、下行波所在的时间(或空间位置)是一致的。再将D、U延拓到Z3，Z3＞ZR，即当延拓深度Z＞ZR以后，不会再出现时间(或深度位置)一致的现象。在上、下行波延拓过程中，若求下行波场D和上行波场U的零移位互相关，在满足时间(或空间位置)一致性条件时，相关值最大，而在其他情况下相关值很小或为零，延拓过程中的相关结果就为叠前偏移成像剖面。

3.4.3.2 叠后波动方程偏移方法

叠后偏移是在叠加剖面的基础上进行偏移处理。叠后波动方程偏移是用某些数学手段求解波动方程，对叠后波场延拓归位，达到偏移的目的。针对求解波动方程的方法，可将波动方程偏移分为三大类主要方法：有限差分法波动方程偏移、F－K域波动方程偏移和克希霍积分法波动方程偏移。

3.4.3.2.1 15°有限差分法波动方程偏移

15°有限差分法波动方程偏移是以地面上获得的水平叠加时间剖面作为边界条件，用差分代替微分，对只包含上行波的近似波动方程求解以得到地下界面的真实图像。这也是一个延拓和成像的过程。

3.4.3.2.1.1 延拓方程的推导

由下述二维波动方程出发。

地震勘探原理、方法及解释

根据爆炸反射面模型，将速度缩小一半，即用V/2代替V，可得

地震勘探原理、方法及解释

此方程有两个解，分别对应于上行波和下行波。但地震记录是上行波记录，故不能用此方程进行延拓，必须将它化为单纯的上行波方程才能利用。通常采用的方法是进行坐标变换后取近似值。第一步是坐标变换，令

地震勘探原理、方法及解释

上式中第二式是把方程中的深度坐标变为时间坐标。第三式是上行波的坐标变换。若称t为老时间，t′为新时间。因为坐标变换不改变实际波场，故原坐标系中波场u(x，z，t)与新坐标系中的波场 (x′，τ，t′)一样，即

地震勘探原理、方法及解释

由复合函数微分法，得

地震勘探原理、方法及解释

将上述二阶偏微分结果代入方程(3.4-2)，整理后得

地震勘探原理、方法及解释

为书写方便，以u、x、t分别代替u′、x′、t′，则(3.4-5)式可写为

地震勘探原理、方法及解释

式中：uxx，uττ，uτt分别表示u的二次导数。注意，此方程仍然包含了上行波和下行波，仍不能用来进行延拓，故还有第二步。

经过了坐标变换，虽然波场不变，但在新坐标系下，上、下行波表现出差异，此差异主要表现为uττ的大小不同。当上行波的传播方向与垂直方向之间的夹角较小时(小于15°)，uττ可以忽略，而对下行波来说，uττ不能忽略。忽略掉uττ项，就得到只包含上行波的近似方程

地震勘探原理、方法及解释

此即15°近似方程(因为它只适用于夹角小于15°的上行波，或者只有倾角小于15°的界面形成的上行波才能满足它)，为常用的延拓方程。

为了求解此方程还必须给出定解条件。由于震源强度有限，可给出如下定解条件

1)测线两端外侧的波场为零，即

u(x，τ，t)≡0 当 x＞ xmax或 x＜xmin时

2)记录最大时间以外的波场为零，即

u(x，τ，t)≡0 当 t＞ tmax时

3)自激自收记录(水平叠加剖面)为给定的边界条件，即时间深度τ=0 处的波场值u(x，0，t)已知。

有了这些定解条件就可对方程(3.4-7)求解得到地下任意深度处的波场值u(x，τ，t)，这是延拓过程。再根据前述成像原理，取(3.4-4)中，第三式的老时间t=0时刻时的波场值，即新时间t=τ时刻的波场值u(x，τ，t)就组成了偏移后的输出剖面。

图3-23 12点差分格式

3.4.3.2.1.2 差分方程

为了求解微分方程(3.4-7)，用差分近似微分，采用如图3-23所示的12点差分格式，将uxx、uτt表示为差分表达式，可得差分方程

地震勘探原理、方法及解释

式中：I和T为向量

I=[0，1，0] T=[-1，2，-1] (3.4-9)

α和β为标量

地震勘探原理、方法及解释

3.4.3.2.1.3 计算步骤和偏移结果

差分方程(3.4-8)形式上是一个隐式方程。即时间深度τ=(j+1)Δτ处的波场值不能单独地用时间深度τ=jΔτ处的波场值组合得到，方程右边仍有τ=(j+1)Δτ 的项。为了求得一排数据u(x，j+1，l)必须用到三排数据u(x，j+1，l+1)，u(x，j，l)和u(x，j，l+1)(图3-24)。

图3-24 有限差分法偏移求解中的一步

①u(x，j，l+1)；②u(x，j，l)；③u(x，j+1，l+1)；④u(x，j+1，l)

利用第二个定解条件，在计算新的深度τ=(j+1)Δτ处波场值时，由最大时间开始，首先计算t=tmax的那一排值。因u(i，j+1，tmax+Δt)≡0和u(i，j，tmax+Δt)≡0，有

地震勘探原理、方法及解释

计算u(i，j+1，tmax)只用到已知的u(i，j，tmax)值，十分容易。然后再利用(3.4-8)式递推地求τ=(j+1)Δτ深度处任何时刻的波场值就没有任何困难了。

具体计算时由地面向下延拓，计算深度Δτ处的波场值。首先计算此深度处在t=tmax时的波场，然后向t减小的方向进行。一个深度计算结束，再向下延拓一个步长Δτ继续计算。依此类推，可以得到地下所有点在不同时刻的波场值。

如前所述，在新时间t=τ时刻的波场值正是所欲求的“像”。因此，每次递推计算某一深度τ处的波场值时，由t=tmax向t减小的方向计算至t=τ时就可以结束。不同深处的“像”u(x，τ，t)组成偏移后的输出剖面。

图3-25 画出了偏移时的计算关系及结果取值位置。A 表示地面观测到的叠加剖面。由A计算下一个深度Δτ处的波场值 B，计算 B 时先算第1′排的数值(只用到A中第1排数值)，再算第2′排数值(要用A 中第1、2 排和B 中第1′排数值)，依此类推，直到 t=τ 为止。再由 B算下一个深度2Δτ处波场值C，……在二维空间(x，t=τ)上呈现出需要的结果剖面信息。

图3-25 偏移结果取值位置图

当延拓计算步长Δτ与地震记录的采样间隔Δt一样时，由图3-25 的几何关系可以看到，偏移剖面是该图中45°对角线上的值。实际工作中 Δτ 不一定要与Δt相等，可根据界面倾角大小确定Δτ，倾角较大时应取较小的Δτ，倾角较小时Δτ可取的较大些，以减少计算工作量。中间值可用插值求得。

与其他波动方程偏移方法相比，有限差分法有能适应横向速度变化，偏移噪声小，在剖面信噪比低的情况下也能很好地工作等优点。但15°有限差分法对倾角太大的情况不能得到好的偏移效果。因此，相继又研究发展了45°、60°有限差分偏移方法和适应更大倾角的高阶有限差分分裂算法。

3.4.3.2.2 频率波数域波动方程偏移

有限差分偏移方法是在时间空间域中进行的。利用傅里叶变换也可使偏移在频率波数域中实现。

与有限差分法偏移思想完全一样，认为水平叠加剖面是由界面上无数震源同时向上发出的上行波在地面处的波场值u(x，0，t)，用它反求地下任一点的波场值u(x，z，t)，这是延拓；据成像原理，取其在t=0时刻的值u(x，z，0)，组成偏移后的输出剖面。

仍由速度减半后的波动方程(3.4-2)出发，对方程两边做关于x和t的二维傅里叶变换，得到一个常微分方程

地震勘探原理、方法及解释

式中： = (kx，z，ω)为波场函数u(x，z，t)的二维傅里叶变换，ω=2πƒ为圆频率，kx为x方向上的空间波数。

式(3.4-11)是常微分方程，其解有两个，分别对应于上行波和下行波。偏移研究的是上行波的向下延拓问题，故只考虑上行波解

地震勘探原理、方法及解释

其中U(kx，0，ω)为解的初值，即上行波在z=0处的记录的傅里叶变换。因此，式(3.4-12)表示由z=0处波场的傅里叶变换求出任何深度处波场傅里叶变换的过程，是频率波数域中的波场延拓方程。

通过傅里叶反变换可由 (kx，z，ω)求出地下任何深度处的波场值

地震勘探原理、方法及解释

根据成像原理，偏移结果应是这些点处t=0时刻的波场值

地震勘探原理、方法及解释

这就是频率波数域偏移的数学模型。由于该式不是傅里叶变换公式，为了能利用快速傅里叶变换求解，经变量置换后，上式可变为一个傅里叶反变换公式。

3.4.3.2.3 克希霍夫积分偏移

克希霍夫积分偏移是一种基于波动方程克希霍夫积分解的偏移方法。

三维纵波波动方程的克希霍夫积分解(可见原理部分)为

地震勘探原理、方法及解释

式中：Q为包围点(x，y，z)的闭曲面，n为Q的外法线，r为由(x，y，z)点至Q面上各点的距离，[ ]表示延迟位，[u]=u(t-r/V)。

此解的实质是由已知的闭曲面Q上各点波场值计算面内任一点处的波场值。它正是惠更斯原理的严格数学形式。

选择闭曲面Q由一个无限大的平面Q0和一个无限大的半球面Q1所组成。Q1面上各点波场值的面积分对面内一点波场函数的贡献为零。因此，仅由平地面Q0上各点的波场值计算地下各点的波场值

地震勘探原理、方法及解释

此时，原公式中的 项消失，积分号前的负号也因z轴正向与n相反而变为正。

以上是正问题的克希霍夫积分计算公式。偏移处理的是反问题，是将反射界面的各点看作为同时激发上行波的源点，将地面接收点看作为二次震源，将时间“倒退”到t=0时刻，寻找反射界面的源波场函数，从而确定反射界面。反问题也能用上式求解，差别仅在于[ ]不再是延迟位而是超前位， 。根据这种理解，克希霍夫积分延拓公式应为

地震勘探原理、方法及解释

按照成像原理，此时t=0时刻的波场值即为偏移结果。只考虑二维偏移，忽略掉y坐标，将空间深度z转换为时间深度t0=2z/V，得到克希霍夫积分偏移公式

地震勘探原理、方法及解释

式中：τ= ]1/2，xl为地面记录道横坐标，x为偏移后剖面道横坐标，r=[z2+(x-xl)2]1/2(见图3-26)。

由 =-cosθ，得

地震勘探原理、方法及解释

由此可见，克希霍夫积分偏移与绕射扫描叠加十分相似，都是按双曲线取值叠加后放在双曲线顶点处。不同之处在于：①不仅要取各道的幅值，还要取各道的幅值对时间的导数值 参加叠加；②各道相应幅值叠加时不是简单相加，而是按(3.4-18)式的加权叠加。

正因如此，所以虽然形式上克希霍夫积分法与绕射扫描叠加类似，但二者有着本质区别。前者的基础是波动方程，可保留波的动力学特性，后者属几何地震学范畴，只保留波的运动学特征。

图3-26 克希霍夫偏移公式中各量示意图

与其他波动方程偏移法相比，克希霍夫积分法具有容易理解，能适应大倾角地层等优点。但它在速度横向变化较大的地区难以使用，且偏移噪声较大。

3.4.3.3 叠前波动方程偏移简介

叠后偏移需经过水平叠加处理才能进行，水平叠加本身是以射线理论为基础的近似处理方法，随着构造的复杂程度以及波场的复杂程度增加而误差越来越大，叠后偏移效果也随构造的复杂度而降低。叠前偏移是直接对野外接收的波场偏移归位，不受动校叠加的影响，理论和实践均证明其偏移效果明显优于叠后偏移。叠前偏移是偏移成像领域的发展方向。叠前偏移有二维或三维偏移，三维偏移可实现三维空间归位成像，成像质量优于二维。实现叠前偏移的方法同样有差分法、F—K法和积分法以及混合方法。下面以相移法三维叠前深度偏移为例，讨论叠前偏移的原理及实现方法。

由三维纵波方程

地震勘探原理、方法及解释

设 (kx，ky，z，ω)为p(x，y，z，t)的三维傅里叶变换，对(3.4-19)式作三维傅里叶变换，可求解得

地震勘探原理、方法及解释

式中：

式(3.4-20)称为相移延拓公式，仅适应V为常数的情况。

设地下为水平层状界面，在某一深度Z处ΔZ厚度层内的层速度为常数的条件下，该层的延拓公式为

地震勘探原理、方法及解释

该式为适应纵向变速V=V(z)的相移延拓公式。

设 (kx，ky，0，ω)为震源函数S(x，y，0，t)的三维傅里叶变换， (kx，ky，0，ω)为地面接收的地震记录R(x，y，0，t)的三维傅氏变换，则将震源函数在(k，ω)域正向延拓z的公式为

地震勘探原理、方法及解释

将记录R反向延拓Z的公式为

地震勘探原理、方法及解释

对(3.4-22)、(3.4-23)式作三维反傅里叶变换，并根据时间或深度一致性成像原理，求两波场在(x，y，z)点的互相关为

地震勘探原理、方法及解释

当相关延迟时间τ=0时，即得成像结果

地震勘探原理、方法及解释

该式也可以在(x，y，z，ω)域计算。

对横向变速介质，当V=V(x，y，z)时，(3.4-20)式中的kz应为

地震勘探原理、方法及解释

该式可写成

地震勘探原理、方法及解释

式中Vα为在Z深度平面的平均速度，则三维相移因子为

地震勘探原理、方法及解释

为满足相移公式条件，先用水平面平均速度Vα做纵向延拓，设延拓后的波场为

地震勘探原理、方法及解释

则横向变速的结果为：

地震勘探原理、方法及解释

在式(3.4-30)中的指数部分用二项式展开并略去高次项，得

地震勘探原理、方法及解释

该式即是相移法延拓后适应速度横向变化的校正因子。根据不同的精度要求保留相应高次项，可分别作一阶、二阶或三阶校正。校正可在F-K域进行，也可在F-x域进行，若在F-x域用差分法进行校正，则称为混合法波动方程叠前深度偏移。以上叠前深度偏移方法的实现过程是对共炮集三维观测记录分别偏移成像，然后按空间位置叠加。

写入方法如下，Offs指令：目销亩标点位置偏移，常用于帆拿机器人码垛、物态斗搭料仓储等摆放有规则的模块中。

格式：Offs(p1,100,50,0)

p1：目标点位置

100：距离p1点X轴的偏差量为100mm

50：距离p1点Y轴的偏差量为50mm

0：距离p1点Z轴的偏差量为0mm

---