基于零偏移距VSP的时变子波反褶积方法

赵 岩，毛宁波

（1.油气资源与勘探技术教育部重点实验室（长江大学），武汉 430100；2.长江大学地球物理与石油资源学院，武汉 430100）

0 引言

由于地层的吸收衰减作用，地震波在传播过程中会逐渐衰减，地面地震数据通常具有有限带宽，这可能会导致无法成像或无法识别细微的地层特征［1-3］。影响地面地震数据分辨率的主要因素之一是近地表松散风化带对地震波的强烈吸收，然而，与地面地震勘探不同，VSP 是在井中接收地震波，其数据受近地表的影响较小，因此，VSP 数据通常具有比地面地震数据更宽的频带以及更高的信噪比和分辨率［4-6］。目前有4 种方法可以通过VSP 数据来提高地面地震数据的分辨率:①从VSP 数据中提取品质因子Q［7-8］，然后将反Q 滤波方法应用于叠后地表地震记录［9-11］；②高频恢复方法［12-13］，从不同深度的VSP 数据中提取滤波算子，然后对叠后地面地震数据进行高频恢复滤波；③VSP 小波匹配方法，分别从叠后地面地震数据和VSP 数据中提取子波，然后将二者的子波进行匹配［14-15］，以改善叠后地面地震数据的分辨率；④VSP 子波反褶积方法，从VSP 数据中提取子波，然后计算反褶积算子，并对叠后地表地震数据进行反褶积［14，16］。

目前的VSP 子波反褶积方法仅从VSP 地震记录中提取一个时不变的反褶积算子，但由于地下介质的黏弹性，VSP 记录的子波是逐渐衰减的，因为它忽略了子波的时变特性，所以并不准确。为了改进该方法，本次研究将零偏移距VSP 数据的下行直达波作为单程地震子波，并将其振幅谱转换为双程地震子波振幅谱，在地震子波是最小相位的假设下，从振幅谱中获得时变反褶积算子，以期其可用于对叠后地面地震数据进行时变子波反褶积。

1 方法原理

本次研究提出的时变子波反褶积方法，首先使用零偏移距VSP 资料的下行直达波来获得时变子波振幅谱，以准确反映单程地震子波的振幅谱随旅行时间的变化；然后利用直达波自身的衰减规律将单程地震子波的振幅谱转换为双程地震子波的振幅谱。该方法的步骤为:①利用零偏移距VSP 数据估算单程地震子波振幅谱；②基于单程地震子波振幅谱和地层衰减函数，估算叠后地面地震数据的双程时变子波振幅谱；③在最小相位假设下计算时变反褶积算子；④通过地面地震记录的Gabor 谱和反褶积算子的乘积得到反射系数的Gabor 谱；⑤通过反Gabor 变换计算反射系数（图1）。

图1 基于零偏移距VSP 的时变子波反褶积方法流程Fig.1 Flowchart of time-varying wavelet deconvolution method based on zero-offset VSP data

1.1 时变子波估算

式中:t为时间，s；f为频率，Hz；Q为品质因子。

定义c(f)=-πf/Q，Δt=tn+1-tn，则式（1）变为

式中:cn(f)表示tn和tn+1间的衰减函数。

对于任意时刻t(tn＜t＜tn+1)，VSP 单程地震子波的振幅谱Wvsp(t,f)可写为

将c1(f)作为0 到t1时刻的衰减函数，就可以计算0 到tn+1任意时刻的单程地震子波振幅谱。

由于VSP 地震子波是单程衰减的振幅谱，因此需要将VSP 单程地震子波的振幅谱转化为与地面地震相对应的双程地震子波的振幅谱。对于任意时刻t，当它满足tn＜t＜tn+1时，则有Wsurface(2t,f)=Wvsp(t,f)exp ⋅

式中:Wsurface(2t,f)为叠后地面地震双程地震子波的频谱。

1.2 时变地震子波反褶积

信号u(t)的Gabor 变换定义为［18］

式中:U(t,f)为信号u(τ)的Gabo r谱；g(τ)为Gabor分析时窗；t为分析时窗的中心位置，s。

地震记录的Gabor 谱可以近似写为3 个部分的乘积［19-21］，表达式为

式中:U(t,f)是地震记录u(t)的Gabor 谱；W(f)是子波w(t)的傅里叶变换；R(t,f)是反射系数r(t)的Gabor 变换。

一旦计算得到时变地震子波谱Wsurface(t,f)，在最小相位假设下，反褶积算子可写为［19］

式中:H[]表示Hilbert 变换。

根据式（6），反射系数Rg(t,f)的Gabor 谱可以由下式计算得到

再利用Gabor 反变换即可计算得到反射系数。

2 理论模型数据试算

图2（a）为层状模型和VSP 观测系统。在该模型中，品质因子Q 的数值是逐渐增加的。在100 m到2 080 m 深度内布置100 个检波器，检波器间距为20 m。震源与井口之间的水平距离为20 m，地震记录长度为1.2 s，采样间隔为0.002 s。由图2（b）可以看到，该记录中直达波的振幅和主频随深度减小，表现出单程VSP 子波的时变特性。

图2 层状模型和VSP 观测系统（a）及零偏移距VSP 下行直达波记录（b）Fig.2 A layered model and VSP geometry（a）and synthetic down-going direct wave record of zero-offset VSP（b）

从图3 可看出，直达波的振幅和主频都随着旅行时的增加而减小。利用式（1）—（3）可计算出每个时间采样点的VSP 单程子波振幅谱［图3（b）］，利用式（4），可得到双程地震子波的振幅谱［图3（c）］。比较图3（b）和图3（c）可知，双程地震子波振幅谱的衰减明显强于单程地震子波的振幅谱。为了验证从VSP 单程地震子波振幅谱计算的双程地震子波振幅谱的正确性，使用图2（a）的层状模型直接计算理论双程地震子波的振幅谱［图3（d）］。比较图3（c）和图3（d）可以看到，通过计算得到的双程地震子波的振幅谱与理论子波振幅谱几乎一致。

图3 利用VSP 下行直达波的振幅谱估算双程地震子波的振幅谱Fig.3 Estimated two-way seismic wavelet amplitude spectrum by using the amplitude spectrum of VSP down-going direct wave

图4 理论合成衰减地震记录的时变子波反褶积结果Fig.4 Result of time-varying wavelet deconvolution to synthetic attenuation seismic record

在估算得到双程地震子波振幅谱后，可以对理论合成地震记录进行时变子波反褶积（图4）。图4（b）中的理论合成衰减地震记录是由图3（d）中的时变子波和图4（a）中的理论随机反射系数褶积产生。利用式（7）和式（8）来进行时变子波反褶积，其中反褶积算子由图3（c）中所示的时变子波振幅谱计算得到。由时变子波反褶积后的地震记录［图4（c）］可以看到，时变子波反褶积后的记录更接近理论反射系数。

3 实际地震数据应用

由图5（a）可以看到，下行直达波的振幅能量很强，而上行波的振幅能量较弱，随着传播距离的增加，下行直达波的振幅在逐渐衰减。对零偏移距VSP 地震记录进行波场分离，并将下行直达波截断为单程VSP 地震子波，地震记录如图5（b）所示。

图6 为利用实际VSP 数据下行直达波的振幅谱计算得到双程地震子波的振幅谱。由图6（a），可以看到，随着深度和旅行时的增加，振幅减小，频带变窄。图6（b）和图6（c）分别为针对每个采样点计算得到的VSP 单程子波的幅度谱和双程子波的幅度谱。将二者进行对比可知，双程子波的振幅谱具有比单程子波的振幅谱更弱的振幅和更窄的频带。

在获得双程时变子波振幅谱后，在最小相位假设下计算反褶积算子。由图7（a）可以看到，随着深度的增加，地震记录的振幅能量减弱，同相轴变粗，分辨率降低。利用上述得到的反褶积算子对其进行时变子波反褶积处理［图8（b）］可以看到，时变反褶积后，中深层地震记录的能量得到有效补偿，且同相轴被压缩，波组信息更加丰富，分辨率明显提高。

图5 实际零偏VSP 地震记录和VSP 下行直达波地震记录Fig.5 Real zero-offset VSP data and VSP downgoing direct wave record

图6 利用图6（b）中实际VSP 数据下行直达波的振幅谱估算双程地震子波的振幅谱Fig.6 Estimated two-way seismic wavelet amplitude spectrum by using the amplitude spectrum of real VSP down-going direct wave of Fig.5（b）

图7 对叠后地面地震数据进行时变子波反褶积结果Fig.7 Results of time-varying wavelet deconvolution to real post-stack surface seismic data

4 结论

（1）提出一种基于零偏移距VSP 的时变子波反褶积方法，该方法利用零偏移距VSP 数据的直达波估算双程时变子波的振幅谱，然后利用双程时变子波的振幅谱估算得到时变子波的反褶积算子，将其用于叠后地面地震记录，其中时变地震子波的估算精度会影响最终的反褶积效果。

（2）理论合成数据和实际地震资料的处理结果证明本文方法可以有效恢复地震记录的高频衰减能量，并可提高地震记录的分辨率。