低信噪比地震數(shù)據圖像結構引導去噪方法與應用

鄭啟明,李 琦,都小芳,吳高奎

(1.北京派特杰奧科技有限公司,北京100044;2.中國地質大學(北京)海洋學院,北京100083;3.中石油東方地球物理公司研究院,河北涿州072751;4.中石化石油勘探開發(fā)研究院,北京102249)

基于地震資料刻畫斷層、裂縫等結構特征,對于發(fā)現(xiàn)和描述油藏具有重要意義,可支撐決策最佳開發(fā)方案、減少探井風險,提高石油勘探成功率。然而,在地震勘探區(qū)域,由于地質結構比較復雜,斷裂帶的地層形態(tài)多變,因而常導致采集的地震數(shù)據信噪比較低,使得地震數(shù)據處理、解釋十分困難。高品質的地震資料是高質量解釋、油氣儲層精準預測的重要依據和保證。因此,在地震數(shù)據處理時需要采用有效的去噪技術來提高地震資料的品質和可靠性。地震資料處理和解釋人員探索了一系列的地震數(shù)據濾波技術以提高地震資料信噪比,如頻率-波數(shù)域濾波、時域中值濾波等[1],這些技術可以有效去除地震數(shù)據中的部分噪聲,從而提高地震資料解釋和油氣儲層預測的準確性。相比通用的圖像處理算法,地震數(shù)據處理更關注利用反射同相軸的橫向一致性,處理過程中可能會壓制傾斜斷層陰影部分的幅值、破壞小尺度斷層和裂縫等非連續(xù)性的結構信息,導致錯誤匹配斷層兩側的同相軸。因此,進一步探索研究可以保護斷裂、裂縫、尖滅等地質信息的疊后地震保邊去噪技術對地震解釋工作尤為重要。

為更好地利用地震反射的局部信息,LUO等[2-3]首次將保邊光滑濾波的思想引入地震資料解釋性處理領域。采用Kuwahara多窗分析技術得出,此類方法可以考慮局部結構信息設計自適應光滑濾波器[4-5],在提高信噪比的同時,保護地震圖像中的有效地質信息。此后,國內外學者引入基于各向異性擴散濾波的構造導向濾波技術[6-10],通過結構張量評估局部結構走向、計算各向異性算子,設計增強地震資料沿同相軸方向的信噪比、抑制斷點位置光滑強度的自適應濾波器。在此基礎上,國內外學者優(yōu)化了相關算法,提高了結構張量評估局部層位走向的有效性[11-16],通過偏微分方程的近似解法提高了計算效率[17],并將相關技術擴展應用到三維地震數(shù)據的處理[18-20]。考慮到擴散濾波的有效性依賴于結構張量的準確估計,HALE[21]提出了引入高度自適應的雙邊濾波技術來實現(xiàn)保邊光滑處理的方法,但該方法的計算代價相對較高。國內學者也探索了在變換域[22-24],如曲波域、變分模態(tài)域或多道奇異譜上的數(shù)據濾波去噪處理技術。TRINH等[25]提出了利用Bessel函數(shù)逼近各向異性Laplace濾波的結構導向光滑算法;在圖像處理領域,HE等[26]提出了利用引導圖像濾波進行自適應保邊光滑處理的方法,計算效率顯著高于雙邊濾波方法,該方法尚未在地震數(shù)據處理中廣泛應用。此外,結構導向光滑濾波技術在全波形反演中也具有較明顯的應用效果[27-28]。近年來人工智能技術飛速發(fā)展,高性能濾波算法[29]、神經網絡地震數(shù)據去噪方法[30-31]和智能化斷層表征技術[32-33]也得到了深入的探索。

在前人研究的基礎上,將“二維圖像濾波”首次擴展應用到二維或三維地震數(shù)據處理中,該計算過程有別于傳統(tǒng)的二維圖像濾波。為抑制隨機與相干噪聲、增強地震反射同相軸的一致性,建立了雙邊顯式加權濾波器的框架,引入了分時窗的引導濾波核函數(shù)與正規(guī)化算子,建立輸入與輸出的目標方程,利用線性回歸顯式求解方程系數(shù),迭代計算直至收斂。該算法采用高效率的保邊光滑處理,使得處理后的地震資料能更加準確地刻畫地質結構,從而更好地服務于地質勘探開發(fā)。

1 圖像結構引導濾波原理

針對地震數(shù)據去噪開發(fā)的濾波器可分為非平均濾波器和加權平均濾波器,前者以中值濾波為代表,其濾波過程一般為非線性過程,后者可進一步分為顯式濾波器和隱式濾波器。顯式濾波器通過設計濾波核函數(shù),然后根據核函數(shù)對待處理的數(shù)據進行加權平均實現(xiàn)特定濾波的目的,其中較典型的有均值濾波、高斯濾波。在地震數(shù)據處理中,以隱式加權平均的擴散濾波方法較為常見,其將包含在矩陣中的濾波器作用于待處理數(shù)據,可以較為直接地在矩陣中引入各向異性參數(shù)并求解隱式矩陣來實現(xiàn)結構導向濾波。需要注意的是,各向同性擴散濾波與高斯濾波是等價濾波器。

定義問題Iout=Iin+n,其中,Iin,Iout,n分別為含噪聲數(shù)據、去噪后的數(shù)據和噪聲。參考顯式加權平均濾波器的設計思路,可利用線性加權平均定義濾波核函數(shù):

(1)

(2)

(3)

其中,k為窗wk的中心坐標,(ak,bk)是表征該線性關系的參數(shù),假定其在窗wk中為常量。若輸入數(shù)據Iin和輸出數(shù)據Iout的差異為隨機或相干噪聲,則確定參數(shù)(ak,bk)的最小化目標函數(shù)為:

(4)

其中,ξ是正則化算子系數(shù),用以確保ak的解在可接受的數(shù)值范圍內。利用線性回歸算法,可顯式地求解(ak,bk)的值,即:

(5)

(6)

分時窗處理過程中,由于在時窗交疊區(qū)的輸出數(shù)據Iout具有非唯一性,因此,在同一點處多個交疊窗的平均輸出表達式為:

(7)

(8)

傳統(tǒng)的雙邊濾波器可以利用引導數(shù)據通過顯式加權平均進行保邊結構光滑,其核心在于使用線性加權平均濾波器,可以引入濾波核函數(shù)。引導濾波器采用同樣思想,但通過構建最優(yōu)化問題引入了隱式的加權平均,通過線性最優(yōu)化解獲得了引導濾波器算子,同時在優(yōu)化問題中引入正則化算子形成超參數(shù),可控制對于邊界結構的敏感程度。

以上推導過程適用于二維或三維地震數(shù)據體。在具體應用過程中,需要準備引導數(shù)據G(原始輸入數(shù)據或濾波處理后的窄帶高信噪比數(shù)據),試驗確定半窗寬度r和正則化因子ξ兩個參數(shù)。當r和ξ變大時,處理后的數(shù)據具有更強的光滑效果,可能會損傷地震反射的斷裂結構信息,實際計算中需要根據保邊處理的目的與需求,綜合考慮確定合適的處理參數(shù)。由公式(5)至公式(7)可知,該算法可多次迭代應用方框濾波器實現(xiàn)地震數(shù)據的去噪處理,自動具備O(n)的計算效率。

具體應用于實際地震數(shù)據處理時,采用以下技術流程:

1) 準備待處理的地震數(shù)據Iin(二維或三維地震數(shù)據)、結構引導數(shù)據G(可以與Iin相同)、設置算法參數(shù)半窗寬度r和正則化因子ξ;

2) 計算f(Iin),f(G),f(Iin·Iin)與f(Iin·G),其中,f為均值濾波器;

3) 計算參數(shù)a=[f(Iin·G)]/[f(Iin·Iin)-f(Iin)f(Iin)+ξ];

4) 計算參數(shù)b=f(G)-af(Iin);

5) 計算濾波后的數(shù)據Iout=f(a)Iin+f(b)。

此流程主要計算6次均值濾波,由于均值濾波器可以使用方框濾波器的方法構建,因此本流程具有O(n)復雜度的計算效率。

該流程中結構信息的引入依賴于結構引導數(shù)據G的選擇。引導濾波應用過程中利用均值濾波器提取G的結構信息。假如G由多個均質的區(qū)塊構成,那么引導濾波會在各個區(qū)塊內分別進行操作,而不會在光滑區(qū)塊的邊界進行均值濾波,從而有效保留邊緣信息,實現(xiàn)保邊濾波的效果。若G是輸入數(shù)據本身,則引導濾波器會強化數(shù)據本身的走向或者傾斜結構,同時具有保留邊界、突出主要結構的功能。導向濾波器中的引導數(shù)據G的選取可根據具體應用調整,本次研究選擇輸入數(shù)據作為引導數(shù)據,是較為直接和有效的方法。

2 模型數(shù)據測試

為驗證圖像結構引導濾波算法的穩(wěn)定性與結果的可靠性,對比該算法與傳統(tǒng)濾波器及主流濾波器的優(yōu)劣,基于Marmousi模型通過卷積15Hz主頻的Ricker子波和反射系數(shù),在設計的二維合成地震剖面(圖1a)中加入一定噪聲得到信噪比(signal to noise ratio,SNR)約8dB的含噪數(shù)據(圖1b),以此評估后續(xù)各類濾波算法的去噪效果。其中,SNR定義如下:

(9)

圖1 基于Marmousi模型的合成地震剖面去噪測試對比a 無噪剖面; b 含噪剖面; c 圖像結構引導濾波后的結果; d 構造導向濾波后的結果; e 中值濾波后的結果; f 均值濾波后的結果; g 高斯濾波后的結果

其中,Isignal表示不含噪聲的有效信號,Isignal+noise表示含噪聲的信號,SNR表示計算得到的信噪比。

為方便對比不同方法的去噪效果,試驗中選用了相同的半窗寬度,其中圖像結構引導濾波與構造導向濾波皆以原始數(shù)據本身作為引導數(shù)據提取結構張量進行處理,含噪地震數(shù)據的自適應引導圖像濾波、構造導向濾波、中值濾波、均值濾波和高斯濾波的處理結果分別如圖1c至圖1g所示。

由模型去噪結果可知,圖像結構引導濾波(圖1c)、構造導向濾波(圖1d)、中值濾波(圖1e)、均值濾波(圖1f)、高斯濾波(圖1g)均有效消除了部分噪聲,提高了模型數(shù)據的信噪比。前兩種方法的信噪比由8.00dB分別提高至18.31,14.77dB,與構造導向濾波相比,圖像結構引導濾波的信噪比高約4dB,對地層邊界具有更好的保護作用;而后3種方法的信噪比由8.00dB分別提高至9.19,9.58,10.09dB,信噪比雖有所提高,但在一定程度上模糊了地層與斷裂的反射邊界。

為便于清晰展示上述去噪效果,圖2更詳細地對比了不同方法分離的噪聲剖面。圖1與圖2的截斷幅值不同。在層位相對平緩的區(qū)域,如距離為8km,時間為1s處,圖像結構引導濾波和構造導向濾波均較好地保持了地震數(shù)據的結構信息(圖2b、圖2c),而中值濾波、均值濾波和高斯濾波因過度光滑破壞了有效地震反射的真實結構信息(圖2d、圖2e、圖2f)。針對本合成數(shù)據中其它斜率較大的結構,從噪聲剖面可知,構造導向濾波(圖2c)保留邊界的能力相對中值濾波(圖2d)、均值濾波(圖2e)和高斯濾波(圖2f)更強,但是圖像結構引導濾波(圖2b)更好地保留了地震反射細微的結構。總之,相比其它處理方法,圖像結構引導濾波不但能有效抑制噪聲,而且在處理方法與計算效率上也更加簡單高效。

圖2 不同方法分離的噪聲剖面對比a 原始噪聲; b 圖像結構引導濾波; c 構造導向濾波; d 中值濾波; e 均值濾波; f 高斯濾波

3 實際應用效果

以信噪比較低、含有斷層的三維地震數(shù)據為例,利用含噪數(shù)據為引導,基于核濾波函數(shù)自適應地對其進行濾波處理。圖3和圖4分別為自適應濾波前、后不同頻帶三維地震剖面和沿層振幅切片。由圖3和圖4可見,自適應濾波處理前,低頻地震數(shù)據的地震剖面和沿層振幅切片中存在明顯低頻噪聲(圖3a、圖4a);受噪聲影響,中頻信號與小尺度斷層反射信息相混疊,導致屬性邊界模糊不清(圖3c、圖4c);高頻地震數(shù)據的地震剖面和沿層振幅切片中的高頻有效信號受噪聲影響更加嚴重,信噪比極低,同相軸連續(xù)性極差(圖3e、圖4e)。

圖3 自適應濾波前、后不同頻帶三維地震剖面a 低頻含噪信號; b 低頻圖像結構引導濾波結果; c 中頻含噪信號; d 中頻圖像結構引導濾波結果; e 高頻含噪信號; f 高頻圖像結構引導濾波結果

圖4 濾波前、后不同頻帶沿層振幅切片屬性a 低頻含噪信號; b 低頻圖像結構引導濾波結果; c 中頻含噪信號; d 中頻圖像結構引導濾波結果; e 高頻含噪信號; f 高頻圖像結構引導濾波結果

與圖3a、圖4a低頻信息相比,圖像結構引導濾波很好地壓制了低頻數(shù)據中的散射噪聲,同時很好地保持大尺度斷層的錯斷邊界(圖3b、圖4b);其次,對比圖3c、圖3d和圖4c、圖4d可見,本文濾波處理方法可有效保邊光滑去除噪聲、提高信噪比,同時保持小尺度斷層、沉積體邊界的有效反射信息(圖3d、圖4d 標識區(qū)域);對比圖3e、圖3f與圖4e、圖4f可見,本文方法可自適應利用高信噪比優(yōu)勢頻帶信號作為引導圖像數(shù)據,有效壓制高頻噪聲、提高地震同相軸連續(xù)性與信噪比(圖3f與圖4f方框區(qū)域),有助于解釋人員提高小斷層和地震屬性的解釋與描述精度。

構造導向與圖像結構引導濾波的全頻帶處理結果如圖5所示。在含噪地震數(shù)據(圖5a)上提取結構張量的傾斜構造導向濾波結果如圖5b與圖5c所示;以含噪數(shù)據自身為引導的圖像結構引導濾波處理結果如圖5d 所示。由圖5可見,采用構造導向濾波與自適應圖像結構引導濾波處理后的地震剖面信噪比顯著提高,除大斷層特征更加明顯之外,小斷裂的細節(jié)更加清晰。對比圖5b、圖5c可見,當構造導向濾波的相關半徑較大時,會導致計算效率降低、橫向保邊結果略差;當相關半徑較小時,局部產生相干噪聲異常(圖5c中紅圈),散射噪聲壓制不徹底(圖5b、圖5c黃框)。

圖5 原始地震數(shù)據、構造導向濾波與自適應圖像結構引導濾波結果a 原始地震數(shù)據; b 小半徑構造導向濾波結果; c 大半徑構造導向濾波結果; d 自適應圖像結構引導濾波結果

自適應圖像結構引導濾波可較為徹底地壓制相干噪聲、保護有效反射信號(圖5d)。在算法上,其相比于構造導向濾波,可自適應地從原始含噪地震數(shù)據中優(yōu)選引導圖像數(shù)據與處理參數(shù),無需象構造導向濾波那樣為獲取傾角信息而多次調整參數(shù)來計算結構張量并重復進行引導濾波;在計算效率上,其利用的線性函數(shù)可有效降低計算復雜度、提高運算效率,計算代價僅相當于數(shù)次均值濾波。

對原始地震數(shù)據、構造導向濾波和自適應圖像結構引導濾波結果分別利用相同的特征值參數(shù)提取相干體剖面與沿層相干屬性信息進一步評價處理效果(圖6、圖7)。

圖6 原始地震數(shù)據、構造導向濾波與自適應引導圖像濾波后相干體剖面a 原始地震數(shù)據相干體剖面; b 小半徑構造導向濾波后相干體剖面; c 大半徑構造導向濾波后相干體剖面; d 自適應圖像結構引導濾波后相干體剖面

圖7 原始地震數(shù)據(a)與自適應圖像結構引導濾波后(b)沿層相干屬性

在斷裂位置附近,對比原始數(shù)據與不同構造導向半徑獲取的相干結果(圖6a、圖6b與圖6c),因受散射噪聲影響,導致其對斷層的分辨能力大大降低,給后續(xù)的解釋工作造成困難。此類污染主要由于斷裂帶附近的地震反射資料信噪比低、同相軸的連續(xù)性不夠明顯,導致計算的相干體屬性難以區(qū)分地震層位與斷裂結構特征;與其它方法相比,本文方法可有效壓制地震數(shù)據中的噪聲,因此可獲得信噪比更高的相干體,從而能夠更加清晰地描述復雜的斷裂結構(圖6d);同時,沿層相干屬性也能夠很好地描述細微斷裂結構,多級斷裂在橫向上的層次關系也更加清晰(圖7)。該結果有助于后續(xù)的層位與斷裂解釋。

4 討論

圖像結構引導濾波的引導數(shù)據G的設計具有關鍵作用。傳統(tǒng)的構造導向濾波需要利用結構張量等方法,顯式地提取結構體在某個區(qū)域內的傾角和方位角,而本文方法的引導圖像數(shù)據不需要顯式表征結構傾角和方位角,只需要其梯度具備與期望輸出圖像Iout一致的邊界信息。同時由于其具有滑動窗實現(xiàn)的特點,以含噪聲的圖像作為引導圖像,可在有效保留結構邊界的前提下壓制隨機噪聲。然而,在實際地震數(shù)據處理中,相干噪聲、散射噪聲是影響地震圖像質量的關鍵,僅以含噪圖像作為引導可能會增強相關噪聲的強度。但由于引導圖像本身設計的靈活性,可采用以下步驟壓制相干噪聲和散射噪聲:①先采用自引導的方式壓制隨機噪聲;②再利用均值、中值、高斯或構造導向濾波處理剖面,獲得邊界受到一定破壞但相干、散射噪聲得到壓制的剖面;③再以步驟②的輸出結果作為引導圖像處理步驟①的輸出結果。但是,此方法仍具有一定局限性,其有效范圍僅限于相干或散射噪聲相對較弱的狀態(tài),否則其它濾波器在壓制相干或散射噪聲時會過度傷害主要結構信息。需要在未來工作中進一步改進圖像結構引導濾波,形成自適應壓制相干或散射噪聲的方法。

5 結論

在進行地震資料的斷層、裂縫解釋時,噪聲的存在會明顯影響解釋精度。通過引入圖像結構引導濾波方法,能自適應地提取有效地震結構引導信息,快速、有效地對地震數(shù)據進行保邊、保幅去噪處理。圖像結構引導濾波基于理論推導可通過方窗濾波器實現(xiàn),對于大尺度地震數(shù)據不存在顯著的計算瓶頸,具有O(n)的計算效率,相比工業(yè)界更常用的構造導向濾波器也有一定優(yōu)勢。含噪合成模型數(shù)據的處理分析結果證實了圖像結構引導濾波可保留反射結構邊界、有效壓制噪聲,結果顯著優(yōu)于常規(guī)非結構性光滑濾波器。針對實際地震數(shù)據的處理,利用引導圖像可靈活設計引導特征,在一定程度上壓制相干噪聲與隨機噪聲,提高地震數(shù)據信噪比,增強地震反射同相軸的一致性、保留并突出斷裂與地層反射的非連續(xù)特征。計算的相干屬性較其它方法斷裂結構特征更加清晰,地質規(guī)律更加明顯,進一步證明了方法的有效性與實用性。