丁 丹，吳 雙，昌彥君，呂淑娟，張 偉

(1.安徽省地質(zhì)礦產(chǎn)勘查局 321地質(zhì)隊(duì),安徽 銅陵 244033;2.中國(guó)地質(zhì)大學(xué) 地球物理與空間信息學(xué)院,湖北 武漢 430074;3.貴州水利水電職業(yè)技術(shù)學(xué)院 土木工程分院,貴州 清鎮(zhèn) 551416;4.武漢捷探科技有限公司,湖北 武漢 430011)

1 引言

電法勘探作為地球物理勘探方法之一,具有原理簡(jiǎn)單,勘探結(jié)果準(zhǔn)確、直觀等優(yōu)點(diǎn)[1]。近年來,國(guó)內(nèi)外實(shí)驗(yàn)及實(shí)踐研究表明,水域電法具有較好的勘探效果,水的導(dǎo)電率能夠傳導(dǎo)較大的電流,對(duì)于解決水域范圍內(nèi)工程地質(zhì)問題[2,3]、探明水底資源分布[4]及水下構(gòu)造勘探[5]等領(lǐng)域具有較好的勘探效果。進(jìn)行水域電法勘探,需要通過電極向水中注入電流,通過接收水中電流引起的電場(chǎng)變化來進(jìn)行分析研究,該方法具有實(shí)現(xiàn)橫向與縱向分辨率高,電極排列布置與移動(dòng)簡(jiǎn)單,施工效率高且勘探速度快等優(yōu)點(diǎn)。但在實(shí)際測(cè)量過程中,儀器供電往往使用50Hz交流電進(jìn)行整流后輸出,因此會(huì)導(dǎo)致采集到的信號(hào)耦合工頻噪聲,對(duì)探勘結(jié)果有較大的影響[6,7]。由于工頻噪聲無處不在且噪聲強(qiáng)度大,在進(jìn)行微弱信號(hào)采集中往往會(huì)淹沒有效信號(hào)[8]。對(duì)工頻噪聲的特點(diǎn)進(jìn)行分析,可以發(fā)現(xiàn)其能量主要集中在基頻及其奇次諧波上[9,10]。為了抑制工頻干擾,在電法勘探領(lǐng)域,對(duì)工頻噪聲的抑制有多種方法,如自適應(yīng)濾波[11]、噪聲估計(jì)法[12,13],這些工頻噪聲抑制方法都需要明確的工頻諧波干擾參數(shù),以獲得可直接消除的噪聲干擾,并且算法實(shí)現(xiàn)復(fù)雜,可能對(duì)水域電法勘探的實(shí)時(shí)性造成影響。

為了解決水域電法勘探存在工頻噪聲干擾的問題,本文結(jié)合諧波濾波器和小波去噪方法對(duì)水域電法采集過程中存在的工頻干擾進(jìn)行處理,提高探測(cè)數(shù)據(jù)的質(zhì)量。之后,通過水槽模擬實(shí)驗(yàn)驗(yàn)證諧波濾波器及小波去噪在水域電法數(shù)據(jù)預(yù)處理中的可行性和有效性。

2 原理

2.1 諧波濾波器

諧波濾波器除了結(jié)構(gòu)簡(jiǎn)單、適應(yīng)性強(qiáng)以外,同時(shí)通過對(duì)參數(shù)的調(diào)整能夠使工頻噪聲及其諧波有較大的衰減,同時(shí)有一定的帶寬能夠滿足諧波干擾的基頻擾動(dòng)。本文的諧波濾波器的設(shè)計(jì)方法是基于梳狀濾波器的原理。梳狀濾波器是指按照一定頻率間隔排列的通帶或者阻帶,只能允許特定頻率的信號(hào)通過,而梳狀濾波器一般可以視為不同中心頻率窄帶帶通濾波器的疊加,根據(jù)文獻(xiàn)[14]中提出的窄帶帶通濾波器傳遞函數(shù):

(1)

式(1)中,N為濾波器的極點(diǎn);j為虛數(shù)單位;ω0=2πf0為濾波器的截止頻率;f0為工頻噪聲的中心頻率,在研究工頻噪聲時(shí)取值為50 Hz。

由于工頻干擾及其諧波具有倍頻關(guān)系,根據(jù)式(1)使得工頻干擾奇次諧波中心頻率fa均設(shè)置為濾波器零點(diǎn),通過零極點(diǎn)的設(shè)計(jì)將多個(gè)帶通濾波器互相疊加,有

fa=(2i-1)ω0

(2)

式(2)中,i=1,2,...,n/2;n=f0/fs,fs為采樣頻率。

通過疊加后的梳狀濾波器傳遞函數(shù)為式(3),即

(3)

將零點(diǎn)代入式(3)可以得到:

1+z-N=0

(4)

令

z=e-j2πn

(5)

有

1+e-j2πNn=0

(6)

e-j2πNn=-1

(7)

2πNn=π

(8)

所以有

圖1 小波去噪流程Fig.1 Flow chart of wavelet denoising

(9)

同時(shí),可以疊加窗函數(shù)對(duì)通帶和阻帶的紋波壓制和衰減進(jìn)行控制。常用的窗函數(shù)包括漢寧窗、海明窗和高斯窗。因?yàn)楦咚勾澳軌蛘{(diào)節(jié)帶寬參數(shù)實(shí)現(xiàn)衰減頻率的調(diào)節(jié),因此被廣泛使用[15]。

加入經(jīng)過歸一化的高斯窗參數(shù)a,取值范圍為(0,1),該參數(shù)能夠改善濾波器的性能,取值增大的同時(shí),濾波器通帶越平緩,對(duì)其他信號(hào)影響較小,但是衰減減小,因此高斯窗參數(shù)a的取值也要不同問題具體分析。通過引入高斯窗函數(shù)后,得到的梳狀帶通濾波器傳遞函數(shù)為,

(10)

根據(jù)式(11),使用全通濾波器減去梳妝帶通濾波器得到諧波濾波器:

(11)

該濾波器的性能取決于高斯窗中高斯窗參數(shù)a的取值。

2.2 小波去噪

小波變換自提出就應(yīng)用于噪聲信號(hào)處理中,2011年,小波變換方法應(yīng)用在高密度電法資料去噪研究中[16]。但是小波變化需要選擇合適小波基函數(shù),不同的小波基函數(shù)有著不同的分離噪聲的效果,如果母小波選擇不正確,容易導(dǎo)致信號(hào)發(fā)生畸變。

利用小波變換的水域直流電探測(cè)信號(hào)降噪包括三個(gè)步驟:首先,提取出有用的水域直流電探測(cè)信號(hào)特征;然后,在小波域進(jìn)行濾波;最后,進(jìn)行水域直流電探測(cè)信號(hào)重構(gòu)就可以達(dá)到消噪的目的(圖1)。直流電測(cè)深測(cè)量得到的電信號(hào)Ua為

(12)

水域直流電探測(cè)信號(hào)特征提取的核心是在小波域中對(duì)探測(cè)信號(hào)進(jìn)行分解。直流電探測(cè)信號(hào)通過選擇一個(gè)小波基函數(shù)并給定分解的層次,對(duì)含噪電信號(hào)進(jìn)行小波變換分解。每一層分解的結(jié)果是上次分解得到的低頻信號(hào)再分解成低頻和高頻兩個(gè)部分。

在小波分解過程中,將某一小波基函數(shù)做位移后,通過式(13)計(jì)算不同尺度下的基函數(shù)與原始水域直流電探測(cè)信號(hào)的內(nèi)積,這樣信號(hào)就被分解成一系列信號(hào)的疊加。

(13)

對(duì)于同樣的探測(cè)信號(hào),不同的基函數(shù)會(huì)得到不同的結(jié)果[17]。目前常用的基函數(shù)可以分成三類:Daubechies小波基,Symlets小波基和Coiflet小波基。Daubechies小波基具有較好的正則性,計(jì)算量小且去噪效果好,可以用于暫態(tài)噪聲和瞬態(tài)噪聲去除,但是隨著小波階數(shù)增加,Daubechies小波基計(jì)算量會(huì)增加,同時(shí)實(shí)時(shí)性變差;Symlets小波基是從Daubechies小波基演化而來,在連續(xù)性、濾波器階數(shù)等方面與Daubechies小波基一致,但Symlets小波基有更好的對(duì)稱性,能夠降低在信號(hào)重構(gòu)時(shí)對(duì)信號(hào)的失真;Coiflet小波基雖然是從Daubechies小波基演化而來,但犧牲了支撐長(zhǎng)度來?yè)Q取一定的對(duì)稱性,一般不適用于實(shí)時(shí)信號(hào)處理[18,19]。

經(jīng)過I層分解后,水域直流電探測(cè)信號(hào)Ua被分解為

(14)

式中,j=1, 2, ...,I;Yj為每層分解得到的估計(jì)小波系數(shù);Mj(k)為第j層分解得到的低頻水域電法探測(cè)信號(hào),以測(cè)量信號(hào)為主(即有效電信號(hào));Nj(k)為第j層分解得到的高頻水域電法探測(cè)信號(hào),以噪聲為主。

隨后,對(duì)分解得到估計(jì)小波系數(shù)進(jìn)行閾值處理,利用硬閾值法來得到小波系數(shù),通過一個(gè)濾波閾值實(shí)現(xiàn)噪聲和有效信號(hào)的分離。閾值的選擇將對(duì)去噪效果的好壞產(chǎn)生決定性的作用。人們依據(jù)水域電法電探測(cè)數(shù)據(jù)采樣點(diǎn)數(shù)進(jìn)行閾值選取,即長(zhǎng)度對(duì)數(shù)閾值(式15)。

(15)

式中,λ為閾值;σ為首層小波分解系數(shù)絕對(duì)值的中間值。

硬閾值法(式16)可以盡可能地保留信號(hào)邊緣等局部特征,使重構(gòu)的信號(hào)具有更好的逼近性。若估計(jì)小波系數(shù)Yj小于λ,認(rèn)為該系數(shù)由背景噪聲引起,將其置零;若估計(jì)小波系數(shù)Yj大于等于λ,認(rèn)為該系數(shù)由有效信號(hào)引起,保留該系數(shù),并用于小波逆變換[20]。

(16)

在小波域完成電法探測(cè)信號(hào)與水域背景噪聲的濾波分離后,就可以將分解后的信號(hào)進(jìn)行重構(gòu),將新產(chǎn)生的Yj重新代入公式(14),所得到的Ua即為去噪后的有效信號(hào)。

將信號(hào)進(jìn)行諧波濾波器和小波去噪處理后,利用公式(17)來計(jì)算信噪比評(píng)價(jià)去噪信號(hào)的優(yōu)劣[21],

(17)

式中,x(n)為原始信號(hào);y(n)為去噪后信號(hào)。

3 實(shí)際應(yīng)用

3.1 實(shí)驗(yàn)平臺(tái)

圖2 勘探模擬實(shí)驗(yàn)平臺(tái)示意圖Fig.2 Schematic diagram of exploration simulation experiment platform

為了更好地驗(yàn)證水域電法勘探效果,搭建了一個(gè)物理模擬實(shí)驗(yàn)平臺(tái)如圖2所示,用于后續(xù)實(shí)驗(yàn)的展開。水域電法勘探物理模擬實(shí)驗(yàn)平臺(tái)能夠模擬海洋環(huán)境下拖曳式直流電法實(shí)際勘探情況。實(shí)驗(yàn)平臺(tái)主要由四大部分組成,分別為上位機(jī)(LabVIEW編程)、四軸三維運(yùn)動(dòng)系統(tǒng)、電源模塊和電壓采集模塊,如圖2所示,為模擬勘探平臺(tái)框架圖。實(shí)驗(yàn)用長(zhǎng)方體水槽規(guī)格為2 m(長(zhǎng))×1 m(寬)×0.5 m(高)。為了直觀展示測(cè)量結(jié)果,實(shí)驗(yàn)采用對(duì)稱四極裝置對(duì)石墨塊(約10 cm)進(jìn)行測(cè)量。

為了更加直觀地展示去噪效果,利用圖2所示平臺(tái),獲取一組水域直流電測(cè)量數(shù)據(jù),并進(jìn)行頻譜分析,得到圖3的電壓曲線,對(duì)其進(jìn)行頻譜分析,結(jié)果如圖4所示。可以發(fā)現(xiàn),直接采集的信號(hào)噪點(diǎn)比較多,當(dāng)電極運(yùn)動(dòng)到在石墨塊上方時(shí),電壓出現(xiàn)下降,但下降幅度為0.6 V,噪聲幅度為0.3 V,由于噪聲的耦合導(dǎo)致后半段變化不明顯,利用式(17)可以求得信號(hào)信噪比為25.514 9。分析頻譜圖可以發(fā)現(xiàn),噪聲信號(hào)在50 Hz附近有一個(gè)較大的峰值,同時(shí)100 Hz附近有較小的波峰,是由工頻噪聲及其諧波導(dǎo)致的。因此如果在實(shí)際勘探中直接采集直流電壓信號(hào),信號(hào)過于微弱,有可能導(dǎo)致異常點(diǎn)淹沒在噪聲中,因此濾波和去噪是有必要的。

圖3 模擬實(shí)驗(yàn)直接采集電壓信號(hào)Fig.3 Direct acquisition of voltage signal in simulation experiment

圖4 模擬實(shí)驗(yàn)直接采集電壓信號(hào)頻譜Fig.4 Spectrum of voltage signal directly acquired in simulation experiment

3.2 高斯窗參數(shù)

為了明確不同的高斯窗參數(shù)a對(duì)諧波濾波器去噪效果的影響,分別選取了a為0.1,0.5和0.8時(shí)的諧波濾波器幅頻曲線。得到結(jié)果如圖5。

圖5 不同a值的去噪效果Fig. 5 Noise removal effect of different a values

通過對(duì)比發(fā)現(xiàn),當(dāng)a增大的時(shí)候,諧波濾波器頻率響應(yīng)曲線會(huì)變得平坦,對(duì)其他頻率信號(hào)影響較小,但濾波效果會(huì)變差,衰減會(huì)減小,當(dāng)a=1時(shí)會(huì)變成全通濾波器。當(dāng)a減小的時(shí)候,可以發(fā)現(xiàn)諧波濾波器的頻率響應(yīng)變得凸出,但是信號(hào)在中心頻率處衰減迅速,濾波效果較好,同時(shí)對(duì)其他頻率的信號(hào)影響較大。因此,對(duì)高斯窗參數(shù)a的選擇為0.5比較合適。對(duì)設(shè)計(jì)好的濾波器進(jìn)行實(shí)際信號(hào)的濾波,得到如圖6的去噪信號(hào),其頻譜如圖7所示,對(duì)比結(jié)果可以發(fā)現(xiàn)該濾波器對(duì)工頻噪聲及其諧波有一定的抑制作用。

圖6 諧波濾波器去噪信號(hào)Fig.6 Harmonic filter denoising signal

圖7 諧波濾波器去噪信號(hào)頻譜Fig.7 Spectrum of harmonic filter denoising signal

3.3 小波基選取

利用圖2所示的測(cè)試平臺(tái),獲取了一組水槽直流電探測(cè)數(shù)據(jù)。將Coiflet小波基、Symlet小波基和Daubechies小波基分別進(jìn)行噪聲處理測(cè)試。對(duì)比結(jié)果(圖8),為了更清晰展示細(xì)節(jié),將圖8紅框中進(jìn)行放大得到圖9,結(jié)果表明Coiflet小波基的效果是不滿意的,無法有效去噪;Daubechies小波基的效果比Coiflet小波基的有所改善,但是仍然存在明顯的噪聲波動(dòng)且信號(hào)不夠平滑;Symlets小波基有較好的效果,噪聲波動(dòng)被明顯減少,同時(shí)水槽實(shí)驗(yàn)中采集的直流電壓信號(hào)變得平滑。因此,這個(gè)結(jié)果驗(yàn)證了Symlets小波基是比較適合水域直流電探測(cè)。通過計(jì)算,相較于直接采集的信號(hào),經(jīng)過處理后的信號(hào)信噪比為33.587 1。

圖8 不同小波基去噪效果Fig.8 Effect diagram of noise removal with different wavelet bases

圖9 不同小波基去噪效果圖局部細(xì)節(jié)Fig.9 Partial details of different wavelet basis denoising effect

4 探測(cè)結(jié)果

在水槽中,使用自動(dòng)跑極裝置進(jìn)行三次測(cè)試,電極排列裝置為對(duì)稱四極裝置,供電電極AB間隔為10 cm,測(cè)量電極MN間隔為2 cm。第一次采集不放置石墨塊,用于采集水槽無目標(biāo)狀態(tài)下的探測(cè)效果;第二次實(shí)驗(yàn)將石墨塊用繩子吊起來,放置于水槽中間位置;第三次實(shí)驗(yàn)將石墨塊放置于水槽底部,在不進(jìn)行任何去噪處理的情況下,得到圖10的采集電壓數(shù)據(jù),此時(shí),三次采集信噪比分別為23.782 7、25.514 9和24.623 1,可以發(fā)現(xiàn):采集電壓信號(hào)有較大的波動(dòng),尤其當(dāng)探測(cè)目標(biāo)放置在水槽底部時(shí),采集電壓較小,采集的信號(hào)被噪聲淹沒,對(duì)后續(xù)視電阻率的計(jì)算存在困難。對(duì)其進(jìn)行諧波濾波器和小波變換去噪后得到結(jié)果如圖11,可以發(fā)現(xiàn)信號(hào)變化更加明顯,三次采集信噪比提高到了30.365 2、33.587 1和31.745 3,這將有利于視電阻率的計(jì)算。

對(duì)水槽模擬實(shí)驗(yàn)采集的電壓數(shù)據(jù)進(jìn)行諧波濾波器和小波去噪處理后,利用式(18)計(jì)算視電阻率:

(18)

圖11 去噪處理的探測(cè)數(shù)據(jù)Fig.11 Detected data after de-noise processing

最后得到結(jié)果如圖12所示,通過對(duì)比實(shí)驗(yàn)分析得到:由于水槽使用淡水,在無目標(biāo)的情況下測(cè)得水的電阻率為21.8 Ω·m,符合理論實(shí)際。對(duì)比石墨塊放置于半水深和水底測(cè)得的視電阻率,可以發(fā)現(xiàn)當(dāng)電極陣列支架越靠近測(cè)量被測(cè)石墨塊時(shí),視電阻率變化越大,勘探效果越好。

圖12 勘探實(shí)驗(yàn)視電阻率計(jì)算結(jié)果Fig.12 Calculation results of apparent resistivity of exploration experiment

5 結(jié)論

針對(duì)目前水域電法勘探存在工頻噪聲及其諧波干擾導(dǎo)致采集數(shù)據(jù)精度不高的問題,提出了一種基于諧波濾波器和小波去噪的方法來對(duì)采樣數(shù)據(jù)進(jìn)行去噪處理,結(jié)合水域電法原理進(jìn)行視電阻率計(jì)算。首先確定了諧波濾波器中高斯窗因子a的取值為0.5時(shí)以及小波去噪中小波基采用Symlets小波基時(shí),有較好的去噪效果。而后基于搭建的實(shí)驗(yàn)平臺(tái)進(jìn)行數(shù)據(jù)采集,對(duì)水槽模擬實(shí)驗(yàn)采集電壓數(shù)據(jù)中的工頻噪聲及其諧波進(jìn)行抑制,使用小波方法對(duì)采集的電壓信號(hào)再次去噪,在水槽模擬實(shí)驗(yàn)中,對(duì)于異常最大的采集信號(hào)的信噪比從25.514 9提高到了33.587 1,而無異常情況下信噪比從23.782 7提高到了30.365 2。最后,計(jì)算得到無目標(biāo)時(shí),水底視電阻率計(jì)算結(jié)果為21.8 Ω·m,放置石墨塊時(shí),上方的水下視電阻率降為10.1 Ω·m,符合預(yù)期結(jié)果,驗(yàn)證了該去噪方法的可行性,具有深入研究的價(jià)值。