黃梅紅, 李月

吉林大學通信工程學院， 長春　130012

基于方向可控濾波的地震勘探隨機噪聲壓制

黃梅紅, 李月*

吉林大學通信工程學院， 長春130012

摘要針對地震勘探隨機噪聲的壓制，本文應用拉伸厄米特高斯函數(shù)設計出方向可控濾波器.根據(jù)時空域上隨機噪聲的無方向性與有效信號的有向性的區(qū)別，通過局部數(shù)字特征，對數(shù)據(jù)進行選擇后重組信號.方向選擇性的增加，使得濾波過程能與不同方向的軸進行匹配，噪聲被壓制的同時保持信號的幅度；方向可調性，使得計算效率提高，且所需存儲空間減少.仿真實驗表明，采用此方法，信號保幅性和去噪效果均比傳統(tǒng)的小波算法以及Curvelet變換好，在-5 db信噪比下，本文方法保幅度為92.99%，信噪比提升221.774%，在實際地震信號處理中有明顯的抑制噪聲、保持有用信號的效果.

關鍵詞拉伸厄米特高斯函數(shù)； 方向可控濾波器； 地震信號； 隨機噪聲壓制

1引言

近年來，小波憑借其時頻域上良好的局部化特性，在地震數(shù)據(jù)隨機噪聲壓制方面應用非常廣泛.常用的二維小波是由一維小波張量積構成的可分離小波，只能對水平、垂直、對角三個方向進行濾波，其基無方向性，對各個方向進行相同程度濾波，濾波系數(shù)難以稀疏表示同相軸的取向信息，對方向邊緣細節(jié)具有模糊作用，一維小波精確的時頻定位性在二維中難以體現(xiàn)(Zhang and Yin，2004；Yu and Whitcombe,2008).為了克服小波無方向性的不足，引進了Ridgelet變換(Candès，1998)，先進行Radon變換，把面上的直線映射到面上的點，再進行小波變換，一定程度提高小波對線的檢測.然而實際中，同相軸曲線居多，不能理想的轉化為集中的點，在地震信號中較少使用.后來，引進了塊剖分的構造思想，在局部中，曲率足夠小，曲線可近似作直線來處理，即Curvelet變換.Curvelet解決了小波方向性不足的缺點，進行了子帶劃分，具有明確方向的基，能夠稀疏的表示圖像平滑區(qū)域，也能較為稀疏的表示邊緣部分，有較好的適用性.但Curvelet的方向性是靠波數(shù)增加、方向濾波器個數(shù)增多的方式來增強，冗余度增加，數(shù)據(jù)存儲空間龐大，同時分塊過程必然引入邊緣效應(Candès and Donoho，1999；Candès et al.，2006；林春和王緒本，2009；Oliveira et al.,2012；袁艷華等，2013).Contourlet變換在低冗余度的情況下，對Curvelet進行離散化實現(xiàn)，采用拉普拉斯金字塔實現(xiàn)多尺度、多方向的分析，雖然在對地震資料處理中有很好的效果，但是依舊無法克服Curvelet變換本身的不足(Do and Vetterli，2002；Li et al.,2008).

我們考慮是否有一種濾波器，在方向上連續(xù)，能夠不進行插值計算得到任意方向的濾波器.Freeman和Adelson(1991)等提出了方向可控濾波器的理論框架，其基本思想是沿任意方向的濾波器都可以表示成一組基濾波器函數(shù)的線性組合形式，這樣一來好似濾波器可以沿任意方向轉動，因此得名.按照方向可控濾波器的思想，輸入圖像與固定的基濾波器卷積即可高效的估計出一副圖像與任一旋轉角度濾波器的卷積結果，極大的減少了計算量；在對方向特征(邊緣、脊等)分析上，方向可控濾波器有很多應用.結合Hilbert變換構成正交濾波器組，鄭林等(2002)提出基于方向可控濾波器的圖像融合新方法；施俊等(2009)提出一種封閉幾何圖形的自動檢測算法；張大明等(2005)則將方向可控濾波器與小波分析結合，也達到了很好的圖像融合的目的，并證明了該方法效果優(yōu)于傳統(tǒng)的小波融合圖像方法.根據(jù)定向濾波器的思想，方向可控濾波器經常應用于金字塔算法中，對圖像進行多尺度、多方向分析.董立娟和練秋生(2008)根據(jù)此思路提出了一種改進的視覺信息處理模型，有效的克服了人臉識別中的一個重要難題；林春等(2011，2012)將方向可控金字塔用于對三維地震切片的分解中，不僅能夠分析不同尺度方向上斷層的特征及走向，還能獲取斷層的信息.通過方向可控濾波器獲取紋理走向，被應用于掌紋(王艷霞和阮秋琦，2008)、虹膜(李志慧，2007)等模式識別.可見，方向可控濾波器主要集中于圖像處理方面，如果將其成功引入到地震信號保幅去噪中，將非常有意義.本文提出了一類新型方向可控濾波器用于地震勘探信號去噪處理，其高方向選擇性使得可以通過分析輸入圖像每個像素的局部模式角度信號來分離不同方向特征以及噪聲，從而達到軸的增強，噪聲的壓制.

Freeman和Adelson(1991)采用的濾波核為厄米特-高斯函數(shù)(HGFs：Hermite-Gauss Functions)，由厄米特多項式和無向性的高斯函數(shù)之積來表示.HGFs的一個固有局限是由于高斯函數(shù)的無方向性，濾波時對各方向都進行了平滑作用，細節(jié)保護不好.于是在高斯諧波函數(shù)中引入一個參數(shù)對HGFs進行拉伸，形成如圖1d的拉伸厄米特高斯函數(shù)(E-HGFs：Elongated Hermite-Gauss Functions).E-HGFs具有明確方向性的特征，對軸敏感度增強的同時，對梯度方向的平滑作用也減弱(Perona，1995).圖1給出了HGFs和E-HGFs的簡單對比，從圖1e和圖1f的濾波結果來看， E-HGFs濾波去噪的同時能更好定位軸的位置.E-HGFs有個重大的缺陷，即不具備可控性，需要設計出可控的E-HGFs的擬合函數(shù)，使其具有方向可控的特性，在2.3部分將給出擬合過程.利用E-HGFs的高方向選擇性和高效性，在局部進行濾波處理，根據(jù)噪聲無方向性，考慮構造和層序的方向，結合統(tǒng)計特性，我們就可以沿著地震數(shù)據(jù)同相軸方向進行濾波.在第3部給出的理論模型和共炮點資料處理結果表明，本文算法可以很好的達到壓制隨機噪聲、保持幅度的目的.

2方向可控濾波器

2.1方向可控濾波器理論

如果一組濾波器f(r,λ)(r為笛卡爾坐標向量，λ為這組濾波器所依賴的參數(shù))能夠展開為一組基濾波器Vn(r)(n=1,2,…N)的線性組合，且線性系數(shù)只與λ有關，那么就說該濾波器具有可控性.可用公式(1)表示為

(1)

其中a(λ)為λ的函數(shù)，V(r)稱為基濾波器.

如果我們用這組濾波器卷積處理一副圖像I(r)時，可用公式(2)表達為

(2)

由(2)式可看出只要I(r)*Vn(r)能夠得到， a(λ)確定，那么任意λ值的濾波結果就能被估計出來.通過線性卷積，對λ無需離散化取值，可以認為Y(r,λ)是關于λ的連續(xù)函數(shù)，具有精確可導可積的性質，克服了傳統(tǒng)通過離散化參數(shù)λ得到濾波器產生的離散插值誤差，同時能減少大量的卷積操作.

現(xiàn)在，我們把研究重點放在方向上，濾波器組在二維中表現(xiàn)為角度旋轉方式，即λ=θ，旋轉矩陣為

(3)

如果一個濾波核T(r)可以展開成一系列適當?shù)腃HF之和，公式為

(4)

那么對于T(r)的旋轉變換可化為

(5)

對比(1)式和(5)式，可以得到濾波核T(r)具有方向可控性，我們稱之為方向可控濾波器.處理圖像結構(Perona，1995)如圖2所示.

圖2　對圖像任意方向處理結構Fig.2　Structure of processing an image in any orientation

2.2二維E-HGFs及其擬合可控函數(shù)

若二維函數(shù)Hk,l(x,y;μ)滿足式(6)，Hk,l(x,y;μ)即為拉伸厄米特高斯函數(shù)(E-HGFs)，公式為

(6)

其中Hp(x)=(-1)pex2(d/dx)pe-x2是p階Hermite多項式，通過參數(shù)μ來控制拉伸程度.當μ=1時，高斯項退化為無方向性，即經典的HGFs.HGFs函數(shù)簡單，用起來較為方便，然而由于其無方向性的局限，濾波時對各個方向均進行平滑作用，降低了檢測精度.E-HGFs憑借其高方向性，一方面克服了折中檢測效果和定位精度的Canny(1986)局限，另一方面可以通過分析不同方向的響應來區(qū)分噪聲和感興趣的特征，應用前景更為廣泛.圖3為圓盤的濾波效果(此圖引用Papari(2012)文獻中Fig.2)，(c)中響應只在垂直方向的中間位置有較大值，很好地展示了E-HGFs對方向定位的高精確度.

然而，E-HGFs是不可控函數(shù)，即在式(5)中N取值無限，因此必須找出可以擬合E-HGFs的可控函數(shù).在方向可控濾波器發(fā)展史上，主要有基于李群理論的方法以及基于奇異值分解的方法.李群法結果為頻譜連續(xù)的函數(shù)，無法進行有限因子分解.奇異值分解法得到的是一系列正交可控性函數(shù)，解決了李群法的無限問題，但結果經常是無法解析的，難以解出數(shù)值解.為克服以上方法的不足，本文采用Giuseppe Papari(2012)等提出的在奇異值分解理論框架上的通過引入特定矩陣算子來構建緊湊可控的小波函數(shù)的方法，以下為其推導：

(7)

Uk,l(ρ,μ)=e-αρ2Hk(μρR)Hl(μ-1ρJ)?(βρ2)，

(9)

根據(jù)(9)式以及厄米特多項式與高斯函數(shù)的關系，可以得到Uk,l(ρ,μ)用貝塞爾函數(shù)展開的式子，Papari(2012)等給出了幾個詳細的結果.本文采用k=0,l=0的濾波器對地震記錄進行方向濾波，其表達式為

U2n(ρ;μ)=e-α ρ2?n(βρ2)，U2n+1(ρ;μ)=0.

(10)

3理論模型和共炮點資料處理

本文采用E-HGFs對地震數(shù)據(jù)進行方向濾波，首先根據(jù)每一點各個方向的濾波響應方差初步區(qū)分噪聲點與有效信號點，然后再與局部小塊均值進行對比，進一步確認，最后按照方差大小比例對數(shù)據(jù)進行非線性濾波，這樣處理的實質是對信號引入了方向維度，濾波器沿著地震數(shù)據(jù)同相軸方向進行濾波.算法流程總結如下：

(1) 采用Papari等(2012)提出的方法擬合出方向可控的E-HGFs濾波器.

(2) 輸入共炮點記錄，進行方向濾波，構建“時間-空間-方向”三維數(shù)據(jù).

圖3　HGFs和E-HGFs一階導對合成圓盤的濾波響應(a) 合成圓盤圖片; (b)(c) HGFs、E-HGFs一階導分別與(a)的卷積結果(Papari et al.，2012).Fig.3　Filter response to the first-order derivative of HGFs and E-HGFs for synthetic image(a) Synthetic image; (b) and (c) Results of convolution with first-order derivative of HGFs and E-HGFs with (a)(Papari et al.，2012).

圖4　k=0,l=0的基濾波器(a,b,c,d,e,f分別對應Fig.

(3) 根據(jù)噪聲的無向性，在方向維度進行方差計算，方差小于一定門限時，初步判定為噪聲點.

(4) 計算鄰域局部均值，結合第(3)步的方差，再次判定，防止誤判.

(5) 對判定的噪聲點進行壓制，其余點取其最大的方向濾波器結果，最后用帶寬足夠大的低通濾波器去除毛刺.

圖7a繪制了在-15db到5db噪聲下，數(shù)據(jù)間隔為1db，采用三種方法對含噪記錄進行處理的結果曲線，從圖中可以看出，在不同信噪比記錄中，本文方法較傳統(tǒng)小波變換和Curvelet去噪后的記錄的信噪比均有所提高.圖7b給出了三種方法在信噪比分別為-15db到5db的處理結果的平均幅值保持，從圖中可以看出，本文方法保幅性明顯存在極大優(yōu)勢.從圖7b中可以看出，噪聲對保幅性影響不大，算法決定了幅值保持度.

通過仿真實驗驗證了本文方法的有效性及可行性，最后是對實際地震數(shù)據(jù)的處理，選取一個共炮點記錄進行實驗.實際數(shù)據(jù)一共168道，每一道的采樣點數(shù)為2501，采樣間隔為1ms.圖8a是原始記錄，

圖5　(a)(b)為原始信號和信噪比為-5dB含噪信號; (c)(d)(e)分別為小波變換、Curvelet變換以及本文方法的濾波結果; (f)(g)(h)分別為小波變換、Curvelet變換以及本文方法的差圖像Fig.5　(a) Raw seismic record with two reflection events; (b) Record with the SNR of -5 dB; (c),(d)and(e) Results after the conventional Wavelet transform, Curvelet transform and steerable filters; (f),(g) and (h) Difference diagrams of conventional wavelet transform, curvelet transform and steerable filters

圖6　(a)單道時域波形對比(b)單道頻譜對比Fig.6　(a) Waveforms of signal in the 102nd channels. (b) Spectrum of signal in 102 channel

圖7　(a)不同噪聲背景下的去噪信噪比(b)不同噪聲背景下的平均幅值保持Fig.7　(a) SNR after processing by different methods under different background noise， (b) Keeping average amplitudes after processing by different methods under different background noise

圖8　(a)原始地震記錄(b)本文方法處理結果Fig.8　(a) Original seismic records， (b) Results of processing by proposed method

圖9　(a) 局部原始地震記錄; (b) 本文方法處理結果Fig.9　(a) Local original seismic records; (b) Results of processing by proposed method

圖10　實際頻譜與處理后頻譜對比Fig.10　Comparison of spectra of original seismic records and processed signals.

圖8b是本文處理結果，可以看到，隨機噪聲已經明顯得到壓制(圖中黑色代表噪聲)，第二幅圖像黑色明顯去除，可以看到同相軸信息且更加連續(xù)，一些被強噪聲湮沒的軸經過處理也恢復出來，如圖方框標出的區(qū)域.圖9為圖8方框區(qū)域內局部第15道到第50道上4000采樣點到5500采樣點的放大圖，圖9a原始記錄中，噪聲湮沒了軸的信息，而經過濾波之后的圖9b噪聲受到壓制，能夠直觀看到同相軸.圖10是截取方框內第22道記錄中3600～5000點的數(shù)據(jù)做的頻譜圖.從頻譜圖中可以看出，頻率在5～40Hz范圍內的信號得到保持，還有一定程度被增強，而頻率超過40Hz之后受到壓制.在實際中，地震信號頻率主要分布在低頻段，噪聲分布在高頻段，結合圖10的頻譜圖，可以看出，本文方法能夠有效的壓制噪聲，對信號起到很好的保持作用.

4結論

本文首次將拉伸厄米特高斯核函數(shù)擬合的方向可控濾波器用于地震信號去噪處理，這種算法出發(fā)點在于根據(jù)同相軸具有空間幾何特性，利用拉伸厄米特高斯核函數(shù)具有明確的方向性，使得構造出來的方向可控濾波器具有高方向分辨率并且能夠準確定位軸的位置，在減少計算量的同時沿著軸向進行濾波處理，有效地壓制隨機噪聲，保持信號的幅值.本文算法使處理后的地震記錄信噪比明顯提升，地震同相軸更清晰、連續(xù)性更好.方向濾波后的數(shù)據(jù)選擇是根據(jù)統(tǒng)計特征進行局部閾值處理，今后，將進一步研究如何對數(shù)據(jù)進行自適應選擇，取得更好的濾波效果.

致謝感謝實驗室?guī)熜謳熃愕臒嵝膸椭绕涓兄x李佳洧師姐的引導，感謝朋友們的鼎力支持,感謝所有給予我無私幫助的人們.

References

Candès E J. 1998. Ridgelets: Theory and applications[Ph.D. Thesis]. USA: Stanford University. Candès E J, Donoho L. 1999. Curvelets. USA: Stanford University. Candès E J, Demanet L, Donoho D, et al. 2006. Fast discrete curvelet transforms.MultiscaleModelingandSimulation, 5(3): 861-899, doi: 10.1137/05064182X.

Canny J. 1986. A computational approach to edge detection. IEEE Transactions on Pattern Analysis and Machine Intelligence, PAMI-8(6): 679-698, doi: 10.1109/TPAMI.1986.4767851.

Do M N, Vetterli M. 2002. Contourlets//Stoeckler J, V-Welland G. Beyond Wavelets. New York: Academic Press. Dong L J, Lian Q S. 2008. Face feature extraction based on steerable filters.Microprocessors(in Chinese), 29(5): 132-134. Freeman W T, Adelson E H. 1991. The design and use of steerable filters.IEEETransactionsonPatternAnalysisandMachineIntelligence, 13(9): 891-906, doi: 10.1109/34.93808.

Li K, Gao J H, Wei W, et al. 2008. Local angle extraction and noise attenuation for seismic image using Contourlet transform. ∥ ETT and GRS 2008. Shanghai: IEEE, 349-352.

Li Z H. 2007. Methodology on Feature Extraction and Descriptors in Iris Recognition[Ph. D. thesis] (in Chinese). Changchun: Jilin University.

Lin C, Wang X B. 2009. Application of second generation Curvelet transform in random noise decaying of seismic data.ComputerEngineeringandApplications(in Chinese), 45(25): 222-224.

Lin C, Wang X B, Liu L H. 2012. Noise decaying of seismic data by Steerable Pyramid decomposition based on Laplace prior.ComputerEngineeringandApplications(in Chinese), 48(2): 222-226.

Lin C, Wang X B, Liu S B, et al. 2011. Steerable Pyramid decomposition method and its application in cracks imaging.JournalofLanzhouUniversity(NaturalSciences) (in Chinese), 47(1): 25-30.Oliveira M S, Henriques M V C, Leite F E A, et al. 2012. Seismic denoising using curvelet analysis.PhysicaA:StatisticalMechanicsanditsApplications, 391(5): 2106-2110, doi: 10.1016/j.physa.2011.04.009.

Papari G, Campisi P, Petkov N. 2012. New families of Fourier eigenfunctions for steerable filtering.IEEETransactionsonImageProcessing, 21(6): 2931-2943, doi: 10.1109/TIP.2011.2179060.

Perona P. 1995. Deformable kernels for early vision.IEEETransactionsonPatternAnalysisandMachineIntelligence, 17(5): 488-499, doi: 10.1109/34.391394.

Shi J, Xiao Z H, Chang Q. 2009. An algorithm for recognizing geometrical shapes automatically based on tunable filter.JournalofNorthUniversityofChina(NaturalScienceEdition) (in Chinese), 30(5): 467-471. Teo P C, Hel-Or Y. 1998. Lie generators for computing steerable functions.PatternRecognitionLetters, 19(1): 7-17, doi: 10.1016/S0167-8655(97)00156-6.

Wang Y X, Ruan Q Q. 2008. A method of line structure features representation and matching of Palmprint.JournalofElectronics

&InformationTechnology(in Chinese), 30(6): 1281-1285. Yu Z, Whitcombe D. 2008. Seismic noise attenuation using 2D complex wavelet transforms. ∥ 70th EAGE Conference & Exhibition. SPE, EAGE.Yuan Y H, Wang Y B, Liu Y K, et al. 2013. Non-dyadic Curvelet transform and its application in seismic noise elimination.ChineseJournalofGeophysics(in Chinese), 56(3): 1023-1032, doi: 10.6038/cjg20130330.

Zhang D M, Hu M L, Zhang C Y. 2005. Image fusion approach based on steerable filters and wavelet analysis.MicrocomputerDevelopment(in Chinese), 15(7): 27-29, 33.

Zhang G Z, Yin X Y. 2004. Noise suppressing by wavelet transform in seismic data processing. ∥ Proceedings of the 2004 7th International Conference on Signal Processing. Beijing: IEEE, 2553-2555.

Zheng L, Han C Z, Zuo D G, et al. 2002. An algorithm for noise image merging based on steerable filters.JournalofXi'anJiaotongUniversity(in Chinese), 36(12): 1236-1239.

附中文參考文獻

董立娟, 練秋生. 2008. 基于可調濾波器金字塔算法的人臉特征提取. 微處理機, 29(5): 132-134.

李志慧. 2007. 虹膜識別特征提取與表達的理論和方法研究[博士論文]. 長春: 吉林大學.

林春, 王緒本. 2009. 第二代Curvelet變換在地震隨機噪聲衰減中的應用. 計算機工程與應用, 45(25): 222-224.

林春, 王緒本, 劉力輝. 2012. Steerable Pyramid分解地震隨機噪聲衰減-基于局部Laplace先驗概率密度模型. 計算機工程與應用, 48(2): 222-226.

林春, 王緒本, 劉四兵等. 2011. Steerable Pyramid分解方法及其在斷層檢測中的應用. 蘭州大學學報 (自然科學版), 47(1): 25-30. 施俊, 肖至恒, 常謙. 2009. 基于可調濾波器的幾何圖形自動識別算法研究. 中北大學學報(自然科學版), 30(5): 467-471.

王艷霞, 阮秋琦. 2008. 一種掌紋紋線結構特征的描述和匹配方法. 電子與信息學報, 30(6): 1281-1285.

袁艷華, 王一博, 劉伊克等. 2013. 非二次冪Curvelet變換及其在地震噪聲壓制中的應用. 地球物理學報, 56(3): 1023-1032, doi: 10.6038/cjg20130330.

張大明, 胡茂林, 張長耀. 2005. 基于方向可調濾波器和小波分析的圖像融合. 微機發(fā)展, 15(7): 27-29, 33.

鄭林, 韓崇昭, 左東廣等. 2002. 基于方向可調濾波器的含噪圖像融合算法. 西安交通大學學報, 36(12): 1236-1239.

(本文編輯張正峰)

Suppression of seismic random noise based on steerable filters

HUANG Mei-Hong, LI Yue*

CollegeofTelecommunicationEngineering，JilinUniversity，Changchun130012,China

AbstractFor seismic random noise suppression, this work designs a steerable filter by taking advantage of elongated Hermite-Gauss functions. According to the different directional responses between valid signal and random noise, we can reconstruct signal by the local characteristics of selected data. With the added directional selectivity, the filtering process can be matched with different directions axes, which makes sure that noise is suppressed without reducing the signal fidelity. The property of directional steerability makes computation more efficient and requires less storage space. Simulation results show that we can get a better signal amplitude and denoising effects than traditional wavelet transform and Curvelet transform algorithm by using this method. In -5db SNR, this method can ensure that the average amplitude reaches 92.99% and SNR enhances 221.774%, which can significantly suppress noise as well as keep the useful signal in processing of real seismic signals.

KeywordsElongated Hermite-Gauss functions； Steerable filters； Seismic signal； Random noise suppression

基金項目國家自然科學基金重點項目(41130421)資助.

作者簡介黃梅紅，女，1991年生，在讀碩士，主要從事信號與信息處理和地震勘探噪聲壓制研究.E-mail:huangmeihong11@sina.com *通訊作者李月，女，1958年生，教授，博士生導師，主要從事信號與信息處理和地震勘探噪聲壓制研究. E-mail: liyue@jlu.edu.cn

doi:10.6038/cjg20160524 中圖分類號P631

收稿日期2015-07-16，2016-02-16收修定稿

黃梅紅, 李月. 2016. 基于方向可控濾波的地震勘探隨機噪聲壓制.地球物理學報，59(5):1815-1823,doi:10.6038/cjg20160524.Huang M H, Li Y. 2016. Suppression of seismic random noise based on steerable filters.ChineseJ.Geophys. (in Chinese),59(5):1815-1823,doi:10.6038/cjg20160524.