黃 躍，許 多，文雪康

（1.中國(guó)石油化工股份有限公司多波地震技術(shù)重點(diǎn)實(shí)驗(yàn)室，四川德陽 618000；2.中國(guó)石油化工股份有限公司西南油氣分公司，四川成都 610016；3.中國(guó)石油化工股份有限公司西南油氣分公司勘探開發(fā)研究院德陽分院，四川德陽 618000；4.中國(guó)石油化工股份有限公司西南油氣分公司監(jiān)理中心，四川德陽 618000）

地震道多子波分解（簡(jiǎn)稱多子波分解）就是將地震道分解為各種不同主頻、不同寬度以及不同時(shí)間位置的子波。將這些子波簡(jiǎn)單疊加在一起即可精確地恢復(fù)原有的地震道。針對(duì)目標(biāo)層位精細(xì)構(gòu)造解釋或油氣藏特征描述的需要，可以有選擇地對(duì)分解后的子波進(jìn)行疊加，重構(gòu)新的地震道，以突出相應(yīng)的波場(chǎng)和波組特征。Castagna等［1］將多子波分解所獲得的瞬時(shí)譜的低頻陰影用于含氣性檢測(cè)取得了較好的效果；趙爽等［2］應(yīng)用該方法進(jìn)行含氣性檢測(cè)；王昕等［3］將該方法應(yīng)用于儲(chǔ)層預(yù)測(cè)。

目前多子波分解方法主要采用確定性子波進(jìn)行地震道分解，人們對(duì)于子波的選擇以及不同子波類型對(duì)地震重構(gòu)信號(hào)分辨率、連續(xù)性等方面的影響尚不夠重視。事實(shí)上，在實(shí)際地震資料的多子波分解與重構(gòu)過程中，選擇不同的子波能夠突出不同的地震屬性特征。以復(fù)數(shù)Paul子波和復(fù)數(shù)Morlet子波為例，理論上復(fù)數(shù)Paul子波具有更好的時(shí)間分辨率，適合處理信噪比高但相位變化快的信號(hào)；而復(fù)數(shù)Morlet子波具有更好的頻率分辨率，時(shí)間分辨率相對(duì)較差，適合地震波吸收衰減研究。有鑒于此，根據(jù)不同的地質(zhì)任務(wù)和研究目的合理進(jìn)行子波選擇無疑是十分必要的。我們從多子波分解常用的匹配追蹤算法出發(fā)，通過理論計(jì)算和實(shí)際資料處理，對(duì)比分析Paul子波和Morlet子波對(duì)重構(gòu)信號(hào)的分辨率及連續(xù)性等方面的影響，探討子波選擇的重要性，以期為多子波分解方法的應(yīng)用，尤其是子波選取，提供指導(dǎo)和借鑒。

1 基本原理

地震道多子波分解中最常用的算法是匹配追蹤算法，它通過迭代將地震道分解成一系列與其時(shí)頻特征相匹配的子波［4－8］，是從一系列連續(xù)小波集合（子波字典）基函數(shù)中選取出與地震信號(hào)相匹配的一系列子波（也稱為子波原子）的運(yùn)算過程。在此過程中，每一次迭代均可獲得一個(gè)最優(yōu)的子波，經(jīng)過N次迭代之后，地震道S（t）被分解為N個(gè)最優(yōu)子波的代數(shù)和［9－14］，即

式中：an是第n個(gè)小波grn的振幅值；R（N）為殘差；S（0）為初始地震道，R（0）＝S（0）。

對(duì)于每一個(gè)小波，可以用設(shè)定的4個(gè)參數(shù)γ＝｛u，σ，ω，φ｝來刻畫，其中，u是時(shí)移，σ是時(shí)間軸尺度，ω是中心頻率，φ是相移。具體匹配追蹤過程如下。

1）對(duì)實(shí)際地震道進(jìn)行Hilbert變換后，構(gòu)建復(fù)地震道，將復(fù)地震道最大振幅包絡(luò)處所對(duì)應(yīng)的時(shí)間un當(dāng)作時(shí)移，用其相應(yīng)時(shí)刻的瞬時(shí)頻率和瞬時(shí)相位作為中心頻率ωn和相移φn。地震道經(jīng)Hilbert變換后構(gòu)建的解析信號(hào)為

其振幅包絡(luò)為

2）在相應(yīng)的時(shí)間范圍內(nèi)利用

求取尺度控制參數(shù)σn。其中，D＝｛gr（t）｝γ∈Γ是一系列連續(xù)小波的集合（子波字典）；〈s，g〉表示函數(shù)s和g的內(nèi)積；為子波grn的標(biāo)準(zhǔn)化。在參數(shù)un，ωn，φn確定后，在一定范圍內(nèi)選取σn進(jìn)行掃描，當(dāng)滿足公式（4）最大值條件時(shí)停止，此時(shí)的σn就是最佳值。當(dāng)｛u，σ，ω，φ｝確定后，對(duì)每一個(gè)參數(shù)給定一較小的擾動(dòng)范圍及固定的變化步長(zhǎng)，再次利用公式（4）對(duì)其進(jìn)行優(yōu)化。

3）通過下式來估算最佳小波grn的振幅值：

對(duì)于當(dāng)前迭代分解和正交投影而言，R（n）＝angrn＋R（n＋1），其 中，grn與R（n＋1）是正交函數(shù)。選取使最大，即最小時(shí)的grn。在經(jīng)過全部的N次迭代之后，最終所得的最小殘差R（N）（t）被當(dāng)作噪聲。

由于復(fù)數(shù)Paul子波和復(fù)數(shù)Morlet子波的速度頻散和振幅衰減屬性與地震子波相似，因此，通常選擇這兩種子波來進(jìn)行多子波分解與重構(gòu)［1－3］。復(fù)數(shù)Paul子波和復(fù)數(shù)Morlet子波的表達(dá)式 為［1，8］

式中：m為控制Paul子波寬度的參數(shù)；k為控制Morlet子波寬度的參數(shù)；fm為主頻。

在實(shí)際的匹配追蹤算法中，復(fù)數(shù)Paul子波和復(fù)數(shù)Morlet子波的選擇對(duì)算法本身無任何影響，但在迭代匹配的過程中，選擇不同類型的子波對(duì)最終結(jié)果有一定影響——每次迭代分離出的子波在形狀上稍有差異，因此，利用這些子波重構(gòu)的剖面也有所不同。

2 不同子波分辨率差異性理論分析

為了說明不同子波類型對(duì)地震重構(gòu)信號(hào)分辨率影響的差異性，我們利用公式（6）和公式（7），通過理論計(jì)算來分析復(fù)數(shù)Paul子波和復(fù)數(shù)Morlet子波對(duì)重構(gòu)信號(hào)頻率分辨率的影響。

所謂頻率分辨率是指將兩個(gè)相鄰譜峰分開的能力。對(duì)單個(gè)子波而言，可以用子波譜的寬度來描述。復(fù)數(shù)Paul子波和復(fù)數(shù)Morlet子波的表達(dá)式中都有一個(gè)控制子波寬度的參數(shù)，子波的寬度與分辨率有關(guān)，寬度越小說明頻率分辨率越高。

2.1 復(fù)數(shù)Paul子波

m為控制Paul子波頻率分辨率的參數(shù)，因此，通過改變m值計(jì)算具有不同m值的Paul子波頻譜，利用子波譜的寬度來討論它的頻率域分辨率。圖1分別給出了m值分別為2，4，5，10四種情形下復(fù)數(shù)Paul子波的三維頻譜。為便于子波頻譜形態(tài)和寬度的對(duì)比分析，部分圖的頻率軸從30Hz開始顯示。

從圖1可見，當(dāng)m較小時(shí)（圖1a，m＝2），盡管Paul子波頻譜在頻率域的寬度相對(duì)較小，也就是說頻率域分辨率相對(duì)較高，但因其頻譜寬度很窄，在低頻固定頻率切片上接近一個(gè)脈沖形狀。由公式（6）可知，m太小時(shí)，子波寬度幾乎為零，不能形成完整的振幅譜曲線，因此，其三維頻譜上顯示不完整，有數(shù)據(jù)的缺失。當(dāng)m較大時(shí)（圖1c，圖1d；m＝5，10），子波寬度很大，顯示出較低的頻率分辨率，且子波譜形狀會(huì)發(fā)生畸變。特別是圖1d中所顯示的子波譜形狀，與理想的光滑對(duì)稱子波譜相去甚遠(yuǎn)，存在嚴(yán)重扭曲現(xiàn)象，若利用這種子波進(jìn)行地震道多子波分解與重構(gòu)，必將導(dǎo)致地震信號(hào)的嚴(yán)重失真。

通過m的改變可以看出，隨著m的增加，子波從窄變化到寬，但m太大時(shí)子波形狀會(huì)發(fā)生畸變。因此，在實(shí)際應(yīng)用中，m值的選擇必須適中。通過對(duì)m值的連續(xù)掃描分析，理論上我們認(rèn)為，m＝4是Paul子波頻率分辨率的最佳控制參數(shù)（圖1b）。

2.2 復(fù)數(shù)Morlet子波

k為控制Morlet子波頻率分辨率的參數(shù)，通過改變k值可以計(jì)算不同k值的Morlet子波頻譜（圖2）。

圖2顯示了k值分別為0.3，0.5，0.7和3.0時(shí)，復(fù)數(shù)Morlet子波的三維頻譜。從圖2可以看出，隨著k值增大，子波的寬度越來越小，說明復(fù)數(shù)Morlet子波在頻率域的分辨率越來越高，時(shí)間分辨率越來越低，用復(fù)數(shù)Morlet子波對(duì)地震信號(hào)分析時(shí)既要考慮頻率分辨率也要考慮時(shí)間分辨率。通過對(duì)圖2 進(jìn)行對(duì)比分析，理論上我們認(rèn)為，k＝0.5是Morlet子波頻率分辨率的最佳控制參數(shù)。

3 實(shí)際應(yīng)用效果對(duì)比分析

為了更好地揭示子波選擇對(duì)地震道分解與重構(gòu)的重要性，采用匹配追蹤算法，選擇復(fù)數(shù)Paul子波、復(fù)數(shù)Morlet子波分別對(duì)四川新場(chǎng)工區(qū)實(shí)測(cè)地震資料進(jìn)行了多子波分解重構(gòu)處理和對(duì)比分析。

圖3為不同m值復(fù)數(shù)Paul子波分解所得主頻為45～55 Hz的高頻子波重構(gòu)地震剖面（m＝2，4，5，10）。通過對(duì)比分析可知，當(dāng)m＝2時(shí)，剖面連續(xù)性很好，低頻成分明顯，視覺上與子波主頻不一致。其主要原因是m＝2時(shí)子波的低頻能量缺失，導(dǎo)致最終結(jié)果低頻能量明顯增強(qiáng)。當(dāng)m＝10時(shí)，重構(gòu)剖面主頻似乎很高，但信噪比極差，幾乎看不到有效信號(hào)。m＝5時(shí)與m＝4重構(gòu)信號(hào)有些相似，僅連續(xù)性較差。m＝4時(shí)，子波譜形狀的對(duì)稱性和完整性都很好，并且重構(gòu)剖面也比較符合實(shí)際；m＝5時(shí)子波形狀有一些扭變。通過綜合分析認(rèn)為，m＝4為最佳參數(shù)，這與理論分析的結(jié)果（圖1）一致。

圖4 給出了不同k值復(fù)數(shù)Morlet子波分解所得主頻為45～55 Hz的高頻子波重構(gòu)地震剖面（k＝0.3，0.5，0.7，3.0）。對(duì)比可見，當(dāng)k＝3.0時(shí)，由于重構(gòu)子波頻譜寬度很窄（圖2d），導(dǎo)致重構(gòu)地震信號(hào)的信噪比和連續(xù)性較差，一些強(qiáng)同相軸明顯變?nèi)酰▓D4d）。而當(dāng)k為0.7時(shí)，重構(gòu)信號(hào)雖然在連續(xù)性上要好于k＝3.0，但與k＝0.5（圖4b）相比仍然有信號(hào)失真現(xiàn)象（圖4c）。當(dāng)k＝0.3時(shí)（圖4a），雖然重構(gòu)信號(hào)的連續(xù)性和信噪比都要高于k＝0.5，但考慮到子波主頻為45～55 Hz，有部分低頻信號(hào)的缺失才是合理的，再對(duì)比復(fù)數(shù)Morlet子波頻譜，分析認(rèn)為k＝0.5 為復(fù)數(shù)Morlet子波的最佳分辨率參數(shù)。

為了進(jìn)一步對(duì)比分析基于復(fù)數(shù)Morlet子波（k＝0.5）和復(fù)數(shù)Paul子波（m＝4）分解重構(gòu)結(jié)果的差異性，圖5 和圖6 分別給出了高、低頻復(fù)數(shù)Morlet子波和復(fù)數(shù)Paul子波的重構(gòu)剖面及其頻譜。通過對(duì)比不難看出，復(fù)數(shù)Paul子波重構(gòu)剖面的連續(xù)性和分辨率均優(yōu)于復(fù)數(shù)Morlet子波重構(gòu)剖面。低主頻情形下復(fù)數(shù)Paul子波重構(gòu)剖面包含了更多的高頻信息（圖6a）；高主頻情形下復(fù)數(shù)Paul子波重構(gòu)剖面具有相對(duì)較好的抗噪性。究其原因，可能是由于復(fù)數(shù)Paul子波重構(gòu)對(duì)高頻噪聲的放大作用較弱所致（圖6b）。

4 結(jié)論

通過理論和實(shí)際應(yīng)用效果對(duì)比分析，得到如下認(rèn)識(shí)：①m＝4是Paul子波最佳頻率域分辨率控制參數(shù)，k＝0.5為Morlet子波最佳頻率域分辨率控制參數(shù)；②復(fù)數(shù)Paul子波能使重構(gòu)剖面保留更多的低頻能量，同相軸更連續(xù)，構(gòu)造更清晰，而復(fù)數(shù)Morlet子波高低頻頻譜差異大，重構(gòu)剖面同相軸連續(xù)性相對(duì)較差。因此，基于復(fù)數(shù)Paul子波的重構(gòu)剖面適合反射同相軸連續(xù)的地下構(gòu)造研究，基于復(fù)數(shù)Morlet子波的重構(gòu)剖面適合地震波吸收衰減研究。

［1］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

［2］趙爽，李仲東，許紅梅，等.多子波分解技術(shù)檢測(cè)含煤砂巖儲(chǔ)層［J］.天然氣工業(yè)，2007，27（9）：44－47 Zhao S，Li Z D，Xu H M，et al.Using multiple wavelet decomposition technique to detect the sandstone reservoir with coal layer［J］.Natural Gas Industry，2007，27（9）：44－47

［3］王昕，魏剛，劉歆，等.多子波分解與重構(gòu)技術(shù)在復(fù)雜儲(chǔ)層預(yù)測(cè)中的應(yīng)用［J］.斷塊油氣田，2010，17（3）：289－292 Wang X，Wen G，Liu X，et al.Application of seismic trace－wavelet decomposition and reconfiguration technology in complex reservoir prediction［J］.Fault－Block Oil ＆Gas Field，2010，17（3）：289－292

［4］Partyka G，Gridley J，Lopez J.Interpretational applications of spectral decomposition in reservoir characterization［J］.The Leading Edge，1999，18（3）：353－360

［5］Liu J，Wu Y F，Han D H，et al.Time－frequency decomposition based on Ricker wavelet［J］.Expanded Abstracts of 74thAnnual Internat SEG Mtg，2004，1937－1940

［6］Liu J，Marfure K J.Matching pursuit decomposition using Morlet wavelets［J］.Expanded Abstracts of 75thAnnual Internat SEG Mtg，2005，786－789

［7］Wang Y.Multichannel matching pursuit for seismic trace decomposition［J］.Geophysics，2010，75（4）：61－66

［8］Wang Y.Seismic time－frequency spectral decomposition by matching pursuit［J］.Geophysics，2007，72（1）：V13－V20

［9］徐天吉，沈忠民，文雪康.多子波分解與重構(gòu)技術(shù)應(yīng)用研究［J］.成都理工 大 學(xué)學(xué)報(bào)，2010，37（6）：660－666 Xu T J，Shen Z M，Wen X K.Research and application of multi－wavelet decomposition and reconstructing technology［J］.Journal of Chengdu University of Technology，2010，37（6）：660－666

［10］陳發(fā)宇，楊長(zhǎng)春.對(duì)依據(jù)頻率的匹配追蹤快速算法的改進(jìn)［J］.石油物探，2009，48（1）：80－83 Chen F Y，Yang C C.An improvement for frequency－dominated fast matching pursuits algorithm［J］.Geophysical Prospecting for Petroleum，2009，48（1）：80－83

［11］陳發(fā)宇，楊長(zhǎng)春.基于MP方法的地震信號(hào)快速分解算法［J］.地球物理學(xué)進(jìn)展，2007，22（6）：1692－1697 Chen F Y，Yang C C.Seismic signal＇s decomposition based on matching pursuits method［J］.Progress in Geophysics，2007，22（6）：1692－1697

［12］馬朋善，王繼強(qiáng)，劉來祥，等.Morlet小波分頻處理在提高地震資料分辨率中的應(yīng)用［J］.石油物探，2007，46（3）：283－287 Ma P S，Wang J Q，Liu L X，et al.The application of frequency decomposition technique of morlet wavelet in enhancing seismic data resolution［J］.Geophysical Prospecting for Petroleum，2007，46（3）：283－287

［13］高靜懷，吳茜，陳文超，等.小波變換域地震資料瞬時(shí)頻率分析方法［J］.石油物探，2007，46（6）：534－540 Gao J K，Wu Q，Cheng W C，et al.Instantaneous frequency analysis of seismic data in wavelet transform domain［J］.Geophysical Prospecting for Petroleum，2007，46（6）：534－540

［14］鄒問，陳愛萍，賀振華，等.基于S變換的地震相分析技術(shù)［J］.石油物探，2006，45（1）：48－51 Zhou W，Chen A P，He Z H，et al.Seismic faces analysis based on S transform［J］.Geophysical Prospecting for Petroleum，2006，45（1）：48－51