魯光銀， 羅 帥， 朱自強， 趙云威， 席飛雁

(中南大學 地球科學與信息物理學院，長沙 410083)

0 引言

探地雷達數(shù)據(jù)中隨機噪聲通常為高頻非平穩(wěn)信號，且分布于整個數(shù)據(jù)剖面，特別是在環(huán)境噪聲嚴重、淺層干擾大或探測深度較大時，剖面數(shù)據(jù)往往不能用于直觀解譯[1]。實際工作中雖然可以在一定程度上壓制隨機干擾，但是要求數(shù)據(jù)達到解譯的標準，就需要對其進行深入的去噪處理。

目前對隨機噪聲的壓制方法主要有頻率域濾波、時空域濾波、中值濾波、奇異值分解法等[2-5]，但這些方法大都是利用隨機噪聲的空間不相干特點對其進行處理，壓制效果往往不佳。小波分析在時頻兩域都具有表征信號局部特征的能力，具有自適應符合低頻信號變化緩慢而高頻信號變化迅速的特點[6]。因此，不少學者將小波變換引入地球物理數(shù)據(jù)的噪聲壓制中，劉秀娟等[7]將一維小波變換用于探地雷達數(shù)據(jù)的噪聲壓制中，并研究了小波基函數(shù)，和分解層數(shù)對壓制效果的影響；陳臻懿等[8]研究了自適應分層閾值法對雷達信號噪聲壓制效果，但其研究的信號是一維的；王超等[9]研究了小波變換在探地雷達弱信號去噪中的應用，并提出通過對比不同尺度下小波系數(shù)的疊加結(jié)果來選取小波基函數(shù)的方法；楊立強等[10]將二維小波變換應用于地震數(shù)據(jù)的隨機噪聲壓制中，并驗證了其實用性；Ray Abma等[11]進行了利用復小波變換壓制震源干擾的研究，并達到了在不影響帶寬的條件下提高信噪比以及在不影響信噪比的前提下保留低頻和近偏移距信號的目的；Julio Cesar S.O.Lyrio[12]提出一維自適應小波去噪方法，并將其成功應用于重力梯度數(shù)據(jù)的噪聲壓制過程中。筆者基于二維小波變換，采用自適應分層閾值法對含一定量高斯白噪聲的探地雷達正演模擬數(shù)據(jù)進行噪聲壓制，驗證了該方法的有效性；對比分析了三種不同去噪方法對探地雷達隨機噪聲的壓制效果。其中，中值濾波去噪效果最差，自適應分層閾值法去噪效果最好，全局閾值法去噪效果介于二者之間。最后，將自適應分層閾值法應用于探地雷達實測數(shù)據(jù)處理中，應用效果良好。

1 方法原理

1.1 二維小波分解基本理論[13-14]

連續(xù)二維小波變換公式為式(1)。

(1)

其逆變換表示為：

(2)

其中：ψa,bx,by(x,y)為二維小波基；a為尺度函數(shù)；bx、by分別為x、y方向的平移函數(shù)；f(x,y)為一二維函數(shù)。

為方便計算將a、b離散，則得到二維離散小波序列：

(3)

其中:k=1、2、3，分別表示可以在不同尺度和不同方向情況下提取二維數(shù)據(jù)的不同成分;ψ1(x,y)主要提取水平低頻，垂直高頻分量；ψ2(x,y)主要提取水平高頻，垂直低頻分量；ψ3(x,y)主要提取水平高頻，垂直高頻分量。

1.2 隨機噪聲壓制方法

將探地雷達數(shù)據(jù)看作離散二維含噪信號[15]，則有：

S(i,j)=E(i,j)+N(i,j)

(i=1,2…N,j=1,2…M)

(4)

其中:S(i,j)、E(i,j)、N(i,j)分別為含噪信號、有效信號、噪聲信號；i、j分別為時間采樣點和空間采樣點坐標。

實際情況下，有效信號常以中、低頻成分為主，而隨機噪聲頻帶寬且通常以高頻成分為主[10]。因此，可將通過二維小波分解將信號轉(zhuǎn)化到二維小波域[16]，并根據(jù)隨機噪聲和有效信號在二維小波域中的差異對小波系數(shù)進行處理，其主要有三種方式：①強制去噪，即將所有高頻系數(shù)置零，僅利用低頻系數(shù)進行小波重構(gòu)，該方式最大的一個缺點就是可能將有效信號視為噪聲，從而導致有效信號丟失；②默認閾值去噪，也叫全局閾值去噪，其處理方式是含噪聲信號采用統(tǒng)一的一個閾值進行去噪處理，這種方式對噪聲的壓制往往不夠徹底；③獨立閾值去噪，即在去噪過程中根據(jù)需要在每層選取不同的閾值進行去噪。該方法不僅能夠較為徹底的壓制噪聲，而且可最大程度的保留有效信號[15]。

筆者采用上述第三種去噪方式，具體的方法如下：

1)選取適當?shù)男〔▽π盘栠M行n(n為正整數(shù))層二維小波分解。

2)將每層視為相互獨立，對每層分解計算出相應的閾值，并采用軟閾值函數(shù)[17]對該層高頻系數(shù)進行處理。

3)利用低頻系數(shù)和閾值量化處理后的高頻系數(shù)進行二維小波重構(gòu)，獲得去噪重構(gòu)信號。

2 模型驗證

此次研究采用Daubechies(dbN)小波(N=1…10),當N=0時為haar小波。在給定消失矩階數(shù)時，Daubechies小波是所有正交小波中具有最短支集的小波。小波函數(shù)ψ的消失矩為N，尺度函數(shù)φ和小波函數(shù)ψ的有效支撐長度為2N-1[18]。圖1為db9的小波函數(shù)和尺度函數(shù)圖像。

2.1 斷層模型

圖2(a)、圖2(b)分別為一斷層模型及其探地雷達正演數(shù)據(jù)，正演計算區(qū)域大小為30 m×10 m，模型背景介質(zhì)為花崗巖，斷層內(nèi)含淡水，斷層寬度為1 m，中心點坐標為(5 m,15 m)，傾角為45°，正演計算采用的中心頻率為250 MHz，時窗長度為500 ns，共400道，每道2 000個采樣點。

圖2(c) 為加入一定量高斯白噪聲后得到的含噪聲記錄，其信噪比為1.0；圖2(d)、圖2(e)、圖2(f)分別為中值濾波法去噪、采用db9小波對數(shù)據(jù)進行3層分解后的全局閾值法去噪和自適應分層閾值法去噪后的結(jié)果，其信噪比分別為8.7、16.0、17.5。

由圖2(d)可以看出，中值濾波法對信噪比的改善效果較差，全局閾值法和自適應分層閾值法的去噪效果明顯。但是對比圖2 (e) 、 圖2(f)中的信號可以看出，自適應分層閾值法去噪后的信號對原始信號的細節(jié)保留較全局閾值法更加完整。

圖1 db9小波圖像Fig.1 db9 wavelet image(a)小波函數(shù)；(b)尺度函數(shù)

圖2 含水斷層模型Fig.2 Aquifer fault model(a)斷層模型；(b)正演數(shù)據(jù)；(c)含噪數(shù)據(jù)；(d)中值濾波去噪；(e)全局閾值去噪；(f)自適應分層閾值去噪

2.2 空洞模型

圖3(a)、圖3 (b)分別為一正方形含水空洞模型及其探地雷達正演數(shù)據(jù)，模型邊長為1 m，中心點坐標為(5 m，15 m)，其余參數(shù)與斷層模型相同。

圖4(a) 為加入一定量高斯白噪聲后得到的含噪聲記錄，其信噪比為1.0；圖4(b)、圖4 (c)、圖4 (d)分別為中值濾波法去噪、采用db9小波對正演數(shù)據(jù)進行3層分解后的全局閾值法去噪和自適應分層閾值法去噪后的結(jié)果，其信噪比分別為8.7、15.1、18.4。對比中值濾波法與全局閾值法去噪發(fā)現(xiàn)，采用自適應分層閾值法的去噪效果較另外兩種方法更加明顯。

圖5是信噪比等于0.1時的結(jié)果，中值濾波法、采用db9小波對正演數(shù)據(jù)進行3層分解后的全局閾值法和自適應分層閾值法去噪后的信噪比分別為7.8、14.4、17.3。從圖5(d)中可以看出，自適應分層閾值法去噪后的信號在信噪比的改善和原始信號特征保留方面效果更佳。

圖6是信噪比等于10時的去噪結(jié)果，采用中值濾波法、db9小波對正演數(shù)據(jù)進行3層分解后的全局閾值法和自適應分層閾值法去噪后的信噪比分別為17.4、22.5、25.3。對比分析圖6(c)、6(d)中的弱反射信號發(fā)現(xiàn)，采用自適應分層閾值法去噪后的數(shù)據(jù)對弱信號的特征保留更加完整。

圖3 含水空洞模型及正演數(shù)據(jù)Fig.3 Water-bearing hollow model and forward data(a)含水空洞模型；(b)正演數(shù)據(jù)

圖4 信噪比為1.0時的去噪結(jié)果(含水空洞模型)Fig.4 Denoising results(SNR=1.0)(a)含噪數(shù)據(jù)；(b)中值濾波去噪；(c)全局閾值去噪；(d)自適應分層閾值去噪

圖5 信噪比為0.1時的去噪結(jié)果(含水空洞模型)Fig.5 Denoising results(SNR=0.1)(a)含噪信號；(b)中值濾波去噪；(c)全局閾值去噪;(d)自適應分層閾值去噪

通過含水斷層和含水空洞兩個模型的含噪聲正演模擬數(shù)據(jù)的驗證可知，自適應分層閾值法壓制隨機噪聲適用于不同的模型；通過與中值濾波法和全局閾值濾波法對隨機噪聲的壓制效果的對比，可以知道自適應分層閾值法對隨機噪聲的壓制更為徹底，對信噪比的改善效果更好，同時在對信號的細節(jié)保留方面較全局閾值法更加完整；三種不同信噪比下數(shù)據(jù)的去噪結(jié)果表明自適應分層閾值法對低、中、高信噪比數(shù)據(jù)均適用。

圖6 信噪比為10時的去噪結(jié)果(含水空洞模型)Fig.6 Denoising results(SNR=10)(a)含噪信號；(b)中值濾波去噪；(c)全局閾值去噪；(d)自適應分層閾值去噪

圖7 實際數(shù)據(jù)驗證Fig.7 Actual data validation(a)預處理數(shù)據(jù)；(b)中值濾波去噪；(c)全局閾值去噪；(d)自適應分層閾值去噪

3 實際數(shù)據(jù)驗證

圖7(a)為一經(jīng)過預處理的實測探地雷達剖面，可以看出隨機噪聲存在于整個數(shù)據(jù)剖面；圖7(b)為采用中值濾波法去噪后的結(jié)果，從中可以看出，中值濾波法對數(shù)據(jù)質(zhì)量有一定程度地改善，但剖面中仍存在大量隨機噪聲。

圖7(c)為采用db9小波進行3層小波分解后的全局閾值法去噪后的結(jié)果。對比中值濾波法去噪，可知全局閾值法去噪對數(shù)據(jù)質(zhì)量有更為明顯的改善，但在深度約12 m到18 m區(qū)間仍然存在少量隨機噪聲。

圖7(d)為自適應分層閾值法去噪后的剖面。對比中值濾波法和全局閾值濾波法的結(jié)果發(fā)現(xiàn)，自適應閾值分層法處理過后的剖面對隨機噪聲的壓制更為徹底，信噪比明顯提高，細節(jié)顯示更加清晰，利于地質(zhì)解譯。深度約12 m到15 m區(qū)間可以明顯看出反射波同相軸錯斷、雜亂，且能量較強 ,推測為含水較為嚴重的破碎帶，最后實際開挖結(jié)果驗證了這一推斷。

4 結(jié)論

通過不同模型、不同信噪比下的正演數(shù)據(jù)和實際數(shù)據(jù)的驗證以及與常規(guī)隨機噪聲壓制方法的比較分析，得出以下結(jié)論：

1)相對于常規(guī)濾波方法，小波變換能同時在時頻兩域?qū)π盘栠M行分析，從而可以更加有效地區(qū)分有效信號與噪聲,達到更好的去噪效果。

2)采用基于二維小波變換的自適應分層閾值法對隨機噪聲的壓制更為徹底，同時也能較好地保留信號的原始特征，且對不同的模型數(shù)據(jù)以及不同信噪比下的資料均適用，說明其適用范圍廣，實用性強。

3)在選取閾值對小波系數(shù)進行處理時需要充分考慮到噪聲與有效信號的細節(jié)差異，避免有效信號丟失。