利用區域自適應極化濾波壓制多分量面波干擾

王志農 孫成禹* 伍敦仕 閆月鋒

(①中國石油大學(華東)地球科學與技術學院，山東青島 266580； ②中國石油勘探開發研究院西北分院，甘肅蘭州 730020)

0 引言

在反射波法地震勘探中，面波是主要的規則干擾波，具有低速、低頻、強振幅和橢圓極化等特征。隨著多分量地震勘探技術的發展，現場數據采集過程中通常會采用單點數字檢波器接收方式，舍棄了傳統檢波器組合接收方式[1]。單點接收雖然會帶來一些優勢，但采集到的地震數據往往信噪比低、面波干擾嚴重。當采用寬方位或全方位觀測系統時，面波和部分線性干擾在遠離炮點的排列上，其視速度與轉換波的傳播速度接近；另外，轉換波與面波在頻率上有較大重疊范圍，常規去噪方法難以在不損害有效波的情況下去除面波干擾[2]。陸上三分量地震記錄中的面波干擾往往會掩蓋較弱的反射信號，給后續處理與反演解釋帶來困難[3]。

在地震數據處理中，傳統的面波壓制方法主要有高通濾波、FK濾波、Radon變換等[4-5]。區域濾波是在面波分布區內對地震記錄進行高通濾波以達到壓制面波的目的。此類方法簡單實用，只對面波區域內信號進行處理，對區域外的資料無影響[6]。但由于其本質仍是高通濾波，因此難免會對低頻信號造成一定程度的損害。1965年Flinn[7]最早將極化濾波應用于天然地震中不同波型的分離，提出了協方差矩陣的統計方法。隨著多分量地震勘探技術的發展，極化濾波也開始應用于地震記錄的去噪。Perelberg等[8]對Flinn的經典極化濾波法做了改進，并引入基于橢圓率和方向性的權重函數。

Deighan等[9]率先將一維小波變換應用于面波壓制，隨后以小波變換為代表的時頻分析法逐漸被廣泛地應用于地震數據處理[10-14]。Meersman等[15]利用SVD極化濾波成功地壓制了低頻地震數據中的面波； Diallo等[16]利用極化濾波在小波域實現了三分量地震數據波型分離； Chen等[17]提出一種自適應三分量地震記錄極化濾波法，取得了較好去噪效果； 馬見青等[18]構建了時頻域自適應協方差矩陣，提出一種時頻域極化濾波壓制面波的方法。

在利用極化濾波壓制面波時，極化分析的時窗選擇對面波壓制效果有很大影響，使用固定/單一時窗會帶來很大局限性。實際應用中，時窗大小若能隨面波主頻對應的波長變化，則會在很大程度上改善極化濾波效果。

本文針對三分量地震數據特點，綜合考慮面波低頻、低速和橢圓極化的特征，提出一種基于區域自適應極化濾波的面波壓制方法： 首先利用S變換進行區域濾波，但去除的面波干擾中不可避免地還殘留一部分低頻有效信號； 再通過HHT計算面波干擾的瞬時頻率從而設計自適應極化濾波窗口，使窗口大小可隨面波波長的變化而改變； 最后，利用自適應極化濾波所得低頻有效信號對區域濾波結果進行補償，得到面波壓制后的最終地震記錄。模型測試和實際數據處理結果表明，本文方法能在有效壓制三分量面波干擾的同時，避免有效信號遭受損害。

1 方法原理

1.1 基于S變換的區域濾波

區域濾波是一種針對目標區域的帶通濾波。在面波壓制的過程中，利用面波的低頻特性對面波區域進行帶通濾波。它與常規帶通濾波的區別在于：帶通濾波是針對整個數據體，而區域濾波則需根據面波的視速度設計合適的濾波窗口，從而減少對區域外信號的影響，最大限度地保留有效信息[19]。由于面波速度較低(通常為100～1000m/s)，在地震記錄中往往呈雙曲線型或扇型分布，因此宜選擇雙曲線型或線型濾波窗口。面波線型濾波窗口定義為

(1)

式中：t(x)表示炮檢距為x的地震道窗口起始時間；v為面波視速度；t0為零炮檢距處的窗口起始時間。在實際地震數據處理中，由于面波存在頻散特性，其視速度會隨著炮檢距的增大而變化。為了得到更好的處理效果，可采用分段線型窗口(圖1)。

S變換可看作是小波變換的一種延伸，其主要優點是分辨率可隨頻率的變化而改變，且在時頻表示中各頻率分量的相位譜與原始信號保持直接聯系[20]。由于S變換具有可變的時窗寬度，能有效克服常規傅里葉變換的時頻域分辨率差的缺點。因此，選用S變換做區域濾波，并在時頻域設計濾波器，能提高區域濾波的精度。圖2b為利用S變換對圖1地震記錄中第400道數據做時頻分析的結果，其中紅色框標示面波干擾的主要分布區域。

圖1 線型窗口示意圖

1.2 自適應極化濾波基本原理

面波除了具有低速、低頻特征外，另一重要特征是橢圓極化。在傳播過程中，體波和面波質點的極化方式不同，可利用其不同的極化特征，采用基于極化分析技術的濾波方法壓制面波干擾。在陸上三分量地震勘探中，采用單點數字檢波器能記錄到完整的波場矢量信息，這為利用極化濾波法壓制面波創造了有利條件。在利用極化濾波壓制面波過程中，極化分析的時窗大小是固定的，這種固定時窗在濾波時具有很大局限性。本文利用HHT計算地震記錄中面波的瞬時頻率，并根據面波主頻求取對應波長，從而設計自適應極化濾波時窗。

HHT主要分為以下兩步： 首先對信號做經驗模態分解(EMD)，得到若干固有模態函數(IMF)，并假設這些固有模態函數均為窄帶信號； 再對這些信號進行Hilbert變換，從而得到具有明確物理意義的瞬時參數。固有模態函數需同時滿足兩個條件： ①該信號極值點的數目相等或至多相差一個； ②在任何時間點上，由該信號的局部極大值點形成的包絡線和由局部極小值點形成的包絡線的平均值為零。

假設單道地震信號為x(t)，找出x(t)的所有極值點，將所有的局部極大值點和極小值點利用三次

圖2 S變換時頻域中的面波干擾(a)圖1中第400道地震數據波形； (b)S變換時頻譜

樣條函數進行擬合，得到上、下包絡線u(t)和v(t)，計算上、下包絡的平均值得到原信號的平均包絡線

(2)

將原信號減去平均包絡，可得不含低頻成分的新信號

h(t)=x(t)-m(t)

(3)

若h(t)滿足IMF條件，則為第一固有模態函數； 否則，循環k次，使h1k(t)滿足IMF條件，記為

c1=h1k(t)=h1(k-1)(t)-m1k(t)

(4)

從x(t)中分離出c1，得到去高頻成分殘余信號

r1=x(t)-c1

(5)

將r1作為原始信號重復以上步驟，得到第2，3，…，n個IMF分量，直至rn成為一個單調函數，不能再從中提取IMF分量時，循環停止； 此時，原信號可表示為

(6)

則原始地震信號x(t)被分解為n個IMF分量和一個殘余量rn(t)。對IMF分量做Hilbert變換，得到具有物理意義的瞬時參數

(7)

(8)

則瞬時振幅、瞬時相位及瞬時頻率可表示為

(9)

自適應極化濾波的時窗Tn(n為樣點數)定義為

(10)

M=[x,y,z]

(11)

式中：x,y,z分別為列向量，若每列都含有n個采樣點，則M就是一個n×3矩陣。令協方差矩陣Mc=MTM，“T”表示矩陣的轉置，那么

(12)

由于MT是3×n的矩陣，則Mc為一個3×3的矩陣。極化分析通過對Mc進行分解，得到其特征值和特征向量，分別用于判別能量大小及傳播方向。奇異值分解可用下式表示

M=UMVT

(13)

式中：V表示Mc的特征向量，為一個3×3的矩陣；W是一個對角矩陣，其對角元素是矩陣M的奇異值，每一個奇異值是Mc對應特征值的正平方根。濾波器的設計主要是依據對橢圓率及傳播方向主軸的估計，因此分別定義橢圓權重函數G1(t)和方向權重函數G2(t)為

(14)

式中：ε是反橢圓率； Δθ(t)是期望入射角與實際入射角的夾角。圖3為G1(t)、G2(t)示意圖。

圖3 權重函數示意圖(a)橢圓方向(0

引入高斯權重函數代替上式中的線性和余弦函數，可將式(14)改寫為

(15)

質點的振動方向是協方差矩陣主特征值所對應的特征向量的方向。濾波過程可表示為

R(t)=G1(t)G2(t)S(t)

(16)

式中：S(t)表示輸入的含有面波的地震記錄；R(t)表示自適應極化濾波的結果。

1.3 基于區域自適應極化濾波的面波壓制方法

由前文的分析可知，基于S變換的區域濾波考慮了面波的低速、低頻特征，但是會損害面波區域內與面波頻帶相同的低頻有效信號。極化濾波考慮到面波與體波極化特征的不同， 可據此將體波和面波進行分離。兩種方法針對面波的不同特征，在面波壓制的過程中具有一定的互補效果，因此本文提出了一種基于區域自適應極化濾波的面波壓制方法，該方法通過自適應極化濾波提取的低頻有效信號對區域濾波結果進行低頻補償，能夠取得較好的面波壓制效果。該方法的流程如圖4所示。

圖4 三分量面波干擾壓制方法流程圖

2 模型測試

2.1 一維合成信號測試

為了測試自適應極化濾波法的有效性，利用正弦掃描信號和指數函數模擬得到面波干擾(圖5a)，與通過褶積合成的有效信號(圖5b)疊加，就得到圖5c所示的合成信號。利用HHT求出該信號的瞬時頻率(圖5d，綠色線段表示自適應時窗大小)，可見低頻時取較大時窗，而在高頻時取較小時窗。

圖6為利用固定時窗和自適應時窗做極化濾波得到的信號分離結果。可見雖然兩種方法都可分離出有效信號，但固定時窗極化濾波法對有效信號造成了損害，有效信號發生畸變(圖6b)； 相比之下，自適應極化濾波能更徹底地分離有效信號，且有效信號未遭受明顯損害(圖6d)。

2.2 正演模型測試

為驗證本文方法的正確性和有效性，設計了一個水平層狀模型，共包含三層，具體參數見表1。針對該模型進行有限差分正演模擬，得到一個含有面波的二分量單炮記錄(圖7)。圖7b中箭頭處為與直達S波混疊在一起的面波干擾，可見其頻散現象明顯，嚴重干擾了有效信號的識別。

圖5 一維合成信號(a)面波干擾； (b)有效信號； (c)合成測試信號； (d)測試信號對應的瞬時頻率

圖6 固定時窗(a、b)與自適應時窗(c、d)極化濾波分離得到的面波干擾(a、c)及有效信號(b、d)

圖7 含面波干擾的X(a)、Z(b)二分量正演記錄

表1 地層參數模型

利用本文方法對正演記錄做面波壓制處理，得到如圖8所示的地震記錄。從去噪后記錄可見，箭頭所指處(圖8b)面波得到了很好的壓制，面波殘留能量較少，直達波和面波被較徹底分離開，有效信號基本未受影響，易于識別。

圖9為分離出的面波成分，圖中除了極少量的反射S波殘留，大部分為面波噪聲。從而驗證了本文方法的有效性。

圖8 使用本文方法壓制面波后的X(a)、Z(b)二分量記錄

圖9 去除的面波的X(a)、Z(b)二分量成分

3 實際數據處理

為驗證本文方法的實際應用效果，針對實際地震資料進行了測試。圖10中三分量地震記錄采集于勝利油田M工區，共有125道，記錄長度為2.5s，時間采樣率為0.5ms。從圖10c的Z分量黃色圈中可見明顯的面波噪聲，它基本呈扇形分布，頻率分布范圍是5～15Hz，嚴重干擾了有效信號的識別。

圖11為利用區域濾波處理得到的結果，可見面波干擾雖被去除，但低頻有效信號也受到了損害。對分離出的面波干擾應用HHT計算，得到三個分量地震數據的瞬時頻率(圖12a～圖12c)及其平均值(圖12d)。

根據瞬時頻率做自適應極化濾波，利用分離出的低頻有效信號對區域濾波結果進行補償(圖13)。

圖10 X(a)、Y(b)、Z(c)三分量原始單炮記錄

圖11 單獨使用區域濾波壓制面波后的X(a)、Y(b)、Z(c)三分量記錄

圖12 利用HHT計算的X(a)、Y(b)、Z(c)三分量地震數據的瞬時頻率及其平均值(d)

對比去噪前、后地震記錄，可見原始記錄中的面波成分已經得到較好壓制，如圖13c的Z分量黃色圈中，被面波掩蓋的有效反射波能清晰地顯現出來，且同相軸連續性較好。

圖14為濾除的面波噪聲，在其中未見到有效波同相軸，說明有效信號受到了較好保護。

圖15為單獨使用固定時窗極化濾波得到的面波壓制記錄。從圖中容易看到，雖然大部分的面波能量已被濾除，但記錄中仍殘留了部分面波，不利于對有效信號的識別。

圖13 本文方法壓制面波后的X(a)、Y(b)、Z(c)三分量記錄

圖14 本文方法分離出的X(a)、Y(b)、Z(c)三分量面波

利用本文方法壓制面波前(黑線)、后(藍線)第20道地震記錄的波形及其放大顯示的對比如圖16所示。可見其中的低頻面波被較徹底地壓制，有效信號得到了較好保護。尤其通過去噪前、后的波形放大對比，發現面波區域的噪聲得到壓制，其他區域的有效信號并未受到損害，故具有較好保幅性。

將本文方法與區域濾波法和極化濾波法壓制面波前、后第20道地震記錄的頻譜(圖17～圖19)進行對比，可以看到：區域濾波法壓制面波后的頻譜(圖17)中，低頻信號明顯被損害； 采用極化濾波法處理所得頻譜(圖18)中，低頻部分仍有少量面波殘留，且其他頻段的有效信號也受到一定程度的影響； 采用本文方法在有效壓制低頻面波的同時，能夠保護其他頻段的有效信號(圖19)。 這就再次驗證了本文方法的有效性和保幅性。

圖16 本文方法去噪前、后第20道X(a)、Y(c)、Z(e)三分量地震記錄波形及其對應的局部放大(b、d、f)

圖17 區域濾波去噪前、后第20道的X(a)、Y(b)、Z(c)三分量地震記錄頻譜

圖18 極化濾波去噪前、后第20道的X(a)、Y(b)、Z(c)三分量地震記錄頻譜

圖19 本文方法去噪前、后第20道的X(a)、Y(b)、Z(c)三分量地震記錄頻譜

4 結束語

全面考慮了面波具有的低速、低頻和橢圓極化的特征，并綜合區域濾波與極化濾波的優勢，本文提出了一種基于區域自適應極化濾波的面波壓制方法： 利用S變換對區域濾波法進行改進，提高了其時頻域的分辨率，摒棄單獨使用區域濾波法，以免對低頻有效信號造成損害； 通過HHT變換計算面波干擾的瞬時頻率，設計自適應極化濾波窗口，克服固定窗口在極化濾波過程中的局限性； 利用自適應極化濾波的多分量去噪優勢對區域濾波結果進行低頻補償，以取得更佳的面波壓制效果。

模型測試和實際數據處理結果表明，本文方法能在有效壓制地震記錄中面波干擾的同時，較好地保護有效信號,效果明顯，實用性強。