韓海英，王志章，王宗俊,張雨晴

(1.中國石油大學(北京)地球科學學院，北京102249;2.中國石油天然氣股份有限公司勘探開發(fā)研究院，北京100083;3.中海油研究總院，北京100027)

近年來，信號時頻分解技術(shù)在地震資料處理和解釋中得到了廣泛的應(yīng)用，如高分辨率處理[1]、烴類檢測[2-3]等。常規(guī)的時頻分解(如窗口傅里葉變換和小波變換等)，時頻窗受不確定準則約束，時頻分辨率難以同時達到最佳。為了克服常規(guī)方法的缺點，Mallat等[4]和Qian等[5]提出了匹配追蹤分解算法(Matching Pursuit，MP)。MP算法基于最佳匹配原則，通過不斷迭代，從原子庫中選取與殘余信號結(jié)構(gòu)最佳匹配的原子，將原始信號分解為一系列原子。由于在迭代尋優(yōu)過程中，時窗長度均由信號局部特征決定，因此這些最優(yōu)原子反映了信號的局部特征，能較好地描述時變地震信號的時頻分布特性。

經(jīng)典的MP算法基于Gabor函數(shù)構(gòu)建超完備子波庫，為了適應(yīng)信號的不同結(jié)構(gòu)特征，Liu等[6-7]先后提出了基于Ricker子波和Morlet子波的匹配追蹤算法。對研究薄互層儲層特征而言，Morlet小波分析薄層的效果欠佳;宋新武等[8]認為Ricker子波波形簡單，延遲時間短，收斂較快，基于Ricker子波匹配追蹤算法的頻譜分解技術(shù)在薄互層分析中具有較高的時頻分辨率，能較好地刻畫地質(zhì)體形態(tài)。原子庫的合理選取，對于信號的匹配追蹤分解效果具有重要作用。

Ricker類地震子波是由Ricker子波演變得到的零相位理論子波，具有3個控制參量，旁瓣能量的收斂速度可變，延遲時間可控，子波庫波形更豐富，與地震子波的匹配較為靈活[9]。我們在研究子波波形隨3個控制參量變化關(guān)系的基礎(chǔ)上，對該子波作適當?shù)暮喕?，建立參量與波形特征之間的半定量關(guān)系。最后，以簡化的子波構(gòu)建原子庫，通過三步法原理[10]實現(xiàn)信號的匹配追蹤分解。

1 基本原理

1.1 Ricker類地震子波

Ricker類地震子波頻域解析式為[9]

(1)

對于確定的k，m越大，高頻能量相對越強;對于確定的m，k越大，高頻能量相對越弱，同時頻譜的有效頻帶變窄，導致子波信號的旁瓣變化復雜。文獻[9]建議k不宜過大。事實上，參量k和m都是控制波譜高低頻相對強弱的參數(shù)，一般情況下，相同的增(減)量，m對頻譜形態(tài)的改變較為平緩。

(2)

這里記C(f)為二參量Ricker類地震子波，簡稱C子波。其中，fm和c為待定參量;A為歸一化系數(shù)。通過簡化和改寫，將子波譜表述成Ricker譜(主頻為fm)的c次方(c=1時對應(yīng)Ricker子波)。fm為波譜主頻，控制波譜的等比例伸縮;c控制波譜主頻兩側(cè)能量的衰減速度。

下面考察參量c與子波波形的關(guān)系。首先，針對連續(xù)相位數(shù)小于1的對稱波形(圖1中藍線)，我們定義主、旁瓣極值比PR和主、旁瓣寬度比WR兩個參量，用于描述子波的波形特征。

主旁瓣極值比(PR)：旁瓣極值與主瓣極值之比的絕對值;

主旁瓣寬度比(WR)：旁瓣寬度(極小值點間距)與主瓣寬度(零值點間距)之比。

圖1 對稱C子波示意圖解

圖2為C子波波形特征值(紅色為極值比PR線，藍色為寬度比WR線)隨參量c的變化曲線。隨著c的增加，主旁瓣極值比增加，旁瓣能量越強;主旁瓣寬度比增加，說明在相同主瓣寬度的情況下，旁瓣能量達到極值需要的時間變長。因此，對C子波有如下結(jié)論：參數(shù)fm控制子波的橫向伸縮，fm越大，子波越短;參數(shù)c控制子波的波形形態(tài)，c越大，子波的主旁瓣極值比、寬度比越大，連續(xù)相位數(shù)增加。當c→0時，C子波趨近于脈沖信號，PR→0，WR→1；當c→∞時，C子波趨近于單頻信號，PR→1，WR→2。

圖2 主旁瓣寬度比(藍)和極值比(紅)曲線

1.2 基于C子波的匹配追蹤

匹配追蹤是基于原子庫掃描的信號自適應(yīng)分解。記D為進行信號分解的超完備子波庫，可由Hilbert空間的窗口函數(shù)通過擴張、調(diào)制和平移生成，定義D={wγ:γ=1,2,3,…}，wγ為單個原子，滿足歸一化條件‖wγ‖=1。記待分解信號為s(t)，經(jīng)N次迭代分解如下：

(3)

式中：wγ,n為第n次迭代得到的匹配子波;an為wγ,n對應(yīng)的振幅；RNs為迭代N次后的殘差。

匹配追蹤算法中子波庫的波形與原信號結(jié)構(gòu)越相似，分解效果越好，這里選擇C子波來構(gòu)建匹配子波庫。匹配子波wγ,n由參量γn={tn,fm,n,cn,φn}來控制，其中，tn，φn為子波的中心時間和相位;fm,n，cn為波譜控制參量。匹配追蹤算法基于不斷迭代的貪婪算法，常規(guī)算法的每一次迭代都要從子波庫中尋找最佳子波。目前，通常采用“三步法”[10]原理來提高分解速度。

(4)

(5)

3) 估算最優(yōu)子波的幅值。按(6)式得到最優(yōu)子波的幅值an：

(6)

每一次迭代，都將生成匹配子波anwγ,n和殘差信號Rn+1s。分解完成后，即可用提取的匹配子波近似表示原始信號，即

(7)

1.3 時頻表征

常規(guī)的匹配追蹤時頻表征通過各匹配子波的Wigner-Ville分布疊加實現(xiàn)，但只能給出振幅分布，無法得到相位信息。這里，采用(8)式所示的時頻表征方法[12]，對每個匹配子波的復譜求和，獲取信號的時頻譜：

(8)

式中：Wγ,n(f)是匹配子波wγ,n(t)的頻譜;env[wγ,n(t)]是wγ,n(t)的瞬時包絡(luò)。

圖3是單個C原子(主頻fm取50Hz，c取0.7)的時頻譜。該時頻表征法具有較高的時頻分辨率，能同時獲取振幅和相位等信息的時頻分布，并且避免了窗口截斷效應(yīng)。

圖3 C原子時頻譜

2 模型試算

設(shè)計兩個合成信號進行基于C子波庫的匹配追蹤算法測試，時間采樣率均為1ms。

圖4為信號1基于C子波的匹配追蹤重構(gòu)結(jié)果，其中,圖4a為原始信號，圖4b為重構(gòu)信號，圖4c 為重構(gòu)相對誤差，重構(gòu)誤差保持在2%以內(nèi)，對信號的損傷較小(重構(gòu)精度與算法終止條件有關(guān)，精度越高，運行時間越長)。信號組成方式如圖5a 所示，其中，第1道為原始信號(圖4a)，為其余4道的疊加;圖5b為信號1的匹配追蹤時頻譜。圖6a為信號1的短時Fourier變換(STFT)譜，圖6b 為S變換譜。對比發(fā)現(xiàn)，匹配追蹤與常規(guī)時頻分解方法在分辨率上存在較大的差異。在300ms和400ms處，3個時頻譜上都有單個能量團出現(xiàn)，無論是時間尺度還是頻率尺度，匹配追蹤的分辨率都是最高的;在500ms處，3種方法都沒能較好地把兩個頻率分辨出來，但匹配追蹤的時間分辨率還是最高的;在600ms處，S變換和STFT譜的時間分辨率不夠，匹配追蹤則較好地識別出2個反射;在700ms和800ms處，匹配追蹤譜的時頻分辨能力明顯高于STFT譜和S變換譜，能反映更多的地層薄層信息。因此，匹配追蹤時頻譜具有較高的時頻聚焦性。

圖4 信號1重構(gòu)圖解a 原始信號; b 重構(gòu)信號; c相對誤差

信號2(圖7)由Ricker類理論子波與隨機反射系數(shù)褶積生成。為符合實際情況，由淺到深理論子波的主頻和帶寬逐漸減小，且主頻與帶寬的比值隨機改變。圖7a為隨機反射系數(shù)，圖7b為合成記錄，從上往下信號的分辨率降低，圖7c為匹配追蹤結(jié)果，圖7d為相對誤差百分比。為作對比，對信號2同時進行基于Ricker子波的匹配追蹤(算法的終止條件不變)，時頻譜如圖8所示。其中，圖8a為圖7b基于C子波的匹配分解結(jié)果，圖8b為基于Ricker子波的匹配分解結(jié)果。對比發(fā)現(xiàn)，兩個分解結(jié)果存在一定的差異：C子波匹配追蹤譜主頻隨時間呈近似下降的趨勢，符合模型設(shè)計，而Ricker子波匹配譜在600ms后主頻呈遞增的趨勢，與實際不符;C子波匹配譜的強能量團與信號的峰值分布符合度較高，如信號2在680～740ms有兩個強峰值，C子波匹配譜上有2個強能量團，而Ricker子波匹配譜上則有3個，在800ms附近存在強峰值，C子波匹配譜上存在強能量顯示，Ricker子波匹配譜上能量不明顯。

圖5 信號1的分解結(jié)果a 原始信號組成示意; b 基于C子波的匹配追蹤時頻譜

圖6 信號1的短時Fourier變換(a)和S變換(b)時頻譜

圖7 信號2a 反射系數(shù); b 合成記錄; c 匹配結(jié)果; d 相對誤差

圖8 信號2匹配追蹤時頻譜a 基于C子波; b 基于Ricker子波

對比結(jié)果表明，匹配追蹤的分解效果依賴于原子庫與信號基本結(jié)構(gòu)的匹配度。C原子庫的參量調(diào)控靈活，波形較Ricker子波庫豐富，對實際地震信號基本結(jié)構(gòu)(時變地震子波)的適應(yīng)性更強。因此，對于子波特征值(主旁瓣極值比和寬度比)時變的地震信號，基于C子波的匹配效果優(yōu)于基于Ricker子波的匹配效果。

3 實際資料處理

圖9為陸上某探區(qū)的過井測線地震剖面。該區(qū)地下構(gòu)造簡單，地層平緩，斷裂較少，淺表氣遍布全區(qū)，第四系沉積以砂泥巖互層為主，成巖性較差，特殊的地質(zhì)條件造就了該區(qū)地震資料的特殊性。經(jīng)勘探發(fā)現(xiàn)該井存在工業(yè)氣流，綜合研究后認為該地區(qū)氣藏為鼻狀構(gòu)造背景上的巖性氣藏。根據(jù)以往經(jīng)驗，含油氣地區(qū)一般會有明顯地震異常。從圖9 可以看出，在井附近同相軸的連續(xù)性變差，有下拉趨勢，主頻變低，在主測線和聯(lián)絡(luò)線上地震相外形呈塔狀。這種地震異常雖然存在，但并非特別明顯。針對該地區(qū)的氣藏特征，應(yīng)用基于C子波的匹配追蹤算法，并對資料進行分頻處理，為該地區(qū)開展巖性氣藏研究提供依據(jù)。

圖9 某工區(qū)過井測線地震剖面

提取該過井剖面井旁道目的層段(1000～2000ms)地震記錄，分別用Ricker原子和C原子對該記錄進行匹配追蹤，結(jié)果如圖10所示。對比發(fā)現(xiàn)，兩種匹配結(jié)果都能實現(xiàn)地震信號的高精度匹配分解，較好地反映信號的時頻特征;在局部細節(jié)上，C原子匹配的時間分辨率優(yōu)于Ricker原子匹配，如：在1900ms附近(黑框標注部分)，C原子匹配較好地將20Hz主頻的兩個能量團區(qū)分開。因此，利用C原子對地震記錄進行匹配追蹤，更有助于薄層的精細刻畫。

圖11是對過井地震剖面進行分頻處理得到的單頻調(diào)諧體剖面，不同頻率剖面顯示的振幅相對強弱各不相同。通過頻率切片，可以清晰地觀察到各層位能量不同的展布和變換趨勢。圖11中黑色橫線為氣層所處深度(1650ms附近)，觀察井(紅線)附近目的層的能量，在12Hz左右的分頻剖面上達到最大，之后隨著頻率的增加而降低，這是因為薄層的調(diào)諧頻率在12Hz左右。在紅線與黑線相交位置附近，右側(cè)的能量在12Hz之后逐漸減小，而左側(cè)在25Hz剖面上仍有較強的能量，這是由于右側(cè)含氣層的影響。匹配追蹤分解結(jié)果為目的層油氣藏的預測和描述提供了可靠的依據(jù)。

圖10 過井地震剖面目的層段井旁單道匹配結(jié)果a 地震記錄; b Ricker原子匹配時頻譜; c C原子匹配時頻譜

圖11 過井地震剖面分頻處理得到的單頻調(diào)諧體剖面

4 結(jié)束語

Ricker類子波是一類新型的理論子波，波譜主頻以及高低頻的衰減速度由三參量聯(lián)合控制。根據(jù)實際需要，將該類子波簡化成二參量C子波，波譜主頻、帶寬由二參量獨立控制。針對子波波形，提出以主旁瓣極值比PR和主旁瓣寬度比WR來定量表征，并闡述了C子波波形比Ricker子波、寬帶Ricker子波更豐富，更能匹配時變地震子波。

根據(jù)C子波建立原子庫，用于信號的匹配追蹤算法。模型試算表明，與短時Fourier變換和S變換相比，匹配追蹤得到的時頻譜圖具有較好的時頻聚焦性;對于時變地震記錄，C子波與信號基本結(jié)構(gòu)的匹配度較Ricker子波高，匹配追蹤效果更好。最后，將基于C子波的匹配追蹤應(yīng)用于實際資料的分頻處理，取得了較好的解釋效果。

但同時需要指出的是，基于C子波的匹配追蹤算法的運算量較大，運行速度的提升是需要進一步改進優(yōu)化的方向。

參 考 文 獻

[1] 熊曉軍，賀振華，黃德濟，等．廣義S變換在地震高分辨處理中的應(yīng)用[J]．勘探地球物理進展，2006，29(6)：415-418

Xiong X J,He Z H,Huang D J,et al.Application of generalized S transform in seismic high resolution processing[J].Progress in Exploration Geophysics,2006,29(6):415-418

[2] Castagna J P,Sun S,Siegfried R W.Instantaneous spectral analysis:detection of low-frequency shadows associated with hydrocarbons[J].The Leading Edge,2003,22(2):120-127

[3] Sinha S,Routh P S,Anno P D,et al.Spectral decomposition of seismic data with continuous-wavelet transforms[J].Geophysics,2005,70(6):19-25

[4] Mallat S,Zhang Z.Matching pursuits with time-frequency dictionaries[J].IEEE Transactions on Signal Processing,1993,41(12):3397-3415

[5] Qian S,Chen D.Signal representation via adaptive normalized Gaussian function[J].Signal Processing,1994,36(1):1-11

[6] Liu J,Wu Y,Han D,et al.Time-frequency decomposition based on Ricker wavelet[J].Expanded Abstracts of 74thSEG Annual Internat Mtg,2004,1937-1940

[7] Liu J,Marfurt K J.Matching pursuit decomposition using Morlet wavelets[J].Expanded Abstracts of 75thSEG Annual Internat Mtg,2005,786-789

[8] 宋新武，鄭浚茂，范興燕，等．基于Ricker子波匹配追蹤算法在薄互層砂體儲層預測中的應(yīng)用[J]．吉林大學學報(地球科學版)，2011，41(增刊1)：387-392

Song X W,Zheng J M,Fan X Y,et al.Application of the thin-interbedded reservoir prediction based on Ricker wavelet match tracing algorithm[J].Journal of Jilin University(Earth Science Edition),2011,41(S1):387-392

[9] 梅金順，王潤秋．Ricker類地震子波[J]．石油地球物理勘探，2012，47(增刊1)：8-14

Mei J S,Wang R Q.Ricker-like seismic wavelets[J].Oil Geophysical Prospecting,2012,47(S1):8-14

[10] Wang Y.Seismic time-frequency spectral decomposition by matching pursuit[J].Geophysics,2007,72(1):13-20

[11] Barnes A E.Instantaneous frequency and amplitude at the envelope peak of a constant-phase wavelet[J].Geophysics,1991,56(7):1058-1060

[12] 張繁昌，李傳輝．非平穩(wěn)地震信號匹配追蹤時頻分析[J]．物探與化探，2011，35(4):546-552

Zhang F C,Li C H.Matching pursuit time-frequency analysis of non-stationary seismic signals[J].Geophysical and Geochemical Exploration,2011,35(4):546-552

[13] 張繁昌,李傳輝.匹配追蹤保幅地震AVF剖面及同頻率剖面的構(gòu)建[J].石油物探,2013,52(3):234-239

Zhang F C,Li C H.Amplitude-preserved AVF and isofrequency section construction based on seismic matching pursuit[J].Geophysical Prospecting for Petroleum,2013,52(3):234-239

[14] 李傳輝,張繁昌.地震信號可變分辨率匹配追蹤頻譜成像方法[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

[15] 張顯文,韓立國,王宇,等.地震信號譜分解匹配追蹤快速算法及其應(yīng)用[J].石油物探,2010,49(1):1-6

Zhang X W,Han L G,Wang Y,et al.Seismic spectral decomposition fast matching persuit algorithm and its application[J].Geophysical Prospecting for Petroleum,2010,49(1):1-6

[16] 陳發(fā)宇,楊長春.對依據(jù)頻率的匹配追蹤快速算法的改進[J].石油物探,2009,48(1):80-83

Chen F Y,Yang C C.An improvement for frequency-domainated fast matching pursuits algorithm[J].Geophysical Prospecting for Petroleum,2009,48(1):80-83

[17] 趙天姿,宋煒,王尚旭.基于匹配追蹤算法的時頻濾波去噪方法[J].石油物探,2008,47(4):367-371

Zhao T Z,Song W,Wang S X.Time-frequency filtering de-noise method based on matching pursuit algorithm[J].Geophysical Prospecting for Petroleum,2008,47(4):367-371