常 青， 徐秀剛， 孫道朋

(海底科學(xué)與探測技術(shù)教育部重點(diǎn)實(shí)驗(yàn)室，中國海洋大學(xué)海洋地球科學(xué)學(xué)院，山東 青島 266100)

多次波是一種規(guī)則干擾波，在海洋地震勘探的資料處理中普遍存在多次波的問題[1]。多次波處理不當(dāng)會混淆、干擾有效波成像，出現(xiàn)多次波同相軸掩蓋有效波同相軸的情況，在很大程度上干擾地震資料的偏移成像效果[2]，甚至?xí)?dǎo)致對地震資料解釋和地下構(gòu)造的錯誤認(rèn)識，影響著探井井位的部署和勘探的成功率。

多次波壓制方法分為基于有效波和多次波的特征差異的濾波方法和基于波動方程的預(yù)測減去法[3]。當(dāng)有效波和多次波時差較大時，濾波方法能夠很好地壓制多次波，但是當(dāng)有效波和多次波的特征差異較小時，該方法則會損傷有效波?；诓▌臃匠痰念A(yù)測減去法卻能很好的處理這種情況，該方法是數(shù)據(jù)驅(qū)動的，在不需要速度模型等先驗(yàn)信息下，就能夠處理地下介質(zhì)較復(fù)雜時的多次波與有效波的分離[4-5]。

最初，Anstey等[6]利用單道的自相關(guān)壓制自由表面多次波。Kennett[7]提出一維表面多次波的正演和反演模擬方法，但是由于該方法對于數(shù)據(jù)需要做許多假設(shè)，不能成功應(yīng)用于實(shí)際資料的多次波壓制。后來，荷蘭Delft大學(xué)的Berkhout[8]提出反饋模型和數(shù)據(jù)矩陣的概念，奠定了反饋迭代法預(yù)測多次波的基礎(chǔ)。由于該方法是是利用疊前地震數(shù)據(jù)預(yù)測自由表面多次波，無需地下介質(zhì)的任何信息，受到眾多地球物理學(xué)家們的關(guān)注。

常規(guī)的多次波壓制技術(shù)是基于輸出信號能量最小準(zhǔn)則進(jìn)行匹配相減的[9]，這一準(zhǔn)則不能很好的處理有效波和多次波非正交的情況下的多次波壓制，而實(shí)際地震數(shù)據(jù)基本不滿足多次波和有效波是正交的(兩個波形不重疊且其相位差為90°)。近年來，許多學(xué)者將ICA(Independent Component Analysis，獨(dú)立分量分析法)算法運(yùn)用到地球物理中，實(shí)現(xiàn)對于非正交數(shù)據(jù)的處理。劉喜武等[10]提出利用ICA方法實(shí)現(xiàn)地震轉(zhuǎn)換波與多次反射波的分離，陸文凱等將ICA方法用于多次波的壓制中并提出“幾何獨(dú)立分量分析方法”[11-12]，Kaplan等[13]也將ICA思想應(yīng)用于多次波的分離中并取得理想的壓制結(jié)果。

為避免壓制效果受速度模型、地層信息等先驗(yàn)信息的影響，本文借鑒了預(yù)測減去法的思路，采用基于波動方程的反饋迭代法實(shí)現(xiàn)對于自由表面多次波的預(yù)測，然后基于高階統(tǒng)計(jì)量的快速ICA算法實(shí)現(xiàn)多次波與有效波的匹配相減，完成最終的自由表面多次波壓制。本文壓制方法完全是數(shù)據(jù)驅(qū)動的，而且克服了常規(guī)匹配算法對于有效波和多次波非正交情況不能很好壓制多次波的缺陷，實(shí)現(xiàn)了對于自由表面多次波的良好壓制。資料試算的結(jié)果，驗(yàn)證了本文方法壓制多次波的正確性，有效性和實(shí)用性。

1 基于反饋迭代和獨(dú)立分量分析的多次波壓制

1.1 反饋迭代預(yù)測多次波的原理

基于波動方程的反饋迭代預(yù)測多次波方法不需要地下地層的介質(zhì)信息，可以直接由地震數(shù)據(jù)本身預(yù)測自由表面多次波。地震數(shù)據(jù)本身被用作波場延拓算子，通過對其進(jìn)行多維褶積運(yùn)算，有效波變成一階多次波，一階多次波變成更高階多次波，該方法通過數(shù)據(jù)驅(qū)動，可以直接利用疊前地震數(shù)據(jù)構(gòu)建表面多次波，無需地下介質(zhì)的信息。

在多次波預(yù)測中，引入了波場延拓的思想，波場延拓過程中，格林函數(shù)表述了延拓路徑。在表面多次波預(yù)測中，一次波的脈沖響應(yīng)函數(shù)X0將作為格林函數(shù)預(yù)測所有可能的多次波。在二維多次波預(yù)測中，采用Rayleigh積分的非穩(wěn)態(tài)空間褶積表示形式，即表面多次波M0(xr,xs,f)的預(yù)測寫成公式(1)所示的空間褶積的積分形式：

M0(xr,xs,f)=

(1)

式中：xs為炮點(diǎn)位置；xr為檢波點(diǎn)位置；xk為求和所進(jìn)行的橫向坐標(biāo)；X0(xr,xk,f)為無表面多次波的脈沖響應(yīng)。

寫成離散求和形式如下：

M0(xr,xs,f)=

(2)

如果對所有的炮點(diǎn)和檢波點(diǎn)組合重復(fù)該過程，則多次波的預(yù)測可表示成矩陣乘法形式：

M0=-X0P。

(3)

M0(xr,xs,f)=

(4)

寫成矩陣乘法形式：

(5)

多次波迭代預(yù)測多次波方法可表述為：

(6)

(7)

1.2 基于獨(dú)立分量分析方法的多次波剔除

在多次波預(yù)測減去法中，多次波相減是最為關(guān)鍵的一環(huán)，利用不同的數(shù)學(xué)方法來獲得自適應(yīng)參數(shù)，使得預(yù)測的多次波在匹配處理后與實(shí)際地震數(shù)據(jù)中的多次波盡可能的一致，然后再從實(shí)際的地震數(shù)據(jù)中將多次波減去，從而得到有效波的信息。

獨(dú)立分量分析方法是近二十年發(fā)展起來的一種高效盲信號分離方法。它的基本思想是從一組觀測到的混合信號中恢復(fù)無法直接觀測到的獨(dú)立原始信號，其用一些基函數(shù)來表示一系列的隨機(jī)變量并且假設(shè)各成分之間是統(tǒng)計(jì)獨(dú)立的或者盡可能的獨(dú)立。其最簡單的形式是：給定N個觀測信號x1,x2,…,xn,假設(shè)它們?yōu)閚個未知的統(tǒng)計(jì)獨(dú)立的源信號s1,s2,…,sn的線性組合：

xj=aj1s1+aj2s2+…+ajnsn,j=1,2,…,n。

(8)

式中aj1,aj2,…,ajn為輸入觀測量的線性組合系數(shù)。

寫成矩陣形式如下：

X=AS。

(9)

設(shè)觀測到的地震信號為X，由一次波P和多次波M組成，則有：

X=P+M。

(10)

設(shè)通過多次波預(yù)測技術(shù)得預(yù)測的多次波為Mp，為了消除Mp和實(shí)際地震數(shù)據(jù)中的多次波M之間的不同子波因素，采用了基于最小能量的子波估計(jì)方法[14]，用如下的子波校正濾波器進(jìn)行子波的校正：

(11)

式中：W(ω)是一個一維窗函數(shù)；上標(biāo)*表示復(fù)數(shù)共扼；?表示褶積運(yùn)算。

利用該子波校正濾波器，得到子波校正后的預(yù)測多次波：

(12)

設(shè)子波因素消除后的近似多次波為：

Mp=αM。

(13)

式中：α為尺度系數(shù)。

壓制地震信號X中的多次波所用的常規(guī)的基于二階統(tǒng)計(jì)量能量函數(shù)的方法是最小化如下的函數(shù)：

(14)

上式對應(yīng)的最小二乘解為：

(15)

考慮到峰度對信號中的大值的敏感性較高，采用了如下的對數(shù)范數(shù)度量[9]：

圖1 傾斜模型速度場Fig.1 Velocity field of the tilt model

(16)

上式中：N是信號的采樣長度；qi是對應(yīng)信號y的振幅歸一化信號，即：

(17)

ICA就是通過最大化高階統(tǒng)計(jì)量的能量函數(shù)：

(18)

來實(shí)現(xiàn)多次波的自適應(yīng)相減。

2 數(shù)值模型試算

為驗(yàn)證本文方法的多次波壓制效果，本文通過傾斜層模型以及Sigsbee2B模型分析等實(shí)現(xiàn)了對于該方法技術(shù)的試算，通過資料試算，驗(yàn)證了方法的正確性、適應(yīng)性和有效性。

本文采用的第一個模型是傾斜模型，大小為2 000 m×1 200 m，水平界面深度300 m，傾斜界面左端點(diǎn)深度960 m，右側(cè)430 m，如圖1所示，表層速度是漸變的，第一層為1 500 m/s，第二層為2 000 m/s，第三層為2 500 m/s。炮記錄采用左邊放炮的激發(fā)方式，共計(jì)96炮，每炮96道接收，每道1 000個采樣點(diǎn)，采樣率為2 ms，子波主頻40 Hz。采樣反饋迭代法和ICA匹配相結(jié)合的方法對該炮記錄進(jìn)行試算，第1炮單炮記錄多次波壓制結(jié)果如圖2所示。

圖2(a)為原始地震記錄(箭頭所指部分為多次波)，圖2(b)是利用基于波動方程的反饋迭代方法預(yù)測的多次波。圖2(c)是利用基于二階統(tǒng)計(jì)量的常規(guī)匹配方法壓制多次波后的地震記錄，可以看到在有效波和多次波非正交部分，多次波仍然有部分殘留；圖2(d)是利用基于高階統(tǒng)計(jì)量的ICA方法壓制多次波后的地震記錄，對比圖2(c)可知，利用基于高階統(tǒng)計(jì)量的ICA方法對于有效波和多次波非正交部分壓制效果要明顯優(yōu)于二階統(tǒng)計(jì)量的常規(guī)方法。

圖2 原始地震記錄與不同方法下的多次波壓制結(jié)果對比Fig.2 Comparison between original record and suppressed result of multiples based on different ways

為了測試方法對于復(fù)雜構(gòu)造的適應(yīng)性，本文還采用國際上常用的檢測多次波壓制效果的Sigsbee2B模型對算法進(jìn)行了測試。該模型的參數(shù)為：縱向共有1 201個采樣點(diǎn)，采樣間隔為7.62 m；縱向共有2 133個采樣點(diǎn)，采樣間隔為11.43 m。該模型的正演炮記錄共有500炮，為了計(jì)算起來更加方便，本文僅抽取其中的348炮、每炮174道、每道1 500個采樣點(diǎn)、8 ms采樣進(jìn)行多次波壓制的實(shí)驗(yàn)。Sigsbee2B模型速度場如圖3所示，該模型表層模擬海水層，中間是由許多正斷層和逆沖斷層組成的沉積地層以及埋藏在沉積地層中的一個高速度的巖丘體，可以看出該地質(zhì)構(gòu)造非常復(fù)雜。

該模型的正演炮記錄總共有12 s，為了使得顯示效果清晰便于對比，本文選取其中的第51炮并截取其中的3.2～12 s進(jìn)行顯示對比。

圖3 Sigsbee2B模型速度場Fig.3 Velocity field of the Sigsbee2B model

圖4(a)中箭頭所指部分為多次波，圖4(b)為利用基于波動方程的反饋迭代方法預(yù)測的多次波，比較兩圖可以看出，原始數(shù)據(jù)中的多次波被很好地預(yù)測出來。圖5(a)為利用基于二階統(tǒng)計(jì)量的常規(guī)匹配方法壓制多次波后的地震記錄，圖5(b)為利用基于高階統(tǒng)計(jì)量的獨(dú)立分量分析自適應(yīng)匹配方法法壓制多次波后的地震記錄。對比兩圖中方框部分可以看出，對于有效波和多次波非正交部分以及深層多次波，后者的壓制效果要明顯好于前者。

圖4 原始單炮記錄和反饋迭代法預(yù)測的多次波Fig.4 Original record and predicted multiples based on feedback iteration method

圖5 基于最小二乘和ICA方法匹配后的一次波記錄Fig.5 Matched primary wave record based on least-square method and ICA method

圖6 多次波壓制前后的偏移距剖面Fig.6 The offset section before and after multiples being suppressed

3 結(jié)論與認(rèn)識

本文實(shí)現(xiàn)了基于反饋迭代技術(shù)和基于高階統(tǒng)計(jì)量的獨(dú)立分量分析技術(shù)的結(jié)合，更好的對于自由表面多次波進(jìn)行了壓制，取得的結(jié)論與認(rèn)識是：

(1)反饋迭代法是基于數(shù)據(jù)驅(qū)動的多次波壓制方法，可以直接利用疊前地震數(shù)據(jù)預(yù)測表面多次波。由于該方法不需要地層的信息，因此可對復(fù)雜地下的表面多次波進(jìn)行預(yù)測。

(2)將獨(dú)立分量分析技術(shù)應(yīng)用于多次波的匹配壓制，克服了常規(guī)匹配算法對于有效波和多次波非正交情況不能很好壓制多次波的缺陷，實(shí)現(xiàn)了對于自由表面多次波的有效壓制。傾斜模型和Sigsbee2B模型等模型資料的處理分析表明，本文方法既能夠有效保持地震記錄中的有效波能量，還能夠?qū)崿F(xiàn)對于自由表面多次波，特別是中、深層多次波的有效壓制。