組合廣義形態(tài)濾波在鹿井地區(qū)音頻大地電磁測(cè)量數(shù)據(jù)去噪中的應(yīng)用

李于波，劉祜，段書新，汪碩

（核工業(yè)北京地質(zhì)研究院 中核集團(tuán)鈾資源勘查與評(píng)價(jià)技術(shù)重點(diǎn)實(shí)驗(yàn)室，北京 100029）

鹿井地區(qū)位于湖南省郴州市汝城縣境內(nèi)，處于湖南、江西、廣東三省交界處，具有鐵、鎢、銅、稀土、鈾等多種礦產(chǎn)資源。為查明鹿井鈾礦田深部成礦要素，開展了音頻大地電磁（Audio Magnetotelluric Method，簡(jiǎn) 稱AMT）測(cè)量。但由于工作區(qū)內(nèi)電磁干擾嚴(yán)重，部分采集到的時(shí)間序列上疊加了與地質(zhì)體無(wú)關(guān)的噪音信號(hào)，需進(jìn)行去噪處理以便反演處理得到地下的真實(shí)電性結(jié)構(gòu)。

針對(duì)音頻大地電磁去噪，前人從時(shí)間域或頻率域的角度先后開展了互功率譜法、小波分析法、遠(yuǎn)參考技術(shù)、小波自適應(yīng)閾值、改進(jìn)閾值的TI 小波、時(shí)間序列依賴關(guān)系等去噪技術(shù)［1-7］。上述系列去噪方法取得了一定的應(yīng)用效果，但亦存在若干問(wèn)題，如：互功率譜法只能去除不相關(guān)噪聲、小波變換參數(shù)選取不當(dāng)會(huì)導(dǎo)致去噪效果適得其反等。為實(shí)現(xiàn)對(duì)大量實(shí)測(cè)數(shù)據(jù)的快速去噪，筆者以鹿井地區(qū)多個(gè)測(cè)點(diǎn)的時(shí)間序列為研究對(duì)象，探討了組合廣義形態(tài)濾波方法在該區(qū)AMT 去噪中的應(yīng)用。

1 鹿井地區(qū)AMT 工程部署及電磁干擾特征

1.1 AMT 工程部署

為查明鹿井鈾礦田深部地質(zhì)情況，在該區(qū)布設(shè)了4 條AMT 探測(cè)剖面（圖1）。野外利用V8（主機(jī)）結(jié)合AMTC30（磁探頭）的方式開展數(shù)據(jù)采集，獲取到TS2（高頻）、TS3（中頻）、TS4（低頻）3 個(gè)時(shí)間序列文件，對(duì)應(yīng)的采樣率分別為24 000、2 400、150 Hz。

圖1 鹿井地區(qū)AMT 測(cè)線布置圖Fig.1 Layout of AMT survey line in Lujing area

1.2 AMT 電磁干擾特征

據(jù)現(xiàn)場(chǎng)實(shí)際踏勘，測(cè)區(qū)內(nèi)有高速、省道、鄉(xiāng)道等公路貫穿，同時(shí)分布了較多的廣播電臺(tái)、高壓線、發(fā)射塔等，預(yù)計(jì)會(huì)對(duì)AMT 實(shí)測(cè)信號(hào)產(chǎn)生較多干擾。對(duì)實(shí)測(cè)的時(shí)間序列進(jìn)行分析，證實(shí)了上述干擾的存在。

圖2a 中，Ex電場(chǎng)信號(hào)由于受到工頻噪聲干擾，其背景值明顯高于其他測(cè)點(diǎn)的背景水平。圖2b 中，Ey電場(chǎng)信 號(hào)在180 和400 采 樣點(diǎn)處出現(xiàn)明顯的尖波，為三角波噪音在時(shí)間序列上的反映。圖2c 中，Hx磁場(chǎng)信號(hào)在30 和330 采樣點(diǎn)處亦出現(xiàn)2 個(gè)格值超過(guò)150 000 的脈沖。圖2d中，Hx磁場(chǎng)信號(hào)出現(xiàn)規(guī)律性的似充放電三角波噪聲干擾。為獲取地下真實(shí)的電性結(jié)構(gòu)，對(duì)上述含噪音數(shù)據(jù)進(jìn)行去噪處理就顯得尤為重要。

圖2 鹿井地區(qū)采集電場(chǎng)信號(hào)受多種噪聲干擾圖Fig.2 The acquired signal interfered by various noises in Lujing area

2 組合廣義形態(tài)濾波原理及其在AMT去噪中的具體實(shí)現(xiàn)

形態(tài)濾波是由法國(guó)數(shù)學(xué)家Matheron G和Serra J于1964年共同創(chuàng)立的信號(hào)分析方法，其核心思想是建立適當(dāng)?shù)慕Y(jié)構(gòu)元素，讓待處理的序列（或集合）依次穿過(guò)結(jié)構(gòu)元素，從而提取到該序列（或集合）的數(shù)學(xué)形態(tài)。在數(shù)據(jù)分析過(guò)程中，若同時(shí)使用不同的結(jié)構(gòu)元素進(jìn)行組合處理，即為組合廣義形態(tài)濾波。

在AMT 時(shí)間序列處理過(guò)程中，為克服時(shí)間序列的基線漂移，采用正負(fù)結(jié)構(gòu)元素串行的方法構(gòu)建組合廣義形態(tài)濾波器［8］，其流程如圖3 所示。以Ex 電場(chǎng)信號(hào)為例，將信號(hào)Input_Ex 依次通過(guò)ΨGOC(GCO)(s1，s2)和ΨGOC(GCO)(-s1，-s2)兩個(gè)基本濾波單元（其中S1，S2 為結(jié)構(gòu)元素），即得到Ex 電場(chǎng)中所包含的噪聲序列Exlb(n)。

圖3 組合廣義形態(tài)濾波流程圖Fig.3 Flow chart of combined generalized morphology filtering

圖（3）中，ΨGOC(GCO)(s1，s2) 濾波單元的定義為

其中GOC(f(n))、GCO(f(n))分別為廣義形態(tài)開-閉濾波器和廣義形態(tài)閉-開濾波器，“°”為形態(tài)開運(yùn)算，“?”為形態(tài)閉運(yùn)算，其定義如式（2）、（3）所示：

式（2）和式（3）中的Θ 為腐蝕運(yùn)算，其目的是減少待處理序列的峰值，即剔除掉不平滑的峰值數(shù)據(jù)，縮小峰域和拓寬谷域，使目標(biāo)序列的值收縮；⊕為膨脹運(yùn)算，其目的是拓寬峰域和縮小谷域，使目標(biāo)序列的值膨脹，兩運(yùn)算的公式分別如式（4）、（5）所示：

在本次AMT 去噪處理中，將S1 和S2 分別設(shè)定為5 點(diǎn)圓盤型結(jié)構(gòu)元素（圖4a）和5 點(diǎn)拋物線型結(jié)構(gòu)元素（圖4b），兩者的函數(shù)表達(dá)式分別如式（6）、（7）所示：

圖4 五點(diǎn)圓盤型（a）和五點(diǎn)拋物線型（b）結(jié)構(gòu)元素示意圖Fig.4 Diagram of structuring elements of disc type（a）and parabolic type（b）

據(jù)圖3 獲取到噪聲序列Exlb(n)后，將其從Ex電場(chǎng)信號(hào)中剔除，即得到重構(gòu)后的時(shí)間序列信號(hào)Exdn(n)。

3 實(shí)測(cè)數(shù)據(jù)的去噪應(yīng)用

為實(shí)現(xiàn)對(duì)大量實(shí)測(cè)數(shù)據(jù)的快速去噪，筆者以Python 編程語(yǔ)言為平臺(tái)，按照?qǐng)D5 所示流程實(shí)現(xiàn)了鹿井地區(qū)AMT 實(shí)測(cè)數(shù)據(jù)的組合廣義形態(tài)濾波去噪。整個(gè)流程的核心是對(duì)4個(gè)電磁場(chǎng)分量做正、負(fù)結(jié)構(gòu)元素的廣義形態(tài)濾波處理，并重構(gòu)TS 文件。其中TS2 時(shí)間序列的結(jié)構(gòu)元素參數(shù)L2=3，k2=0.5；TS3 的結(jié)構(gòu)元素參數(shù)L3=5，k3=0.4；TS4 的結(jié)構(gòu)元素參數(shù)L4=2，k4=0.3。

圖5 數(shù)據(jù)處理及組合廣義形態(tài)濾波去噪流程圖Fig.5 Flow chart of data processing and filtering denoising by generalized morphology combination

以line2-05、line2-07、line2-27 測(cè)點(diǎn)為例，重構(gòu)后的時(shí)間序列（圖6、7、8）分別去除了階躍噪聲、三角波、似充放電三角波的干擾，表明組合廣義形態(tài)濾波對(duì)該類電磁干擾具有較好的去噪效果。

圖6 鹿井地區(qū)采集的含噪時(shí)間序列去噪前、去噪后信號(hào)對(duì)比圖（line2-05）Fig.6 The comparison of original series and denoised series in Lujing area（line2-05）

圖7 鹿井地區(qū)采集的含噪時(shí)間序列去噪前、去噪后信號(hào)對(duì)比圖（line2-07）Fig.7 The comparison of original series and denoised series in Lujing area（line2-07）

圖8 鹿井地區(qū)采集的含噪時(shí)間序列去噪前、去噪后信號(hào)對(duì)比圖（line2-27）Fig.8 The comparison of original series and denoised series in Lujing area（line2-27）

上述去噪效果在各個(gè)測(cè)點(diǎn)的功率譜密度上亦有所表現(xiàn)。以Line2-07 測(cè)點(diǎn)為例，圖9a 展示了TS2 序列中各分量的功率譜密度，曲線在50 Hz 及其倍頻處出現(xiàn)了明顯的峰值情況，表明該點(diǎn)電磁場(chǎng)信號(hào)受到了工業(yè)諧波干擾；同時(shí)，低頻的功率譜密度小、高頻的功率譜密度大，推測(cè)是三角波等大格值噪聲淹沒(méi)了低頻有用信號(hào)所致。對(duì)其進(jìn)行組合廣義形態(tài)濾波處理，得到重構(gòu)后的時(shí)間序列功率譜密度曲線如圖9b 所示。可以看出，去噪后的時(shí)間序列在150、250、350、450 Hz 等頻率處的工業(yè)噪聲得到了壓制，且低頻1～300 Hz 范圍內(nèi)的功率譜密度得到了加強(qiáng)，Ex、Ey、HxHy 四道信號(hào)的功率譜密度均較為穩(wěn)定，分別保持在1 000、100、10 000、1 000 w/Hz 左右。

圖9 Line2-07 測(cè)點(diǎn)TS2、TS3 序列去噪前、后功率譜密度對(duì)比圖Fig.9 The comparison of original and denoised TS2 and TS3 PSD of Line 2-07 point

相比較而言，該測(cè)點(diǎn)處TS3 時(shí)間序列受噪聲影響較小，經(jīng)組合廣義形態(tài)濾波處理前后的功率譜密度曲線變化不大。

圖10 為測(cè)點(diǎn)Line2-05 去噪前后的視電阻率曲線對(duì)比。由圖可知：在小于1 Hz 的頻率范圍，rhoxy 和rhoyx 經(jīng)廣義形態(tài)濾波去噪后變得更加光滑、穩(wěn)定。針對(duì)去噪后rhoxy 曲線在10 Hz 左右出現(xiàn)的跳變，推測(cè)是由于形態(tài)濾波時(shí)將部分有用信號(hào)進(jìn)行了剔除所致，需要在后續(xù)數(shù)據(jù)處理過(guò)程中修改濾波參數(shù)，以減少有效信號(hào)的損失，進(jìn)一步提高廣義形態(tài)濾波效果。

圖10 鹿井Line2-05 測(cè)點(diǎn)去噪前、去噪后視電阻率對(duì)比圖Fig.10 The comparison of original and denoised apparent resistivity of Line 2-05 point

圖11為測(cè)點(diǎn)Line2-07去噪前后的視電阻率曲線對(duì)比。由圖可知：去噪后的視電阻率曲線更平滑穩(wěn)定，該現(xiàn)象在rhoxy曲線10 Hz以低頻率表現(xiàn)的更為明顯。整體來(lái)看，去噪后50 Hz處的干擾依舊存在，后續(xù)可采用FIR數(shù)字濾波、傅里葉逆變換［9］、希爾伯特-黃變換［10］等方法進(jìn)行進(jìn)一步的去除。

圖11 Line2-07 測(cè)點(diǎn)去噪前、去噪后視電阻率對(duì)比Fig.11 The comparison of original and denoised apparent resistivity of point line 2-07

圖12 為測(cè)點(diǎn)Line2-27 去噪前后的視電阻率曲線對(duì)比?？梢钥吹剑家曤娮杪是€在3 Hz 以后出現(xiàn)了上下波動(dòng)的不穩(wěn)定現(xiàn)象，去噪后的視電阻率曲線則變得連續(xù)、穩(wěn)定。

圖12 鹿井Line2-27 測(cè)點(diǎn)去噪前、去噪后視電阻率對(duì)比Fig.12 The comparison of original and denoised apparent resistivity of point Line 2-27

4 結(jié)論

通過(guò)開展組合廣義形態(tài)濾波在鹿井地區(qū)音頻大地電磁去噪中的應(yīng)用研究，取得以下成果：

1）采用圓盤型結(jié)構(gòu)元素和拋物線型結(jié)構(gòu)元素構(gòu)建了組合廣義形態(tài)濾波器，在處理鹿井地區(qū)AMT實(shí)測(cè)數(shù)據(jù)時(shí)分別賦以不同的結(jié)構(gòu)元素參數(shù)K、L，對(duì)不同頻率采集的時(shí)間序列取得了較好的去噪效果。

2）數(shù)學(xué)形態(tài)濾波方法在做濾波處理后能提取到噪聲的波形形態(tài)，對(duì)于三角波、似充放電三角波、方波、階躍噪聲等波形形態(tài)明顯的電磁干擾，具有較好的去噪效果。

3）將組合廣義形態(tài)濾波法應(yīng)用于鹿井地區(qū)AMT實(shí)測(cè)數(shù)據(jù)，能便捷高效地完成時(shí)間域TSn時(shí)間序列的去噪處理，為AMT數(shù)據(jù)去噪提供了一種新的思路。