朱振宇 高佳倫 姜秀娣 孫文博 薛東川 王清振

(中海油研究總院有限責(zé)任公司,北京100028； ②中國石油大學(xué)(北京),北京 102200)

1 引言

近年來，隨著油氣勘探工作的不斷深入，頻譜分解方法已成為油氣儲層預(yù)測中的常用手段。Widess[1]利用楔形模型揭示了薄層反射與薄層厚度的關(guān)系。隨后，Neidell等[2]引入了調(diào)諧厚度的概念，即當(dāng)層厚為地震子波波長的四分之一時，合成地震記錄的振幅最大，此時的地層厚度即為調(diào)諧厚度，這是利用地震反射振幅解釋薄層厚度的基礎(chǔ)。Partyka[3]發(fā)現(xiàn)，薄層反射在短時窗傅里葉變換頻譜上會出現(xiàn)陷頻現(xiàn)象，利用這種現(xiàn)象可以估算薄層厚度；同時，根據(jù)地震薄層反射的調(diào)諧理論，不同的調(diào)諧頻率對應(yīng)不同的薄層厚度[4]，因此通過頻譜分解獲得的單頻體可以反映該頻率成分對應(yīng)的地質(zhì)現(xiàn)象[5]，挖掘更豐富的信息，還可以突出關(guān)鍵的地質(zhì)目標(biāo)，從而提高地震資料的解釋精度[6]。

利用頻譜分解技術(shù)進行地震成像時，需要選擇合適的頻譜分解方法[7]。Morlet提出的小波變換可以很好地刻畫非平穩(wěn)信號的局部特征，長期以來一直是頻譜分解的有力工具[8]。小波變換的效果取決于基小波的選擇，由于Morlet小波具有較高的時間、頻率分辨率[9]，一直被廣泛使用。但是當(dāng)Morlet小波的中心頻率較小時，其修正項就不可忽略。為此，Harrop等[10]提出了一種改進的Morlet小波，該小波在中心頻率較小時仍然可以滿足小波的允許條件。在Harrop等的研究基礎(chǔ)上，高靜懷等[11]提出一種靈活性高、具有三個可調(diào)參數(shù)的三參數(shù)小波，主要研究了σ(調(diào)制頻率)和τ(能量衰減因子)兩個參數(shù)。

本文研究了三參數(shù)小波的每個參數(shù)及其參數(shù)組合，重點分析了參數(shù)β(能量延遲因子)的影響。通過改變參數(shù)組合可以獲得多種形態(tài)的基小波，從而滿足不同的處理需求。首先闡述了頻譜分解方法的基本原理；然后介紹了三參數(shù)小波，重點研究了每個參數(shù)對小波的影響；隨后設(shè)計了正演模型，利用三參數(shù)小波進行時頻分析；最后利用基于三參數(shù)小波的頻譜分解方法預(yù)測儲層，取得了很好的效果。

圖1 楔形地質(zhì)模型及其地震響應(yīng)特征(a)楔形地質(zhì)模型； (b)合成地震記錄；(c)薄層濾波器的振幅譜；(d)調(diào)諧曲線楔形地質(zhì)模型的蓋層和底層具有相同的波阻抗

2 頻譜分解方法

頻譜分解方法的理論基礎(chǔ)來源于薄層反射的調(diào)諧現(xiàn)象。由楔形地質(zhì)模型(圖1a)的正演地震記錄(圖1b)可見，薄層的存在對入射波而言相當(dāng)于一個濾波器，由薄層濾波器響應(yīng)的傅里葉振幅譜可以看到周期性頻陷現(xiàn)象(圖1c)。提取層厚度為25、50ms對應(yīng)的兩條調(diào)諧曲線(圖1d)，曲線的第一個峰值頻率即為該厚度對應(yīng)的調(diào)諧頻率，曲線的陷頻頻率Pf與薄層厚度T(雙程旅行時)存在如下對應(yīng)關(guān)系[3]

(1)

利用式(1)就可以計算薄層厚度。

通過上述分析(圖1)可知，不同的薄層厚度對應(yīng)不同的調(diào)諧頻率，換言之，不同的頻率成分包含的地質(zhì)信息也會有所差異。所以在頻譜分解時可分頻解釋[12](圖2)，在不同的單頻體上突顯關(guān)鍵目標(biāo)。利用頻譜分解可以開展針對性的解釋工作，包括刻畫砂體的展布和河道的邊界等[13-16]，其工作流程如下：

(1)選定目的層段；

(2)選擇恰當(dāng)?shù)臅r頻分析方法進行頻譜分解，借助分頻處理手段獲得離散的單頻體；

(3)在每個單頻體上制作沿層切片，利用頻率切片預(yù)測砂體展布、刻畫河道邊界；

(4)將不同的單頻切片進行RGB混色融合，提高儲層的解釋精度(此步驟可以靈活取舍)。

圖2 分頻解釋

3 三參數(shù)小波變換

在頻譜分解時需要選擇合適的時頻分析方法。目前已發(fā)展的時頻分析方法種類繁多[17]，小波變換是其中的一種方法。20世紀(jì)80年代初期，Grossmann 等[18]發(fā)現(xiàn)傳統(tǒng)的傅里葉變換在分析地震信號的局部特征時遇到了瓶頸，無法刻畫信號在某種頻率出現(xiàn)的時間位置。為此，引入小波變換的概念彌補這一缺陷。定義信號s(t)的小波變換(WT)為

(2)

式中：ψ(t)為基小波，上角“*”表示取共軛；a≠0為尺度因子；b為平移因子。小波變換的特點是既考慮頻率分辨率，也考慮時間分辨率，在分析信號時具有“變焦“功能，可很好地刻畫信號的局部特征[19]。小波變換結(jié)果與基小波的選擇關(guān)系緊密[20]。在實際的應(yīng)用中，Morlet小波具有良好的聯(lián)合時頻分辨率，使用較為廣泛。Morlet小波的時間域表達式為

(3)

式中ω0≥5rad/s。由式(3)可以看出，Morlet小波只有一個可調(diào)參數(shù)，因此在使用時不夠靈活。并且當(dāng)ω0取較小值時，其時域局部化效果不能令人滿意。三參數(shù)小波可以解決上述問題，其表達式為

ψ(t,Γ)= e-τ(t-β)2{p(Γ)[cos(σt)-k(Γ)]+

iq(Γ)sin(σt)}

(4)

其中

(5)

(6)

(7)

Γ=(σ,τ,β)

(8)

式中：σ為小波的調(diào)制頻率；τ為能量衰減因子；β為能量延遲因子[11]。

通過改變?nèi)齾?shù)小波的3個參數(shù)，可以得到不同形態(tài)的基小波。因此三參數(shù)小波的變化形式多樣，選擇靈活[21]，可以滿足不同的處理需求。下面著重討論每個參數(shù)對小波形態(tài)的影響。分析式(4)可以看出，三參數(shù)小波是由衰減函數(shù)和復(fù)變函數(shù)相乘得到的，這個復(fù)變函數(shù)的實部和虛部均為三角函數(shù)，而衰減函數(shù)相當(dāng)于一個窗口。所以三參數(shù)小波使時間軸上無限分布的函數(shù)變?yōu)榫o支撐集函數(shù)，從而符合小波“小”的概念。σ作為調(diào)制頻率，控制三角函數(shù)的頻率，影響小波的震蕩程度，即σ越大，小波的震蕩越劇烈。

圖3為不同σ的三參數(shù)小波。由圖可見，隨著σ的不斷增大，小波旁瓣增多，震蕩加快。τ作為能量衰減因子，控制衰減函數(shù)的衰減速度，即τ越大，衰減越快，窗口越窄，小波的波形就越窄。圖4為不同τ的三參數(shù)小波。由圖可見，隨著τ的不斷增大，小波在時間域的展布變窄，波形變窄。值得注意的是，σ對小波波形的影響程度與τ的取值有關(guān)，當(dāng)τ取值較大(≥3)時，σ對小波的影響變小。這是因為τ較大時，窗口過窄，改變σ，就不易體現(xiàn)小波的周期變化，即震蕩程度變化不明顯。圖5為τ較大時不同σ的三參數(shù)小波，與τ較小時不同σ的三參數(shù)小波(圖3)相比，當(dāng)τ較大時，改變σ，小波震蕩程度變化不明顯。

圖3 不同σ的三參數(shù)小波(a)Γ=(1,1,0)；(b)Γ=(3,1,0)； (c)Γ=(6,1,0)

圖4 不同τ的三參數(shù)小波(a)Γ=(2.5,1,0)；(b)Γ=(2.5,2,0)； (c)Γ=(2.5,3,0)

圖5 τ較大時不同σ的三參數(shù)小波(a)Γ=(1,3,0)；(b)Γ=(3,3,0)； (c)Γ=(6,3,0)

圖6 不同β的三參數(shù)小波(a)Γ=(π,3,0)； (b)Γ=(π,3,0.6)； (c)Γ=(π,3,1)； (d)Γ=(π,3,2)

最后討論β對小波的影響，β作為能量延遲因子對小波形態(tài)的影響較復(fù)雜。如果改變β，會使衰減函數(shù)產(chǎn)生時間延遲，即窗口在時間軸滑動。當(dāng)β取三角函數(shù)周期的整數(shù)倍時，窗口就按周期的整數(shù)倍移動，則小波只會發(fā)生時移；當(dāng)β為周期的非整數(shù)倍時，窗口內(nèi)被改造的函數(shù)成分發(fā)生變化，因此小波不僅發(fā)生時移，還產(chǎn)生變形。由不同β的三參數(shù)小波(圖6)可見：①固定σ為π(周期為2)，當(dāng)β=0.6時，小波不僅發(fā)生時移，還產(chǎn)生變形(圖6b)，這種變形即為相位延遲。可以理解為由于β改變，導(dǎo)致窗口平移(也可以認為窗口不變，三角函數(shù)平移β個單位)，則相位延遲σβ個單位。②當(dāng)β=1.0(為周期的一半)時，小波實部與虛部都發(fā)生反轉(zhuǎn)，且產(chǎn)生1s的時移(圖6c)。③當(dāng)β=2.0(圖6d)時，小波只產(chǎn)生2s的時移，波形并沒有發(fā)生改變。由此可見：三參數(shù)小波在β為周期的整數(shù)倍時，可以匹配零相位子波；在β為周期的非整數(shù)倍時，可以匹配非零相位子波。

4 模型測試及實際地震資料應(yīng)用

4.1 模型測試

通過模型測試三參數(shù)小波和Morlet小波時頻分析的差異，測試模型是反射系數(shù)變化的薄互層地質(zhì)模型(圖7a)。圖7為測試模型及不同小波的時頻譜。由圖可見： ①對比ω0=5rad/s的Morlet小波(圖7c)和Γ=(5,2,0)的三參數(shù)小波時頻譜(圖7d)發(fā)現(xiàn)，前者只有一個強能量團，看不到任何薄層信息，后者顯示高頻能量呈鋸齒狀分布，即利用三參數(shù)小波做基小波可反映部分薄層信息，而Morlet小波不能精確地反映地下的地質(zhì)細節(jié)。②對比ω0=1rad/s的Morlet小波(圖7e)和Γ=(1,2,0)的三參數(shù)小波時頻譜(圖7f)發(fā)現(xiàn)，當(dāng)ω0較小時提高了Morlet小波的時間分辨率，但這是以犧牲頻率分辨率為代價的，其頻率分辨率幾乎為零(圖7e)，也就失去了時頻域聯(lián)合解釋的意義； 在三參數(shù)小波時頻譜(圖7f)中，主頻帶仍然在50Hz附近，且能量團發(fā)生分裂，出現(xiàn)多個錐形譜，可識別薄互層地質(zhì)結(jié)構(gòu)，即三參數(shù)小波的靈活性高，能更好地刻畫復(fù)雜地質(zhì)現(xiàn)象。③對比Γ=(5,2,0)(圖7d)和Γ=(1,2,0)的三參數(shù)小波時頻譜(圖7f)發(fā)現(xiàn)，后者(即當(dāng)σ取較小值時)的時間分辨率更高，刻畫薄層的能力更強。

圖7 測試模型及不同小波的時頻譜(a)反射系數(shù)序列； (b)圖a數(shù)據(jù)與50Hz雷克子波褶積得到的地震記錄； (c)ω0=5rad/s的Morlet小波時頻譜；(d)Γ=(5,2,0)的三參數(shù)小波時頻譜； (e)ω0=1rad/s的Morlet小波時頻譜； (f)Γ=(1,2,0)的三參數(shù)小波時頻譜圖a中的反射系數(shù)正、負相間，其絕對值從0.1增大到0.5，再從0.5減小至0.1，共計10個反射系數(shù)，相鄰兩個反射系數(shù)的時間間隔為10ms

4.2 實際地震資料應(yīng)用

4.2.1 三參數(shù)小波頻譜分解

實際資料來源于渤海某靶區(qū)，目的層砂體薄厚不均，疊置關(guān)系復(fù)雜，砂體橫向分布不穩(wěn)定，砂、泥巖交替頻繁。靶區(qū)內(nèi)河道眾多，多期河道縱、橫分布，橫向變化快。對于河流相儲層而言，刻畫河道的展布及邊界，弄清河道之間的疊置關(guān)系是儲層預(yù)測中非常重要的工作。

三參數(shù)小波頻譜分解的關(guān)鍵在于參數(shù)的選取，結(jié)合前文的研究以及實際地震資料，參數(shù)選取時需要注意以下問題： 第一，σ取較小值更利于刻畫地質(zhì)細節(jié)。第二，當(dāng)τ值較大(≥3)時，由于σ對小波形態(tài)的影響變小，所以不同σ的頻譜分解結(jié)果差異不大(圖5)。如由三參數(shù)小波頻譜分解切片(圖8)可見， 即使σ存在很大差異，處理結(jié)果也沒有太大區(qū)別。第三，前文研究表明，令β=0可匹配零相位地震子波，β≠0雖然可以更好地匹配非零相位子波，但由于β非零造成基小波時移，則頻譜分解后的數(shù)據(jù)體也會產(chǎn)生時移，導(dǎo)致對地質(zhì)體埋深的錯誤解釋。

圖9為不同β的三參數(shù)小波頻譜分解切片。由圖可見：當(dāng)β=0時，在目的層下方38ms的切片(圖9a)上可見河道信息；當(dāng)β=5時，在目的層下方38ms的切片(圖9b)上見不到任何河道信息，而在目的層下方92ms的切片(圖9c)上可見河道信息?？梢?，β非零會造成對河道埋深的錯誤解釋。

圖8 三參數(shù)小波頻譜分解切片(40Hz)(a)Γ=(1,3,0)； (b)Γ=(6,3,0)

圖9 不同β的三參數(shù)小波頻譜分解切片(40Hz)(a)Γ=(2.5,0.5,0)(目的層之下38ms)； (b)Γ=(2.5,0.5,5)(目的層之下38ms)； (c)Γ=(2.5,0.5,5)(目的層之下92ms)

針對該區(qū)目的層段的實際情況提取地震子波，并與不同的三參數(shù)小波做相關(guān)，優(yōu)選相關(guān)性大的參數(shù)組合作為該目的層段的小波[7]。本文選擇Γ=(1,0.5,0)進行頻譜分解，并與Morlet小波的頻譜分解效果進行對比。對于河流相儲層而言，主要從河道展布及邊界刻畫是否清晰作為評價標(biāo)準(zhǔn)。圖10為三參數(shù)小波和Morlet小波的頻譜分解結(jié)果。由圖可見：在20Hz切片上，Morlet小波的信噪比低，對河道邊界(紅色圓圈內(nèi))的刻畫不如三參數(shù)小波清晰；在80Hz切片上，在Morlet小波展示的河道(藍色圓圈內(nèi))展布不連續(xù)，邊界不清晰，而三參數(shù)小波的頻譜分解效果較好。綜上所述，在該區(qū)的河流相儲層預(yù)測中，三參數(shù)小波頻譜分解的效果優(yōu)于Morlet小波。

4.2.2 RGB混色融合

利用頻譜分解預(yù)測儲層時，可采用分頻解釋技術(shù)獲得離散的單頻體，但由于某個地質(zhì)體在時頻域的響應(yīng)不僅僅對應(yīng)一種頻率，所以每個單頻體只反映地質(zhì)體的一個時間厚度信息，為了提高解釋精度，需要綜合這些單頻信息，RGB混色融合技術(shù)的基本思想即源于此[22]。這里的R代表紅色，G代表綠色，B代表藍色，將三種顏色按照不同比例疊加，可以產(chǎn)生不同的顏色。本文用20Hz數(shù)據(jù)體作為紅色分量，40Hz數(shù)據(jù)體作為綠色分量，80Hz數(shù)據(jù)體作為藍色分量，即R、G、B分別代表低、中、高頻。圖11為20、40和80Hz切片的RGB融合顯示。由圖可見，通過RGB融合顯示，河道的能量更強，邊界更清晰，展布更連續(xù)。

圖10 三參數(shù)小波(PWT)和Morlet小波(WT)的頻譜分解結(jié)果

圖11 20Hz、40Hz和80Hz切片的RGB融合顯示

5 結(jié)論

(1)三參數(shù)小波的靈活性高，σ影響小波的震蕩程度。τ影響小波的寬窄，當(dāng)τ較大時，σ對小波的影響減弱。β的影響頗為復(fù)雜，當(dāng)其為周期的整數(shù)倍時，小波只產(chǎn)生時移，可以匹配零相位子波；當(dāng)其為周期的非整數(shù)倍時，小波不僅產(chǎn)生時移，而且產(chǎn)生相位延遲，可以匹配非零相位子波。同時需要注意，β非零造成河道埋深解釋錯誤。

(2)三參數(shù)小波變換較Morlet小波變換具有較好的時頻分辨率，能更好地刻畫薄互層內(nèi)細微的地質(zhì)沉積結(jié)構(gòu)。三參數(shù)小波變換處理結(jié)果的信噪比高，展示的河道邊界清晰、砂體展布連續(xù)。因此基于三參數(shù)小波的頻譜分解方法可以高精度地預(yù)測儲層。

(3)實際處理中可將從目的層段提取的地震子波與三參數(shù)小波做相關(guān)，優(yōu)選相關(guān)性大的參數(shù)組合，從而獲得最佳的三參數(shù)小波。