地震波在地下地层中传播时, 由于介质的黏弹性性质, 地震波会出现明显的吸收衰减现象, 表现为振幅衰减和波形畸变, 降低地震资料的分辨率。自1962年Futterman首次提出吸收衰减理论以及补偿方法^[1]以来, 在此基础上发展了许多补偿算法, 主要有频率域^[2-3]、时间域^[4]和时频域^[5]的Q值求取方法, 但这些方法求取的Q值精度不佳, 为提高叠前叠后衰减补偿的精度, 非线性最优化反演^[6]、正则化反演^[7-9]、反Q滤波^[10-13]、Gabor域振幅平滑^[14]等方法被发展起来。对叠前地震资料来说, 吸收衰减使振幅的横向反射特征发生改变, 很容易影响AVA的正确分析, 针对这一问题, 许多基于NMO拉伸校正^[15]、常Q模型^[16]、射线路径^[17-20]的补偿方法发展了起来。虽然这些方法都能对叠前衰减进行补偿, 但是都有一个共同的问题, 即在补偿过程中高频成分振幅会被过度放大, 产生高频噪音, 降低地震资料信噪比。匹配追踪算法^[21]是一种高精度的地震信号分解与重构算法, 但计算效率一直制约着该算法的应用, 为提高算法计算效率, 对其时频原子库进行了改进^[22-24], 得到了较好的效果。笔者在前人的研究基础上, 结合匹配追踪算法构建随匹配子波中心时间和中心频率乘积值变化的振幅衰减曲线, 对叠前道集进行振幅补偿和频率校正, 将二维信息时间和频率降低为一维信息时间和频率的乘积值, 以增强算法的稳健性。

1 方法原理 1.1 吸收衰减机制

在黏弹性介质中, 波数是一个满足因果关系的与频率相关的复函数。Wang^[25]通过研究发现, 为增加频散校正的准确度, 复波数可表示为

$ k\left( f \right) = \frac{{2{\rm{ \mathsf{ π} }}f}}{v} + {\rm{i}}\gamma . $

(1)

式中, γ和v分别表示与频率有关的衰减因子和相速度; f为频率。在Kolsky^[26]衰减模型中, 假设衰减是一个关于频率的线性函数, 衰减因子表示为

$ \gamma \left( f \right) = \frac{{\rm{ \mathit{ π} }}}{{vQ}}f. $

(2)

相速度表示为

$ v\left( f \right) = {v_0}\left( {{f_0}} \right)|f/{f_0}{|^{\frac{1}{{{\rm{ \mathit{ π} }}Q}}}}. $

(3)

式中, Q为品质因子; f₀为参考频率; v₀为参考频率处的参考相速度。根据Bickel和Natarajan^[27]的理论, 公式(3)可以简写为

$ v\left( f \right) = {v_0}\left( {{f_0}} \right)\left( {1 + \frac{1}{{{\rm{ \mathit{ π} }}Q}}{\rm{ln}}\left| {\frac{f}{{{f_0}}}} \right|} \right). $

(4)

对公式(3)作一阶级数近似可得到慢度P(f):

$ P\left( f \right) = \frac{1}{{v\left( f \right)}} = \frac{1}{{{v_0}\left( {{f_0}} \right)}}\left( {1-\frac{1}{{{\rm{ \mathit{ π} }}Q}}{\rm{ln}}\left| {\frac{f}{{{f_0}}}} \right|} \right). $

(5)

将公式(2)、(5)代入式(1)中, 可得复波数为

$ k\left( f \right) = \frac{{\pi f}}{{{v_0}\left( {{f_0}} \right)}}\left( {1-\frac{1}{{\pi Q}}{\rm{ln}}\left| {\frac{f}{{{f_0}}}} \right|} \right)\left( {2 + \frac{1}{Q}{\rm{i}}} \right). $

(6)

在频率域内求解单程波动方程的平面波, 可得衰减信号X(f)和原信号X₀(f)有如下关系:

$ X\left( f \right) = {X_0}(f){\rm{exp}}[{\rm{i}}k\left( f \right)\xi]. $

(7)

式中, ξ为波前的传播距离。将公式(6)代入式(7)可得

$ X\left( f \right) = {X_0}(f){\rm{exp}}(2{\rm{ \mathsf{ π} }}ft{\rm{i}}){\rm{exp}}\left( { - \frac{{{\rm{ \mathsf{ π} }}ft}}{Q}} \right){\rm{exp}}\left( { - \frac{{2ft}}{Q}{\rm{ln}}\left| {\frac{f}{{{f_0}}}} \right|{\rm{i}}} \right). $

(8)

式中, t为地震波传播时间。定义衰减算子如下:

$ \alpha \left( {f, t} \right) = {\rm{exp}}\left[{-\left( {\frac{{{\rm{ \mathit{ π} }}ft}}{Q} + {\rm{i}}\frac{{2ft}}{Q}{\rm{ln}}\left| {\frac{f}{{{f_0}}}} \right|} \right)} \right]. $

(9)

式中, 实部为振幅衰减项, 虚部为速度频散项。

1.2 零偏移距数据的吸收补偿

对于零偏移距地震道, 其吸收衰减只与深度有关。当只考虑振幅衰减时, 零偏移距衰减算子为

$ \left| {\alpha \left( {f, {t_0}} \right)} \right| = {\rm{exp}}\left( {-\frac{{\rm{ \mathit{ π} }}}{Q}f{t_0}} \right). $

(10)

式中, t₀为零偏移距反射时间。由式(10)不难发现, 衰减算子是关于时间和频率的函数, 且与两者的乘积有关, 因此可以寻找衰减振幅与时间和频率乘积值间的关系进行吸收补偿。

匹配追踪是一种自适应信号分解算法, 可以对信号进行匹配分解。对零偏移距地震信号匹配分解后, 其可表示为一系列匹配子波的线性叠加:

$ s\left( t \right) = \sum\limits_i {{A_i}g\left( {{t_{0i}}, {f_{0i}}, {\varphi _i}} \right)} . $

(11)

式中, A_i、t_0i、φ_i和f_0i分别表示不同匹配子波的振幅、中心时间、相位和中心频率。将各个匹配子波中心时间和中心频率乘积值t_0if_0i和相应匹配子波振幅投影到坐标平面内; 把t_0if_0i分为等间隔的条段并求取各个条段内的子波振幅平均值, 即得到随t_0if_0i变化的振幅衰减曲线a(t_0if_0i); 然后由振幅衰减曲线计算振幅增益曲线:

$ \tilde a\left( {{t_{0i}}{f_{0i}}} \right) = \frac{{a\left( {{t_{0i}}{f_{0i}}} \right)}}{{{a^2}\left( {{t_{0i}}{f_{0i}}} \right) + \alpha {A_{{\rm{max}}}}}}. $

(12)

式中, α为白噪系数; A_max为振幅衰减曲线的最大值。通过式(12)即可求得各个匹配子波衰减振幅的增益系数, 对衰减振幅进行补偿:

$ \overline {{A_i}} = \tilde a\left( {{t_{0i}}{f_{0i}}} \right){A_i}. $

(13)

1.3 非零偏移距数据的吸收补偿

对于非零偏地震道, 其振幅衰减算子可以表示为

$ \left| {\beta \left( {f, {t_f}} \right)} \right| = {\rm{exp}}\left( {-\frac{{\rm{ \mathit{ π} }}}{Q}f{t_f}} \right). $

(14)

式中, t_f为地震波非零偏移距反射时间。分析可知, 对于同一反射层, 零偏移距反射时间t₀与非零偏移距反射时间t_f有如下关系:

$ {t_f}{\rm{ = }}\frac{1}{{\cos \;\theta }}{t_0}. $

(15)

式中, θ为地震波传播方向与垂直方向的夹角, (°)。则式(14)可写为

$ \left| {\beta \left( {f, {t_f}} \right)} \right| = {\rm{exp}}\left( {-\frac{{\rm{ \mathit{ π} }}}{Q}f{t_0}\frac{1}{{\cos \;\theta }}} \right) = {\left| {\alpha \left( {f, {t_0}} \right)} \right|^{\frac{1}{{\cos \;\theta }}}}. $

(16)

由式(16)可知:对于同一反射层, θ越大, 振幅衰减越严重; 并且非零偏移距振幅衰减算子可由零偏移距振幅衰减算子与夹角θ表示, 即其衰减同时与深度和偏移距有关。

结合零偏移距吸收补偿方法, 首先计算与深度有关的匹配子波振幅增益系数, 然后结合夹角θ得到非零偏移距各个匹配子波振幅增益系数, 可表示为

$ b\left( {{t_{0i}}{f_{0i}}} \right) = {\left[{\tilde a\left( {{t_{0i}}{f_{0i}}} \right)} \right]^{\frac{1}{{\cos \;\theta }}}}. $

(17)

式(17)同时补偿了深度衰减和偏移距衰减。当仅补偿与偏移距有关衰减, 即非零偏移距与零偏移距的吸收差异时, 增益系数可写为

$ \tilde b\left( {{t_{0i}}{f_{0i}}} \right) = {\left[{\tilde a\left( {{t_{0i}}{f_{0i}}} \right)} \right]^{\frac{1}{{\cos \;\theta }} -1}}. $

(18)

然后, 用各个匹配子波振幅增益系数对衰减振幅进行补偿, 最后, 用振幅补偿后的各个匹配子波重构地震道。

1.4 频散校正

对于速度频散, 通过匹配追踪算法求取地层品质因子^[28]:

$ {\rm{ln}}a\left( \chi \right) =-\frac{{\rm{ \mathit{ π} }}}{Q}\chi + {\rm{ln}}{a_0}. $

(19)

式中, χ=t_0if_0i为匹配子波中心时间和中心频率的乘积值; a(χ)=a(t_0if_0i)为随χ变化的振幅衰减曲线; a₀为未衰减振幅。这样就可以通过拟合斜率得到地层品质因子值, 然后运用公式(9)的虚部进行频散校正。

2 方法测试与应用效果分析 2.1 理论模型测试

为测试本文吸收补偿方法的有效性, 建立衰减道集模型对方法进行测试。图 1为建立该衰减道集模型所需参数, 包括密度、横波速度、纵波速度以及品质因子模型。根据佐普里兹方程计算得到反射系数, 并用60 Hz雷克子波合成叠前道集; 在此过程中不考虑透射损失和几何扩散, 振幅衰减和速度频散用式(9)给出, 其中速度频散的参考频率为200 Hz。

图 2(a)为没有地层吸收的叠前道集, 图 2(b)为包含地层吸收的叠前道集。通过对比可以发现, 由于受地层吸收的影响, 纵向上, 随着传播时间的增加, 地震波振幅衰减越严重, 分辨率降低; 横向上, 地震反射特征发生明显变化, 改变了地震记录的同相性, 不能反映真实的振幅随入射角变化特征。

图 3(a)是对图 2(b)通过匹配追踪求取品质因子然后进行频散校正后的结果, 由于频散造成的同相轴下拉现象得到了校正, 改善了叠前道集的同相性; 图 3(b)是对图 3(a)的结果进一步作与偏移距相关吸收补偿后的结果, 由于消除了与偏移距相关的吸收差异, 地震记录的横向一致性得到明显的提高, AVA反射特征得到恢复; 图 3(c)是图 3(b)作与深度有关吸收补偿后的结果, 可以看出深部衰减能量得到有效补偿, 并且与没有地层吸收及速度频散效应的图 2(a)完全一致。

图 4是图 2(b)模型数据中480 ms处反射子波与偏移距相关吸收补偿前后对比。图 4(a)是补偿前后反射振幅随偏移距的变化对比, 可以看出:由于吸收的影响, 补偿前反射振幅随偏移距增大而减小, 改变了横向反射特征, 产生错误的AVA信息, 补偿后反射振幅随偏移距增大而增大, 恢复了正确的横向反射特征。图 4(b)是补偿前后反射主频随偏移距的变化对比, 可以看出:由于频散的影响, 反射波主频随偏移距增大逐渐降低, 削弱了叠前道集横向一致性, 补偿后不同偏移距地震反射波主频没有发生变化, 这可以增强叠前道集的聚焦性。

2.2 实际资料处理

图 5(a)为某气井附近的叠前道集, 由于地层吸收的影响, 深部能量衰减严重, 分辨率较低; 横向上反射振幅特征并没有随偏移距的增大而增大, 不能代表真实的AVA特征。图 5(b)是进行频散校正和与深度有关吸收补偿后的道集, 频散引起的同相轴相位畸变问题得到解决, 且深部衰减能量得到有效补偿, 分辨率得到提高; 图 5(c)是对图 5(b)结果进一步作与偏移距有关吸收补偿后的结果, 从红框中对比可以看出, 补偿后反射强度随偏移距的增大而增大, 显现出典型的含气储层AVA反射特征。

图 6是与深度有关吸收衰减不同补偿方法的结果对比, 图 6(c)是反Q滤波方法同时进行振幅补偿和频散校正的结果。可以看出, 由于频散校正是无条件稳定的, 而振幅补偿是不稳定的, 所以该方法在补偿深部信号的同时, 高频噪音也被相应放大, 严重降低了补偿后地震数据的信噪比; 而本文方法的补偿结果如图 6(b)所示, 在补偿衰减能量的同时, 有效抑制了高频噪音的放大, 使补偿后的地震数据具有更高的分辨率和信噪比。

3 结束语

基于匹配追踪算法, 提出一种叠前资料的吸收补偿和频散校正方法。首先, 对叠前道集进行匹配分解, 用零偏移距匹配子波求取随匹配子波中心时间和中心频率乘积变化的匹配子波振幅衰减曲线, 以此获取匹配子波振幅增益系数, 对与深度有关吸收进行补偿; 然后通过零偏移距地震数据匹配子波振幅增益系数和非零偏地震道集角度计算非零偏移距匹配子波振幅增益系数, 对与偏移距有关的吸收衰减进行补偿。本文中方法通过降维处理, 将二维信息时间和频率降低为一维信息时间和频率的乘积, 增加了算法的稳健性。对于速度频散, 同样通过匹配子波求取地层品质因子, 从而进行频散校正。模型数据和实际资料的处理结果表明, 本文中方法有效补偿了深部地层的衰减能量, 极大地改善了地震数据的AVA反射特征, 有利于提取正确的AVA信息, 为叠前反演提供更可靠的道集数据。