基于二維小波變換的隨機(jī)噪聲壓制方法在GPR數(shù)據(jù)中的應(yīng)用

魯光銀， 羅 帥， 朱自強(qiáng)， 趙云威， 席飛雁

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

0 引言

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

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

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ù)。

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

(3)

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

1.2 隨機(jī)噪聲壓制方法

將探地雷達(dá)數(shù)據(jù)看作離散二維含噪信號(hào)[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)分別為含噪信號(hào)、有效信號(hào)、噪聲信號(hào)；i、j分別為時(shí)間采樣點(diǎn)和空間采樣點(diǎn)坐標(biāo)。

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

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

1)選取適當(dāng)?shù)男〔▽?duì)信號(hào)進(jìn)行n(n為正整數(shù))層二維小波分解。

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

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

2 模型驗(yàn)證

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

2.1 斷層模型

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

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

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

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

圖2 含水?dāng)鄬幽Ｐ虵ig.2 Aquifer fault model(a)斷層模型；(b)正演數(shù)據(jù)；(c)含噪數(shù)據(jù)；(d)中值濾波去噪；(e)全局閾值去噪；(f)自適應(yīng)分層閾值去噪

2.2 空洞模型

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

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

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

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

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

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

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

通過含水?dāng)鄬雍秃斩磧蓚€(gè)模型的含噪聲正演模擬數(shù)據(jù)的驗(yàn)證可知，自適應(yīng)分層閾值法壓制隨機(jī)噪聲適用于不同的模型；通過與中值濾波法和全局閾值濾波法對(duì)隨機(jī)噪聲的壓制效果的對(duì)比，可以知道自適應(yīng)分層閾值法對(duì)隨機(jī)噪聲的壓制更為徹底，對(duì)信噪比的改善效果更好，同時(shí)在對(duì)信號(hào)的細(xì)節(jié)保留方面較全局閾值法更加完整；三種不同信噪比下數(shù)據(jù)的去噪結(jié)果表明自適應(yīng)分層閾值法對(duì)低、中、高信噪比數(shù)據(jù)均適用。

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

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

3 實(shí)際數(shù)據(jù)驗(yàn)證

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

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

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

4 結(jié)論

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

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

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

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