汪偉明　賀巍

摘要：本文擬應(yīng)用MATLAB軟件對容易受到外界環(huán)境干擾和由于儀器自身干擾的影響，對地磁信號進行降噪處理。通過小波分析和數(shù)學(xué)形態(tài)學(xué)融合算法，同時將該方法與數(shù)學(xué)形態(tài)學(xué)算法兩者進行對比分析，并且計算出降噪后信號相對于含噪信號的均方誤差的值，然后對其不同算法的降噪信號的處理結(jié)果與未受干擾的原始正常信號進行對比分析，分別求其信噪比的大小，綜合考慮兩種算法對地磁含噪信號處理的優(yōu)勢和劣勢。進而得出小波分析與數(shù)學(xué)形態(tài)學(xué)算法對地磁信號的降噪處理，能最大化的減少信號的失真程度，保留其原始形態(tài)。為以后批量處理地磁觀測的干擾數(shù)據(jù)提供新方法和新思路，同時也能抑制對地磁信號的干擾和提高觀測數(shù)據(jù)的內(nèi)在質(zhì)量。

關(guān)鍵詞：MATLAB軟件；降噪處理；小波分析；數(shù)學(xué)形態(tài)學(xué)；均方誤差；信噪比

中圖分類號：TN911.4 文獻標(biāo)識碼：A 文章編號：1007-9416（2018）02-0122-03

隨著數(shù)學(xué)信號處理技術(shù)的不斷提升，基于MATLAB軟件對地磁信號降噪處理的方法也越來越多，截止目前為止，降噪處理的方法主要有：基于fft硬件濾波法，信號相關(guān)法，自適應(yīng)濾波法，圖譜識別法，有限沖激響應(yīng)濾波法和小波變換法等等。對于先前已有的降噪方法，均有明顯的缺點，主要表現(xiàn)在對降噪信號處理不徹底，或者波形容易失真等。所以開始對地磁信號研究新的降噪方法非常重要。數(shù)學(xué)形態(tài)學(xué)算法其主要應(yīng)用于算法開發(fā)、數(shù)據(jù)可視化、數(shù)據(jù)分析以及數(shù)值計算的高級技術(shù)計算語言和交互式環(huán)境[1]。該算法將一個含噪的信號分解成具有物理意義的各個部分，然后將其背景剝離，從而有效抑制含噪信號的毛刺和縫隙，保留原始信號的主要形態(tài)。應(yīng)用數(shù)學(xué)形態(tài)學(xué)的算法處理含噪信號時，不容易發(fā)現(xiàn)本身微弱信號的尖峰和毛刺，同時由于含噪信號因其有振蕩衰減的特點，所以在含噪信號特別微弱的時候，數(shù)學(xué)形態(tài)學(xué)算法很難在形態(tài)上區(qū)分是振蕩衰減的信號還是信號的尖峰和毛刺。因此單純應(yīng)用數(shù)學(xué)形態(tài)學(xué)算法對地磁信號的降噪處理無法取得理想效果。小波分析對地磁信號的研究較為廣泛，因其本身具有可以調(diào)整時間窗和頻率窗的特點，可以精確找到某一段地磁信號的突變點，對地磁信號的降噪算法研究具有深刻意義。通過小波分析，進行門限閾值量化處理后進行小波重構(gòu)，進而抑制信號的噪聲。結(jié)合小波分析和數(shù)學(xué)形態(tài)學(xué)兩者算法的優(yōu)點，將其兩者融合，來對地磁信號進行降噪處理，進而實現(xiàn)抑制地磁信號的干擾信息，最大化，最有效的恢復(fù)原始的地磁真實信號。

本文利用小波分析與數(shù)學(xué)形態(tài)學(xué)融合算法對地磁的高頻干擾信息進行降噪處理，同時將該方法與數(shù)學(xué)形態(tài)學(xué)算法進行對比分析，綜合分析兩者的濾波效果，達到消除地磁信號因外界環(huán)境受到的高頻干擾信息，兩者通過定量對比分析其均方誤差值和信噪比的大小，分析其優(yōu)勢和劣勢，最大化的提取真實地磁信號的初步分析和探討。

1 基于數(shù)學(xué)形態(tài)學(xué)算法的混合濾波器設(shè)計

1.1 數(shù)學(xué)形態(tài)學(xué)算法介紹

數(shù)學(xué)形態(tài)學(xué)是積分幾何研究成果的基礎(chǔ)上創(chuàng)立的，是基于集合論的數(shù)學(xué)分支，它提供了非常有效的非線性濾波技術(shù)[2]。數(shù)學(xué)形態(tài)學(xué)的運算主要包括腐蝕Θ，膨脹，形態(tài)開OC和形態(tài)閉CO的級聯(lián)組合。對于地磁信號而言，以GM4儀器為例，記錄到的值屬于一維的離散信號f（j）。其主要包括地磁真實信號x（j）和含噪信號y（j），其表達式為：

所以對于地磁信號所有的高頻干擾含噪信號而言，其主要就是平滑和抑制峰值和谷底的噪聲，由此引出來開-閉OC濾波器和閉-開CO濾波器，兩者任意一個濾波器均可以實現(xiàn)對地磁信號高頻干擾的降噪處理。單個的濾波器雖然可以消除高頻干擾信息，但是開-閉OC濾波器和閉-開CO濾波器因分別具有開運算的收縮性和閉運算的擴張性，會使地磁信號的濾波輸出結(jié)果偏大或者偏小，進而影響其濾波效果。因此本文擬采用開-閉OC濾波器閉-開CO濾波器輸出結(jié)果的平均值來組成混合濾波器，通過混合濾波器的設(shè)計，對地磁信號進行降噪處理，可以更好地消除高頻干擾。得到的濾波效果會更接近真實信號，其表達式為：

1.3 結(jié)構(gòu)元素的確定

設(shè)計好混合濾波器的算法后，組成混合濾波器的結(jié)構(gòu)元素直接影響其地磁信號的濾波結(jié)果。常見的結(jié)構(gòu)元素有余弦，直線，三角，半圓等。本文擬采用余弦結(jié)構(gòu)和直線結(jié)構(gòu)分別對高頻干擾信息進行消除。此余弦結(jié)構(gòu)的的形狀由幅值S和寬度W共同決定，而直線結(jié)構(gòu)的形狀僅有寬度W決定。所以應(yīng)用表達式（8），選取不同結(jié)構(gòu)元素的形狀，對于含有地磁行的高頻干擾進行消除，對地磁信號進行降噪處理。分析其不同種類的結(jié)構(gòu)元素，對其濾波的效果。

2 小波分析與數(shù)學(xué)形態(tài)學(xué)算法的融合濾波器設(shè)計

基于數(shù)學(xué)形態(tài)學(xué)算法的混合濾波器設(shè)計后，接著開始設(shè)計小波分析與數(shù)學(xué)形態(tài)學(xué)算法相融合的濾波器，其主要包括以下幾個步驟：

①選擇合適的小波基db和尺度p，確定小波分解層數(shù)N，并將含有高頻信息的地磁信號f（j）分解，分解后提取每一層的小波系數(shù)WT（p1，k1）；

②然后利用已經(jīng)設(shè)計好的混合濾波器對每一層的小波系數(shù)WT（p1，k1）濾波，得到新的小波系數(shù)WT”（pn，kn）；

③確定閾值的λp，對小波系數(shù)WT”（p，k）進行閾值處理。在閾值前乘系數(shù)k（0

λp是第p層小波分解的閾值；mp是第p層的小波分解系數(shù)的中值；np時第p層的小波系數(shù)的個數(shù)。

④重構(gòu)第③步處理后的小波系數(shù)WT”（p，k），通過重構(gòu)小波系數(shù)，還原濾波后的地磁信號。得到的結(jié)果即為小波分析與數(shù)學(xué)形態(tài)學(xué)融合濾波器的濾波結(jié)果。

3 基于兩種不同算法所設(shè)計濾波器濾波效果的衡量標(biāo)準(zhǔn)

對于含噪地磁信號分別通過基于數(shù)學(xué)形態(tài)學(xué)的混合濾波器和小波分析與數(shù)學(xué)形態(tài)學(xué)相融合算法的濾波器，輸出濾波結(jié)果，通過得到濾波后的曲線，計算其地磁信號濾波的均方誤差ESE的值和信噪比SNR的大小，進而綜合分析考慮兩者的濾波效果。其計算公式如下：

4 仿真分析

為驗證兩種算法所設(shè)計的濾波器對地磁高頻干擾信息的濾除效果，擬采用MATLAB軟件進行仿真分析。

4.1 不同結(jié)構(gòu)元素選取對地磁信號降噪處理的仿真分析

如2017年11月1日榆林臺GM4-1儀器記錄到Z分量正常的地磁信號，選取世界時00時00分-16時00分時段的波形曲線，將該曲線加入高斯白噪聲干擾，然后分析余弦結(jié)構(gòu)和直線結(jié)構(gòu)兩種不同結(jié)構(gòu)元素下，所設(shè)計的混合濾波器對榆林臺高斯白噪聲的降噪效果情況。對于兩種不同結(jié)構(gòu)元素，最優(yōu)的S和W組合時的濾波效果如表1可知，由此可知對于不同結(jié)構(gòu)元素，余弦結(jié)構(gòu)降噪后最大信噪比SNR=28.139dB，最小均方誤差ESE=2.523X10-2。而對于直線結(jié)構(gòu)降噪后最大信噪比SNR=10.589dB，最小均方誤差ESE=1.492X10-2。同時由圖1可見對于高斯白噪聲而言，雖然余弦結(jié)構(gòu)元素的混合濾波器能得到最大的信噪比，但是地磁信號嚴(yán)重失真，反而對于直線結(jié)構(gòu)元素，均方誤差值較小，信噪比較大，最大化的還原真實的地磁信號，對高斯白噪聲的降噪效果較好。

4.2 兩種不同算法的濾波器對地磁信號降噪的仿真分析

松江佘山臺FGE[K]儀器記錄到的地磁信號的Z分量的變化，因受到地鐵高頻干擾，地磁信號嚴(yán)重失真。同樣選取2017年11月1日世界時00時00分-16時00分時段的波形曲線，因其榆林臺與松江佘山臺基本在同一緯度上，地磁的日變化曲線基本一致，但存在一定的相位差。如下圖2所示，松江佘山臺在該時間段內(nèi)，Z分量受到高頻干擾，地磁信號曲線嚴(yán)重失真，而榆林臺地磁信號的Z分量未受干擾，能較好的還原地磁信號的真實形態(tài)。所以可以應(yīng)用這兩種不同算法所設(shè)計的濾波器，對該臺的高頻干擾信息降噪處理，降噪處理的結(jié)果同榆林臺同時段的曲線做對比分析，最后通過仿真結(jié)果綜合分析兩種算法所設(shè)計的濾波器對高頻干擾的降噪效果情況。

由表2可知，對于高頻干擾信息的濾除，余弦結(jié)構(gòu)的混合濾波器的信噪比明顯優(yōu)于直線結(jié)構(gòu)的混合濾波器，并且兩種結(jié)構(gòu)的均方誤差比較接近。因此選余弦結(jié)構(gòu)的數(shù)學(xué)形態(tài)學(xué)混合濾波器和小波分析與數(shù)學(xué)形態(tài)學(xué)相融合的濾波器兩者做對比分析，分析其降噪效果情況。

由表3可知，對于兩種不同算法所設(shè)計的濾波器，小波分析與數(shù)學(xué)形態(tài)學(xué)相融合所設(shè)計濾波器的降噪后最大信噪比SNR=250.621dB，最小均方誤差ESE=3.8519X10-2。而對于數(shù)學(xué)形態(tài)學(xué)混合濾波器的降噪后最大信噪比SNR=260.893dB，最小均方誤差ESE=3.8327X10-2。由此可見，小波分析與數(shù)學(xué)形態(tài)學(xué)融合算法的濾波器對高頻干擾信息的處理要優(yōu)于數(shù)學(xué)形態(tài)學(xué)混合濾波器。同時由圖3可見，對于松江佘山臺受高頻干擾的處理結(jié)果而言，小波分析與數(shù)學(xué)形態(tài)學(xué)算法相融合的濾波器明顯最大化的保留地磁的原始形態(tài)，與榆林臺在同一時間段的地磁信號變化曲線基本一致。

5 結(jié)語

通過MATLAB軟件對兩種不同算法所設(shè)計的濾波器進行仿真測試，分析其地磁信號受到高頻干擾的外界環(huán)境干擾時，對其降噪處理，所得到的濾波結(jié)果明顯有差異，經(jīng)分析，得出小波分析與數(shù)學(xué)形態(tài)學(xué)相融合算法的濾波器的明顯優(yōu)于數(shù)學(xué)形態(tài)學(xué)算法的濾波器。小波分析與數(shù)學(xué)形態(tài)學(xué)算法相融合算法設(shè)計的濾波器能保證信號最大化且無失真還原地磁真實信號，其濾波得到的結(jié)果信噪比高，均方誤差值也小。為以后批量處理地磁干擾數(shù)據(jù)提供新思路和新方法，同時也有助于提高地磁相對觀測資料的內(nèi)在質(zhì)量。

參考文獻

[1]許可，劉瑞瑞，等.基于matlab的背景噪聲計算程序的設(shè)計與應(yīng)用[J].山西地震，2015，2（1）：43-45.

[2]魯權(quán)，刑西淳，李西京，等.涇陽臺數(shù)字地磁信號的干擾分析及去噪處理[J].地球物理學(xué)報，2013，33（1）：49-55.

[3]潘泉，孟晉麗，張磊，等.小波濾波方法及應(yīng)用[J].電子與信息學(xué)報，2007，29（1）：236-242.

[4]李成錄，張永仁，基于小波分析與數(shù)學(xué)形態(tài)學(xué)融合的降噪算法研究[J].價值工程，2012，19（2）：14-16.

[5]趙昭，劉利林，張承學(xué)，等.形態(tài)學(xué)濾波器結(jié)構(gòu)元素選取原則研究與分析[J].電力系統(tǒng)保護與控制，2009，37（14）：21-25.

[6]萬永革.數(shù)字信號處理的MATLAB實現(xiàn)[M].科學(xué)出版社，2012.

[7]胡廣書.數(shù)字信號處理—理論、算法與實現(xiàn)[M].北京：清華大學(xué)出版社，1997.