Q值是描述波传播过程中介质吸收特性的重要参数之一。在地震勘探中,Q值可用于分析地下介质的衰减特征,提供能量补偿信息,并辅助识别岩性和孔隙流体。传统的Q值估算方法大多基于零偏移距VSP和叠后地震数据。Tonn^[1]利用合成VSP地震数据对几种Q值估计方法进行了对比分析,认为对数谱比(logarithmic spectral ratio,LSR)方法抗噪能力较弱。为减少噪声对Q值估计的影响,多位学者采用了积分、正则化和滤波等手段提高Q值估计方法的抗噪能力^[2-4]。Cheng等^[5]提出了对数谱面积双重差(logarithmic spectral area double difference,LSADD)方法。该方法通过两次计算对数谱面积差异来估计Q值,提高了方法的抗噪能力。尽管该方法具有较高的准确性和和抗噪能力,但目前仅适用于VSP和叠后地震数据。叠后地震数据是经过动校正、叠加后得到的地震数据,利用叠后数据进行Q值估计时都假设地震波是按垂直方向进行传播的,叠后数据相当于自激自收地震记录。这和实际地震波的传播路径有较大的差异,由于野外实际地震数据接收时都是以一定的炮检距进行接收的,同时叠加的平均效应也会改变波形和波的衰减特征,因而基于叠后数据的Q值估计必然有较大的误差。

相比叠后地震数据,叠前数据保留了原始地震波形,可用于更准确地估算Q值。利用叠前数据进行Q值估计时,应充分考虑不同接收点地震数据的地震波传播路径的差异。叠后Q值估计方法假设地震资料是零炮检距接收的,不适应叠前数据。一些研究人员针对叠前数据开展了Q值估计方法研究,这些方法可以分为三大类。第一类是基于共中心点道集的Q值随炮检距变化(Q versus offset,QVO)方法^[6-8]。Hackert等^[9]借助测井数据校正频谱对QVO方法进行改进,降低了薄层效应对衰减定量计算的干扰。Zheng等^[10]利用四阶累积量匹配从叠前数据中提取瞬时子波谱,并结合QVO方法估算Q值,降低了反射系数对估算结果的影响。Wu等^[11]借助广义S变换,提出了连续谱比斜率法,增强了QVO方法的稳定性。该类QVO方法抗噪能力较弱,且容易受到动校拉伸和射线路径差异的影响。第二类方法是由Zhang等^[12]提出的峰值频率漂移(peak frequency shift,PFS)和层剥离技术,该方法在雷克子波的假设下,推导频率峰值随旅行时变化的解析式,利用层剥离技术从CMP道集估算Q值。Behura等^[13]采用层剥离技术,估算了反射波地震数据的层间各向异性衰减量。然而,层剥离技术在迭代计算Q值时容易产生累积误差,并且此类方法基于直射路径理论计算视Q值,对于大偏移距和深层的Q值估计有较大偏差。第三类方法主要是利用速度信息,基于坐标变换和波场延拓等方法消除射线路径对叠前Q值估计的影响,然后采用常规方法进行Q值估计。Reine等^[7-8]提出基于τ-p变换的叠前数据Q值估计方法,利用τ-p变换将同一水平慢度下的反射投影到t-x域,消除射线路径的影响,并改进LSR方法,利用反演思想估计Q值。在Reine方法的基础上,Beckwith等^[14]提出了一种频率相关的衰减因子估算方法。Oliveira等^[15]基于重建基准面算法,利用速度信息进行波场延拓,校正了射线路径引起旅行时间差异对Q值估算的影响,估计的Q值精度更高。坐标变换类方法能充分消除射线路径对Q值估计的影响,相较于前两类方法精度更高。但该类方法往往需要预先知道较准确的地层速度,速度模型的误差对Q值估计影响较大。

径向道(radial trace,RT)变换由Claerbout提出^[16-17],最初主要用于偏移和成像,Li等^[18]将地震记录扩展到广义法自激自收模型,从而减小射线路径的影响。为了充分考虑射线路径对Q值估计的影响,本文在径向道变换域研究了叠前Q值估计方法。利用RT变换技术将叠前数据从炮检距—旅行时(offset-time,X-T)域变换到视速度—旅行时(apparent velocity-time,R-T)域,即扩展叠前数据到广义自激自收模型,减小上覆地层射线路径差异对Q值估计的影响。在R-T域中推导了Q值随视速度变化关系(Q versus apparent velocity,QVAV)。将QVAV关系与LSADD方法相结合,提出了QVAV_LSADD方法。通过选择不同频带分割点来计算出多个截距值,将多个截距值的均值作为最终估计,从而克服了叠前数据信噪比较低引起的估计值不稳定问题。

波在非完全弹性介质中传播的吸收衰减表达式为:

式中:C表示与频率无关的振幅干扰项,包括透射损失、反射系数和几何扩散因子等其他因素的影响;A₁(f)与A₂(f)是两个不同时间的地震波的振幅谱;f是频率;Δt是旅行时;Q是两个界面之间的品质因子。

对式(1)两边取自然对数:

式(2)即为常用的LSR方法,利用线性拟合求取斜率值即可计算地层Q值。LSR方法应用的主要限制是其抗噪性较弱。

为了准确测量两个反射界面之间的Q值并消除与频率无关的振幅干扰项的影响,Cheng等^[5]提出了LSADD方法。该方法通过计算两个界面之间的对数谱面积差(logarithmic spectral area difference,LSAD),并设置一个分割点将LSAD分为两部分。通过对两部分LSAD的差值运算,可以消除干扰项的影响,从而准确估算地层的Q值。首先,离散化后的LSAD表达式为:

式中:N表示在频带范围[f_min,f_max]内对频率f离散采样后的采样点数;Δf是频率采样间隔。

设置频带分割点f_m将LSAD分为两部分,在频带分割点前一部分有m-1个频率采样点,在分割点后面有N-m个采样点,我们定义前一部分的对数谱面积为:

同理,定义频带分割点后一部分对数谱面积为:

观察式(4)和式(5),对两个LSAD再求差,从而消除了幅度干扰项C的影响,其表达式如下:

借助带宽[f_min,f_max],我们可以简化式(6)的右侧项,从而获得Q值与双差之间的关系:

观察上述推导过程可以很容易发现,通过计算对数谱面积的差值,使用式(7)可以估算出精确的Q值。由于求和运算能够很好地压制随机噪声,因LSADD方法具有很好的抗噪能力。

由于叠后地震数据假设炮检距为零,第L层反射界面的反射波是第L-1层反射界面反射波继续传播经过第L层的结果。因而可以利用第L-1、L层反射界面反射波频谱由式(7)计算出第L层的Q值。然而,由于射线路径差异的影响,对于叠前数据,对于同一道数据相邻两点的LSADD与Q值之间的关系比式(6)要复杂得多。主要原因在于,同一偏移距上两个反射界面的反射角和射线路径不同,无法消除上覆层射线路径的差异,从而导致Q值与LSADD之间呈现非线性且复杂的关系。如果仍然使用叠后方法,估算的Q值将随偏移距而变化,导致显著的估算误差。为了准确估算Q值,本文考虑使用RT变换来消除射线路径的影响。

RT变换R能够将一个X-T域的地震波场s(x,t)映射到R-T域波场$\stackrel{~}{s}$(v_a,t'),而其逆变换R^-1能将$\stackrel{~}{s}$(v_a,t')映射回X-T域,即:

式中:x为炮检距;t为双程旅行时;v_a为视速度;t'为映射后的双程旅行时。RT变换就是一个坐标变换,假设RT变换原点坐标为(x₀,t₀),对于X-T域中的一个点(x,t),变换到R-T域后变为(v_a,t'),v_a及t'利用式(10)进行计算:

这就实现了RT正变换。同样,对于R-T域中的点(v_a,t')通过下式实现RT逆变换:

我们设计了图1所示的地层模型,炮检距范围为0~1 200 m,道间距为10 m,根据Zoeppritz方程计算不同炮检距的反射系数,并利用45 Hz的雷克子波合成CMP道集,如图2所示。在0~3 000 m/s的范围内计算301个视速度,将合成地震数据转换到R-T域,如图2b所示,可看出同相轴在R-T域变为近似垂直的轴。将相同视速度下R-T域的反射映射到X-T域,图2中的A、B、C、D点展示了视速度为1 000 m/s时对应的映射关系。X-T域的径向射线(红色虚线)上的A、B两点映射为R-T域同一视速度(1 000 m/s)上的C、D两点。

图3显示了RT变换域Q值分析的优势。为了估计第二层的Q值,我们要分析地震子波进入第二层前后的频谱。叠后Q值估计方法是分析一道上前后反射波的频谱来估计Q值,也就是利用E点接收的F点反射波和C点反射波进行分析,这和实际地震波的传播过程相差较大。准确的Q值估计应该严格按照地震波的传播路径进行分析,即分析P点接收到的B点反射波和E点接收到的C点的反射波,通过分析这两个反射波的频谱来估计Q值。但是P点和E点之间的相对位置关系要利用地下介质的层速度才能得到。RT域Q值估计是利用Q点接收的G点反射波和E点接收的C点反射波,由图3可以看出,它和准确的Q值估计过程非常接近,它相当于假设地震波是按直线传播的。因而在RT域进行Q值估计可以降低射线路径的影响。

前面使用一个3层的简单模型说明了RT变换在减小射线路径差异对Q值估计影响方面的优势。然而,利用RT变换映射后反射波的射线路径与直射线路径一致。对于多层模型而言,对深层大偏移距的射线路径差异消除不够充分。针对这一劣势,我们利用QVO思想推导了QVAV关系,从而给出了零视速度处的地层Q值估计公式。

我们设计了如图4所示的五层模型,炮检距范围为0~1 200 m,间隔10 m。同样,我们根据Zoeppritz方程计算反射系数,利用45 Hz的Ricker子波合成CMP道集。X-T和R-T域的CMP道集分别如图5和5b所示。

使用LSADD和LSR方法对已映射的地震道进行Q值估计,并与PFS方法进行比较。图6绘制了不同层的Q值估算结果。随着层数的增加,大炮检距处的射线路径影响消除不足。主要原因是RT变换映射后的射线路径类似于直射线路径,与真实的曲射线路径相比,当传播层数和炮检距增加时,会出现更显著的误差。由于PFS方法也是基于直射线理论,所以也会出现类似的现象。这种现象会导致深层和大炮检距处的Q值估计存在偏差,降低Q值估计的精度。因此,受QVO方法的启发,本文提出了QVAV方法来减弱这种影响。

QVO方法是在NMO后的CMP道集上应用的,该方法计算每个炮检距的对数谱比斜率,通过线性拟合求取截距值,进而估算零炮检距Q值。如上文所述,R-T域的射线路径可以近似为具有一定入射角的广义自激自收模型,类似于NMO后的CMP道集。此外,RT变换减小了上覆射线路径对Q值估算的干扰,优于QVO理论的近似假设。因此,我们受到QVO思想的启发,推导了在R-T域CMP道集中Q值随视速度的变化关系。

首先,当RT变换原点为(x₀,t₀)=(0,0)时,炮检距x和视速度v_a之间的关系满足:

根据水平层状介质中t与x的近似关系,有:

式中:v_r为均方根(root-mean-square,RMS)速度;t₀为零炮检距旅行时。

将式(12)代入式(13),有:

式中:t_0,n和t_a,n分别表示零视速度(炮检距)和视速度为v_a时从地表到第n层的旅行时;v_r,n表示从地表到第n层的RMS速度。对式(14)进行移项化简,并级数展开,可得到反射波旅行时t_a,n与视速度v_a的关系如下:

当视速度或者炮检距较小时,关于${v}_{a}^{4}$/${v}_{r,n}^{4}$之后的高阶项可以忽略掉,因此有:

由式(2)知,第n-1层与第n层之间的零视速度(炮检距)处的频谱比斜率值(slope of spectral ratio,SSR)${p}_{0}^{n}$为:

式中:Δtⁿ₀和Qⁿ₀分别为零视速度(炮检距)处第n-1层到第n层之间的旅行时差和层Q值。由于RT变换已经大幅减弱了射线路径的影响,对于视速度v_a处的SSR,满足:

式中:${Q}_{a}^{n}$表示与视速度v_a相关的第n层的拟层Q值。

假设在一个小排列的情况下,不同视速度v_a处${Q}_{a}^{n}$变化不大,可以近似为零视速度处的${Q}_{a}^{n}$值,式(18)可近似为:

由式(19)可知,求取的斜率值将随着视速度的平方呈线性变化,我们称之为Q值随视速度的变化关系(QVAV)。通过将不同视速度v_a下的SSR值做线性回归,我们可以得到截距值,即零视速度的SSR,进而可以求得地层Q值。我们称这一方法为QVAV方法。

基于QVAV方法进行叠前Q值估计的计算步骤如下:

1)在给定的视速度范围内,按照设置的视速度采样间隔选取v_a1,v_a2,…,v_an,对CMP道集进行RT变换,将其转换到R-T域;

2)利用式(10),将一个恒定视速度v_ai下R-T域的反射(v_ai,t'_i)映射到X-T域的(x_i,t_i);

3)选择相同的上覆地层射线路径(恒定v_ai)的两个反射界面对应X-T域的两道;

4)对两道数据进行时频分析得到顶底界面的振幅谱,并计算不同v_ai时的SSR值;

5)将SSR值关于视速度的平方做线性回归,得到截距,即零视速度的SSR,并根据式(16)计算地层Q值。

LSR方法是QVAV方法的基础,其抗噪能力较弱。对于含噪数据而言,不同视速度v_a下的SSR值随噪声波动剧烈,严重影响Q值估计的稳定性。然而,LSADD方法具有较强的抗噪性,因此我们受到LSADD方法的启发,提出了QVAV_LSADD方法,实现了对低信噪比的叠前数据Q值的准确估计。

首先,由式(6)可知,当频带分割点为f_m(m=1,2,…,k)时,第n-1层与第n层之间的零视速度处的对数谱面积双重差${\left.\Delta {{\stackrel{~}{S}}^{n}}_{0}\right|}_{{f}_{m}}$为:

由于RT变换已经大幅减弱了射线路径的影响,对于视速度v_a处的对数谱面积双重差,满足:

同理,式(21)可以被近似为:

由式(22)可知,求取的Δ$\stackrel{~}{S}$将随着视速度的平方呈线性变化。因此,通过将不同视速度v_a下的Δ$\stackrel{~}{S}$做线性回归,并求取截距值,即bⁿ|_fm,得到:

式(23)是当频带分割点为f_m时,求取的截距值,利用其可以估算地层Q值。根据LSADD理论,Δ$\stackrel{~}{S}$是关于分割点f_m的函数,通过滑动选取多个分割点取均值可以提高抗噪性。利用这一优势,分别计算不同f_m(m=1,2,…,k)时的截距bⁿ|_fm,并取其平均作为最终估计,进而增强叠前Q值估算的抗噪性。最终Q值计算公式如下:

我们称式(24)的方法为QVAV_LSADD。

下面,给出基于QVAV_LSADD方法估计Q值的步骤:

1)在给定的视速度范围内,按照设置的视速度采样间隔选取v_a1,v_a2,…,v_an,对CMP道集进行RT变换,将其转换到R-T域;

2)利用式(10),将一个恒定视速度v_ai下R-T域的反射(v_ai,t'_i)映射到X-T域的(x_i,t_i);

3)选择相同的上覆地层射线路径(恒定v_ai)的两个反射界面对应X-T域的两道;

4)对两道数据进行时频分析得到顶底界面的振幅谱;

5)固定一个分割点f_m,计算不同视速度v_ai下的Δ$\stackrel{~}{S}$;

6)将Δ$\stackrel{~}{S}$关于视速度的平方做线性回归,得到截距bⁿ|_fm;

7)滑动选取多个分割点f_m,重复步骤5)~6)计算不同的截距值,利用式(24)得到最终估计的Q值。

我们分别使用QVAV_LSADD、QVAV、QVO和PFS方法对图5中的共中心点(CMP)道集进行Q值估计。图7展示了各种方法Q值估计的结果。

可以看到:QVAV类方法估计的Q值都非常准确,而QVO方法由于没有消除上覆射线路径的影响,估计结果与真实值有一定的偏差。另一方面,由于PFS方法采用直射线理论假设,在深层处与真实射线路径误差较大,进而导致估计的Q值与真值相差较大。同时,PFS方法在估计Q值时采用了逐层剥离的思想,Q值误差的逐层累积也会影响深层的Q值估计。总体而言,本文提出的QVAV类方法在大炮检距和多层地质模型下能够获得较高的估计精度,充分证明了该方法的有效性。

为了验证本文方法在实际数据应用中的可行性和实用性,我们对实际的CMP集进行Q值估计。图8为经过去噪预处理之后的CMP道集,图8b展示了变换到R-T域之后的道集。选择图8b中红线所示的两层层位,取恒定视速度下的不同反射层位点,将其映射到X-T域,即为图8中所示层位,定位旅行时间及反射波位置,从而进行Q值估算。

选择了5个反射界面分别位于3 380、3 780、4 100、4 300和4 780 ms。利用QAVA_LSADD和QVAV方法得到了层Q值和等效Q值曲线,与RMS速度曲线共同展示在图9中。我们可以看到,QVAV_LSADD方法在深层的Q值估计比QVAV方法有更小的扰动。这一结果与我们模型测试实验的结论相似,证明QVAV_LSADD方法在信噪比较差的深层有较好的估计结果。此外,我们还可以看出Q值趋势与速度曲线具有一定的相关性,这进一步证明了QVAV_LSADD方法在应用于实际地震数据时的有效性。

准确的Q值估计对于提高地震资料的分辨率及储层预测都有重要的意义。本文针对利用叠前地震资料提取Q值问题,研究了一种RT变换叠前Q值估计方法。该方法通过RT变换减小射线路径的影响,并结合LSADD方法提高抗噪性。在RT变换域,推导了Q值与视速度之间的变化关系,将QVAV与LSADD方法结合,通过滑动选择多个分割点取平均值来提高抗噪性。理论模型数据测试结果表明该方法精度高、抗噪性强,实际资料验证了该方法的有效性和实用性。