中国西部沙漠地区的油气勘探中, 近地表地层多为疏松的沙砾层, 其对地震波的吸收衰减作用严重, 因此导致地震资料的分辨率低, 频带窄。目前, 提高地震资料分辨率的方法很多, 最常见的方法为地表一致性反褶积^[1]。人们根据不同的估计地震子波的方法, 发展出谱模拟反褶积^[2]、预测反褶积^[3]等多种方法。但实际地震资料处理中, 子波是未知的, 且通常难以满足上述方法对子波估计的假设, 所以反褶积方法对于消除近地表影响的能力有限。

为消除因近地表低速层引起的地震波吸收衰减作用, 研究学者提出了假设条件简单、结果稳定的品质因子Q补偿方法。WANG^[4]提出的稳定补偿方法, 很好地解决了高频端补偿过度的问题。在此基础上, 越来越多的学者重视对地层品质因子Q的估算以及反Q滤波稳定补偿的研究, 并将研究成果有效地应用于实际生产。Q值估算的方法包括上升时间法、振幅衰减法、谱比法、解析信号法、频移法等^[5]。王晓涛等^[6]利用地震波能量的指数衰减速率估算Q值。蒋立等^[7]对比了地表一致性反褶积方法与近地表Q补偿方法的适用条件及稳定性, 实验结果证明戈壁区适合采用地表一致性反褶积方法, 沙漠区适合采用近地表Q补偿方法。石战结等^[8]利用大炮震源微测井数据分析了沙漠区近地表地震波衰减规律。翟桐立等^[9]采用基于井地一体化测量的观测系统来估算Q值, 改变了传统的微测井测量方式, 解决了虚反射、近场效应的问题。LI等^[10]利用井地联合观测系统减小了检波器耦合效应对Q值估算的影响。张固澜等^[11]利用零偏移距VSP数据证实了地震波在传播过程中的速度频散, 并利用反Q滤波相位补偿有效消除了地震子波的速度频散, 使得地面地震资料与零偏移距VSP走廊叠加剖面更加吻合。宋吉杰等^[12]提出一种基于信息融合的高精度Q值估算的方法, 该方法将调查点处的实测Q值与采用地震初至折射波法估算的相对Q值进行融合, 来构建高精度Q场。李伟娜等^[13]提出了一种双线性回归品质因子Q估计算法, 并在实际应用中获得了与近地表速度分层一致且较为稳定的、地质意义明确的近地表Q估计值。郭平^[14]采用峰值频移法解决了地震数据频散和子波空间不一致的问题, 并在实际生产中得到了良好的应用。苏勤等^[15]将质心频移法和峰值频移法的优点相结合构建了近地表Q场, 然后进行Q值空变补偿, 并且得到了良好的实际应用效果。

本文主要针对中国西部沙漠区的地震数据开展Q补偿应用研究。为了提高沙漠区地震数据成像精度, 首先选取谱比法和峰值频移法进行对比分析, 确定适合的Q值估计方法; 然后对比相位补偿、振幅补偿和相位振幅同时补偿3种方法的应用效果, 明确适用于研究区地震数据的Q补偿方法; 最后总结出一套近地表Q补偿的技术流程, 并将其应用于实际地震数据, 以提高地震数据的分辨率并改善地震数据成像效果。

1 方法技术

品质因子Q可以定量描述地下介质对地震波的吸收衰减作用, 是地层的属性之一, 反映了地层岩性的物理特性^[5]:

$ Q = 2{\rm{ \mathit{ π} }}\frac{E}{{\Delta E}} $

(1)

式中: ΔE表示在一个波长距离内地震波损耗的能量; E为总能量。由(1)式可知, Q值越大, 地层对地震波的吸收作用越小, 反之Q值越小, 地层对地震波的吸收作用越大。我们分别采用谱比法和峰值频率偏移法计算地层的Q值, 并优选出最适合研究区的Q值估算方法和稳定Q补偿方法, 用于实际地震数据处理。

1.1 谱比法

谱比法应用范围广泛, 其原理是利用子波衰减前后的振幅谱比率来估算Q值。以单层模型为例, t₁时刻子波的振幅谱B(f, t₁)可表示为:

$ B(f, {t_1}) = A({t_1}){B_0}(f)\exp ( - \frac{{{\rm{ \mathit{ π} }}f{t_1}}}{Q}) $

(2)

式中: B₀(f)表示初始时刻地震波振幅谱; A(t₁)表示与地层吸收作用无关的振幅因子; f表示频率。因地震子波穿过单层模型时受吸收衰减作用的影响, 故t₂时刻地震子波相对于t₁时刻地震子波能量衰减, t₂时刻振幅谱表示为:

$ B(f, {t_2}) = A({t_2}){B_0}(f)\exp ( - \frac{{{\rm{ \mathit{ π} }}f{t_2}}}{Q}) $

(3)

式中: A(t₂)表示与地层吸收作用无关的振幅因子。地震子波衰减前、后振幅谱的比率为:

$ \frac{{B(f, {t_2})}}{{B(f, {t_1})}} = \frac{{A({t_2})\exp ( - \frac{{{\rm{ \mathit{ π} }}f{t_2}}}{Q})}}{{A({t_1})\exp ( - \frac{{{\rm{ \mathit{ π} }}f{t_1}}}{Q})}} $

(4)

可简写为:

$ R = {S^ * }\exp ( - \frac{{{\rm{ \mathit{ π} }}f\Delta t}}{Q}) $

(5)

式中: Δt=t₂－t₁; R=B(f, t₂)/B(f, t₁), 表示相对衰减系数; S=A(t₂)/A(t₁), 表示与地层吸收衰减作用无关的量。

将(4)式两边同时取自然对数:

$ \ln [\frac{{B(f, {t_2})}}{{B(f, {t_1})}}] = \ln [\frac{{A({t_2})}}{{A({t_1})}}] - \frac{{{\rm{ \mathit{ π} }}f\Delta t}}{Q} $

(6)

将ln[A(t₂)/A(t₁)]看作常数, 则(6)式可看作ln[B(f, t₂)/B(f, t₁)]与频率f的一次函数, 进行线性拟合得到与地层品质因子Q有关的斜率k, 可通过斜率k推导出品质因子Q的表达式为:

$ Q = - \frac{{{\rm{ \mathit{ π} }}\Delta t}}{k} $

(7)

利用(7)式可快速估算地层的Q值。通过该方法计算出的Q值, 为等效Q值, 只与地层介质的吸收性质有关, 排除了非固有衰减对Q值估计的影响。

1.2 峰值频移法

对大多数介质而言, 地震子波衰减随频率的增加而增加, 且高频成分比低频成分衰减得更快, 因此, 信号的频率中心在传播过程中向低频端偏移^[14](图 1)。图 1中红线表示初始状态子波频谱, 蓝线表示经过地层吸收衰减后的子波频谱, 衰减后的子波频谱向低频端偏移, 且振幅变小, 其中f_m表示衰减前峰值频率, f_p表示衰减后峰值频率。

ZHANG等^[16]利用地震波在频率域衰减后主频向低频端偏移的特点, 推导了Q值与峰值频率的关系:

$ Q = {\rm{ \mathit{ π} }}t\frac{{{f_{\rm{p}}}f_{\rm{m}}^2}}{{2(f_{\rm{m}}^2 - f_{\rm{p}}^2)}} $

(8)

如果f_m已知, 则可以根据(8)式计算出品质因子Q, 但实际地震数据中, 因衰减前频谱峰值频率f_m是未知的, 故可以利用质心频率f_c根据如下公式换算得到。

$ {f_{\rm{c}}}{\rm{ = }}\frac{{\int_0^\infty {fA(f){\rm{d}}f} }}{{\int_0^\infty {A(f){\rm{d}}f} }} $

(9)

$ {f_{\rm{m}}}{\rm{ = }}\frac{{\sqrt {\rm{ \mathit{ π} }} {f_{\rm{c}}}}}{2} $

(10)

1.3 稳定Q补偿方法

由于振幅补偿量随频率和时间的增加而增加, 一般情况下, 补偿算法主要对地震数据的高频端进行补偿, 这使得其中的噪声被过量补偿, 导致地震数据失真。因此, 本文使用WANG^[4]提出的稳定Q补偿方法, 对完全吸收的频率成分不再进行补偿, 这样可以避免对高频噪声的过量补偿。该方法同时对振幅进行补偿和相位调整, 可解决由地层吸收造成的地震波能量损失和相位一致性问题。稳定反Q滤波如下:

$ \begin{array}{c} U(\tau + \Delta \tau , \omega ) = U(\tau , \omega )\exp [{(\frac{\omega }{{{\omega _{\rm{h}}}}})^{ - \gamma }}\frac{{\omega \Delta \tau }}{{2Q}}] \cdot \\ \exp [{\rm{i}}{(\frac{\omega }{{{\omega _{\rm{h}}}}})^{ - \gamma }}\omega \Delta \tau ] \end{array} $

(11)

式中: U(τ, ω)表示未经补偿的频率域地震数据, 即对原始地震数据进行傅里叶变换的结果; U(τ+Δτ, ω)表示经过振幅和相位补偿后的频率域地震数据; τ为旅行时, Δτ为旅行时增量; ω为原始地震数据角频率; ω_h为中心频率, 作为调整参数与地震波频带的高频率相关; γ=(1/π)Q^-1。

在实际情况中即使输入地震数据是无噪声的, 随着Q值的降低和成像时间的增加, 全反Q滤波器将导致不稳定, 并在地震数据中产生较强的伪影。因此利用如下的稳定性条件公式解决Q补偿的不稳定性:

$ \mathit{\Lambda} (\omega ) = \frac{{\beta (\omega ) + {\sigma ^2}}}{{{\beta ^2}(\omega ) + {\sigma ^2}}} $

(12)

式中: σ²为稳定因子; Λ(ω)=exp[(ω/ω_h)^－γωΔτ/2Q], Λ(ω)表示稳定的振幅补偿量, 如果Λ(ω)=1, 则表示仅对相位的反Q滤波器补偿是无条件稳定的。其中, 参数β(w)和σ²的表达式为:

$ \beta (\omega ) = \exp [ - {(\frac{\omega }{{{\omega _{\rm{h}}}}})^{ - \gamma }}\frac{\omega }{{2Q}}t] $

(13)

$ {\sigma ^2} = \exp [ - (0.23{G_{{\rm{im}}}} + 1.63)] $

(14)

式中: G_im是增益限制, 为可调节参数。

1.4 近地表Q补偿处理流程

近地表Q值估算方法众多, 但实际应用中很难找到一种适用于所有情况的方法。本文基于叠前地震数据, 利用谱比法和峰值频移法进行Q值估算, 最终选取了适用于研究区的Q值估算方法。处理流程如下: ①建立近地表速度模型, 由静校正量和高速层顶高程换算近地表旅行时; ②拾取叠前地震数据的主频和均方根振幅值, 得到相对衰减系数R值; ③利用主频和近地表旅行时采用频移法估算Q值, 利用相对衰减系数R和近地表旅行时采用谱比法估算Q值, 并与实测Q值标定; ④选取近地表Q补偿最佳参数。利用Q补偿参数与补偿方法实现近地表Q补偿。具体流程如图 2所示。

2 应用实例 2.1 工区概况

工区位于中国西部塔里木盆地, 地表较为平坦, 如图 3a所示, 工区内最大高程差约为27 m, 地形呈西高东低、北高南低之势, 为典型的沙漠环境。近地表结构横向上复杂多变, 主要有红柳植被区、沙丘浮土区、水域沼泽区和硬板地区。图 3b为低降速带厚度分布情况, 红柳植被区的低降速带厚度最大, 约为80~88 m。水域沼泽区的低降速带厚度最小, 约为33~50 m。沙丘浮土区a与硬板地区的低降速带厚度约为54~70 m。沙丘浮土区b的低降速带厚度约75~85 m。工区地层纵向结构简单, 低速层速度通常为400 m/s, 高速层速度通常为1 600~1 900 m/s, 高速层顶界面随地势由东向西、由南向北缓慢抬升。不同近地表地层类型的疏松情况不同, 因此地层对地震波的吸收衰减作用不同。当近地表介质对地震波的吸收衰减作用积累到深层时, 出现地震子波变形, 地震剖面分辨率降低的情况, 从而影响地震成像质量。

2.2 Q值计算

根据图 2中的近地表Q补偿处理流程, 首先建立近地表速度模型, 求得静校正量和高速层顶高程, 进而换算得到近地表旅行时, 如图 4a所示, 近地表旅行时与厚度之间几乎呈正相关关系, 厚度越大, 旅行时越大, 厚度越小, 旅行时越小, 说明低降速带的横向速度变化不大。然后利用叠前地震数据计算主频(图 4b), 采用峰值频移法估算近地表Q值; 最后利用均方根振幅, 换算相对衰减系数R(图 4c), 采用谱比法估算近地表Q值。

对比分析图 4中3幅图的关系, 结果如表 1所示。假设不考虑其它因素的影响, 则旅行时越大, 地层对地震波的吸收衰减作用就越强, 对应的主频和相对衰减系数越小。表 1中对比组1中大体存在此趋势, 只有水域沼泽区不符合此趋势, 其近地表旅行时高于红柳植被区, 低于硬板地区, 但其主频最小、相对衰减系数最小, 衰减最严重。我们大胆推测此现象可能与水域沼泽区近地表地层含水有关, 地层岩石孔隙中的流体改变了地层对地震波的吸收衰减作用, 岩石的含水量增大, 增加了岩石的非弹性, 进而增强了地震波的衰减^[17]。对比组2中, 地表结构类型虽然相同, 但由于沙丘浮土区b位于水域沼泽的下游, 也受地层含水的影响, 所以沙丘浮土区b旅行时高于沙丘浮土区a, 其主频和相对衰减系数反而高于沙丘浮土区a。将不受水域影响的区域进行对比, 在对比组3中, 我们发现旅行时越大, 主频和相对衰减系数越小, 地层对地震波的吸收衰减作用越强, 这与理论规律相符。

根据(7)式和(8)式分别求得近地表Q值, 结果如图 5所示, 整个研究区近地表Q值范围为1~4。从平面分布结果上看, 在水域沼泽区利用峰值频移法和谱比法分别求得的Q值均为最小, 吸收衰减作用最严重; 在红柳植被区, 利用峰值频移法求得近地表Q值小于利用谱比法求得的近地表Q值。红柳植被区含水量比水域沼泽区低, 但高于沙丘浮土区b, 所以近地表Q值应大于水域沼泽区, 小于沙丘浮土区b。因此, 我们认为峰值频移法求得的近地表Q值适用于本研究区。

从红柳植被区任意抽取一炮, 利用峰值频移法和谱比法分别求得近地表Q值, 进行相同的稳定近地表Q补偿, 结果如图 6所示。将两种近地表Q补偿的结果与补偿前单炮记录对比, 发现利用峰值频移法估算的近地表Q补偿后的地震数据分辨率更高, 同相轴更清晰, 同时信噪比也得到了一定程度的提高。选取红框中的目的层反射同相轴进行频谱分析, 结果如图 7所示, 可以看出, 相较谱比法, 峰值频移法对高频端的拓频更为有效。

2.3 稳定Q补偿方法应用效果

图 8为不同Q补偿方法得到的地震剖面与原始剖面, 图 9为在图 8中目的层附近提取的地震子波。图 8a为补偿前地震剖面, 图 8b只进行了相位校正, 红色箭头所指的同相轴连续性增强; 从图 9中黄色曲线提取的子波可以看到, 地震子波未被压缩。图 8c通过振幅补偿, 得到的断点更加清晰, 由图 9中绿色曲线可以看出, 地震子波被压缩, 但是蓝色箭头所指的同相轴分辨率变低。图 8d综合了两种方法的优势, 既压缩了地震子波(图 9中蓝色曲线), 又实现了相位校正, 地震剖面上的断点更加清楚, 地震剖面的分辨率也更高。数值计算结果表明采用振幅相位同时补偿的方法更合理。

2.4 成像效果对比

图 10为西部沙漠区地震数据近地表Q补偿前、后的时间偏移剖面, 可以看出, 地震数据分辨率明显提高, 断点清晰, 强反射轴下方的波谷得到了有效刻画, 达到了既定的处理目标, 这对后续的解释工作具有重要意义。图 11为图 10目的层补偿前、后频谱对比结果, 同样表明本文方法可以有效拓宽地震数据的高频分量, 同时保持低频分量。

近地表Q补偿前、后地震数据相干切片如图 12所示, 近地表Q补偿后的相干切片整体分辨率得到了提高, 特别是红圈标识的位置, 清晰地刻画了断层, 突破了前期地震数据处理的效果, 对储层预测具有重要意义。

图 13为近地表Q补偿前、后联井纵波阻抗剖面, 可以看出, 近地表Q补偿前井数据T1层位置与波阻抗界面吻合, 但T2、T3、T4层未出现与之对应的波阻抗界面; 近地表Q补偿后地震剖面分辨率明显提高, T1与T4之间的薄层清晰可见, 与井资料高度吻合, 层位识别精度明显优于近地表Q补偿前的联井剖面。利用本文方法可以对复杂地下介质地层连续变化的区域进行较为精确的刻画, 且在纵向分辨率提升的基础上保持良好的横向连续性。

3 结论

本文对比分析了谱比法和峰值频移法两种Q值估算方法的原理, 总结了一套近地表Q补偿流程, 并结合实际工区数据, 分析了两种方法的应用结果。

1) 对于本工区受地层含水影响的区块, 其近地表旅行时越小, 对地震波的吸收衰减作用越强, 证明了近地表低速层中的流体对地震波的吸收衰减作用有增强效果。

2) 将近地表Q补偿应用于本文工区的地震单炮记录, 补偿结果表明峰值频移法分辨率更高, 高频端展宽效果更好, 更适用于本文工区。相位振幅同时补偿的方式结合了相位补偿和振幅补偿的优势, 成像效果更好。

3) 将近地表Q补偿后的叠前地震数据应用于后续的叠前时间偏移、相干分析等, 其结果较前期地震处理结果断点更清晰, 分辨率更高, 目的层得到了精确成像, 从频谱分析结果也可以看出, 高频端得到了展宽。补偿后的相干切片分辨率明显提高, 补偿前无法刻画的断层在补偿后的相干切片上得到了清晰的显示。

本文通过频谱图、成像剖面以及相干切片定性地评价了近地表Q补偿后的地震成像效果, 在下一步研究中, 还需进一步讨论如何定量地评价近地表Q补偿结果。