萬玲， 張揚， 林君， 蔣川東， 林婷婷

吉林大學(xué)儀器科學(xué)與電氣工程學(xué)院/地球信息探測儀器教育部重點實驗室， 長春　130026

?

基于能量運算的磁共振信號尖峰噪聲抑制方法

萬玲， 張揚， 林君， 蔣川東， 林婷婷*

吉林大學(xué)儀器科學(xué)與電氣工程學(xué)院/地球信息探測儀器教育部重點實驗室， 長春130026

摘要磁共振探測信號微弱，使用高靈敏度的核磁共振地下水探測儀，極易受到環(huán)境噪聲干擾.其中，工頻諧波噪聲和尖峰噪聲，是影響信號質(zhì)量最嚴重的兩類噪聲.國內(nèi)外研究表明，通過參考線圈的布設(shè)，依據(jù)探測線圈和參考線圈中噪聲相關(guān)性，利用自適應(yīng)參考對消算法，能夠有效抑制工頻諧波噪聲.然而，尖峰噪聲的存在嚴重影響了主通道與參考道的數(shù)據(jù)相關(guān)性，成為了參考對消算法應(yīng)用的難題與障礙.為解決這一問題，本文提出磁共振信號中尖峰噪聲的抑制方法，推導(dǎo)了能量域磁共振信號表達式，通過計算信號能量，可有效檢測尖峰噪聲并突出不易識別的小幅度尖峰噪聲，采用基于中位數(shù)的絕對偏差法確定閾值，進而剔除尖峰噪聲.為了驗證去噪效果，與應(yīng)用較廣的統(tǒng)計疊加法進行對比研究.仿真結(jié)果表明，對于干擾幅度較大、持續(xù)時間較長的尖峰噪聲，能量運算法和統(tǒng)計疊加法均能識別并剔除，且不損失有效的磁共振信號，標準差偏差控制在0.3%以內(nèi)，可以滿足實際應(yīng)用的要求；對于小于信號幅度1.5倍的尖峰噪聲，能量運算法可有效識別并剔除，而統(tǒng)計疊加法無法實現(xiàn).針對多通道探測系統(tǒng)，使用能量運算法剔除尖峰噪聲后，可明顯提高主通道和參考道的數(shù)據(jù)相關(guān)性，為后續(xù)自適應(yīng)參考對消算法的應(yīng)用奠定了基礎(chǔ).實測數(shù)據(jù)處理結(jié)果進一步證明了本文方法的實用性.

關(guān)鍵詞磁共振信號； 參考對消算法； 能量運算； 尖峰噪聲； 噪聲相關(guān)性

In the present paper, we provide a detailed insight into the technique of spike noise cancellation based on transferring the signal from time domain to energy domain and using the median absolute deviation method (MAD). First, we calculate the energy of MRS signal. After doing this, the signal is emphasized which is much more clear than in time domain. Then, a threshold principle based on MAD is provided to cancel all the noise above it.

In addition, this paper contains a comparison of the spike noise cancellation effects of Statistical Stacking and Energy cancellation. It is found that the two methods provide identical noise cancellation performance when the spike noise is strong and long enough. Statistical Stacking method is limited when the spike noise get weaker, especially weaker than 1.5 times of the signal. But Energy cancellation we proposed can still remove the spike noise and keep the MRS signal effectively. Moreover, from the noise properties analysis, correlation between the detection loop and the reference loop can be improved after spike noise cancellation. We first apply it on the synthetic data to see the correlation improvement and find that the correlation is improved from 0.1369 to 0.4941 after using the Energy cancellation.Then we obtain the true field data using the multi-channel MRS instrument in Changchun suburb, and apply the Energy cancellation method on the field data. The correlation is improved from 0.1575 to 0.2481. We also find that the harmonic noise is much more evident after spike noise cancellation, which is advantageous to use the adaptive noise cancellation later.

We anticipate that better noise cancelling results can be obtained with the Energy cancellation and adaptive noise cancellation methods used together and therefore a wider application of multi-channel MRS instrument can be made in future.

1引言

地下水中的氫原子被電磁場激發(fā)，會形成宏觀磁矩，在地面鋪設(shè)線圈可拾取宏觀磁矩進動產(chǎn)生的核磁共振信號(Legchenko and Valla, 2002; Legchenko et al.,2002；Lubczynski and Roy,2003；Roy and Lubczynski,2003；Schirov et al.,1991).信號的強弱與含水量大小有關(guān)，通過觀察拾取信號的強弱可判斷地下水含量.另外，不同強度的電磁場可激發(fā)地下不同深度的含水體.通過改變激發(fā)電磁場的強度，可獲得地下不同深度的含水量信息.與常規(guī)的物探找水方法相比，利用核磁共振原理探尋地下水，不打鉆就能確定地下含水層的深度和含水量大小等信息(林君等，2010；潘玉玲和張昌達，2000).因此，核磁共振測深(Magnetic Resonance Sounding, MRS)方法，作為一種直接探測地下水的地球物理方法(林君，2010；林君等，2012；陸其鵠等，2009；張榮等，2006)，已經(jīng)得到了世界上越來越多的關(guān)注.

但是，由于MRS信號非常微弱，高靈敏度的接收裝置極易受到周圍環(huán)境噪聲的干擾(林君等，2010)，MRS信號常常被淹沒在噪聲中，主要包括工頻諧波噪聲、尖峰噪聲等.隨著國際多通道磁共振儀器的誕生，工頻諧波噪聲可以通過參考通道的布設(shè)，利用其與主通道的相關(guān)性，通過自適應(yīng)參考對消算法消除(Fabian,2010；Müller-Petke and Yaramanci,2010,2011a，2011b；Radic，2006；Walsh，2008；田寶鳳等，2012).然而，強尖峰噪聲干擾的存在，將影響自適應(yīng)參考對消算法對工頻噪聲的剔除效果，進而阻礙MRS信號的有效提取，降低反演結(jié)果對水文地質(zhì)參數(shù)解釋的準確度.

尖峰噪聲主要是由太陽磁暴、雷暴或者任何物體突然放電等引起的.這種尖峰噪聲的主要特征是，時域上的分布和幅度具有偶然性和隨機性，并且幅度大于或遠大于信號幅度；頻域上頻率不固定，且頻譜分布范圍廣，通常其覆蓋從幾赫茲到100 MHz以上，與MRS信號頻率混疊在一起(林君等，2010).針對尖峰噪聲對MRS信號干擾問題，國內(nèi)外專家和學(xué)者進行了相應(yīng)的研究.王中興等(2009)提出差閾值替代疊加方法削弱和抑制奇異噪聲，提高信號參數(shù)的提取精度；蔣川東等(2012)采用統(tǒng)計疊加方法去除尖峰噪聲干擾；Strehl(2006), Strehl等(2006)根據(jù)尖峰噪聲對MRS信號的干擾特點，提出了基于小波硬閾值變換方法去除MRS信號中的尖峰噪聲.差閾值和統(tǒng)計疊加的方法能夠判別出幅值較大的尖峰噪聲，然而當尖峰噪聲幅度較小時，難以判別噪聲從而不能獲得好的消噪效果；采用小波變換的方法進行噪聲去除，其去噪效果受到分解級數(shù)和閾值策略的限制，此外，當尖峰噪聲持續(xù)時間較長時，容易損失有效磁共振信號.綜上，磁共振尖峰噪聲處理中需要一種對其進行有效識別與去除的技術(shù).

本文借鑒了尖峰噪聲剔除方法在醫(yī)學(xué)領(lǐng)域(Mukhopadhyay and Ray,1998)的成功應(yīng)用，提出基于能量運算的磁共振信號尖峰噪聲抑制方法.通過計算信號能量，突出不易識別的小幅度尖峰噪聲，選取適當?shù)拈撝颠M行尖峰識別并剔除，進而有效提高主通道和參考道的數(shù)據(jù)相關(guān)性，為自適應(yīng)參考對消算法的廣泛應(yīng)用奠定基礎(chǔ).

2核磁共振探測噪聲分析及對策

2.1核磁共振地下水探測

將發(fā)射線圈鋪在地表，供入頻率為拉莫爾頻率的脈沖電流，形成激發(fā)磁場，當脈沖終止后，利用接收線圈接收磁共振全波信號E(t)，表達式為：

(1)

圖1　發(fā)射的電流脈沖與接收的磁共振信號示意圖Fig.1　Excited pulse and received MRS signal

磁共振信號微弱，通常為nV級.為提高MRS信號有效獲取能力，可應(yīng)用參考對消算法的多通道磁共振探測系統(tǒng)被廣泛應(yīng)用.其實際工作時在距離探測線圈一定距離處，鋪設(shè)一個或多個參考線圈，滿足探測線圈接收含噪的磁共振信號，參考線圈只接收噪聲的要求.利用探測線圈和參考線圈中噪聲的相關(guān)性，采用自適應(yīng)參考對消算法可有效實現(xiàn)噪聲壓制(田寶鳳等，2012).圖2為帶一個參考線圈的多通道磁共振探測系統(tǒng)的布線示意圖.

圖2　帶一個參考線圈的多通道磁共振探測系統(tǒng)示意圖Fig.2　Multi-channel MRS instrumentation with a reference loop

2.2MRS信號中噪聲分析及消噪對策

分析發(fā)現(xiàn)，影響MRS信號質(zhì)量的噪聲主要包括工頻諧波噪聲、尖峰噪聲和自然噪聲(林君等，2010).其中，自然噪聲類似于均勻白噪聲或高斯白噪聲，通過多次測量，采用疊加算法即可消除，而工頻諧波噪聲和尖峰噪聲則需要選取一定的消噪對策進行剔除.

多通道核磁共振地下水探測系統(tǒng)可以依據(jù)噪聲相關(guān)性，采用自適應(yīng)參考對消方法消噪，消噪效果的好壞主要受噪聲相關(guān)程度的影響(Müller-Petke and Yaramanci,2011b).研究發(fā)現(xiàn)，探測線圈和參考線圈中的工頻諧波噪聲相關(guān)性較好，而尖峰噪聲相關(guān)性很差，導(dǎo)致自適應(yīng)參考對消方法難以有效發(fā)揮作用.因此綜合考慮，MRS信號的消噪對策應(yīng)首先剔除尖峰噪聲，再通過自適應(yīng)參考對消方法壓制工頻諧波干擾，最后采用疊加算法剔除自然噪聲影響，圖3給出了MRS信號的處理流程.

圖3　MRS信號處理流程圖Fig.3　Flow chart of MRS signal processing

3基于能量運算的去尖峰算法

3.1信號能量計算表達式的推導(dǎo)與定義

用xn代表采樣后離散的余弦信號，用式(2)表示：

(2)

式中，Ω=2πf/fs，fs是采樣頻率；A和φ分別是信號的振幅和初始相位.假定A、Ω和φ三個參數(shù)恒定，并選擇與xn相鄰的xn-1和xn+1組成方程組，得到表達式(3)：

(3)

根據(jù)三角函數(shù)的積化和差公式對表達式(3)進行變換，可以得出：

(4)

再根據(jù)三角函數(shù)的二倍角公式對表達式(4)進行變換，可以得出：

xn+1xn-1=A2cos2(Ω n+φ)-A2sin2Ω，

(5)

(6)

整理得：

(7)

當Ω<π/4，即f/fs<1/8，在誤差允許的范圍內(nèi)，

(8)

由于余弦信號的能量正比于振幅的平方與頻率的平方的乘積(Miller,1937)，因此可用表達式(8)表示余弦信號的能量，得到：

(9)

式中，En即為余弦信號xn的能量.

MRS信號是呈指數(shù)規(guī)律衰減的余弦信號，因此MRS信號的振幅A是隨時間按指數(shù)規(guī)律變化的.假設(shè)A=A0·e-a n，Ω=Ω0，φ=φ0，可得到MRS信號為：

(10)

利用信號能量表達式(9)計算表達式(10)所代表的MRS信號的能量，可以得到

(11)

同樣，當采樣頻率大于信號頻率8倍時，滿足

(12)

3.2尖峰噪聲的檢測原理

當有尖峰噪聲干擾時，假設(shè)xn是接收系統(tǒng)采集的信號，x1n是不含尖峰噪聲的MRS信號，x2n是尖峰噪聲.可知，xn是x1n和x2n的線性組合，滿足如下關(guān)系：

(13)

假設(shè)使

(14)

成立，可以得出：

(15)如果x1n和x2n不相關(guān)，用數(shù)學(xué)期望值求信號的能量，滿足如下關(guān)系：

(16)

引入延遲時間為2τ的x(t)的自相關(guān)函數(shù)，表達式為

(17)

式中，Rx(t+τ,t-τ)的傅里葉變換為

(18)

式中，W(t,ω)就是該信號的譜估計.由于W(t,ω)是角頻率w的偶函數(shù)

(19)

式中，B是通頻帶寬.因此，利用公式(15)至公式(19)對表達式(14)進行變換，可以得出：

=Rx(n,n)-Rx(n+1,n-1)

(20)

式中，(1-cos2ωT)起高通濾波器的作用.假設(shè)

(21)

可知，W(n,ω)經(jīng)過高通濾波器后是W′(n,ω).由式(20)可以得到

(22)其中

(23)

(24)

經(jīng)過上述能量運算后，尖峰噪聲被突出，此時需要選取適當?shù)拈撝颠M行尖峰識別并剔除.

圖4a為一組受尖峰噪聲干擾的MRS信號，可以看出MRS信號是按指數(shù)規(guī)律變化的衰減信號，而尖峰噪聲在數(shù)據(jù)列中的表現(xiàn)形式為異常值，其數(shù)值多大于信號值.經(jīng)過能量運算后信號如圖4b所示，可見其依然按指數(shù)規(guī)律衰減.信號和尖峰噪聲的這一特點，使得選取的閾值既不能把信號最大值標記成異常值，又必須識別不明顯的異常值，即幅度值小的尖峰噪聲.因此，常規(guī)采用基于數(shù)據(jù)列標準差的閾值確定方法(吳石林和張玘，2010)，當數(shù)據(jù)列中存在數(shù)值較大的異常值時，閾值受其影響會相對較高，進而降低了對不明顯異常的檢測能力.

基于中位數(shù)的絕對偏差法(Median Absolute Deviation, MAD)(Hoaglin et al.,2000)是一種非常有效的統(tǒng)計學(xué)方法.對于一系列測量值為x1,x2,…,xn的數(shù)據(jù)列，計算MAD，應(yīng)先求某測量值相對于其所在數(shù)據(jù)列中位數(shù)的絕對偏差，再對該絕對偏差求中位數(shù)，定義為：

(25)

圖4　峰值檢測實例(a) 受尖峰噪聲干擾的MRS信號； (b) 信號能量及閾值.Fig.4　Example of spike detection(a) MRS signal disturbed by spike noises； (b) Energy of signal and threshold.

將能量域中的信號與閾值進行比較，超過閾值數(shù)據(jù)視作尖峰噪聲并剔除.由于實際采集的信號包含隨機噪聲，隨機噪聲中的高頻成分在能量域中容易被誤識別為“小尖峰”，導(dǎo)致增加無效的判斷次數(shù)及計算量.因此，可適當?shù)赝ㄟ^低通濾波器，去除隨機噪聲中的高頻成分，以準確識別并剔除尖峰噪聲.

仿真一組含有16000個等精度采樣點、采樣頻率為66667 Hz的四次疊加信號，單次采集的原始信號如圖5a所示，圖中的尖峰噪聲、工頻諧波噪聲和高斯白噪聲的干擾時刻以及干擾強度均是隨機的.使用能量運算法剔除尖峰噪聲主要包括三個步驟：

(1) 計算單次采集信號的能量，將信號由時域轉(zhuǎn)化到能量域.再將計算后的信號的能量通過低通濾波器，如圖5b所示.

(2) 采用MAD方法確定閾值.

(3) 將低通濾波后信號的能量與閾值進行比較，如圖5c所示，超出閾值的部分視作尖峰噪聲，用該時刻無尖峰噪聲的數(shù)據(jù)列的中位數(shù)進行替代，以達到剔除尖峰噪聲的目的.圖5d給出了剔除尖峰噪聲后的信號.

4消噪算法仿真分析

目前在磁共振信號處理中通常采用統(tǒng)計疊加方法進行尖峰噪聲剔除(Jiang et al.,2011)，其基本思想是把尖峰噪聲看作統(tǒng)計學(xué)中的粗大誤差，首先提取一個可疑的測量值，然后按t分布檢驗提取的測量值是否含有粗大誤差，再決定其是否需要剔除，最終剔除后的數(shù)據(jù)參與到疊加運算中.

圖5　基于能量運算的去尖峰算法消噪流程圖(a) 單次采集的原始信號； (b) 計算后的信號能量； (c) 中位數(shù)絕對偏差法確定的閾值； (d) 剔除尖峰噪聲后的信號.Fig.5　The workflow of energy calculation(a) Original signals of each sampled； (b) Energy of signal calculated by energy calculation； (c) Threshold determined by median absolute deviation；(d) MRS signals of each sampled without spikes.

分別采用能量運算法和統(tǒng)計疊加法，對含有不同干擾幅度以及干擾持續(xù)時間尖峰噪聲的信號進行處理，驗證能量運算法的應(yīng)用效果，并觀察多通道系統(tǒng)采集中剔除尖峰噪聲后，主通道和參考道數(shù)據(jù)相關(guān)性的改善情況.

除了能從波形上觀察尖峰噪聲的剔除效果，本文引入信噪比和標準差進行數(shù)值衡量.信噪比的定義為：

(26)

式中，信噪比SNR采用了dB的表示方式.

標準差的定義為：

(27)

4.1干擾幅度不同的尖峰噪聲去噪仿真

圖6　統(tǒng)計疊加法和能量運算法對不同幅度尖峰噪聲的去噪結(jié)果(a) 仿真的含噪信號； (b) 應(yīng)用統(tǒng)計疊加法剔除尖峰噪聲； (c) 應(yīng)用能量運算法剔除尖峰噪聲.虛線是區(qū)分不同幅度尖峰噪聲的標識，幅度值超過該虛線的視作大幅度尖峰噪聲，反之，視作小幅度尖峰噪聲.Fig.6　Statistical stacking approach and energy calculation approach to remove spikes of different amplitude(a) Simulated signals with spikes； (b) Remove spikes by statistical stacking approach； (c) Remove spikes by energy calculation approach. The dash line is the sign of different amplitude，data above the dash line represent big spikes，conversely，represent small spikes.

圖中具有衰減趨勢的虛線，其幅度是信號幅度的1.5倍，將它作為區(qū)別不同幅度尖峰噪聲的標識.幅度值超過該虛線的視作大幅度的尖峰噪聲；反之，視作小幅度的尖峰噪聲.分別采用統(tǒng)計疊加法和能量運算法對圖6a所示的仿真含噪信號進行消噪處理，圖6b和6c為對應(yīng)的消噪結(jié)果.從消噪結(jié)果圖中可以明顯看出，應(yīng)用統(tǒng)計疊加法后，幅值較大的尖峰噪聲被剔除，但幅值較小的尖峰噪聲仍然存在；而經(jīng)過能量運算法后，幅值較大的尖峰噪聲和幅值較小的尖峰噪聲均被剔除.應(yīng)用統(tǒng)計疊加法后，信噪比SNR=6.6617 dB，提高了5.7432 dB.另外，標準差STD=80.1300 nV，與不含尖峰噪聲的MRS信號相比，標準差偏差為3.4610%.而應(yīng)用能量運算法后，信噪比SNR=12.6133 dB，提高了11.6948 dB.標準差STD=77.7195 nV，標準差偏差為0.3486%.

圖7　統(tǒng)計疊加法和能量運算法對不同持續(xù)時間尖峰噪聲的去噪結(jié)果(a) 仿真的含噪信號； (b) 應(yīng)用統(tǒng)計疊加法剔除尖峰噪聲； (c) 應(yīng)用能量運算法剔除尖峰噪聲.Fig.7　Statistical stacking approach and energy calculation approach to remove spikes of different duration time(a) Simulated signals with spikes； (b) Remove spikes by statistical stacking approach；(c) Remove spikes by energy calculation approach.

結(jié)合統(tǒng)計疊加法和能量運算法的消噪原理對其消噪差異進行分析，可知，統(tǒng)計疊加法先剔除一個可疑的測量值，然后檢驗剔除的測量值是否含有粗大誤差.在某些情況下，幅度值較小的干擾被認為不含有粗大誤差，不予剔除.因此，統(tǒng)計疊加法能剔除幅度值較大的尖峰噪聲，卻常常忽略幅度值較小的尖峰噪聲.然而，經(jīng)過能量運算后，原來時域中幅度值較小的尖峰噪聲在能量域中其幅度值會被增加，明顯地超過信號的幅度，很容易被識別并去除.

4.2持續(xù)時間不同的尖峰噪聲去噪仿真

分別采用統(tǒng)計疊加法和能量運算法對圖7a所示的仿真含噪信號進行消噪處理，圖7b和7c為對應(yīng)的消噪結(jié)果.應(yīng)用統(tǒng)計疊加法后，信噪比SNR=18.3899 dB，提高了18.1638 dB.另外，標準差STD=77.6652 nV，與不含尖峰噪聲的MRS信號相比，標準差偏差為0.2785%.而應(yīng)用能量運算法后，信噪比SNR=18.3921 dB，提高了18.1660 dB.標準差STD=77.6744 nV，標準差偏差為0.2904%.從消噪結(jié)果圖以及消噪后的信噪比和標準差中可以明顯看出，統(tǒng)計疊加法和能量運算法均可剔除干擾幅度較大，持續(xù)時間較長的尖峰噪聲，且不損失有效的磁共振信號，標準差偏差控制在0.3%以內(nèi)，滿足實際應(yīng)用的需求.

4.3多通道系統(tǒng)數(shù)據(jù)去噪后提高信號相關(guān)性驗證

模擬探測線圈采集含噪聲的磁共振信號，參考線圈只采集周圍環(huán)境噪聲，二者共同構(gòu)成多通道采集系統(tǒng).參考對消算法與通道間數(shù)據(jù)相關(guān)性有關(guān)，因此需引入相關(guān)性的概念.設(shè)Pxy(f)是信號x和信號y的互相關(guān)功率譜密度，Pxx(f)是信號x的自相關(guān)功率譜密度，Pyy(f)是信號y的自相關(guān)功率譜密度.相關(guān)性的定義為：

圖8　使用統(tǒng)計疊加方法和使用基于能量運算的尖峰噪聲抑制方法剔除主通道中的尖峰噪聲(a) 統(tǒng)計疊加方法； (b) 能量運算法； (c) 頻譜結(jié)果.Fig.8　Remove spikes by statistical stacking approach and energy calculation approach in detection loop(a) Statistical stacking approach； (b) Energy calculation approach； (c) Results in frequency domain.

(28)

式中，r介于0和1之間，當r=1時，信號x和信號y完全相關(guān)；當r=0時，信號x和信號y不相關(guān).

仿真模型參數(shù)如下：MRS信號的中心頻率fL=2326 Hz，工頻諧波干擾頻率f1=2250 Hz，f2=2350 Hz，參考道中工頻干擾的幅度是主通道中工頻干擾幅度的3倍，相位相同.主通道中隨機加入幅度為150 nV到550 nV、持續(xù)時間為4 ms到15 ms的五個尖峰噪聲干擾；同時，參考道中隨機加入幅度為300 nV到1000 nV、持續(xù)時間為5 ms到20 ms的八個尖峰噪聲干擾.使得主通道和參考道的數(shù)據(jù)相關(guān)性r=0.1369.

圖9　使用統(tǒng)計疊加方法和使用基于能量運算的尖峰噪聲抑制方法剔除參考道中的尖峰噪聲(a) 統(tǒng)計疊加方法； (b) 能量運算法； (c) 頻譜結(jié)果.Fig.9　Remove spikes by statistical stacking approach and energy calculation approach in reference loop(a) Statistical stacking approach； (b) Energy calculation approach； (c) Results in frequency domain.

分別使用統(tǒng)計疊加法和能量運算法剔除主通道以及參考道中尖峰噪聲.經(jīng)過計算，應(yīng)用統(tǒng)計疊加法剔除尖峰噪聲后，主通道與參考道信號相關(guān)性r=0.4502，提高了0.3133.應(yīng)用能量運算法剔除尖峰噪聲后，主通道與參考道信號相關(guān)性r=0.4941，提高了0.3572.圖8和圖9分別給出了統(tǒng)計疊加法和能量運算法剔除尖峰噪聲前后主通道以及參考道信號的時域波形及頻譜.消噪前后信號的時域波形進一步驗證了4.1和4.2節(jié)得出的結(jié)論：統(tǒng)計疊加法和能量運算法均可剔除干擾幅度較大，持續(xù)時間較長的尖峰噪聲；另外，能量運算法可以識別并剔除干擾幅度較小的尖峰噪聲.從消噪前后信號頻譜圖以及消噪后的相關(guān)性中可以明顯看出，采用能量運算法剔除尖峰噪聲后，能夠顯著提高主通道與參考道信號相關(guān)性，為后續(xù)應(yīng)用自適應(yīng)噪聲對消算法更好地剔除工頻諧波噪聲奠定了基礎(chǔ).

5實測數(shù)據(jù)處理

為了驗證算法的實用性，在長春市郊區(qū)燒鍋鎮(zhèn)進行了野外試驗.使用核磁共振地下水探測儀(林君等，2010)進行信號采集，當?shù)氐睦獱栴l率fL=2326 Hz，圖10a中藍色線為儀器采集的原始信號時域波形圖.從圖中可以看出：該信號含有較多的尖峰噪聲，并且從信號的時域波形上不能明顯看出MRS信號指數(shù)型的衰減趨勢.經(jīng)計算，原始信號的信噪比SNR=-37.7090 dB.

依次操作3.3節(jié)中使用能量運算法的三個步驟，對圖10a所示實測數(shù)據(jù)進行尖峰噪聲剔除，得到消噪后的時域信號，如圖10a中紅色線所示.可見，原始信號中的尖峰噪聲被明顯剔除，消噪后的信號呈現(xiàn)指數(shù)形式衰減.再從頻域上進行比較，如圖10b所示，消噪后尖峰噪聲的頻率分量被壓制，突出了信號的頻率.經(jīng)計算，消噪后信號的信噪比SNR1=-29.1380 dB，剔除尖峰噪聲后信噪比提高了ΔSNR1=8.5710 dB.

圖10　基于能量運算的尖峰噪聲抑制方法剔除實測數(shù)據(jù)中的尖峰噪聲(a) 剔除尖峰噪聲前后的時域圖； (b) 剔除尖峰噪聲前后的頻域圖.Fig.10　Spikes in MRS experimental data removed by energy calculation approach(a) Remove spikes in time domain；(b) Remove spikes in frequency domain.

為了進一步說明方法的有效性，分別對不同信噪比下的另外兩組實測數(shù)據(jù)進行處理，處理結(jié)果如圖11所示.由圖11消噪前后的時域信號及頻譜可以看出，對于不同信噪比實測數(shù)據(jù)，基于能量運算的尖峰噪聲抑制方法可以獲得較好的消噪效果.

為了驗證算法在多通道系統(tǒng)數(shù)據(jù)處理中的實用性，在長春市凈月潭進行了野外試驗.使用帶參考線圈的核磁共振地下水探測儀進行信號采集，當?shù)氐睦獱栴l率fL=2326 Hz，受到2250 Hz和2350 Hz兩個明顯的工頻噪聲干擾.使用能量運算法剔除主通道和參考道中的尖峰噪聲，結(jié)果如圖12和圖13所示.從圖中可以看出：信號中的尖峰噪聲被剔除，并且去除尖峰噪聲后，工頻噪聲的頻率成分更加明顯.經(jīng)計算，剔除尖峰噪聲前，主通道和參考道的信號相關(guān)性r=0.1575；剔除尖峰噪聲后，主通道和參考道的信號相關(guān)性r=0.2481.可見，對多通道磁共振探測系統(tǒng)，先剔除主通道和參考道中的尖峰噪聲，提高信號相關(guān)性后，再應(yīng)用自適應(yīng)噪聲對消算法去除工頻諧波噪聲將會取得更好的去噪效果.

6結(jié)論

在強尖峰噪聲干擾下如何更好地實現(xiàn)微弱MRS信號可靠提取是核磁共振地下水探測抗干擾技術(shù)中的關(guān)鍵問題，它將直接影響后續(xù)應(yīng)用自適應(yīng)參考對消算法對工頻諧波噪聲的有效濾除.能量運算法可有效檢測尖峰噪聲并突出不易識別的小幅度尖峰噪聲.本文根據(jù)能量運算法的這一特點，針對強尖峰噪聲干擾問題，推導(dǎo)了能量域信號表達式，通過計算信號能量，有效識別尖峰噪聲，然后采用基于中位數(shù)的絕對偏差法確定閾值，實現(xiàn)對尖峰噪聲的有效剔除.與應(yīng)用效果較好的統(tǒng)計疊加法對比研究，通過數(shù)值仿真和實測數(shù)據(jù)處理，得出如下結(jié)論：

圖11　不同信噪比下實測數(shù)據(jù)處理結(jié)果Fig.11　Processing results of MRS signals on different SNRs

圖12　基于能量運算的尖峰噪聲抑制方法剔除主通道中的尖峰噪聲(a) 剔除尖峰噪聲前后的時域圖； (b) 剔除尖峰噪聲前后的頻域圖.Fig.12　Spikes in detection loop removed by energy calculation approach(a) Remove spikes in time domain； (b) Remove spikes in frequency domain.

圖13　基于能量運算的尖峰噪聲抑制方法剔除參考道中的尖峰噪聲(a) 剔除尖峰噪聲前后的時域圖； (b) 剔除尖峰噪聲前后的頻域圖.Fig.13　Spikes in reference loop removed by energy calculation approach(a) Remove spikes in time domain； (b) Remove spikes in frequency domain.

(1) 對于干擾幅度較大，持續(xù)時間較長的尖峰噪聲，能量運算法和統(tǒng)計疊加法均能識別并剔除，且不損失有效的磁共振信號，仿真數(shù)據(jù)去噪后信噪比提高了18 dB，標準差偏差控制在0.3%以內(nèi)，滿足實際應(yīng)用的需求.

(2) 對于小于信號幅度1.5倍的小幅度尖峰噪聲，能量運算法可有效識別并剔除，而統(tǒng)計疊加法無法實現(xiàn).同一組仿真數(shù)據(jù)應(yīng)用統(tǒng)計疊加法剔除尖峰噪聲后，信噪比提高了5.7432 dB；應(yīng)用能量運算法剔除尖峰噪聲后，信噪比提高了11.6948 dB.

(3) 針對多通道系統(tǒng)，仿真數(shù)據(jù)經(jīng)過能量運算法剔除尖峰噪聲后，相關(guān)性提高了0.3572，可見剔除尖峰噪聲后可明顯提高主通道和參考道的數(shù)據(jù)相關(guān)性.

(4) 本文所述能量運算法不受尖峰噪聲干擾幅度以及干擾持續(xù)時間限制.隨機抽取三組不同信噪比的實測數(shù)據(jù)，采用能量運算法剔除尖峰噪聲后，信噪比均得到大幅度提升，獲得了很好的實際應(yīng)用效果.

(5) 針對受工頻噪聲干擾的多通道磁共振探測系統(tǒng)的實測數(shù)據(jù)，使用本文所述能量運算法剔除主通道和參考道的尖峰噪聲后，可提高數(shù)據(jù)相關(guān)性.

在核磁共振地下水實際探測應(yīng)用中，除尖峰噪聲之外，MRS數(shù)據(jù)還會受到工頻諧波噪聲和其他未知噪聲等多種影響，尤其在城市周邊和村莊附近探測時，高強度的工頻噪聲干擾會將MRS信號淹沒.鑒于此，帶有參考線圈的多通道探測系統(tǒng)被推出，自適應(yīng)參考對消算法成為MRS消噪算法的主流方向.在運用本文能量運算法對MRS信號進行尖峰噪聲剔除后，主通道和參考道噪聲相關(guān)性的提高為后續(xù)自適應(yīng)參考對消算法的有效應(yīng)用奠定了基礎(chǔ).下一步將研究基于自適應(yīng)參考對消算法對工頻諧波噪聲進行濾除的問題.

References

Fabian N. 2010. Processing of full time series, multichannel surface NMR signals [Master′s thesis]. Switzerland: ETH Zürich.

Hoaglin D C, Mosteller F, Tukey J W. 2000. Understanding Robust and Exploratory Data Analysis. New York: Wiley.

Jiang C D, Lin J, Duan Q M, et al. 2011. Statistical stacking and adaptive notch filter to remove high-level electromagnetic noise from MRS measurements.NearSurf.Geophys., 9(5): 459-468.

Legchenko A, Valla P. 2002. A review of the basic principles for proton magnetic resonance sounding measurements.J.Appl.Geophys., 50(1-2): 3-19.

Legchenko A, Baltassat J M, Beauce A, et al. 2002. Nuclear magnetic resonance as a geophysical tool for hydrogeologists.J.Appl.Geophys., 50(1-2): 21-46.

Lin J, Duan Q M, Wang Y J, et al. 2010. Theory and Design of Magnetic Resonance Sounding Instrument for Groundwater Detection and Its Applications (in Chinese). Beijing: Science Press.

Lin J. 2010. Situation and progress of nuclear magnetic resonance technique for groundwater investigations.ProgressinGeophysics(in Chinese), 25(2): 681-691, doi: 10.3969/j.issn.1004-2903.2010.02.043. Lin J, Jiang C D, Duan Q M, et al. 2012. The situation and progress of magnetic resonance sounding for groundwater investigations and underground applications.JournalofJilinUniversity(EarthScienceEdition) (in Chinese), 42(5): 1560-1570.

Lu Q H, Wu T B, Lin J. 2009. A research report on development of instrument science for geophysics.ProgressinGeophysics(in Chinese), 24(2): 750-758.

Lubczynski M, Roy J. 2003. Hydrogeological interpretation and potential of the new magnetic resonance sounding (MRS) method.JournalofHydrology, 283(1-4): 19-40.

Mukhopadhyay S, Ray G C. 1998. A new interpretation of nonlinear energy operator and its efficacy in spike detection.IEEET.Biomed.Eng., 45(2): 180-187.

Müller-Petke M, Yaramanci U. 2010. Improving the signal-to-noise ratio of surface-NMR measurements by reference channel based noise cancellation.∥ 16thEAGE European Meeting of Environmental and Engineering Geophysics. EAGE.Müller-Petke M, Yaramanci U. 2011a. Noise Cancellation for surface NMR: Derivation of time and frequency domain approaches.∥ 17thEAGE European Meeting of Environmental and Engineering Geophysics. EAGE.

Müller-Petke M, Yaramanci U. 2011b. Noise cancellation for surface NMR-application of time and frequency domain approaches.∥ 17thEAGE European Meeting of Environmental and Engineering Geophysics. EAGE.

Miller D C. 1937. Sound Waves: Their Shape and Speed. New York: Macmillan Company Press.

Pan Y L, Zhang C D. 2000. Theories and Methods of Surface Nuclear Magnetic Resonance (in Chinese). Beijing: China University of Geosciences Press.Radic T. 2006. Improving the signal-to-noise ratio of surface NMR data due to the remote reference technique.∥ 12thEAGE European Meeting of Environmental and Engineering Geophysics. EAGE.Roy J, Lubczynski M. 2003. The magnetic resonance sounding technique and its use for groundwater investigations.Hydrogeol.J., 11(4): 455-465. Schirov M, Legchenko A, Creer G. 1991. A new direct non-invasive groundwater detection technology for Australia.ExplorationGeophysics, 22(2): 333-338.

Strehl S. 2006. Development of strategies for improved filtering and fitting of SNMR-signals. Germany: Technical University of Berlin, Institute of Applied Geosciences, Department of Applied Geophysics.

Strehl S, Rommel I, Hertrich M, et al. 2006. New strategies for filtering and fitting of MRS signals.∥ Proceedings 3rdInternational MRS Workshop. Madrid, Spain, 65-68. Tian B F, Lin J, Duan Q M, et al. 2012. Variable step adaptive noise cancellation algorithm for magnetic resonance sounding signal with a reference coil.ChineseJ.Geophys. (in Chinese), 55(7): 2462-2472, doi: 10.6038/j.issn.0001-5733.2012.07.030.

Walsh D O. 2008. Multi-channel surface NMR instrumentation and software for 1D/2D groundwater investigations.J.Appl.Geophys., 66(3-4): 140-150. Wang Z X, Rong L L, Lin J. 2009. Spike noise elimination from surface nuclear magnetic resonance signal for underground water.JournalofJilinUniversity(EngineeringandTechnologyEdition) (in Chinese), 39(5): 1282-1287.

Wu S L, Zhang Q. 2010. Error Analysis and Data Processing (in Chinese). Beijing: Tsinghua University Press.

Zhang R, Hu X Y, Yang D K, et al. 2006. Review of development of surface nuclear magnetic resonance.ProgressinGeophysics(in Chinese), 21(1): 284-289.

附中文參考文獻

林君, 段清明, 王應(yīng)吉等. 2010. 核磁共振找水儀原理與應(yīng)用. 北京: 科學(xué)出版社.

林君. 2010. 核磁共振找水技術(shù)的研究現(xiàn)狀與發(fā)展趨勢. 地球物理學(xué)進展, 25(2): 681-691, doi: 10.3969/j.issn.1004-2903.2010.02.043. 林君, 蔣川東, 段清明等. 2012. 復(fù)雜條件下地下水磁共振探測與災(zāi)害水源探查研究進展. 吉林大學(xué)學(xué)報(地球科學(xué)版), 42(5): 1560-1570.

陸其鵠, 吳天彪, 林君. 2009. 地球物理儀器學(xué)科發(fā)展研究報告. 地球物理學(xué)進展, 24(2): 750-758.

潘玉玲, 張昌達. 2000. 地面核磁共振找水理論和方法. 北京: 中國地質(zhì)大學(xué)出版社.

田寶鳳, 林君, 段清明等. 2012. 基于參考線圈和變步長自適應(yīng)的磁共振信號噪聲壓制方法. 地球物理學(xué)報, 55(7): 2462-2472, doi: 10.6038/j.issn.0001-5733.2012.07.030.

王中興, 榮亮亮, 林君. 2009. 地面核磁共振找水信號中的奇異干擾抑制. 吉林大學(xué)學(xué)報(工學(xué)版), 39(5): 1282-1287. 吳石林, 張玘. 2010. 誤差分析與數(shù)據(jù)處理. 北京: 清華大學(xué)出版社. 張榮, 胡祥云, 楊迪琨等. 2006. 地面核磁共振技術(shù)發(fā)展述評. 地球物理學(xué)進展, 21(1): 284-289.

(本文編輯何燕)

基金項目國家重大科學(xué)儀器設(shè)備開發(fā)和應(yīng)用專項(2011YQ030113)，國家自然科學(xué)基金面上項目(41374075)，吉林省科技重點攻關(guān)項目(20150519008JH，20140204022GX),中國博士后科學(xué)基金面上項目(2014M561296)聯(lián)合資助.

作者簡介萬玲，女，1986年生，2013年畢業(yè)于吉林大學(xué)，講師，主要從事磁共振信號處理方法研究.E-mail：wanling@jlu.edu.cn *通訊作者林婷婷，女，1983年生，2010年畢業(yè)于吉林大學(xué)，博士生導(dǎo)師，主要從事磁共振地下水探測方法研究.E-mail：ttlin@jlu.edu.cn

doi:10.6038/cjg20160631 中圖分類號P631

收稿日期2015-02-02，2016-05-03收修定稿

Spikes removal of magnetic resonance sounding data based on energy calculation

WAN Ling， ZHANG Yang， LIN Jun， JIANG Chuan-Dong， LIN Ting-Ting*

CollegeofInstrumentationandElectricalEngineering/KeyLaboratoryofGeo-ExplorationandInstrumentation,MinistryofEducation,JilinUniversity,Changchun130026,China

AbstractMagnetic Resonance Sounding (MRS) signal is extremely easy corrupted by the noise, especially by the harmonic noise and spike noise. Harmonic noise cancellation is often based on remote references and the adaptive noise cancellation(ANC)algorithm to increase the signal-to-noise ratio (SNR). However, ANC algorithm cannot play an effective role when spike noise exists. Because the spike noise often lead to a bad correlation between the detection loop and the remote reference loop, it is necessary to remove the spike noise before using the adaptive noise cancellation algorithm to cancel the harmonic noise.

KeywordsMagnetic Resonance Sounding signal； Noise cancellation algorithm； Energy calculation； Spike noise； Noise correlation

萬玲， 張揚， 林君等. 2016. 基于能量運算的磁共振信號尖峰噪聲抑制方法.地球物理學(xué)報，59(6):2290-2301,doi:10.6038/cjg20160631.

Wan L, Zhang Y, Lin J, et al. 2016. Spikes removal of magnetic resonance sounding data based on energy calculation.ChineseJ.Geophys. (in Chinese),59(6):2290-2301,doi:10.6038/cjg20160631.