基于多窗口自適應雙邊濾波去噪方法研究與應用

彭海龍,鄧 勇,赫建偉,劉 兵,王瑞敏

(中海石油(中國)有限公司湛江分公司,廣東湛江524057)

地震數據中除有效反射波之外的所有信息都被視作為噪聲或者干擾波,根據噪聲信號的不同特征,可以將噪聲分為相干噪聲和隨機噪聲兩類[1]。相干噪聲主要包括多次波、面波、海上鳴震、折射波等;隨機噪聲主要包括各種環境因素和人為因素引起的噪聲[2]。相干噪聲在時空域具有一定的規律性,在頻譜、視速度以及到時等方面與有效信號存在差異;而隨機噪聲則缺乏規律性,在各種域中與有效信號混雜在一起,存在于整個記錄之中,較難去除[3]。因此,要想獲取高信噪比的地震資料,必須進行噪聲壓制。

目前,消除地震記錄中的隨機噪聲的方法主要分為時空域方法和變換域方法。時空域方法主要包括多項式擬合濾波[4-5]、中值濾波[6]、奇異值分解濾波[7]、獨立分量分析濾波[8]、雙曲濾波[9]和時頻分析法[10]等;變換域方法主要包括K-L變換濾波[11]、F-X預測濾波[12]、小波變換濾波[13]和曲波變換濾波[14]等。上述去噪方法都已得到廣泛應用且取得了較好的去噪效果,但受各自假設條件和方法原理的限制,對于地震數據中有效信號不連續的邊緣(斷層)以及突變點(地層或者巖性尖滅點)等位置處的噪聲壓制效果不甚理想。常規的構造導向濾波方法考慮的是單一窗口內的數據信息,且窗口大小很難選定合適的范圍[15],因此,在地層傾角方位角約束下的去噪效果難以滿足需求。而其它常規濾波方法如F-X濾波方法在選取濾波窗口時,沒有很好地考慮地層結構因素,對地層不連續地質體有效信號的損傷較大。雙邊濾波方法由于具有簡單、保護邊緣信息的特點得到了廣泛應用,國內外學者對其進行了大量的研究[16],提出了自適應雙邊濾波方法以及改進型的雙邊濾波方法。然而由于核函數的建立方式對數據的信噪比和邊緣信息較為敏感,濾波窗口的選取也較為困難,過大或者過小都容易損傷有效信號,因而限制其在保持邊緣信息情況下壓制隨機噪聲的效果[17-18],在不連續地質體(斷層、裂縫等)較發育的區域,上述方法很難取得理想的噪聲壓制效果。

針對上述情況,本文采用基于多窗口自適應雙邊濾波方法對地震信號中的復雜斷塊區域噪聲進行壓制。在自適應雙邊濾波方法基礎上,沿局部地層傾角方位角提取濾波窗口,進行多窗口邊緣保持濾波。首先介紹了計算地質體傾角方位角信息的方法原理,然后根據地質信息提取沿層多濾波窗口,再介紹在選定的濾波窗口中執行自適應雙邊濾波的過程。利用建立的三維模型數據驗證了該方法的有效性,最后在南海L盆地某工區進行了實際應用。

1 方法原理

1.1 局部地層傾角方位角計算

基于多窗口自適應雙邊濾波方法(簡稱本文方法),重點在于沿局部地層傾角方位角提取濾波窗口,因此,該方法的第一步在于計算局部地層的法線方向。從梯度結構張量出發[19],分別計算三維地震數據3個方向的梯度,則梯度結構張量T的表達式為:

(1)

其中,

(2)

?u/?x,?u?y,?u/?z分別表示地震數據u在3個方向的梯度,分別記為g1,g2,g3,則平滑張量TS為:

(3)

其中,高斯濾波器G的表達式為:

(4)

式中:σ表示濾波尺度因子。

對平滑張量進行特征值分解:

(5)

式中:Λ表示對角矩陣,對角線上的值為特征值,由大到小分別排列為λ1,λ2,λ3;v表示特征向量矩陣,列向量表示以上特征值對應的特征向量,分別為v1,v2,v3。局部地層的法線方向就是最大特征值λ1對應的特征向量v1。

地層傾角θ和方位角φ表達式如下:

式中:v1x,v1y,v1z分別表示v1在x,y,z方向的分量。地層傾角和方位角計算完畢之后,按公式(8)和公式(9)分別估算x,y方向的延遲時間Tx,Ty。

根據兩個方向的延遲時間來獲取每一個目標點位置處的局部沿層多濾波窗口,然后執行多窗口自適應雙邊濾波。

1.2 沿層提取多濾波窗口

為了進一步壓制地震數據中的隨機噪聲,提高信噪比,LUO等[20]提出一種多窗口保邊濾波方法。該方法在每一個目標點位置處設置一個滑動窗口,以尋找所有包含該目標點的窗口。然后計算每個窗口內所有樣點信息的標準方差,計算方差最小窗口內的所有樣點的平均值,并將該平均值作為目標點的輸出結果[21],能夠取得較好的去噪效果。為利用多窗口保邊濾波的優勢,我們將多窗口濾波的多窗口提取思路引入本文濾波技術流程。具體為:根據上一步計算的局部地層的傾角、方位角,沿層提取目標點周圍多窗口,該多窗口中包含不同空間位置處的數據,多窗口的提取示意如圖1所示。圖1中,紅色區域表示提取的多窗口范圍,紅色區域上下虛線框表示具有傾角θ和方位角φ的地層,V表示地層的法線方向向量。然后計算目標點周圍多窗口的標準方差,找到方差最小的窗口,在該窗口內進行雙邊濾波。

圖1 局部地層濾波窗口提取示意

多窗口的設計參考WANG等[22]的思想,采用窗口多維劃分方法,以二維數據的總窗口5×5為例,劃分為9個小多邊形,分為4個五邊形(藍色虛線框),4個六邊形(綠色實線框)和1個矩形(紅色實線框),具體如圖2所示。

圖2 多窗口劃分方法(二維)示意

1.3 雙邊濾波方法

雙邊濾波方法最初用于圖像信息處理,它結合給定的窗口中數據信息對噪聲進行處理,其公式如下[23]:

(10)

(11)

(12)

(13)

式中:N(x,y)表示去噪后圖像;wp為歸一化參數;I(i,j)表示原始圖像;Ω表示濾波的窗口范圍;ws(i,j)表示高斯距離權值;wr(i,j)表示高斯像素相似度權值;d和δ分別表示兩個像素ξ和c的空間距離差和像素差;σs和σr分別表示相似度方差和空間方差。公式(10)中,主要通過處理目標點位置像素與窗口內其余像素的距離和像素亮度信息作為權重進行加權平均之后,得到目標位置處的像素。

自適應雙邊濾波方法根據圖像的局部特征自適應地設置空間參數,通過計算圖像各個像素點的目標尺度得到目標點像素周圍的平滑區域范圍,用來控制空間距離權重中的相似度權重[24]。在邊緣和不連續信息變化區域,目標尺度小,則相似度權重越小;在平滑連續區域,目標尺度大,則相似度權重越大。相似度方差表達式[25-26]為:

(14)

式中:I(x,y)為目標位置像素值;σu表示圖像梯度分布的統計參數;Bij(Ω)表示濾波窗口的區域范圍。根據前文提取的多窗口,在濾波窗口內執行雙邊濾波。

2 模型測試

為對比常規構造導向濾波方法、F-X濾波方法和本文方法,設計一個三維模型,該模型淺部為連續沉積的傾斜地層,中深部為一個包含斷層的連續沉積波狀地層(圖3a),以此為基礎,添加一定的隨機噪聲,含噪模型如圖3b所示。以含噪模型為處理數據,分別采用上述3種方法進行濾波處理。選取處理后數據的縱橫測線剖面和時間切片(分別見圖3中的虛線、實線位置)進行對比,3種方法的處理結果如圖4所示。

處理結果表明,3種方法都能在一定程度上壓制噪聲,但是不同的去噪方法其壓制噪聲效果不同。采用F-X濾波方法能有效抑制噪聲,但是對斷層邊緣信息有所損傷,使得斷層邊界變得模糊,斷面成像不干脆(圖4a,圖4b);且去除的噪聲中包含明顯的斷層有效信息(圖4c和圖4d)。

圖3 原始無噪三維模型(a)和含噪三維模型(b)

圖4 不同濾波方法的去噪結果以及去除的噪聲a F-X濾波水平切片; b F-X濾波橫測線剖面; c F-X濾波去除的噪聲水平切片; d F-X濾波去除的噪聲橫測線剖面; e 常規構造導向濾波水平切片; f 常規構造導向濾波橫測線剖面; g 常規構造導向濾波去除的噪聲水平切片; h 常規構造導向濾波去除的噪聲橫測線剖面; i 本文濾波結果的水平切片; j 本文濾波結果的橫測線剖面; k 本文方法去除的噪聲水平切片; l 本文方法去除的噪聲橫測線剖面

采用常規構造導向濾波方法的結果中斷層和不整合面的邊緣保持以及噪聲壓制都優于F-X濾波方法,但是斷面邊緣仍然存在一定的模糊(圖4e,圖4f);去除的噪聲中仍然包含一定的斷面和不整合面的有效反射(圖4g,圖4h)。采用本文方法去除噪聲的效果最好,斷面清晰,邊緣保持較好(圖4i,圖4j);去除的噪聲中斷面和不整合面的有效反射信息不明顯(圖4k,圖4l)。

3 實際應用

采用南海L盆地某工區的實際資料進行測試,該地震資料除包含與自由表面相關的多次波和繞射多次波外,在斷層、裂縫等不連續地質體發育的中深層,隨機噪聲較為嚴重[27-28],具體表現為:①該區塊屬于淺海相沉積,水深變化劇烈,水道砂、席狀砂以及各種裂隙、裂縫和斷裂較為發育,斷面反射不清晰,地層連續性較差;②淺層含氣范圍廣,中深層屬于高溫高壓區,波阻抗差異小,各向異性較強,吸收衰減嚴重,模糊區較為發育,地震資料信噪比低。由于該區隨機噪聲嚴重,對地質體邊緣的解釋工作容易形成誤判,而常規去噪方法很難解決復雜構造情況下的噪聲壓制問題,且易傷害地層和斷面的有效反射信號,因而影響后續地震屬性體的品質,不利于圈閉和儲層預測。

分別采用F-X濾波方法、常規構造導向濾波方法和本文方法去除噪聲。采用F-X濾波方法去除噪聲的結果如圖5和圖6所示,從圖5c中可以看出,去除的噪聲中出現有效反射信息(圖中橢圓處),且在去除的噪聲水平切片中也出現了斷面以及地層的邊緣信息(圖6c),表明F-X濾波方法損害了有效反射信息。采用常規構造導向濾波方法去除噪聲的結果如圖7和圖8所示,可以看出,在去除的噪聲剖面中也出現斷層有效反射信息(圖7c),且在去除的噪聲水平切片中也出現了斷面以及地層的邊緣信息(圖8c),表明常規構造導向濾波方法同樣損害了有效反射信息。采用本文方法去除噪聲的結果如圖9和圖10所示。對比圖5c,圖7c和圖9c可以看出,采用本文方法去除噪聲后原始數據中斷面信息基本沒有受到損傷(圖9c橢圓位置);對比圖6c,圖8c和圖10c 可以看出,采用本文方法去除的是不規則無規律的信號,表明濾掉的是隨機噪聲,且沒有損傷地層和斷面的有效反射信號(圖10c)。以上結果充分驗證了本文方法的有效性及合理性。

圖5 采用F-X濾波方法去除噪聲結果及其去除的噪聲剖面a 原始數據; b 濾波結果剖面; c 去除的噪聲剖面

圖6 采用F-X濾波方法去除噪聲結果及其去除的噪聲水平切片a 原始數據; b 濾波結果水平切片; c 去除的噪聲水平切片

圖7 采用常規構造濾波方法去除噪聲結果及其去除的噪聲剖面a 原始數據; b 濾波結果剖面; c 去除的噪聲剖面

圖8 采用常規構造濾波方法去除噪聲結果及其去除的噪聲水平切片a 原始數據; b 濾波結果水平切片; c 去除的噪聲水平切片

圖9 采用本文方法去除噪聲的結果及其去除的噪聲剖面a 原始數據; b 濾波結果剖面; c 去除的噪聲剖面

圖10 采用本文方法去除噪聲的結果及其去除的噪聲水平切片a 原始數據; b 濾波結果水平切片; c 去除的噪聲水平切片

4 結論

本文提出的多窗口自適應雙邊濾波方法,改進了常規雙邊濾波方法的單窗口技術流程,采用多窗口濾波流程,多窗口的選擇根據地質構造信息進行,減小其核函數對地震數據信噪比和地質體邊緣信息的依賴性。與常規構造導向濾波方法和F-X濾波方法相比,本文方法能夠在有效壓制隨機噪聲的基礎上,更好地保持邊緣信息。經模型數據和實際數據驗證,該方法能夠更好地實現邊緣保持,同時有效壓制地震信號中的隨機噪聲,去除噪聲后的資料信噪比得到明顯提高,地震剖面反射同相軸連續,保真性好,具有較強的實用性。

由于該方法是在多窗口中進行雙邊濾波,尋找最小方差的窗口作為濾波窗口進行濾波并輸出結果,與常規濾波方法相比,不足之處在于計算時間較長。因此,下一步工作是對該方法進行優化,在實現保邊和去噪效果的基礎上提高其計算效率。