地球物理学报  2016, Vol. 59 Issue (3): 1095-1101   PDF    
地震复谱分解技术及其在烃类检测中的应用
韩利1, 刘春成1, 张益明1, 韩立国2, 叶云飞1    
1. 中海油研究总院, 北京 100028;
2. 吉林大学, 长春 130026
摘要: 谱分解技术在地震解释领域已得到广泛应用,但常用的谱分解方法存在两方面的不足.一是时间分辨率低,难以对薄层进行刻画;二是在烃类检测中多解性强,难以区分流体类型.为了改善该问题,本文提出一种基于地震复谱分解技术的烃类检测方法.复谱分解是指用一个包含多个不同频率Ricker子波的复子波库对地震道进行分解,从而得到时变子波频率和相位信息的过程.借助稀疏反演技术复谱分解可以获得高分辨率的时频能量谱和时频相位谱.本文首先通过拟合算例验证了复谱分解方法刻画薄层的能力以及求取子波频率和相位的准确性.然后利用基于Kelvin-Voigt模型的黏弹波动方程数值模拟对衰减引起子波相位改变的原因进行了分析.最后通过实际资料应用展示了本文方法在储层预测中的高时间分辨率优势,验证了利用子波相位信息识别气藏的有效性.
关键词: 复谱分解     高分辨率     时频谱     时频相位谱     烃类检测    
Seismic complex spectral decomposition and its application on hydrocarbon detection
HAN Li1, LIU Chun-Cheng1, ZHANG Yi-Ming1, HAN Li-Guo2, YE Yun-Fei1    
1. CNOOC Research Institute, Beijng 100082, China;
2. Jilin University, Changchun 130026, China
Abstract: Spectral decomposition techniques have been widely used in seismic interpretation. However, the traditional methods have two limitations, low time resolution and uncertainty in predicting fluid types. In this paper, a seismic complex spectral decomposition technique is developed to achieve a much higher resolution in time-frequency distributions via inversion strategy and to reduce the uncertainty in hydrocarbon detection with the extracted wavelet phase information. The time-varying wavelet frequency and phase information are obtained by decomposing the seismic trace using a wavelet library composed of a number of complex Ricker wavelets with different dominant frequencies and zero phase. This process is referred to as the seismic complex decomposition. The Synthetic examples show the ability of the proposed method to distinguish closely positioned wavelets and the accuracy to extract the wavelet frequency and phase information from the seismic data. The causes of wavelet phase change due to attenuation are studied by using the viscoelastic wave equation modelling which is based on the Kelvin-Voigt model. The result demonstrates that the phase change not only occurs in the process of wave propagating, but also occurs at the interface where existing impedance or attenuation contrast. The impedance contrast only changes the polarity of the wave, while the attenuation contrast is the key reason that rotates the wavelet phase. The gas reservoir has high attenuation property, so that it could be predicted via the wavelet phase change information.The real data example demonstrates the high time resolution of the proposed method in reservoir prediction and the successful application of the phase information in distinguishing gas saturated layers and water saturated layers.
Key words: Complex spectral decomposition     High resolution     Time-frequency spectrum     Time-frequency phase spectrum     Hydrocarbon detection    
1 引言

地震谱分解技术可以将一个时域信号分解为时间和频率的二维函数,在储层刻画(Partyka et al.,1999)、烃类检测(Castagna et al.,2003; Han et al.,20122014)及频率依赖的AVO分析(Yoo and Gibson,2005)等地震解释领域得到了广泛的应用.

常用的时频谱分解方法有短时Fourier变换、Gabor变换、连续小波变换(CWT)(Mallat,1989; Sinha et al.,2005)、S变换(Pinnegar and Mansinha,2003)等. 这些方法在储层预测和烃类检测领域一直起着非常重要的作用,但同时也存在两方面的不足:一是时间分辨率低,难以对薄层进行精细刻画.当地层比较厚时,地层顶底反射可以通过常规谱分解方法识别,而当地层较薄时,常规谱分解方法很难识别顶底反射.事实上,这些基于变换的方法,不论使用固定时窗还是变尺度时窗,都受测不准原理限制,时间分辨率和频率分辨率此消彼长,即当确定频率分辨率的同时,就决定了最大时间分辨率难以突破(韩利,2013).二是在烃类检测中多解性强.含水储层和含烃类储层常常产生相似的频率响应,仅通过频率信息难以准确识别流体类型.

波在含流体介质中传播时,受到介质衰减和频散的影响(张固澜等,2014),不但会改变子波的频率,而且也会使子波相位发生改变.子波频率的变化特征已被用于储层预测,但子波相位的变化特征一直没有得到应用,其主要原因是现阶段对储层引起的相位响应特征还不够清楚,且缺乏从地震剖面中提取时变子波相位信息的技术手段.近些年,学者们不断提出高分辨率的谱分解方法(Wang,2007朱振宇等2009武国宁等;2012Han and van derBaan,2013; Yuan and Wang,2013; Han et al.,2014),但仍然没有提取和利用相位信息的有效手段.针对这 个问题,本文将基于稀疏反演技术(Bonar and Sacchi,2010)的复谱分解方法用于烃类检测.本文方法具有两方面的优势,一方面提高了谱分解刻画薄储层的时间分辨率,另一方面可以求取用于气藏检测的时变子波相位信息.文章结构如下:首先阐述了复谱分解方法的原理,并用拟合例子展示该方法分辨薄层的能力和子波相位信息求取的准确性;然后通过数 值模拟方法对衰减引起相位改变的原因进行了分析;最后将本文提出的方法应用于实际资料的烃类检测.

2 复谱分解方法原理

传统的地震褶积模型将地震记录 s(t)描述为一个固定频率、相位的地震子波 w(t)和反射系数 r(t)的褶积与随机噪声 n(t)之和:

其中,*表示褶积运算.但实际上地震子波受地下介质衰减和频散的影响,频率和相位都是随时间变化的.为了从地震记录中提取这些与岩性和流体有关的频率和相位信息,我们建立了一个比公式(1)更为复杂的复数褶积模型,

其中,w 不再是一个子波,而是一个由不同频率的零相位Ricker子波构成的复数子波库;K表示参与计算的子波的个数; w k 表示复子波库中第k个复子波,其频率和相位分别为fφ w=0; r k 为对应的携 带有地震记录中子波频率和相位信息的复反射系数.

子波库的构成如图 1所示.图 1a展示了由不同频率的零相位Ricker子波构成的实子波库;通过对实子波库进行Hilbert变换可以得到复子波库,如图 1b所示,实线为复子波库实部,虚线为复子波库虚部.

图 1 子波库示意图
(a) 实子波库,由不同频率的零相位Ricker子波组成; (b) 复子波库,通过对实子波库进行Hilbert变换得到,其中实线为复子波库实部;虚线为复子波库虚部.
Fig. 1 Schematic diagram depicting the construction of complex wavelet library
(a) The real wavelet library composed of many Ricker wavelets with different dominant frequency and zero phase; (b) The complex wavelet library constructed from (a) by Hilbert Transform. The solid lines represent the real components and the dashed lines represent the imaginary components.

将式(2)的褶积形式写成如下线性相乘形式:

其中,W 表示复子波库的褶积矩阵.式(3)中 r 的元素个数多于地震记录的元素个数,是一个欠定问题.文中通过如式(4)所示的L1范数进行正则化(Yuan et al.,2015),并使用快速迭代软阈值算法(FISTA)(Beck and Teboulle,2009)对式(4)进行求解,

式中,‖ r1 为L1范数约束项,λ为稀疏度控制参数.

通过对 r 取模的平方和反正切得到地震记录的时频主频谱 F sdom和时频相位谱 φ s

为了保持复谱分解的高时间分辨率优势,同时避免过高的频率分辨率带来的负面影响,可将时频主频谱乘上Ricker子波库的频谱转置得到时频能量谱,

其中 F Ricker为Ricker子波库频谱,T为转置运算; F s时频能量谱.

3 拟合算例

图 2为一个展示复谱分解方法分辨薄层能力的 算例.图 2a为由不同时间间隔的30 Hz主频的Ricker 子波对组成的拟合信号,用来表征不同厚度地层的顶底反射.图 2b为由CWT方法得到的时频能量谱.当地层较厚时(60 ms间隔),CWT方法可以区分地层顶底反射;随着地层厚度变薄(<30 ms间隔),CWT时频谱不能很好地区分开顶底反射.图 2c和2d分别为复谱分解方法得到的时频主频谱和时频能量谱,可以看出复谱分解方法的时间分辨率很高.由于复谱分解的求解过程是一个迭代收敛的过程,稀疏度参数确定后,区分越薄的地层反演需要的迭代次数越多.在牺牲计算速度的情况下,复谱分解有能力将超薄层(<10 ms间隔)顶底反射识别出来.

图 2 时频谱分解分辨薄层的能力
(a) 红线为拟合信号,是两个黑线信号的组合,用来表征不同厚度地层的顶底反射; (b) 为CWT方法得到的时频谱; (c)和(d)分别为复谱分解方法得到的时频主频谱和时频能量谱.
Fig. 2 The ability of spectral decomposition in distinguishing thin layers
(a) The red signal is a combination of two black signals, representing the top and bottom reflections from strata of different thickness. (b) Time-frequency spectrum generated by the CWT method. (c) and (d) are time-frequency dominant spectrum and time-frequency energy spectrum generated by the proposed method, respectively.

图 3为一个由不同频率和相位的Ricker子波组成的拟合地震信号例子.拟合信号的子波频率从上到下依次为60,40,20,30 Hz和30 Hz,相位依次为0°,-90°,45°,-180°和-180°.图 3b为CWT方法得到的时频能量谱,代表着常规谱分解的分辨率水平.图 3c为CWT方法得到的时频相位谱,从中难以提取与衰减有关的子波相位信息.图 3d3f分别为地震复谱分解方法得到的时频主频谱、时频能量谱和时频相位谱.地震复谱分解得到的三种时 频分布都具有很高的时间分辨率,并且所求结果与 拟合信号真实情况一致.拟合算例表明,相对于传统CWT方法,复谱分解方法的优势不仅在于具有高时频分辨率特征,而且还在于能够准确地计算出时变子波的相位信息.

图 3 拟合算例
(a) 由不同频率和相位的Ricker子波组成的拟合地震信号,子波频率从上到下依次为60, 40, 20, 30 Hz和30 Hz,相位依次为0°,-90°,45°,-180°和-180°; (b)和(c)分别为CWT方法得到的时频能量谱和时频相位谱; (d)—(f)分别为复谱分解方法得到的时频 主频谱、时频能量谱和时频相位谱.
Fig. 3 A synthetic example
(a) The synthetic signal composed of five Ricker wavelets with different dominant frequencies and phase changes; (b) and (c) are the time-frequency energy spectrum and time-frequency phase spectrum obtained with the CWT method, respectively. (d)—(f) are the time-frequency dominant spectrum, time-frequency energy spectrum and time-frequency phase spectrum obtained with the complex spectral decomposition method, respectively.
4 衰减有关的子波相位变化

地震波在能耗介质中传播会产生衰减已经成为一个共识,并且这种衰减会引起子波频率和相位的改变.介质衰减程度越大,引起的频率和相位变化量越大.在勘探中还发现另一个常被忽略的现象,即气藏的顶界面也会引起相位的改变.为此,本文利用基于Kelvin-Voigt模型(Morozov,2011)的黏弹波动方程对波在界面产生相位改变的原因进行了研究.在正演中对界面上下介质阻抗差异的影响和衰减差异的影响进行了独立分析,如图 4所示.层状模型中,界面的左侧介质较右侧介质波阻抗低、衰减程度大.图 4a4b中蓝线表示仅考虑阻抗差异存在时产生的透射波和反射波,可以看到,透射波与入射波一致,仍然是零相位,而反射波相位发生180°改变,说明阻抗差异不会改变透射波的相位,仅引起反射波的极性改变.图 4a中的红线表示仅考虑衰减差异存在时产生的透射波和反射波,可以看出波在衰减介质中传播时相位发生了改变,且衰减程度的差异是引起反射波相位旋转的主要原因.图 4b中黑线表示同时考虑阻抗差异和衰减差异存在时产生的透射波和反射波.从模拟结果中我们可以得出,衰减引起的相位改变不仅发生在波传播过程中,还发生在存在衰减差异的界面处,储层界面反射是阻抗差异和衰减差异共同作用的结果.通常气藏与盖层有明显的阻抗差异和衰减差异,因此,通过相位信息检测气藏存在可能性.

图 4 数值模拟分析引起相位改变的原因. 界面左侧的介质相对于右侧介质波阻抗低、衰减程度高.(a)和(b)中蓝线表示 仅考虑阻抗差异存在时产生的透射波和反射波; (a)中的红线表示仅存考虑 存在衰减差异存在时产生的透射和反射波; (b)中黑线表示同时考虑存在阻抗差异和衰减差异存在时产生的的透射和反射波 Fig. 4 Causes of wavelet phase change due to attenuation. The media on the left side of the interface has lower impedance and higher attenuation than that on the right side. The blue waves in (a) and (b) denote the reflected and transmitted wave in the case that only the impedance contrast is considered. The red waves in (a) denote the case that only the attenuation contrast is considered. The black waves in (b) denote the case that both the impedance contrst and the attenuation contrast are considered
5 实际资料应用

图 5展示了一个实际资料烃类检测的算例.图 5a为目标区过探井的地震剖面,剖面上所投的测井曲线为含水饱和度曲线.根据探井揭示和地质解释的气层和水层位置已在图中标示.图 5b为CWT方法得到的用于烃类检测的频率异常剖面,时间分辨率低于地震分辨率,对薄层层位无法准确刻画;图 5c为使用本文提出的复谱分解方法得到的频率异常剖面,具有很高的时间分辨率,且与测井上揭示的地层位置一致.含气储层在频率异常剖面中有响应,但水层在频率异常剖面上也有响应,所以这会干扰气藏预测,造成解释的多解性.图 5d为复谱分解方法求取的时频相位谱异常剖面,从图中可以看到,气层和水层响应在相位异常剖面上存在一定差异,可用于降低解释的多解性.

图 5 实际资料烃类检测算例. (a) 目标区过探井的地震剖面,所投测井曲线为含水饱和度,钻井揭示和地质解释的气层和水层位置已在图中标示; (b) CWT方法得到的 用于烃类检测的频率异常剖面; (c) 使用本文提出的复谱分解方法得到 的频率异常剖面; (d) 复谱分解方法求取的时频相位谱异常剖面 Fig. 5 A real hydrocarbon detection example. (a) A seismic connecting-wells section. The curve of water saturation curve is posted and the interpreted gas saturated and water saturated layers are marked; (b) The abnormal frequency section generated by CWT method. (c) and (d) are the abnormal frequency and abnormal phase sections generated by the complex decomposition method, respectively
6 结论

本文提出了一种基于地震复谱分解技术的烃类检测方法.借助稀疏反演技术,复谱分解能够得到高分辨率的时频能量谱和时频相位谱.拟合算例证明了复谱分解不仅具有刻画薄层的能力,而且还能够准确求取时变子波相位信息.通过基于Kelvin-Voigt模型的黏弹波动方程数值模拟,证明了与衰减有关的相位变化不仅发生在波传播的过程中,也发生在存在阻抗和衰减差异的界面处.阻抗差异只引起反射子波的极性变化,而衰减差异是引起子波相位旋转的的关键因素.因此,可将子波相位变化用于具有高衰减特征的气藏检测.在实际资料烃类检测应用中,通过复分解方法求取的频率异常剖面具有很高的时间分辨率,并且进一步结合相位谱异常有效地区分了该区含气和含水储层的响应,降低了仅利用频率异常进行烃类检测的多解性.

致谢   感谢萨斯喀彻温省大学博士生邓武兵在黏弹波动方程数值模拟上做的贡献.

参考文献
[1] Beck A, Teboulle M. 2009. A fast iterative shrinkage-thresholding algorithm for linear inverse problems. SIAM Journal on Imaging Sciences, 2(1):183-202.
[2] Bonar D C, Sacchi M D. 2010. Complex spectral decomposition via inversion strategies.//SEG Technical Program Expanded Abstracts. SEG:1408-1412.
[3] Castagna J P, Sun S J, Siegfried R W. 2003. Instantaneous spectral analysis:Detection of low-frequency shadows associated with hydrocarbons. The Leading Edge, 22(2):120-127.
[4] Han J J, van der Baan M. 2013. Empirical mode decomposition for seismic time-frequency analysis. Geophysics, 78(2):O9-O19.
[5] Han L, Han L G, Li Z. 2012. Inverse spectral decomposition with the SPGL1 algorithm. Journal of Geophysics and Engineering, 9(4):423-427.
[6] Han L. 2013. Research on the methods of high-resolution full spectrum decomposition[Ph. D. thesis](in Chinese). Changchun:Jilin University.
[7] Han L, Sacchi M D, Han L G. 2014. Spectral decomposition and de-noising via time-frequency and space-wavenumber reassignment. Geophysical Prospecting, 62(2):244-257.
[8] Mallat S G. 1989. A theory for multiresolution signal decomposition:the wavelet representation. IEEE Transactions on Pattern Analysis and Machine Intelligence, 11(7):674-693.
[9] Morozov I B. 2011. Anelastic acoustic impedance and the correspondence principle. Geophysical Prospecting, 59(1):24-34.
[10] Partyka G, Gridley J, Lopez J. 1999. Interpretational applications of spectral decomposition in reservoir characterization. The Leading Edge, 18(3):353-360.
[11] Pinnegar C R, Mansinha L. 2003. The S-transform with windows of arbitrary and varying shape. Geophysics, 68:381-385.
[12] Sinha S, Routh P S, Anno P D, et al. 2005. Spectral decomposition of seismic data with continuous-wavelet transform. Geophysics, 70(6):P19-P25.
[13] Wang Y H. 2007. Seismic time-frequency spectral decomposition by matching pursuit. Geophysics, 72(1):V13-V20.
[14] Wu G J, Cao S Y, Sun N. 2012. Matching pursuit method based on complex seismic traces and its application of hydrocarbon exploration. Chinese J. Geophys.(in Chinese), 55(6):2027-2034, doi:10.6038/j.issn.0001-5733.2012.06.024.
[15] Yoo S, Gibson R L Jr. 2005. Frequency dependent AVO analysis after target oriented stretch correction.// SEG:293-296SEG Technical Program Expanded Abstracts..
[16] Yuan S Y, Wang S X. 2013. Spectral sparse Bayesian learning reflectivity inversion. Geophysical Prospecting, 61(4):735-746.
[17] Yuan S Y, Wang S X, Luo C M, et al. 2015. Simultaneous multitrace impedance inversion with transform-domain sparsity promotion. Geophysics, 80(2):R71-R80.
[18] Zhang G L, He Z H, Wang X M, et al. 2014. Seismic wave dispersion effects and inverse Q-filter phase compensation. Chinese J. Geophys.(in Chinese), 57(5):1655-1663, doi:10.6038/cjg20140528.
[19] Zhu Z Y, Lü D Y, Sang S Y, et al. 2009. Research of spectrum decomposition method based on physical wavelet transform and its application. Chinese J. Geophys.(in Chinese), 52(8):2152-2157, doi:10.3969/j.issn.0001-5733.2009.08.025.
[20] 韩利. 2013. 高分辨率全谱分解方法研究[博士学位论文]. 长春:吉林大学.
[21] 武国宁, 曹思远, 孙娜. 2012. 基于复数道地震记录的匹配追踪算法及其在储层预测中的应用.地球物理学报, 55(6):2027-2034, doi:10.6038/j.issn.0001-5733.2012.06.024.
[22] 张固澜, 贺振华, 王熙明等. 2014. 地震波频散效应与反Q滤波相位补偿. 地球物理学报, 57(5):1655-1663, doi:10.6038/cjg20140528.
[23] 朱振宇, 吕丁友, 桑淑云等. 2009. 基于物理小波的频谱分解方法及应用研究. 地球物理学报, 52(8):2152-2157, doi:10.3969/j.issn.0001-5733.2009.08.025.