基于三參數(shù)小波的頻譜分解方法

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

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

1 引言

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

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

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

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

2 頻譜分解方法

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

(1)

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

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

(1)選定目的層段；

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

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

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

圖2 分頻解釋

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

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

(2)

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

(3)

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

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

iq(Γ)sin(σt)}

(4)

其中

(5)

(6)

(7)

Γ=(σ,τ,β)

(8)

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

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

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

圖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í)不同σ的三參數(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)

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

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

4.1 模型測(cè)試

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

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

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

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

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

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

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

圖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)

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

4.2.2 RGB混色融合

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

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

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

5 結(jié)論

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

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

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