地震復譜分解技術(shù)及其在烴類檢測中的應用

韓利，劉春成，張益明，韓立國，葉云飛

1 中海油研究總院，北京　100028　2　吉林大學，長春　130026

?

地震復譜分解技術(shù)及其在烴類檢測中的應用

韓利1，劉春成1，張益明1，韓立國2，葉云飛1

1 中海油研究總院，北京1000282吉林大學，長春130026

摘要譜分解技術(shù)在地震解釋領域已得到廣泛應用，但常用的譜分解方法存在兩方面的不足.一是時間分辨率低，難以對薄層進行刻畫；二是在烴類檢測中多解性強，難以區(qū)分流體類型.為了改善該問題，本文提出一種基于地震復譜分解技術(shù)的烴類檢測方法.復譜分解是指用一個包含多個不同頻率Ricker子波的復子波庫對地震道進行分解，從而得到時變子波頻率和相位信息的過程.借助稀疏反演技術(shù)復譜分解可以獲得高分辨率的時頻能量譜和時頻相位譜.本文首先通過擬合算例驗證了復譜分解方法刻畫薄層的能力以及求取子波頻率和相位的準確性.然后利用基于Kelvin-Voigt模型的黏彈波動方程數(shù)值模擬對衰減引起子波相位改變的原因進行了分析.最后通過實際資料應用展示了本文方法在儲層預測中的高時間分辨率優(yōu)勢，驗證了利用子波相位信息識別氣藏的有效性.

關(guān)鍵詞復譜分解； 高分辨率； 時頻譜； 時頻相位譜； 烴類檢測

1引言

地震譜分解技術(shù)可以將一個時域信號分解為時間和頻率的二維函數(shù)，在儲層刻畫(Partyka et al.,1999)、烴類檢測(Castagna et al.,2003; Han et al.,2012,2014)及頻率依賴的AVO分析(Yoo and Gibson,2005)等地震解釋領域得到了廣泛的應用.

常用的時頻譜分解方法有短時Fourier變換、Gabor變換、連續(xù)小波變換(CWT)(Mallat,1989; Sinha et al.,2005)、S變換(Pinnegar and Mansinha,2003)等.這些方法在儲層預測和烴類檢測領域一直起著非常重要的作用，但同時也存在兩方面的不足：一是時間分辨率低，難以對薄層進行精細刻畫.當?shù)貙颖容^厚時，地層頂?shù)追瓷淇梢酝ㄟ^常規(guī)譜分解方法識別，而當?shù)貙虞^薄時，常規(guī)譜分解方法很難識別頂?shù)追瓷?事實上，這些基于變換的方法，不論使用固定時窗還是變尺度時窗，都受測不準原理限制，時間分辨率和頻率分辨率此消彼長，即當確定頻率分辨率的同時，就決定了最大時間分辨率難以突破(韓利，2013).二是在烴類檢測中多解性強.含水儲層和含烴類儲層常常產(chǎn)生相似的頻率響應，僅通過頻率信息難以準確識別流體類型.

波在含流體介質(zhì)中傳播時，受到介質(zhì)衰減和頻散的影響(張固瀾等，2014)，不但會改變子波的頻率，而且也會使子波相位發(fā)生改變.子波頻率的變化特征已被用于儲層預測，但子波相位的變化特征一直沒有得到應用，其主要原因是現(xiàn)階段對儲層引起的相位響應特征還不夠清楚，且缺乏從地震剖面中提取時變子波相位信息的技術(shù)手段.近些年，學者們不斷提出高分辨率的譜分解方法(Wang,2007；朱振宇等2009；武國寧等；2012；Han and van der Baan,2013; Yuan and Wang,2013; Han et al.,2014)，但仍然沒有提取和利用相位信息的有效手段.針對這個問題，本文將基于稀疏反演技術(shù)(Bonar and Sacchi,2010)的復譜分解方法用于烴類檢測.本文方法具有兩方面的優(yōu)勢，一方面提高了譜分解刻畫薄儲層的時間分辨率，另一方面可以求取用于氣藏檢測的時變子波相位信息.文章結(jié)構(gòu)如下：首先闡述了復譜分解方法的原理，并用擬合例子展示該方法分辨薄層的能力和子波相位信息求取的準確性；然后通過數(shù)值模擬方法對衰減引起相位改變的原因進行了分析；最后將本文提出的方法應用于實際資料的烴類檢測.

2復譜分解方法原理

傳統(tǒng)的地震褶積模型將地震記錄s(t)描述為一個固定頻率、相位的地震子波w(t)和反射系數(shù)r(t)的褶積與隨機噪聲n(t)之和：

(1)

其中，*表示褶積運算.但實際上地震子波受地下介質(zhì)衰減和頻散的影響，頻率和相位都是隨時間變化的.為了從地震記錄中提取這些與巖性和流體有關(guān)的頻率和相位信息，我們建立了一個比公式(1)更為復雜的復數(shù)褶積模型，

(2)

其中，w不再是一個子波，而是一個由不同頻率的零相位Ricker子波構(gòu)成的復數(shù)子波庫；K表示參與計算的子波的個數(shù)；wk表示復子波庫中第k個復子波，其頻率和相位分別為f和φw=0；rk為對應的攜帶有地震記錄中子波頻率和相位信息的復反射系數(shù).

圖1　子波庫示意圖(a) 實子波庫，由不同頻率的零相位Ricker子波組成； (b) 復子波庫，通過對實子波庫進行Hilbert變換得到，其中實線為復子波庫實部；虛線為復子波庫虛部.

子波庫的構(gòu)成如圖1所示.圖1a展示了由不同頻率的零相位Ricker子波構(gòu)成的實子波庫；通過對實子波庫進行Hilbert變換可以得到復子波庫，如圖1b所示，實線為復子波庫實部，虛線為復子波庫虛部.

將式(2)的褶積形式寫成如下線性相乘形式：

(3)

其中，W表示復子波庫的褶積矩陣.式(3)中r的元素個數(shù)多于地震記錄的元素個數(shù)，是一個欠定問題.文中通過如式(4)所示的L1范數(shù)進行正則化(Yuan et al.,2015)，并使用快速迭代軟閾值算法(FISTA)(Beck and Teboulle，2009)對式(4)進行求解，

(4)

式中，‖r‖1為L1范數(shù)約束項，λ為稀疏度控制參數(shù).

通過對r取模的平方和反正切得到地震記錄的時頻主頻譜Fsdom和時頻相位譜φs，

(5)

為了保持復譜分解的高時間分辨率優(yōu)勢，同時避免過高的頻率分辨率帶來的負面影響，可將時頻主頻譜乘上Ricker子波庫的頻譜轉(zhuǎn)置得到時頻能量譜，

(6)

圖2　時頻譜分解分辨薄層的能力(a) 紅線為擬合信號，是兩個黑線信號的組合，用來表征不同厚度地層的頂?shù)追瓷洌?(b) 為CWT方法得到的時頻譜； (c)和(d)分別為復譜分解方法得到的時頻主頻譜和時頻能量譜.

其中FRicker為Ricker子波庫頻譜,T為轉(zhuǎn)置運算；Fs時頻能量譜.

3擬合算例

圖2為一個展示復譜分解方法分辨薄層能力的算例.圖2a為由不同時間間隔的30Hz主頻的Ricker子波對組成的擬合信號，用來表征不同厚度地層的頂?shù)追瓷?圖2b為由CWT方法得到的時頻能量譜.當?shù)貙虞^厚時(60ms間隔)，CWT方法可以區(qū)分地層頂?shù)追瓷?；隨著地層厚度變薄(<30ms間隔)，CWT時頻譜不能很好地區(qū)分開頂?shù)追瓷?圖2c和2d分別為復譜分解方法得到的時頻主頻譜和時頻能量譜，可以看出復譜分解方法的時間分辨率很高.由于復譜分解的求解過程是一個迭代收斂的過程，稀疏度參數(shù)確定后，區(qū)分越薄的地層反演需要的迭代次數(shù)越多.在犧牲計算速度的情況下，復譜分解有能力將超薄層(<10ms間隔)頂?shù)追瓷渥R別出來.

圖3為一個由不同頻率和相位的Ricker子波組成的擬合地震信號例子.擬合信號的子波頻率從上到下依次為60,40,20,30Hz和30Hz，相位依次為0°，-90°，45°，-180°和-180°.圖3b為CWT方法得到的時頻能量譜，代表著常規(guī)譜分解的分辨率水平.圖3c為CWT方法得到的時頻相位譜，從中難以提取與衰減有關(guān)的子波相位信息.圖3d—3f分別為地震復譜分解方法得到的時頻主頻譜、時頻能量譜和時頻相位譜.地震復譜分解得到的三種時頻分布都具有很高的時間分辨率，并且所求結(jié)果與 擬合信號真實情況一致.擬合算例表明，相對于傳統(tǒng)CWT方法，復譜分解方法的優(yōu)勢不僅在于具有高時頻分辨率特征，而且還在于能夠準確地計算出時變子波的相位信息.

圖3　擬合算例(a) 由不同頻率和相位的Ricker子波組成的擬合地震信號，子波頻率從上到下依次為60,40,20,30 Hz和30 Hz，相位依次為0°，-90°，45°，-180°和-180°； (b)和(c)分別為CWT方法得到的時頻能量譜和時頻相位譜； (d)—(f)分別為復譜分解方法得到的時頻主頻譜、時頻能量譜和時頻相位譜.

4衰減有關(guān)的子波相位變化

地震波在能耗介質(zhì)中傳播會產(chǎn)生衰減已經(jīng)成為一個共識，并且這種衰減會引起子波頻率和相位的改變.介質(zhì)衰減程度越大，引起的頻率和相位變化量越大.在勘探中還發(fā)現(xiàn)另一個常被忽略的現(xiàn)象，即氣藏的頂界面也會引起相位的改變.為此，本文利用基于Kelvin-Voigt模型(Morozov，2011)的黏彈波動方程對波在界面產(chǎn)生相位改變的原因進行了研究.在正演中對界面上下介質(zhì)阻抗差異的影響和衰減差異的影響進行了獨立分析，如圖4所示.層狀模型中，界面的左側(cè)介質(zhì)較右側(cè)介質(zhì)波阻抗低、衰減程度大.圖4a和4b中藍線表示僅考慮阻抗差異存在時產(chǎn)生的透射波和反射波，可以看到，透射波與入射波一致，仍然是零相位，而反射波相位發(fā)生180°改變，說明阻抗差異不會改變透射波的相位，僅引起反射波的極性改變.圖4a中的紅線表示僅考慮衰減差異存在時產(chǎn)生的透射波和反射波，可以看出波在衰減介質(zhì)中傳播時相位發(fā)生了改變，且衰減程度的差異是引起反射波相位旋轉(zhuǎn)的主要原因.圖4b中黑線表示同時考慮阻抗差異和衰減差異存在時產(chǎn)生的透射波和反射波.從模擬結(jié)果中我們可以得出，衰減引起的相位改變不僅發(fā)生在波傳播過程中，還發(fā)生在存在衰減差異的界面處，儲層界面反射是阻抗差異和衰減差異共同作用的結(jié)果.通常氣藏與蓋層有明顯的阻抗差異和衰減差異，因此，通過相位信息檢測氣藏存在可能性.

5實際資料應用

圖5展示了一個實際資料烴類檢測的算例.圖5a為目標區(qū)過探井的地震剖面，剖面上所投的測井曲線為含水飽和度曲線.根據(jù)探井揭示和地質(zhì)解釋的氣層和水層位置已在圖中標示.圖5b為CWT方法得到的用于烴類檢測的頻率異常剖面，時間分辨率低于地震分辨率，對薄層層位無法準確刻畫；圖5c為使用本文提出的復譜分解方法得到的頻率異常剖面，具有很高的時間分辨率，且與測井上揭示的地層位置一致.含氣儲層在頻率異常剖面中有響應，但水層在頻率異常剖面上也有響應，所以這會干擾氣藏預測，造成解釋的多解性.圖5d為復譜分解方法求取的時頻相位譜異常剖面，從圖中可以看到，氣層和水層響應在相位異常剖面上存在一定差異，可用于降低解釋的多解性.

6結(jié)論

本文提出了一種基于地震復譜分解技術(shù)的烴類檢測方法.借助稀疏反演技術(shù)，復譜分解能夠得到高分辨率的時頻能量譜和時頻相位譜.擬合算例證明了復譜分解不僅具有刻畫薄層的能力，而且還能夠準確求取時變子波相位信息.通過基于Kelvin-Voigt模型的黏彈波動方程數(shù)值模擬，證明了與衰減有關(guān)的相位變化不僅發(fā)生在波傳播的過程中，也發(fā)生在存在阻抗和衰減差異的界面處.阻抗差異只引起反射子波的極性變化，而衰減差異是引起子波相位旋轉(zhuǎn)的的關(guān)鍵因素.因此，可將子波相位變化用于具有高衰減特征的氣藏檢測.在實際資料烴類檢測應用中，通過復分解方法求取的頻率異常剖面具有很高的時間分辨率，并且進一步結(jié)合相位譜異常有效地區(qū)分了該區(qū)含氣和含水儲層的響應，降低了僅利用頻率異常進行烴類檢測的多解性.致謝感謝薩斯喀徹溫省大學博士生鄧武兵在黏彈波動方程數(shù)值模擬上做的貢獻.

圖4　數(shù)值模擬分析引起相位改變的原因.界面左側(cè)的介質(zhì)相對于右側(cè)介質(zhì)波阻抗低、衰減程度高.(a)和(b)中藍線表示僅考慮阻抗差異存在時產(chǎn)生的透射波和反射波； (a)中的紅線表示僅存考慮存在衰減差異存在時產(chǎn)生的透射和反射波； (b)中黑線表示同時考慮存在阻抗差異和衰減差異存在時產(chǎn)生的的透射和反射波

References

Beck A,Teboulle M.2009.A fast iterative shrinkage-thresholding algorithm for linear inverse problems.SIAM Journal on Imaging Sciences,2(1): 183-202.

Bonar D C,Sacchi M D.2010.Complex spectral decomposition via inversion strategies.∥ SEG Technical Program Expanded Abstracts.SEG: 1408-1412.

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.

Han J J,van der Baan M.2013.Empirical mode decomposition for seismic time-frequency analysis.Geophysics,78(2): O9-O19.

Han L,Han L G,Li Z.2012.Inverse spectral decomposition with the SPGL1 algorithm.Journal of Geophysics and Engineering,9(4): 423-427.

Han L.2013.Research on the methods of high-resolution full spectrum decomposition [Ph.D.thesis] (in Chinese).Changchun: Jilin University.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.

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.

Morozov I B.2011.Anelastic acoustic impedance and the correspondence principle.Geophysical Prospecting,59(1): 24-34.

Partyka G,Gridley J,Lopez J.1999.Interpretational applications of spectral decomposition in reservoir characterization.The Leading Edge,18(3): 353-360.

Pinnegar C R,Mansinha L.2003.The S-transform with windows of arbitrary and varying shape.Geophysics,68: 381-385.

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.

Wang Y H.2007.Seismic time-frequency spectral decomposition by matching pursuit.Geophysics,72(1): V13-V20.

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.

Yoo S,Gibson R L Jr.2005.Frequency dependent AVO analysis after target oriented stretch correction.∥ SEG Technical Program Expanded Abstracts.SEG: 293-296.

Yuan S Y,Wang S X.2013.Spectral sparse Bayesian learning reflectivity inversion.Geophysical Prospecting,61(4): 735-746.

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.

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.

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.

附中文參考文獻

韓利.2013.高分辨率全譜分解方法研究[博士學位論文].長春: 吉林大學.

武國寧,曹思遠,孫娜.2012.基于復數(shù)道地震記錄的匹配追蹤算法及其在儲層預測中的應用.地球物理學報,55(6): 2027-2034,doi: 10.6038/j.issn.0001-5733.2012.06.024.

張固瀾,賀振華,王熙明等.2014.地震波頻散效應與反Q濾波相位補償.地球物理學報,57(5): 1655-1663,doi: 10.6038/cjg20140528.

朱振宇,呂丁友,桑淑云等.2009.基于物理小波的頻譜分解方法及應用研究.地球物理學報,52(8): 2152-2157,doi: 10.3969/j.issn.0001-5733.2009.08.025.

(本文編輯胡素芳)

基金項目國家重大專項子課題(2011ZX05025-001-07)資助.

作者簡介韓利，男，1984年生，2013年博士畢業(yè)于吉林大學，目前在中海油研究總院主要從事地震資料解釋和烴源巖預測方面的工作. E-mail：hanli4@cnooc.com.cn

doi:10.6038/cjg20160329 中圖分類號P631

收稿日期2015-04-13，2015-10-26收修定稿

Seismic complex spectral decomposition and its application on hydrocarbon detection

HAN Li1,LIU Chun-Cheng1,ZHANG Yi-Ming1,HAN Li-Guo2,YE Yun-Fei1

1CNOOCResearchInstitute，Beijng100082，China2JilinUniversity，Changchun130026，China

AbstractSpectral 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.

KeywordsComplex spectral decomposition; High resolution; Time-frequency spectrum; Time-frequency phase spectrum; Hydrocarbon detection

韓利，劉春成，張益明等.2016.地震復譜分解技術(shù)及其在烴類檢測中的應用.地球物理學報，59(3):1095-1101,doi:10.6038/cjg20160329.

Han L,Liu C C,Zhang Y M,et al.2016.Seismic complex spectral decomposition and its application on hydrocarbon detection.Chinese J.Geophys.(in Chinese),59(3):1095-1101,doi:10.6038/cjg20160329.