預測反褶積在壓制周期性多次反射波中的應用研究

唐海敏 苗慶庫 王 銀

(鐵道第三勘察設計院集團有限公司，天津 300251)

反褶積也叫反濾波，在反射地震資料處理中已經作為最常用的技術之一。它是通過壓縮基本地震子波以提高地震資料時間分辨率的過程。反褶積一般用于疊前，也廣泛用于疊后資料。反褶積產生更高時間分辨率的剖面，沒有反褶積的疊加剖面上的鳴震性質，即地震波經多次反射形成的干擾，顯著地限制了分辨率。預測褶積的理論用于解決反褶積的問題叫做預測反褶積，它既可以用于壓制海上資料虛反射之類的多次波，又可用于陸上壓制地震信號的層間多次波，現在已經是一種常規的處理步驟。

1969年4月號地球物理雜志發表的Peacock和Treitel的文章中［3］描述了擬定削弱多次波(特別是海上常常遇到的鳴震)的長周期預測反褶積算子的方法。他們描述的預測反褶積算子設計方法致力于保留一次波，而盡力削弱多次波，其程序中的算子是針對削弱一次反射波后一定時間間隔有規則重復出現的周期性多次波而設計的。

本文通過對包含多次反射波的合成地震記錄和Tesseral建模地震記錄做多道預測反褶積，并嘗試使用不同的預測參數(預測步長和預測算子長度)，通過反復對比，取得了較為理想的壓制多次波效果。

1 預測反褶積原理

所謂反褶積是一個濾波過程，在地震勘探中這個濾波過程的作用恰好與大地濾波的作用相反，如圖1所示。

圖1 大地濾波和反褶積的流程對比

地震勘探的基本原理是向地下輸入一個短信號，該信號可從地下兩個地層單元的界面上反射回來。在此類模型中，脈沖函數并不是簡單的速度差而是波阻抗差。由此得出眾所周知的反射系數基本公式

式中，ρ是密度，v是波速。

地震道X(t)是由有效波s(t)和干擾波n(t)疊加組成的，具體公式如下

該有效波是一次反射波，對反射地震勘探而言，除一次反射波以外的波一般可視為干擾波。層狀介質的一次反射波通常用線性褶積模型表示。

褶積模型公式

式中，w(t)為系統子波，r(t)為反射系數函數，n(t)為干擾波。

一道合成地震記錄是震源子波w(t)與一個反射率模型r(t)相褶積的結果。

預測反褶積是用預測的方法，根據地震記錄一次反射和干擾的信息預測出純干擾部分，再從包括一次波和干擾的地震記錄中減去純干擾部分，以壓制一次反射后面的多次波干擾，得到壓制后的一次反射信號。在預測反濾波問題中，設計一個預測因子c(t)，對輸入地震記錄x(t)的過去值x(t－m)，x(t－m+1)，…，x(t－1)和現在值x(t)進行預測，所得到的未來的預測值x(t+a)是海上鳴震或陸地多次波干擾，將它從包括一次反射和干擾的地震記錄x(t+a)中減去，所得到的預測誤差ε(t+a)=x(t+a)－x(t+a)就是壓制干擾后的一次反射信號。

1.1 多道多步預測誤差濾波器

設有K道地震數據，每道時間采樣為M+1，預測步長為L，多步預測濾波器系數個數為N，

則多步向前預測誤差定義為

CN，j是K道多步向前預測誤差濾波器的K*K階(1≤K≤K)矩陣系數。

可以整理為矩陣表達形式

I是K階的特征矩陣。

將預測誤差的平方和最小化得到擴展的多步向前預測誤差濾波常態方程

R0是 K*K 階，RL+i－1，i=1，…，N 是 K*K 階，R－L－i+1

PN是最小向前多步預測誤差協方差矩陣

1.2 預測濾波器的算子長度和預測步長

維納推導出將輸入轉換為任意期望輸出的濾波器。

一個長度為n的矩陣方程的普遍形式是［8］

式中，ri，ai，及 gi(i=0，1，2，3，…，n －1)分別是輸入子波自相關，期望濾波系數及期望輸出與輸入子波的自相關。

最佳維納濾波器ai是最佳的，是指它的實際輸出與期望輸出之間的最小平方誤差最小。當期望輸出是零延遲尖脈沖(1，0，0…0)時，維納濾波就與最小平方濾波等同。換句話說，后者是前者的特例。

給定輸入X(t)，我們需要預測它在某一將來(t+α)的值，這里α是預測步長。維納指出，用于估計X(t+α)的濾波器可用如方程(8)的一個矩陣方程的特殊形式來計算［8］。

方程(8)可一般化為一個長度為n的預測濾波器及一個預測步長為α的情況

現在考慮單位預測步長，α=1，n=6的特殊情況，則方程(9)可改寫如下

式中，L=r0－r1a0－r2a1－r3a2－r4a3－r5a4

尖脈沖反褶積實際上是單位預測步長預測濾波器的特殊情況。

合成資料的例子展示了多道預測反褶積在削弱在資料中周期性出現的多次反射的有效性。

2 合成地震資料的預測反褶積

設計一個合成地震記錄，該記錄有10道，每道有1000個采樣點，采樣時間1 s，共產生三次多重反射波，分別在400 ms，600 ms，800 ms，多次反射波的周期為200 ms，設反射系數c=－0.333，加入隨機噪聲為+/－0.001，

合成地震記錄如圖2所示。

圖2 合成地震記錄

反褶積后如圖3所示。

使用的參數為:反褶積一道所用的道數=5，預測步長=200 ms，多重反射周期=200 ms。

圖3 反褶積后的合成地震記錄

3 Tesseral模型設計及預測反褶積對比研究

新建模型如圖4所示，有兩個水平層，上層為海水，下層為某巖層，橫向寬度為1000m，縱向深度為1000m，炮點放在中間500m處，檢波器范圍從600m到1000m，道間距取4m，炮點距離最近的檢波器長度為100m，地震信號記錄時間從0 s開始到2 s結束，采樣間隔為2 ms，每道有1000個采樣點。

子波設為70 Hz最小相位子波，經計算該子波的波長約為21 ms，試驗中預測步長L不能小于21 ms，否則將造成子波的畸變。

第一層:縱波波速1 500m/s，橫波波速0.25m/s，密度1 g/cm3;

第二層:縱波波速2 500m/s，橫波波速1 450m/s，密度2.2 g/cm3。

根據反射系數計算公式，可得反射系數c=(2 500×2.2－1 500×1)/(2 500×2.2+1 500×1)=0.571 4。

在波阻抗相差較大的情況下，比較容易產生多次反射波。

圖4 多次波形成的地質模型

以下是一個多次反射模型的示意，O為震源，R，M1，M2，M3為檢波器位置，如圖5所示。

圖5 全程多次反射波形成示意

選擇彈性波模擬，生成的原始含有多次反射波的地震單炮記錄，如圖6所示。

圖6 含多次反射波的地震單炮記錄

在程序中涉及到三個參數:

①反褶積一道所使用的相鄰道數nc;

②預測算子長度ncf;

③預測步長L。

多次波周期可以近似為(200×2)/1 500=267 ms，歸一化后相當于P=135。先固定兩個參數nc=5，L=135，取預測算子長度 ncf分別為 30 ms，50 ms，80 ms，100 ms，150 ms，200 ms，可對比得出 150 ms 和 200 ms的預測算子長度后的結果都很好。為了提高計算效率，選定ncf=150 ms為好。

固定nc=5，預測算子長度選擇ncf=150 ms，因為再增大算子長度，反褶積的效果保持穩定不變，而大的算子長度會占用比較多的計算時間，影響運行效率。綜合考慮算子長度選擇150 ms，取預測步長分別為L=8 ms，50 ms，100 ms，110 ms，120 ms，130 ms，160 ms，200 ms，280 ms。

步長的選擇出現兩端差中間好的現象，即預測步長太小會出現很多噪聲和不規則的擾動，步長太大則壓制多次波的效果不好。經比較，預測步長選擇120比較適宜。

單獨考慮nc的選擇，即反褶積一道所用的相鄰道數，固定好其他兩個參數，預測算子長度ncf=150，預測步長固定為L=120，取nc=5和10做對比，二者幾乎沒有任何效果上的差別，如果選nc=10，計算時間過長，效率低，綜合考慮應該選取nc=5。

綜上可知，對于本次試驗，多道預測反褶積壓制周期性地震多次波干擾的三個參數最佳值應是:

反褶積一道選取相鄰道數nc=5;預測算子長度ncf=150 ms;預測步長L=120 ms

4 結論

關于預測反褶積的參數選擇，預測步長以取最遠道的自相關極值為佳，并通過試驗確定最佳預測步長。預測算子長度應該大于多次波的周期長度，如果算子長度取小了，就會使預測結果帶進假的能量，使預測反褶積波形尾部出現波動，但算子長度也不能太大，長度太大不僅增加計算量，而且會增加相鄰反射及噪聲的影響。預測反褶積程序的執行效率較高，宜通過對比試驗綜合確定最佳參數。

［1］［美］渥．伊爾馬滋．地震數據處理［M］．黃緒德，譯．北京:石油工業出版社，1994:71-135．

［2］大港油田科技叢書編委會．地震勘探資料處理和解釋技術［M］．北京:石油工業出版社，1999

［3］Peacock K L， Treitel S．Predictive deconvolution Theory and practices．Geophysics，1969，34(2):155-169

［4］陸邦干，范偉粹．地球物理資料數字處理［M］．北京:石油工業出版社，1987:80-102

［5］黃緒德．反褶積與地震道反演［M］．北京:石油工業出版社，1992:1-10，45-69

［6］N．S．奈德爾，等．褶積模型［M］．牛毓荃，譯．北京:石油工業出版社，1986．:51-54，126-165

［7］牟永光．地震數據處理方法［M］．北京:石油工業出版社，2007:80-90

［8］Robinson，E．A．，and Treitel，S．．Geophysical signal and analysis．Prentice-Hall，Inc．1980