基于曲波變換的面波去噪方法研究

姚恒星，謝凱 （長江大學(xué)電子信息學(xué)院，湖北 荊州434023）

在地震勘探領(lǐng)域中，噪聲嚴(yán)重影響了地震資料的后續(xù)處理，因此如何有效地提高地震資料信噪比是地震資料處理的首要任務(wù)。面波［1］作為一種常見的規(guī)則干擾，在疊前數(shù)據(jù)上呈掃帚形狀，具有低頻能量強(qiáng)等特點(diǎn)，嚴(yán)重影響了地震有效反射波，而且還降低了地震數(shù)據(jù)的信噪比。常見的去除面波的方法有低通或高通濾波、F－K頻率波數(shù)譜濾波［2］等。每種方法有著一定的效果，但這些方法都只側(cè)重考慮了面波單一的特性，有著較大的局限性。二維小波變換去噪［3］雖然可以壓制噪聲，但它是以點(diǎn)為元素來描述信號圖像特征，無法有效地表達(dá)邊緣信息，并且還會損傷有效信號。為此Candès等提出了曲波變換［4，5］（Curvelet），即基于脊波變換的多尺度多分辨率的幾何分析方法［6］。該變換稱為后小波變換，不僅彌補(bǔ)了小波變換的不足，還可以更加有效的表示二維圖像，具有各向異性、方向性和局部性，可以稀疏表達(dá)圖像的平滑區(qū)域和邊緣區(qū)域［7，8］。把地震數(shù)據(jù)看成是二維圖像，可以利用面波的方向性［9，10］，用曲波變換去噪的方法去除面波的同時(shí)能夠有效地保護(hù)有效反射信號。為此，筆者提出基于曲波變換的面波去噪方法研究。

1 曲波變換的基本原理

第1代曲波變換的數(shù)字實(shí)現(xiàn)很復(fù)雜，需要子帶分解、正規(guī)化、平滑分塊和脊波分析等許多步驟，并且在進(jìn)行曲波金字塔分解時(shí)帶來了巨大的數(shù)據(jù)冗余。因此，Strack等在第1代曲波變換的基礎(chǔ)上提出了第2代曲波變換的新框架體系。第2代曲波變換在構(gòu)造上已經(jīng)完全不同于第1代曲波變換，實(shí)現(xiàn)過程不需要用到脊波變換，直接通過構(gòu)造曲波窗函數(shù)來實(shí)現(xiàn)曲波分解，不僅變換意義明確，而且實(shí)現(xiàn)起來更加快速和方便。

曲波變換和小波變換屬于稀疏理論的范疇，都是利用基函數(shù)與信號作內(nèi)積來實(shí)現(xiàn)信號的稀疏表示。曲波變換可表示為：

式中，φj，l，k表示曲波函數(shù)，j，l，k分別表示尺度方向和位置參數(shù)。

下面介紹第2代曲波變換的基本原理。曲波變換在頻域內(nèi)的實(shí)現(xiàn)是采用頻域中的窗函數(shù)來表示的，首先構(gòu)造徑向窗和角度窗W（r），r∈和V（t），t∈［－1，1］，那么對所有尺度j，傅里葉頻率窗Uj定義為：

式中，［］表示的整數(shù)部分；Uj是受W和V支撐區(qū)間的限制而獲得的楔形區(qū)域，如圖1所示的陰影區(qū)域。

定義在尺度j，方向θl，位置參數(shù)k＝（k1，k2）處的連續(xù)曲波變換為：

式中，（ω）是二維有效信號的頻域表示表示對φj，l，k（ω）取共軛；Rθl由θl旋轉(zhuǎn)獲得；＝·2－j，k2·2－j2）。

和小波基礎(chǔ)理論一樣，曲波變換包括粗尺度和精細(xì)尺度。粗尺度的曲波變換不具有方向性，所以整個(gè)曲波變換是由粗尺度下各向同性的小波和精細(xì)尺度下的方向元素組成的。

把笛卡爾坐標(biāo)系下的二維函數(shù)f（t）作為有效信號，離散曲波變換定義：

式中，cD（j，l，k）是曲波變換系數(shù)的離散形式（t）是曲波函數(shù)的離散形式。

第2代曲波變換的快速算法有2種：USFFT算法和Wrapping算法。筆者采用基于Wrapping的快速離散曲波算法，其核心思想是圍繞原點(diǎn)Wrapping，對任意區(qū)域周期化，再一一映射到原點(diǎn)仿射區(qū)域。過程如下：

1）對給定的二維函數(shù)f［t1，t2］進(jìn)行二維 FFT，得到頻域表示［n1，n2］，－≤n1，n2≤。

2）在頻域，對每個(gè)尺度和角度組（j，l），重采樣［n1，n2］得到［n1，n2－n1tanθl］，（n1，n2）∈Pj。

3）內(nèi)插后的與窗函數(shù)相乘得到［n1，n2］＝［n1，n2－n1tanθl］＾Uj［n1，n2］。

4）圍繞原點(diǎn) Wrapping局部化［n1，n2］。

5）對，l進(jìn)行二維IFFT（FFT的逆變換）得到離散曲波系數(shù)集合CD（j，l，k）。

Wrapping算法基本思路如下：首先變換到頻域，再在頻域中局部化，最后采用二維IFFT得到曲波系數(shù)。此外局部化和二維IFFT可合二為一，即用局部化窗口來乘局部傅里葉基。該算法采用3個(gè)參數(shù)，使得理解更容易，運(yùn)算操作更簡便，冗余度更低。

圖1 曲波變換的時(shí)域和頻域

2 面波去噪方法

對于地震數(shù)據(jù)的噪聲去除問題，通常采用的方法是通過變換將信號去噪問題從時(shí)域轉(zhuǎn)換到頻域加以解決。基于曲波變換的面波去噪方法如下：首先對經(jīng)過預(yù)處理的含噪信號進(jìn)行多尺度分解，然后在各尺度下盡可能提取出有效信號的曲波系數(shù)，同時(shí)去除噪聲的曲波系數(shù)，最后用曲波逆變換重構(gòu)出地震信號，從而達(dá)到去噪目的。

設(shè)地震信號為s，有效信號為d，噪聲為n，則含噪聲地震數(shù)據(jù)可表示為：

有效信號可用下述方法估算，即：

式中，C表示曲波變換；C－1表示曲波逆變換；Cs表示對地震信號s作曲波變換后的系數(shù)；F表示閾值函數(shù)，定義為：

式中，m為大小與尺度有關(guān)的常數(shù)；σc為曲波域中噪聲標(biāo)準(zhǔn)差；σ為噪聲標(biāo)準(zhǔn)差。

該方法所述的地震信號去除面波方法的流程包括4個(gè)步驟：

1）通過對大量含面波的實(shí)際疊前地震資料進(jìn)行分析，確定地震信號中包含的噪聲模型及其相應(yīng)的各項(xiàng)參數(shù)，從許多具有代表性的面波里取其平均作為面波噪聲模型。

2）因?yàn)榈卣鹦盘柺嵌S的，在空間域中分析存在很多局限性。由于曲波變換具有多尺度和多方向性的特點(diǎn)，能夠稀疏地表示二維信號，使得在曲波域中能更精細(xì)地分離出面波噪聲和有用地震地震反射信號。該方法對二維地震數(shù)據(jù)進(jìn)行曲波分解，選擇合適的尺度，將空域地震信號變換到曲波域中，得到地震信號的各方向各尺度的曲波系數(shù)。同時(shí)對面波模型也進(jìn)行曲波分解，得到面波曲波系數(shù)，將面波數(shù)據(jù)變換到曲波域。

3）根據(jù)面波噪聲模型的曲波系數(shù)在不同尺度和方向上的分布特點(diǎn)來設(shè)置去除面波的閾值。接著采用閾值去噪法，低于閾值的系數(shù)可以認(rèn)為是面波的曲波系數(shù)，從而將其置零去掉；大于閾值的系數(shù)認(rèn)為是有效反射信號的曲波系數(shù)，將其保留。

4）對濾掉面波后的曲波系數(shù)進(jìn)行曲波反變換，重構(gòu)得到去噪后的地震信號。

3 試驗(yàn)結(jié)果與分析

筆者分別對合成地震資料和疊后實(shí)際資料進(jìn)行了相應(yīng)處理。其測試結(jié)果如圖2所示。

圖2（a）是理論合成數(shù)據(jù)，該模型有2條主頻為35Hz和45Hz的反射同相軸和1條主頻為15Hz的面波同相軸，2條斜線代表面波，2條曲線代表有效反射波。該模型數(shù)據(jù)由100道組成，每道包括600個(gè)采樣點(diǎn)。圖2（b）是采用筆者所描述的方法去噪后結(jié)果，圖2（c）是去除掉的面波。由圖2（c）可見去除的大部分是面波。

圖2 模型地震數(shù)據(jù)的去噪結(jié)果對比

圖3（a）是含面波實(shí)際地震數(shù)據(jù)，共33道，每道1501個(gè)采樣點(diǎn)。從圖3（a）中可以看出，原始實(shí)際地震數(shù)據(jù)受面波干擾嚴(yán)重，許多有效波同相軸無法識別。筆者設(shè)計(jì)了包含面波的噪聲模型，將地震數(shù)據(jù)和面波噪聲模型都進(jìn)行曲波分解，得到其在曲波域中的表示。根據(jù)面波噪聲模型在曲波域中分布的尺度和方向，確定了濾波閾值，然后根據(jù)此閾值對地震信號的曲波系數(shù)進(jìn)行處理。圖3（b）是濾波后的重構(gòu)圖像顯示效果，可以看出有效波的水平同相軸變得更加清晰，部分在處理前無法識別的反射層在處理后顯現(xiàn)出來。圖3（c）是通過該方法處理所去除掉的面波噪聲。

表1是實(shí)際地震數(shù)據(jù)去除面波前后的信噪比對比，抽取的是地震數(shù)據(jù)的前10道數(shù)據(jù)。由表1中可以看出在經(jīng)過曲波變換去面波處理后，地震數(shù)據(jù)的信噪比有了明顯的提升，由此可見曲波去面波方法具有較強(qiáng)的實(shí)用性。

4 結(jié)語

通過以上試驗(yàn)可知，將該方法運(yùn)用于地震信號的去噪處理，能夠充分利用面波在曲波域中的分布特點(diǎn)，有效地壓制面波干擾并保護(hù)反射地震信號，重構(gòu)后的地震有效反射信號同相軸變得更加清晰，不僅提高了信噪比，還提高了成像質(zhì)量。

表1 實(shí)際地震數(shù)據(jù)去噪信噪比對比

圖3 實(shí)際地震數(shù)據(jù)的去噪結(jié)果對比

［1］李晶 .面波在地震波場中的特性研究及其應(yīng)用 ［D］.成都：成都理工大學(xué)，2006.

［2］李彩芹，張華 .小波變換與F－K聯(lián)合濾波在面波分離中的應(yīng)用 ［J］.中國煤田地質(zhì)，2007，19（4）：60～61，84.

［3］林椹尠，宋國鄉(xiāng)，薛文 .圖像的幾種小波去噪方法的比較與改進(jìn) ［J］.西安電子科技大學(xué)學(xué)報(bào)（自然科學(xué)版），2004，31（4）：626～629.

［4］董烈乾，李振春，王德營，等 .第2代Curvelet變換壓制面波方法 ［J］.石油地球物理勘探，2011，46（6）：897～904.

［5］蔡炳煌 .基于曲波分析的圖像處理與應(yīng)用 ［D］.汕頭：汕頭大學(xué)，2007.

［6］才溪 .多尺度圖像融合理論與方法 ［M］.北京：電子工業(yè)出版社，2014.

［7］Starck J L，Candes E，Donoho D L.The Curvelet Transform for Image Denoising ［J］.IEEE Transactions on Image Processing，2002，11（6）：670～684.

［8］Kristof De Meersman.Ground Roll polarization filtering with spatial smoothness constraints ［J］.SEG Technical Program Expanded Abstracts，2008（27）：413～416.

［9］代虎 .地震面波資料處理的基本方法 ［J］.黑龍江水利科技，2012，41（9）：28～30.

［10］張恒磊，劉天佑.Curvelet域蒙特卡羅估計(jì)的噪聲衰減 ［J］.西南石油大學(xué)學(xué)報(bào)（自然科學(xué)版），2011，33（4）：64～68.