在地震资料处理中, 动校正的作用是消除炮检距对反射波传播时间的影响, 将共中心点道集(CMP道集)中的反射波同相轴拉平, 其精确性直接影响干扰波的压制效果以及后续偏移成像的质量。采用常规动校正方法存在同相轴无法拉平和拉伸畸变等问题, 利用分偏移距动校正方法^[1-2]或消除动校正剩余时差的方法^[3-6]能够使同相轴拉平, 但不能消除拉伸畸变。研究表明, 动校正拉伸畸变通常发生在浅层和大炮检距处^[7-11], 具体表现为波形被拉长, 频率被降低。常规地震数据处理中拉伸超过一定比例的数据被切除^[10-11], 这样浅层地震数据会因拉伸严重几乎全部或大部分被切除, 导致覆盖次数降低。为了解决常规动校正的拉伸畸变问题, 尽可能多地保留远道地震信息, 施剑等^[8]、LICHMAN^[12]和崔宝文等^[13]提出基于频谱替换的无拉伸动校正方法, 使远道信息得到了保留。该方法理论上可以消除任何时差, 并且由于采用快速傅里叶变换算法, 具有很高的计算效率。但由于不同的同相轴在进行相位替换后会相互影响, 因此实际数据处理中该方法常用于剩余时差校正^[14]。RUPERT等^[15]通过数据块整体移动并对重叠部分求和的方法进行动校正, 在一个时窗内仅使用一个动校正量, 避免了因动校正量不同而导致的拉伸畸变。但该方法在炮检距较大时容易出现同相轴的过校正和欠校正问题, 在相邻时窗重叠处数据相加再除以叠加次数的方法也导致动校正效果变差。此外, 该方法对弱反射层的动校正会受到相邻强反射层的影响。PERROUD等^[16]、孟庆生^[17]和HASSAN等^[18]通过修改动校正速度使得同一个波形的各时间点具有相同的动校正量, 从而消除了动校正拉伸畸变, 但该方法需要找准每个同相轴的初至时间, 实际资料由于波长从浅至深不断变化等因素, 很难选取合适的时窗和时移因子。同时, 速度修改后两个同相轴之间的速度变化加快, 如果某个反射同相轴被遗漏, 会比常规动校正更容易出现拉伸畸变。李录明^[19]提出先自动检测反射波垂直反射时间和速度、后计算待校正炮检距的反射时间并将该时间处波形整体移动的方法, 较好地消除了拉伸畸变, 但该方法对速度分析的精度要求非常高。赵波等^[20]采用动校正过程中的速度信息进行正演, 再采用滤波方式消除动校正拉伸畸变, 但该方法对速度分析的精度要求较高且计算量大。利用抛物Radon变换进行叠加的方法能够消除拉伸畸变^[21-22], 但必须在Radon域寻找一个合适的路径, 或自动拾取能量团的聚焦点(存在噪声时可靠性较差), 因此该方法在工业上的应用受到一定的限制。基于匹配追踪算法的无拉伸动校正方法^[23]由于子波库的冗余性以及方法本身的问题会出现非正交投影和过匹配现象, 计算效率和动校正效果会受到影响。BIONDI等^[24]提出多次迭代的方法, 每次迭代时进行部分动校正并利用反褶积来修正拉伸的影响, 减小拉伸畸变, 但每次迭代过程中采用的是常规动校正方法, 动校正结果仍然会引入拉伸畸变。KAZEMI等^[25]研究了局部无拉伸的动校正方法, 将地震数据从浅到深划分多个时间窗口, 重新计算出每个窗口内每个采样点的时距曲线后再进行动校正, 但该方法要求多个反射波在零炮检距处不相互干涉, 且大炮检距数据整体在时间方向上存在不连续的情况。

本文分析了常规动校正拉伸畸变的原因及动校拉伸的范围, 给出一种无拉伸畸变的动校正方法, 模型数据和实际资料处理测试证明了方法的可行性。

1 方法原理 1.1 常规动校正拉伸分析

常规动校正方法中, 地震反射波的时距方程为^[7]：

$t=\sqrt{t_{0}^{2}+{{\{x/[v({{t}_{0}})]\}}^{2}}}$

(1)

式中:t是地震波传播时间; t₀是零炮检距反射波双程旅行时间; x是炮检距; vt₀是地震波速度, 通常为均方根速度。若地下介质是均匀介质, 则地震波速度v为常数。如图 1所示, 设t₀₁和t₀₂是零炮检距反射波双程旅行时间, 且t₀₁对应子波的初至, 相应的t₁和t₂分别是炮检距为x的地震波传播时间：

${{t}_{1}}=f({{t}_{01}})=\sqrt{t_{01}^{2}+{{\left( x/v \right)}^{2}}}$

(2)

${{t}_{2}}=f({{t}_{02}})=\sqrt{t_{02}^{2}+{{\left( x/v \right)}^{2}}}$

(3)

根据微分中值定理：

${{t}_{2}}-{{t}_{1}}=f({{t}_{02}})-f({{t}_{01}})={{f}^{\prime }}{{t}_{0\varepsilon }}({{t}_{02}}-{{t}_{01}})$

(4)

式中:t_0ε为t₀₁和t₀₂之间的一点,

$\frac{{{t}_{2}}-{{t}_{1}}}{{{t}_{02}}-{{t}_{01}}}={{f}^{\prime }}\left( {{t}_{0\varepsilon }} \right)=\frac{{{t}_{0\varepsilon }}}{\sqrt{t_{0\varepsilon }^{2}+{{\left( \frac{x}{v} \right)}^{2}}}}＜1$

(5)

即t₂-t₁＜t₀₂-t₀₁, 亦即计算的非零炮检距的子波长度与零炮检距的子波长度不相等。类似地, 非均匀介质也有相同的性质^[7]。所以, 对地震数据进行常规动校正时, 除炮检距为零外, 信号总是要在时间上被拉长, 波形会发生不同程度的畸变, 这也是动校正中客观存在的现象。由于f′(t_0ε)=t_0ε/$\sqrt{t_{0\varepsilon }^{2}+{{\left( x/v \right)}^{2}}}$介于0和1之间, 它随着t_0ε的减小而减小, 随着x的增大而减小, 根据公式(4)可知, 在浅层和大炮检距处的拉伸畸变较大, 深层和近炮检距处的拉伸畸变较小。

1.2 拉伸畸变的消除方法

由上述分析可见, 产生动校正拉伸畸变的根本原因, 是同一子波采用了不同的动校正量。子波内部不同样点间的动校正量变化越大, 拉伸越严重。为有效压制这种拉伸畸变, 我们对子波内各样点的动校正量进行动态调整, 按照某种规则找出最佳动校正量。首先计算动校正量可能存在的修正范围, 然后在修正范围内利用互相关算法筛选出无拉伸畸变的动校正量, 最后使用该动校正量对道集数据进行“逐点搬家”, 实现无拉伸畸变的动校正。

1.2.1 动校正量的修正范围

如图 1所示, 假设t₀₁为子波初至时间, t₀₂是该子波延续时间内的某一时间点, t₀₁和t₀₂之间的时间差为τ=t₀₂-t₀₁。设对应该反射子波的速度v为常数, 根据公式(1)计算的零炮检距时间等于t₀₁和t₀₂的反射波时距曲线为1和2, 可以看出, 当炮检距较大时, 根据公式(1)计算的子波两端分别在t₁和t₂处, t₁和t₂之间的时间长度为τ′=t₂-t₁＜t₀₂-t₀₁, 动校正后其时间位置分别对应于t′₁和t′₂, 因而出现了子波被拉伸的现象(图 1)。

物理上子波的长度应保持不变, t₀₂对应的时距曲线应该是在曲线1上加上t₀₂-t₀₁得到的结果, 对应图 1中的曲线3。所以曲线3的时距曲线方程是：

${{t}_{3}}=({{t}_{02}}-{{t}_{01}})+{{t}_{1}}=\tau +\sqrt{t_{01}^{2}+{{\left( \frac{x}{v} \right)}^{2}}}$

(6)

要使动校正后的子波不出现拉伸畸变, 应计算波形初至时间的动校正量, 且整个波形的校正都使用与之相同的动校正量。然而, 实际动校正过程中每一个波形的初至时间很难获取, 其速度也不是固定不变, 而是随时间变化的函数, 因而t₃也不容易求取。

可以认为, 曲线3是曲线2加一个修正量得到的结果, 当t₀₂=t₀₁时, 修正量为0; 当t₀₁和t₀₂恰好在波形两端, 即子波长度λ=t₀₂-t₀₁时, 修正量最大。考虑到速度随时间变化, 炮检距为x时t₀时间的动校正量最大修正范围是：

$\Delta t_{max}^{\prime }=\lambda +t({{t}_{0}}-\lambda )-t({{t}_{0}})$

(7)

其中,

$\begin{align} & t({{t}_{0}}-\lambda )=\sqrt{{{({{t}_{0}}-\lambda )}^{2}}+{{\left[ \frac{x}{v({{t}_{0}}-\lambda )} \right]}^{2}}} \\ & t\left( {{t}_{0}} \right)=\sqrt{t_{0}^{2}+{{\left[ \frac{x}{v\left( {{t}_{0}} \right)} \right]}^{2}}} \\ \end{align}$

(8)

所以, 理想的动校正量应该是由常规动校正量Δt=$\sqrt{t_{0}^{2}+{{\left[ \frac{x}{v\left( {{t}_{0}} \right)} \right]}^{2}}}-{{t}_{0}}$加上一定修正量得到的结果, 修正量范围为0~Δt_max′。Δt_max′的大小与常规动校正拉伸畸变的程度相对应, 浅层远道处较大, 深层近道处较小。因此, 本文方法在浅层远道的计算量较大, 在深层近道的计算量与常规方法相近, 因而具有较强的实用性。

1.2.2 修正量的筛选

在动校正量的修正范围内, 存在一个最佳修正量, 用它对常规动校正量进行修正, 能够有效消除拉伸畸变。从动校正量的修正范围内筛选最佳修正量的步骤如下。

1) 建立标准道。为减小动校正拉伸畸变对标准道的影响, 可采用常规动校正后近道数据的叠加作为标准道^[5]。

2) 扫描修正量。对于每一个待校正的t₀时间数据, 每个炮检距都需要计算常规动校正量Δt和最大修正量Δt_max′。在0~Δt_max′之间按一定时间间隔对修正量Δt′进行扫描, 并将动校正量更新为Δt+Δt′, 对应时间为t=t₀+Δt+Δt′, 即将t时间的数据移动到t₀时间。

3) 选取数据时窗。对于标准道, 以t₀为中心选取一定长度的数据时窗；对于待校正道, 由于同一个波形应具有相同的动校正量, 故以t为中心确定与标准道相同长度的数据时窗。因为不同波形的动校正量不同, 所以数据时窗的长度不宜太长。

4) 计算互相关量。对待校正道时窗内数据和标准道时窗内数据进行归一化互相关计算, 将互相关值作为评价动校正效果的准则。归一化互相关公式如下：

${{r}_{cor}}=\frac{\sum\limits_{j=1}^{N}{m({{t}_{j}})w({{t}_{j}})}}{\sqrt{\sum\limits_{j=1}^{N}{{{m}^{2}}({{t}_{j}})}\sum\limits_{j=1}^{N}{{{w}^{2}}({{t}_{j}})}}}$

(9)

式中:r_cor表示互相关系数; N为时窗长度; w(t_j)为待校正道时窗内的数据; m(t_j)表示标准道时窗内的数据。

动校正的过程就是将动校正前在时间t处的数据“搬家”, 移到t₀处, 因此, 待校正道时窗内的数据与标准道时窗内的数据互相关, 相当于对待校正的数据进行动校正后再与标准道的数据互相关, 其动校正量是常规动校正量加修正量的结果, 即Δt+Δt′。对于0~Δt_max′之间的每一个修正量, 都会有一个待校正数据的时窗(相当于进行动校正), 如果某个时窗内的数据与标准道时窗内数据的互相关系数最大且为正值, 则说明动校正效果最好, 采用此时的更新动校正量Δt+Δt′对数据进行动校正效果最好。选定此时的动校正量为最终动校正量, 并记录其动校正前的时间t=t₀+Δt+Δt′, 如果最大的互相关系数小于或等于0, 说明标准道时窗内数据全为0或不需要动校正, 此时记录的动校正前的时间为0。

重复步骤2)至步骤4), 确定所有t₀时间和所有道的动校正量, 即可依次完成动校正。

1.2.3 无拉伸畸变的动校正

相邻同相轴之间可能存在无波形信息的空白区域, 该区域的大小随炮检距的不同而变化, 常规动校正对同相轴和同相轴之间的空白区域都存在拉伸。无拉伸畸变的动校正通过对动校正量进行修正, 使得同一反射波组的波形不被拉伸, 但由于同一反射波组的t₀时间唯一, 波组之间的空白区域势必会被拉伸。因此, 无拉伸畸变的动校正相当于把常规动校正对反射波组的拉伸进行修正, 并把拉伸量转嫁给了空白区域。

如果地震数据含有一定的随机噪声, 对于存在反射波组的区域, 采用互相关方法能够计算最佳的动校正量, 实现无拉伸畸变的动校正。但在两个同相轴之间的空白区域, 由于随机噪声的影响, 动校正后数据的相对位置关系可能会发生改变, 需要对此进行修正。无拉伸畸变动校正与常规动校正相比, 各样点动校正后的新位置可能存在差异, 但同一道内各样点的时间先后关系应该保持不变。根据这一原则, 我们利用常规动校正方法中各样点的时间相对关系进行约束, 方法是：对于两个t₀时间t₀₁和t₀₂(t₀₂＞t₀₁), 在炮检距为x的地震道上, 设按照本文方法计算得到的动校正前时间分别为t₁和t₂, 而利用常规方法计算得到的动校正前时间分别为t₁′和t₂′, 若计算结果t₁′＜t₂′, t₁＜t₂, 则认为是合理的, 若t′＜t₂′, 而t₁≥t₂, 则需要逐渐减小t₁₁直到满足条件为止。对所有时间样点完成修正后, 将待校正数据中时间为t处的数据“搬家”, 移到t₀处, 就完成了无拉伸畸变的动校正。

2 模型数据测试 2.1 无同相轴交叉的模型

建立5层介质模型一, 各层速度分别为2000, 2100, 3000, 3500和3700m/s, 反射层深度分别为300, 700, 1200, 1800m。通过正演获得图 2a所示的模型一数据, 共180道, 道间距是10m, 最小炮检距为0, 模型中4个同相轴不出现交叉。图 2b是对模型一数据进行常规动校正的结果, 可以看出同相轴①和②在大炮检距处出现了拉伸畸变, 没有被拉平, 其余两个同相轴的拉伸较小, 基本看不出畸变, 说明常规动校正方法容易在浅层和远道处出现拉伸畸变。图 2c是模型一数据采用本文方法的动校正结果, 可以看出4个同相轴均已拉平且没有出现畸变, 说明本文方法是有效可行的。图 2d是两种方法动校正结果的叠加波形对比, 为清楚地展示动校正效果, 本文仅截取了0.25~0.45s之间同相轴①的数据进行显示。从图 2d可以看出, 常规动校正结果由于远道处出现拉伸, 叠加波形的振幅偏小, 且波形被拉伸; 本文动校正方法由于使同相轴得到拉平, 实现了各道数据的同相叠加, 因此叠加后波形特征得到了保持, 振幅比常规动校正方法大。波形的拉伸会造成频率的降低, 且随着炮检距的增大, 常规动校正的拉伸会越来越严重, 主频向低频方向的移动也越来越多。图 2e是对常规动校正结果(图 2b)逐道进行傅里叶分析得到的振幅谱, 从第80道开始明显出现频率范围向低频方向移动的现象, 这与图 2b中子波的拉伸现象一致, 图 2e中的频率向低频方向的移动越明显, 图 2b中地震道的拉伸畸变越严重。图 2f是对本文方法动校正结果(图 2c)逐道进行傅里叶分析得到的振幅谱, 可以看出, 各道地震数据的主频基本一致, 没有出现频带移动的现象。

2.2 存在同相轴交叉的模型

在地层倾斜等复杂构造条件下, 有可能出现深层界面的反射波旅行时小于浅层的情形, 即出现同相轴交叉的现象。在地层厚度很小时, 反射波同相轴也可能存在一定的干涉。不失一般性, 我们构造图 3a所示的模型二数据, 分析存在同相轴交叉及干涉情形的动校正拉伸。该模型数据共180道, 最小炮检距为0, 道间距为10m, 同相轴①和同相轴②在140道附近的A点交叉, 同相轴③和同相轴④存在干涉。图 3b是模型二数据常规动校正结果, 由于常规动校正根据速度计算动校正量, 在同相轴交叉后, 同相轴②出现的时间比同相轴①出现的时间早, 导致同相轴①上方出现同相轴②的信息^[23-24]。由于速度逐渐增大, 同相轴①在交叉点附近难以显示整个波形, 出现被“劈开”的现象, 在同相轴②下方会出现同相轴①的信息; 同相轴②在与同相轴①发生交叉之后的部分会出现两次过校正, 且这两次过校正在A点相交。因为常规动校正在“搬家”完成后并不抹去原先的信息, 所以同相轴①和同相轴②可能出现多次。图 3c是模型二数据使用本文方法动校正的结果, 图中两个同相轴均已拉平, 虽然在交叉点A附近的校正效果与其它区域相比要略差一点, 但依然没有出现拉伸畸变。与常规动校正方法不同, 由于本文方法只选择互相关值大于0时对应的动校正量, 因此同相轴①的信息不会在不同区域多次出现。图 3d是常规方法和本文方法动校正结果的叠加波形对比, 为清楚地展示动校正效果, 本文仅截取0.80~0.88s之间同相轴①的数据进行显示。由图 3d可见, 常规动校正结果被同相轴②的信息所干扰, 叠加后波形发生畸变; 本文方法由于使同相轴得到了拉平且不受同相轴②的干扰, 各道数据实现了同相叠加, 叠加后波形特征得到保持, 叠加后振幅比常规动校正结果强。图 3e是对常规动校正结果(图 3b)逐道进行傅里叶分析得到的振幅谱, 由于同相轴①在不同区域出现两次, 同相轴②在不同区域出现3次, 导致各道振幅谱之间变化较大, 且在交叉点A处有一条贯穿低频到高频的线, 在50道之后出现频率范围明显向低频方向移动的现象。图 3f是对本文方法动校正结果(图 3c)逐道进行傅里叶分析得到的振幅谱, 可以看出, 各道数据没有发生频率范围移动的现象, 虽然在交叉点A附近依然存在频率范围贯穿低频到高频的线, 但已经得到了很大改善。这条线是由于两个交叉同相轴的波形未能完全分离导致的。

为分析反射波同相轴存在干涉时的动校正效果, 对同相轴③和④的分炮检距叠加波形进行放大显示。图 4a为常规动校正方法和本文动校正方法在0~500m炮检距范围的叠加结果对比, 可见两种方法的叠加波形几乎完全一致, 说明常规方法在近炮检距的拉伸量很小。为了分析两种方法对远炮检距数据的动校正效果, 对常规方法动校正后不同炮检距的叠加波形进行放大显示, 如图 4b所示, 可见炮检距越大, 波形拉伸越严重。图 4c是本文方法动校正后不同炮检距的叠加波形对比, 可见各炮检距的波形都未出现拉伸, 但由于大炮检距比近炮检距的波形干涉严重, 因此, 随着炮检距增大, 叠加波形的振幅略微减小。

2.3 抗噪性分析

为了测试本文方法对随机噪声的敏感性, 在图 2a所示的模型数据中加入随机噪声, 使其信噪比分别为5∶1, 1∶1和0.5∶1.0, 进而分析本文方法对高信噪比、中等信噪比和低信噪比数据的处理效果(图 5)。

图 5a, 图 5d, 图 5g分别是信噪比为5∶1, 1∶1和0.5∶1.0的模型数据, 随着信噪比的降低, 同相轴逐渐被淹没在噪声中。图 5b, 图 5e, 图 5h是对应的常规动校正结果, 同相轴的拉伸畸变几乎不随信噪比的降低而变化。图 5c, 图 5f, 图 5i是本文方法的动校正结果, 对于不同信噪比数据, 同相轴都能够被拉平且无畸变。图 6对比了图 5中动校正结果的叠加波形(为了清楚地展示出动校正效果, 仅显示了第一个同相轴), 可见常规动校正方法由于在大炮检距处存在拉伸畸变, 导致非同相叠加, 因此叠加波形的振幅小于本文方法, 且波形被拉伸, 拉伸程度不随信噪比的降低而改变, 说明常规动校正拉伸畸变对信噪比变化不敏感。本文动校正方法则由于控制了拉伸畸变, 叠加波形未被拉伸, 然而在同相轴以外的区域, 噪声能量随着信噪比的降低而增强, 这是由于选择动校正量时, 与标准道波形相似的噪声被作为有效波进行了处理, 导致相干噪声被增强。因此, 使用本文方法处理低信噪比数据时, 需先进行叠前去噪。

3 实际地震数据试处理

采用我国东部某地区道间距为50m, 信噪比较高的地震数据对本文方法的有效性与实用性进行了测试(图 7)。图 7a为实际CMP道集, 图 7b是常规动校正结果(仅截取0~1.2s之间的数据进行了显示)。可以看出, 经过常规动校正后同相轴基本被校平, 但大炮检距的同相轴比近炮检距的同相轴要宽, 说明波形被拉伸。图 7c为采用本文方法动校正的结果, 同相轴均已拉平且未出现拉伸现象。图 7d是两种方法动校正结果的叠加数据对比(为清楚地展示动校正效果, 只显示了0.9~1.2s之间的数据), 可以看出, 由于本文方法动校正的波形没有被拉伸, 叠加后波形的振幅大于常规方法。图 7e是对常规方法动校正结果(图 7b)逐道进行傅里叶分析得到的振幅谱, 可以看到前20道数据的频率范围基本没变, 但20道之后的频率逐渐向低频方向移动, 这是由波形拉伸所致。图 7f是对本文方法动校正结果(图 7c)逐道进行傅里叶分析得到的振幅谱, 由于波形没有被拉伸, 各道数据的频率范围基本不变, 说明本文方法能够对实际资料有效进行无拉伸畸变的动校正。

图 8比较了实际叠加剖面, 道间距为25m。图 8a为常规动校正后的叠加剖面, 由于存在拉伸畸变, 导致非同相叠加, 因此部分同相轴连续性较差。图 8c为本文方法动校正后的叠加剖面, 由于动校正过程中同相轴均已拉平且不出现拉伸, 叠加后的波形能量得到加强, 图 8a中部分能量较弱的同相轴也被加强, 连续性得到较大的改善, 剖面整体的分辨率得到提高。图 8b和图 8d分别是图 8a和图 8c各道地震数据的振幅谱, 常规方法的主频在20Hz附近, 而本文方法的主频约为35Hz, 比常规方法提高了15Hz。

4 结论

本文分析了地震资料处理中常规动校正产生拉伸畸变的原因, 研究了修改动校正量以减少拉伸畸变的方法, 模型数据测试和实际资料处理验证了方法的正确性与有效性, 得出如下结论：

1) 将常规动校正后近炮检距的叠加道作为标准道, 计算与标准道的互相关, 在一定范围内对常规动校正量进行修正, 可以消除拉伸畸变。

2) 对于含有一定噪声的数据, 本文方法仍然可以实现无拉伸畸变的动校正, 但随着信噪比的降低, 叠加波形中同相轴之外的相干噪声能量增强, 因此, 将本文方法应用于低信噪比数据时, 需先进行叠前去噪处理。

致谢: 本文研究过程中, 中国石油大学(华东)地球科学与技术学院伍敦仕博士与中国石油化工股份有限公司胜利油田分公司物探研究院赵爱国博士参与了讨论并给予了宝贵的帮助, 在此表示感谢。