張繁昌，李傳輝

（中國(guó)石油大學(xué)（華東）地球科學(xué)與技術(shù)學(xué)院，山東青島 266580）

地震資料分頻解釋技術(shù)利用地震信號(hào)豐富的頻率信息，減少了常規(guī)解釋的不確定性，在薄儲(chǔ)層及小尺度地質(zhì)目標(biāo)的識(shí)別中起到了重要作用。唐湘蓉等利用地震波的高頻信息進(jìn)行了薄砂體預(yù)測(cè)［1］；胥德平等探討了基于廣義S變換分頻技術(shù)的儲(chǔ)層識(shí)別方法［2］；張志讓等將頻譜成像技術(shù)應(yīng)用于特殊地質(zhì)體的儲(chǔ)層預(yù)測(cè)［3］。此外，不同頻率下的振幅響應(yīng)（AVF）剖面還攜帶有儲(chǔ)層物性信息，Haitao等根據(jù)Biot飽和流體孔隙介質(zhì)地震波的傳播理論，設(shè)計(jì)了具有不同巖石物性參數(shù)（孔隙度、含油氣飽和度、滲透率等）的模型，并研究?jī)?chǔ)層地震響應(yīng)在低頻、中頻和高頻的變化規(guī)律與機(jī)理［4］；并參照AVO 類型的劃分，將AVF類型劃分為3類［5］，為利用AVF 規(guī)律進(jìn)行儲(chǔ)層物性解釋提供了理論指導(dǎo)。

目前地震資料處理中通常采用兩種分頻技術(shù)，一種是帶通濾波分頻，一種是利用Morlet小波變換的多分辨特性進(jìn)行分頻。由于Morlet小波變換具有優(yōu)良的時(shí)頻局部化性質(zhì)，可以將地震信號(hào)分解為一系列具有中心頻率的窄帶信號(hào)，較好地實(shí)現(xiàn)不同尺度地震信號(hào)的分離，近年來(lái)被廣泛應(yīng)用于地震分頻處理之中［6－9］。

雖然小波分頻技術(shù)在地震資料處理中得到廣泛應(yīng)用，但是Morlet小波變換方法得到的AVF剖面從低頻到高頻存在同相軸樹(shù)形分叉、上下漂移等現(xiàn)象，難以追蹤某一同相軸的振幅隨頻率變化規(guī)律。此外，Morlet小波變換得到的同頻率剖面存在調(diào)諧效應(yīng)，產(chǎn)生平行同相軸等假象。

Mallatt等提出的匹配追蹤方法［10－11］具有很多優(yōu)勢(shì)，由匹配追蹤得到的地震信號(hào)瞬時(shí)譜聚焦性最高［12－13］，具有比Morlet小波變換更高的時(shí)間或頻率分辨率。地震信號(hào)經(jīng)匹配追蹤分解后，可表示為一系列時(shí)頻原子的組合，所以由匹配追蹤時(shí)頻原子完全重構(gòu)地震信號(hào)非常簡(jiǎn)單，但利用時(shí)頻原子構(gòu)建AVF剖面和給定頻率的同頻率剖面卻是研究的難點(diǎn)。

1 方法原理

匹配追蹤分解算法通過(guò)創(chuàng)建超完備時(shí)頻原子庫(kù)，根據(jù)信號(hào)自身的特點(diǎn)將信號(hào)在時(shí)頻原子庫(kù)中展開(kāi)，以實(shí)現(xiàn)信號(hào)的自適應(yīng)分解。

設(shè)地震信號(hào)為s（t），D為進(jìn)行信號(hào)分解的超完備時(shí)頻原子庫(kù)，D中的每個(gè)時(shí)頻原子gγ均滿足有限支撐性質(zhì)。匹配追蹤算法通過(guò)一步步重復(fù)迭代，將信號(hào)s（t）垂直投影到D的時(shí)頻原子上。經(jīng)過(guò)匹配追蹤分解后，地震信號(hào)表示為不同振幅、中心時(shí)間、主頻和相位的時(shí)頻原子的線性組合［10］：

其中，an，tn，fn，φn分別表示第n個(gè)時(shí)頻原子的幅度、中心時(shí)間、主頻和相位。將第n個(gè)時(shí)頻原子用gγn（t）表示，即grn（t）＝gr（t－tn，fn，φn），則（1）式表示為

公式（2）說(shuō)明地震信號(hào)經(jīng)過(guò)匹配追蹤分解后，就可以用一系列時(shí)頻原子的線性組合來(lái)進(jìn)行重構(gòu)，但該式只能進(jìn)行地震信號(hào)的完全重構(gòu)，而不能重構(gòu)給定頻率的同頻率信號(hào)。要利用時(shí)頻原子構(gòu)建同頻率信號(hào)，需借助匹配追蹤瞬時(shí)譜來(lái)實(shí)現(xiàn)。

目前，匹配追蹤瞬時(shí)譜的計(jì)算普遍用各時(shí)頻原子的Wigner－Ville分布表示［10，14］，但這種表示方式只能提供地震信號(hào)的振幅分布。要利用匹配追蹤時(shí)頻原子構(gòu)建同頻率剖面，實(shí)現(xiàn)地震數(shù)據(jù)的分頻處理，需要構(gòu)建包含相位信息的瞬時(shí)譜形式，而不僅僅是振幅信息。

由于時(shí)頻原子具有良好的有限支撐性質(zhì)，其能量集中在以中心時(shí)間和主頻為中心的時(shí)頻點(diǎn)附近。設(shè)時(shí)頻原子gγn（t）的頻譜為Gγn（f），振幅包絡(luò)為env［gγn（t）］，將振幅包絡(luò)在頻率方向按照其頻譜Gγn（f）加權(quán)，就得到該時(shí)頻原子的瞬時(shí)譜

通過(guò)匹配追蹤分解，地震信號(hào)被分解成一系列時(shí)頻原子gγn（t），所有原子瞬時(shí)譜的疊加即為地震信號(hào)的匹配追蹤瞬時(shí)譜：

回溯古時(shí)，有莊周不斷捫心自問(wèn)，在人生路途中尋覓到自己人生的意義。在他人都痛拍欄桿，嗟嘆世道不公，意欲吞吐天地時(shí)，他卻疑惑、思索。充耳為蝸角虛名，滿目為勾心斗角，他毅然轉(zhuǎn)身，洞察了超然物我之外的“蝴蝶”，發(fā)現(xiàn)了摶扶搖而上的“大鵬”，了悟了生死榮辱之外的至理。生而為莊周，他在一路探尋中摸索，在摸索中得到答案?！兜赖陆?jīng)》有言：“夫唯不爭(zhēng)，故天下莫能與之爭(zhēng)。”莊子以其“不爭(zhēng)”之人生，活出了獨(dú)特的風(fēng)骨，發(fā)現(xiàn)了自己的人生意義。

與Wigner－Ville分布計(jì)算匹配追蹤瞬時(shí)譜的方法相比，（4）式表達(dá)的匹配追蹤瞬時(shí)譜不但計(jì)算簡(jiǎn)潔，而且同時(shí)包含了振幅及相位信息。

同理，將所有時(shí)頻原子gγn（t）的波形在頻率方向按照Gγn（f）加權(quán)，即可以得到地震信號(hào)的AVF剖面：

該式反映了地震信號(hào)的振幅隨頻率變化規(guī)律。將（5）式沿頻率方向積分，得：

AVF剖面沿頻率方向的積分應(yīng)當(dāng)為原地震信號(hào)，但比較（6）式和（2）式發(fā)現(xiàn)，二者并不相等，（6）式多了這一項(xiàng)，說(shuō)明由（6）式重構(gòu)的地震信號(hào)的振幅大小與原始地震信號(hào)并不相同，換句話說(shuō)，（6）式不能實(shí)現(xiàn)地震信號(hào)的保幅重構(gòu)。

為了實(shí)現(xiàn)地震信號(hào)的保幅重構(gòu)，將（5）式改為

即在計(jì)算AVF 剖面時(shí)，首先對(duì)Gγn（f）進(jìn)行歸一化。這樣，將（7）式對(duì)頻率積分，得

對(duì)比（8）式和（2）式，可見(jiàn)按照（8）式重構(gòu)的信號(hào)就是原來(lái)的地震信號(hào)，實(shí)現(xiàn)了地震信號(hào)的保幅重構(gòu)。

由（7）式得到振幅保持的地震AVF 剖面后，給定某一頻率fj，就可以提取該頻率的同頻率信號(hào)DF（t，fj）：

與（2）式不同的是，（9）式表示在重構(gòu)過(guò)程中，雖然所有時(shí)頻原子均參與計(jì)算，但每個(gè)時(shí)頻原子被賦予不同的權(quán)重，其大小由頻率fj位置處各自的頻譜值決定。（9）式所表示的同頻率信號(hào)構(gòu)建方法由于包含了所有時(shí)頻原子而使得頻帶較寬，同時(shí)，不同時(shí)頻原子根據(jù)其位于fj頻率處的譜值大小，對(duì)同頻率信號(hào)的貢獻(xiàn)有主次之分，使地震信號(hào)的主要能量集中在頻率fj附近。地震剖面在匹配追蹤分解后，每一道都按（9）式進(jìn)行計(jì)算，就得到頻率為fj時(shí)的同頻率剖面。

2 方法測(cè)試

匹配追蹤將地震信號(hào)分解成眾多時(shí)頻原子，各時(shí)頻原子均有不同的主頻、中心時(shí)間、振幅和相位，這些時(shí)頻原子按照（2）式就可以完全重構(gòu)原來(lái)的地震信號(hào)。那么，要利用匹配追蹤結(jié)果對(duì)地震信號(hào)進(jìn)行分頻處理，即有選擇地構(gòu)建給定頻率的信號(hào)，是否可以通過(guò)將此頻率范圍內(nèi)的時(shí)頻原子也按（2）式實(shí)現(xiàn)呢？根據(jù)這個(gè)設(shè)想，對(duì)圖1a所示的實(shí)際地震剖面先進(jìn)行匹配追蹤分解，然后利用主頻為30 Hz的時(shí)頻原子構(gòu)建了30 Hz同頻率剖面（圖1b）。分析發(fā)現(xiàn)，圖1b的同頻率剖面出現(xiàn)很多空白區(qū)，同相軸時(shí)連時(shí)斷，不能很好地保持原地震剖面的地質(zhì)特征，說(shuō)明僅由給定頻率的時(shí)頻原子并不能合理地構(gòu)建同頻率剖面。

圖1c為利用圖1a地震剖面匹配追蹤分解的時(shí)頻原子按（9）式得到的30 Hz同頻率剖面，通過(guò)與圖1b的對(duì)比可見(jiàn)，由（9）式計(jì)算的同頻率剖面不存在空白區(qū)，很好地保持了原地震剖面的特征。

圖1 地震剖面及不同方法得到的30 Hz同頻率剖面

3 結(jié)果對(duì)比

3.1 瞬時(shí)譜的對(duì)比

從圖1a的地震剖面中任取一道地震信號(hào)，例如第10道，顯示于圖2a左側(cè)，對(duì)其進(jìn)行匹配追蹤分解，然后利用（4）式計(jì)算此地震道的瞬時(shí)譜，其振幅分布如圖2a右側(cè)所示。為了對(duì)比，將該地震信號(hào)匹配追蹤結(jié)果利用Wigner－Ville分布方法計(jì)算其振幅分布（圖2b）。對(duì)比圖2a和圖2b可見(jiàn)，利用（4）式得到的瞬時(shí)譜，無(wú)論是在時(shí)頻分辨率還是在能量聚集性上與Wigner－Ville分布方法的計(jì)算結(jié)果都幾乎完全相同，說(shuō)明了（4）式的正確性。

3.2 AVF剖面的對(duì)比

仍以圖2a左側(cè)的地震信號(hào)為例。利用（7）式計(jì)算得到該地震信號(hào)的AVF 剖面（圖3a右側(cè)）。圖3b為同一地震信號(hào)利用Morlet小波變換得到的AVF剖面。由圖3 可以看出，兩種分頻結(jié)果的主要能量分布一致，但是沿頻率方向，圖3b的AVF剖面表現(xiàn)出“樹(shù)形分叉”現(xiàn)象，同相軸彎曲、分叉，難以追蹤某一同相軸的振幅隨頻率變化規(guī)律；而在圖3a中的匹配追蹤AVF剖面上，同相軸平直，沒(méi)有上下漂移現(xiàn)象，可以方便地進(jìn)行AVF分析。

圖2 不同方法得到的瞬時(shí)譜對(duì)比

圖3 不同方法得到的AVF剖面對(duì)比

3.3 同頻率剖面的對(duì)比

對(duì)圖1a所示的地震剖面進(jìn)行匹配追蹤，利用（9）式的同頻率信號(hào)構(gòu)建方法對(duì)地震剖面逐道分頻，得到每個(gè)地震道的AVF剖面，圖4a為不同地震道的AVF剖面。經(jīng)過(guò)匹配追蹤分頻處理后，得到三維頻率數(shù)據(jù)體，再沿頻率方向逐道取出特定頻率的信號(hào)，便得到同頻率剖面，圖4b為不同頻率的剖面排列在一起的三維顯示。

圖4 匹配追蹤分頻處理獲得的三維頻率數(shù)據(jù)體

圖5為從圖4b中取出的15，30和50 Hz的匹配追蹤同頻率剖面。為了對(duì)比，圖6 給出了Morlet小波變換得到的對(duì)應(yīng)頻率剖面?？梢钥闯?，匹配追蹤分頻方法與Morlet小波變換一樣具有多分辨率特性，能夠分離出不同尺度的地震信號(hào)，利用不同頻率的剖面可以揭示不同厚度、不同規(guī)模的地層反射特征。

進(jìn)一步對(duì)比發(fā)現(xiàn)，圖6所示的Morlet小波變換同頻率剖面中存在平行同相軸假象，這是由于濾波造成的調(diào)諧效應(yīng)。而在圖5 的匹配追蹤同頻率剖面上則不存在這種調(diào)諧效應(yīng)，不同頻率的剖面很好地保持了圖1a原地震剖面的反射特征。

圖5 匹配追蹤同頻率剖面

圖6 Morlet小波變換同頻率剖面

觀察圖5和圖6還可以看出，15Hz低頻剖面表現(xiàn)為淺層較弱、深層較強(qiáng)；而50 Hz高頻剖面表現(xiàn)為深層較弱、淺層較強(qiáng)；只有30 Hz剖面深淺層反射能量比較均衡。這是由于地震波在地下傳播過(guò)程中，高頻成分衰減較快，而地震資料處理通常是依據(jù)地震波的優(yōu)勢(shì)頻帶進(jìn)行，不能很好地兼顧高、低頻成分。將圖5a和圖5c分別進(jìn)行Q掃描［14－15］等處理，得到如圖7所示的結(jié)果。可以看出，不同頻率成分的剖面上，淺、中、深層能量都保持均衡。

圖7 對(duì)圖5數(shù)據(jù)進(jìn)行剩余補(bǔ)償后的同頻率剖面

4 結(jié)束語(yǔ)

匹配追蹤分解方法得到的瞬時(shí)譜具有極高的聚焦性，如何利用高聚焦性的時(shí)頻原子構(gòu)建AVF剖面及給定頻率的剖面卻是研究的難點(diǎn)。我們?cè)谄ヅ渥粉櫵惴ǖ幕A(chǔ)上，提出了一種保幅AVF剖面構(gòu)建方法。該方法將所有時(shí)頻原子在頻率方向上按其歸一化頻譜加權(quán)，即可得到地震信號(hào)的AVF 剖面。同Morlet小波變換得到的AVF剖面相比，本方法得到的AVF 剖面上沒(méi)有出現(xiàn)同相軸彎曲和分叉現(xiàn)象。

在取得AVF 剖面的基礎(chǔ)上，進(jìn)一步構(gòu)建了匹配追蹤同頻率剖面，即利用匹配追蹤時(shí)頻原子對(duì)給定頻率進(jìn)行地震同頻率信號(hào)的構(gòu)建，構(gòu)建過(guò)程中所有時(shí)頻原子均參與計(jì)算，其貢獻(xiàn)大小由該頻率位置處各時(shí)頻原子的譜值決定。利用該方法提取的同頻率剖面，有效避免了以濾波機(jī)制為基礎(chǔ)的分頻方法造成的平行同相軸假象，分頻效果更理想。由此解決了匹配追蹤分頻重構(gòu)技術(shù)難題。

［1］唐湘蓉，蔡涵鵬，賀振華.地震波高頻信息在薄層砂體預(yù)測(cè)中的應(yīng)用［J］.石油物探，2012，51（3）：244－250 Tang X R，Cai H P，He Z H.Thin－bed sand body prediction based on seismic wave high－frequency information［J］.Geophysical Prospecting for Petroleum，2012，51（3）：244－250

［2］胥德平，郭科，文曉濤，等.基于廣義S變換和JADE算法的儲(chǔ)層識(shí)別［J］.石油物探，2011，50（4）：319－323 Xu D P，Guo K，Wen X T，et al.Reservoir identification based on GST and JADE algorithm［J］.Geophysical Prospecting for Petroleum，2011，50（4）：319－323

［3］張志讓，楊建勛，李梅，等.頻譜成像技術(shù)在特殊地質(zhì)體儲(chǔ)層預(yù)測(cè)中的應(yīng)用研究［J］.石油物探，2011，50（1）：69－75 Zhang Z R，Yang J X，Li M，et al.Application of frequency spectrum imaging technique for reservoir prediction of special geological bodies［J］.Geophysical Prospecting for Petroleum，2011，50（1）：69－75

［4］Haitao R，Gennady G，F(xiàn)red H.Amplitude versus frequency variations in thinly layered porous rocks［J］.Expanded Abstracts of 78thAnnual Internat SEG Mtg，2008，1744－1748

［5］Haitao R，Gennady G，F(xiàn)red H.Poroelastic analysis of amplitude－versus－frequency variations［J］.Geophysics，2009，74（6）：41－48

［6］Avijit C，David O.Frequency－time decomposition of seismic data using wavelet－based methods［J］.Geophysics，1995，60（6）：1906－1916

［7］龔洪林，王振卿，李錄明，等.應(yīng)用地震分頻技術(shù)預(yù)測(cè)碳酸鹽巖儲(chǔ)層［J］.地球物理學(xué)進(jìn)展，2008，23（1）：129－135 Gong H L，Wang Z Q，Li L M，et al.Predicting carbonate reservoir by applying seismic spectral decomposition technique［J］.Progress in Geophysics，2008，23（1）：129－135

［8］魏文，王小杰，李紅梅.基于疊前道集小波域Q 值求取方法研究［J］.石油物探，2011，50（4）：355－360 Wei W，Wang X J，Li H M.Study on extraction method for Q based on pre－stack gather in wavelet domain［J］.Geophysical Prospecting for Petroleum，2011，50（4）：355－360

［9］Marcilio C M，Kurt J M，Paulo R S，et al.Wavelet transform Teager－Kaiser energy applied to a carbonate field in Brazil［J］.The Leading Edge，2009，28（3）：708－713

［10］Mallat S，Zhang Z.Matching pursuit with time－frequency dictionaries［J］.IEEE transactions on signal processing，1993，41（12）：3397－3415

［11］張繁昌，李傳輝.基于正交時(shí)頻原子的地震信號(hào)快速匹配追蹤［J］.地球物理學(xué)報(bào)，2012，55（1）：277－283 Zhang F C，Li C H.Orthogonal time－frequency atom based fast matching pursuit for seismic signals［J］.Chinese Journal of Geophysics（in Chinese），2012，55（1）：277－283

［12］Castagna J P，Sun S J，Siegfried R W.Instantaneous spectral analysis：Detection of low－frequency shadows associated with hydrocarbons［J］.The Leading Edge，2003，22（2）：120－127

［13］李傳輝，張繁昌.地震信號(hào)可變分辨率匹配追蹤頻譜成像方法［J］.石油物探，2012，51（3）：213－218 Li C H，Zhang F C.Variable resolution matching pursuit spectrum imaging of seismic signals［J］.Geophysical Prospecting for Petroleum，2012，51（3）：213－218

［14］Wang Y H.Seismic time－frequency spectral decomposition by matching pursuit［J］.Geophysics，2007，72（1）：13－20

［15］魏文，王興謀，李紅梅，等.基于地震波衰減的特征屬性重構(gòu)方法［J］.石油物探，2012，51（3）：219－225 Wei W，Wang X M，Li H M，et al.Method of characteristics attributes reconstruction based on seismic wave attenuation［J］.Geophysical Prospecting for Petroleum，2012，51（3）：219－225