本發(fā)明涉及地面核磁共振地下水探測(cè)信號(hào)噪聲濾除和參數(shù)提取技術(shù)領(lǐng)域，具體是利用基于工頻諧波建模和自相關(guān)的地面核磁共振信號(hào)提取方法。

背景技術(shù)：

地面核磁共振(Magnetic Resonance Sounding,MRS)作為一種地球物理新技術(shù)，近年來(lái)被廣泛應(yīng)用到地下水探測(cè)中。與傳統(tǒng)的地球物理間接找水方法相比，地面核磁共振方法具有高分辨率，高效率，信息豐富和對(duì)水文參數(shù)具有唯一解釋的特點(diǎn)。地面核磁共振方法的基本原理是通過(guò)探測(cè)地下水中氫質(zhì)子共振躍遷產(chǎn)生的MRS信號(hào)，提取其特征參數(shù)，再通過(guò)反演解釋來(lái)獲取地下水的水質(zhì)信息。從一些國(guó)內(nèi)外研究和應(yīng)用中可以看出，MRS方法不僅在尋找地下水和地下水資源評(píng)價(jià)等方面取得了良好的應(yīng)用效果，而且在獲取水文地質(zhì)參數(shù)和定量化解釋方面也取得進(jìn)展。

MRS方法雖然比其他傳統(tǒng)方法具有明顯優(yōu)勢(shì)，但也存在一些不足。實(shí)際探測(cè)到的MRS信號(hào)是非常微弱的納伏級(jí)信號(hào)，信號(hào)中噪聲干擾非常大，包括尖峰噪聲、隨機(jī)噪聲和工頻諧波噪聲等，其中由電力線干擾引起的工頻諧波噪聲對(duì)其影響最大。目前，針對(duì)濾除核磁共振信號(hào)中噪聲的理論和方法有很多。訾彥勇在著作《基于EMD的磁共振探測(cè)信號(hào)噪聲抑制方法》，提出一種新型的信號(hào)時(shí)頻處理方法，可以在沒(méi)有輸入信號(hào)任何先驗(yàn)知識(shí)的情況下，自適應(yīng)地將信號(hào)分解成若干個(gè)固有模態(tài)函數(shù)(Intrinsic Mode Functions,IMF)，實(shí)現(xiàn)信號(hào)趨勢(shì)的有效提取，但存在端點(diǎn)效應(yīng)和模態(tài)混疊問(wèn)題；田寶鳳等人在著作《基于參考線圈和變步長(zhǎng)自適應(yīng)的磁共振信號(hào)噪聲壓制方法》中直接利用觀測(cè)數(shù)據(jù)不斷遞歸更新處理參數(shù)，實(shí)現(xiàn)不同信噪比和信號(hào)強(qiáng)度下MRS信號(hào)噪聲的去除，但該方法需要遠(yuǎn)端線圈同時(shí)采集噪聲數(shù)據(jù)做參考，對(duì)于單通道儀器無(wú)法實(shí)現(xiàn)，而且只能對(duì)相關(guān)噪聲處理，不能保證將多個(gè)工頻諧波噪聲全部消除。

CN104459809A公開(kāi)了“基于獨(dú)立成分分析的全波核磁共振信號(hào)噪聲濾除方法”，該算法能夠有效地實(shí)現(xiàn)全波MRS信號(hào)的信噪分離，且數(shù)據(jù)擬合后初始振幅和弛豫時(shí)間的相對(duì)誤差小于±5.00％，相比于其他經(jīng)典算法該算法更具有消噪性能的優(yōu)越性，但該方法只針對(duì)去除磁共振信號(hào)中的工頻諧波噪聲。CN104898172A公開(kāi)了“一種基于互相關(guān)的核磁共振全波信號(hào)噪聲濾除方法”，該發(fā)明利用噪聲與拉莫爾頻率的正弦信號(hào)不相關(guān)，而MRS幅度衰減正弦信號(hào)與拉莫爾頻率的正弦信號(hào)具有相關(guān)性的特點(diǎn)，通過(guò)互相關(guān)運(yùn)算濾除噪聲，然后擬合互相關(guān)波形的包絡(luò)，并重構(gòu)不含噪聲的互相關(guān)波形，最后利用解卷積算法提取核磁共振全波數(shù)據(jù)中的MRS信號(hào)。該方法運(yùn)算數(shù)據(jù)計(jì)算量小，可以同時(shí)壓制工頻及其諧波噪聲、隨機(jī)噪聲和尖峰噪聲，明顯提高核磁共振全波數(shù)據(jù)信噪比。但是在地磁場(chǎng)不穩(wěn)定時(shí)或信噪比較低時(shí)，拉莫爾頻率不能準(zhǔn)確獲得，利用存在頻率偏差的信號(hào)進(jìn)行互相關(guān)計(jì)算后，信號(hào)參數(shù)提取誤差大。CN104636609A中公開(kāi)了“一種基于經(jīng)驗(yàn)?zāi)B(tài)分解與小波分析的信號(hào)聯(lián)合去噪方法”，根據(jù)信號(hào)的自相關(guān)性，對(duì)混有高斯白噪聲的信號(hào)進(jìn)行EMD分解，由于EMD分解的性質(zhì)，高斯白噪聲已不再是真正的白噪聲，但白噪聲的統(tǒng)計(jì)特性近似存在，即所述混有高斯白噪聲的信號(hào)的自相關(guān)函數(shù)在零點(diǎn)取得最大值，幅值隨著時(shí)間差的變化而變化，但其隨著時(shí)間的衰減很快。利用這種差異可以選取出噪聲起主導(dǎo)作用的IMF分量有效降低噪聲對(duì)信號(hào)的影響。但該方法主要是根據(jù)自相關(guān)特性對(duì)EMD分解后每個(gè)模態(tài)的屬性進(jìn)行判斷，主要用于低信噪比下去除高斯白噪聲，所以其應(yīng)用存在局限性。

技術(shù)實(shí)現(xiàn)要素：

本發(fā)明所要解決的技術(shù)問(wèn)題在于提供一種地面核磁共振信號(hào)提取方法，解決了磁共振測(cè)深找水工作中由于強(qiáng)工頻諧波干擾和隨機(jī)噪聲造成的MRS信號(hào)有效提取的問(wèn)題。

本發(fā)明是這樣實(shí)現(xiàn)的，一種地面核磁共振信號(hào)提取方法，該方法包括以下步驟：

步驟(1):利用地面核磁共振地下水探測(cè)儀器采集到一組觀測(cè)MRS含噪數(shù)據(jù)；

步驟(2):利用統(tǒng)計(jì)方法判斷是否存在尖峰噪聲，如果存在，去除尖峰噪聲并用插值結(jié)果代替，如果不存在，則測(cè)量數(shù)據(jù)保持不變；

步驟(3):將去除尖峰噪聲的數(shù)據(jù)，利用諧波建模的方法去除工頻諧波噪聲；

步驟(4):對(duì)步驟(3)的結(jié)果進(jìn)行自相關(guān)和疊加處理減小隨機(jī)噪聲；

步驟(5):對(duì)步驟(4)后的結(jié)果進(jìn)行MRS信號(hào)參數(shù)提取。

進(jìn)一步地，步驟(3)中的建模方法的具體步驟為：

步驟3a：設(shè)定工頻諧波基頻的搜索范圍以及步長(zhǎng)進(jìn)行粗掃；

步驟3b：將工頻諧波模型轉(zhuǎn)化為Ax＝b的形式，解包含測(cè)量數(shù)據(jù)的線性方程組，得到工頻各諧波頻點(diǎn)的系數(shù)矩陣；

步驟3c：計(jì)算工頻諧波基頻搜索范圍內(nèi)不同掃描值的模型估計(jì)值；

步驟3d：計(jì)算模型估計(jì)值與測(cè)量數(shù)據(jù)差值的2范數(shù)；

步驟3e：確定使誤差2范數(shù)取得最小值的頻點(diǎn)，從而獲得兩側(cè)相鄰的頻點(diǎn)；

步驟3f：將步驟3e中得到的兩側(cè)相鄰頻點(diǎn)作為下一次掃描范圍，選取適當(dāng)步長(zhǎng)再次掃描；

步驟3g：判斷是否掃描設(shè)定的次數(shù)，如果達(dá)到則根據(jù)誤差2范數(shù)最小的頻點(diǎn)進(jìn)行工頻諧波建模；否則重復(fù)步驟3c～3g；

步驟3h:從步驟(2)的結(jié)果中減去工頻諧波模型，去除工頻諧波噪聲。

進(jìn)一步地，步驟3f：步長(zhǎng)為(f_m2-f_m1)/M，其中f_m1，f_m2分別為最小值頻點(diǎn)的兩側(cè)相鄰頻點(diǎn)，其中M＝3～7。

進(jìn)一步地，步驟3g掃描設(shè)定的掃描次數(shù)為3-7次。

進(jìn)一步地，所述的步驟(4)中對(duì)步驟(3)得到的結(jié)果進(jìn)行自相關(guān)和疊加處理減小隨機(jī)噪聲的具體步驟為：

步驟4a：將MRS信號(hào)代入自相關(guān)函數(shù)公式，修改積分上下限，求得MRS信號(hào)的自相關(guān)表達(dá)式；

步驟4b：判斷采集時(shí)間是否大于1s，如果大于1s截取自相關(guān)結(jié)果的前半段用近似，其中e₀為初始振幅、為弛豫時(shí)間、f為拉莫爾頻率和τ為自相關(guān)時(shí)間間隔；如果小于1s自相0關(guān)結(jié)果用近似，其中t_max為時(shí)間最大值；

步驟4c：基于步驟4b的判斷后將多次測(cè)量數(shù)據(jù)的自相關(guān)結(jié)果進(jìn)行疊加，進(jìn)一步減小隨機(jī)噪聲。

進(jìn)一步地，所述的步驟(5)中對(duì)步驟(4)得到的結(jié)果進(jìn)行參數(shù)提取的具體步驟為：

步驟5a:對(duì)步驟(4)中疊加的數(shù)據(jù)進(jìn)行希爾伯特變換，再通過(guò)低通濾波，轉(zhuǎn)化為兩個(gè)正交分量；

步驟5b:利用非線性擬合方法求解MRS信號(hào)的特征參數(shù)A₀和弛豫時(shí)間其中再計(jì)算初始振幅e₀結(jié)果，參數(shù)提取結(jié)果為初始振幅e₀和弛豫時(shí)間

本發(fā)明與現(xiàn)有技術(shù)相比，有益效果在于：本發(fā)明提出了基于工頻諧波建模和自相關(guān)的地面核磁共振信號(hào)提取方法，針對(duì)單通道采集的全波MRS數(shù)據(jù)，不僅可以一次性去除所有諧波干擾，而且還結(jié)合了自相關(guān)方法壓制隨機(jī)噪聲，最后通過(guò)非線性擬合實(shí)現(xiàn)了有效MRS信號(hào)特征參數(shù)的提取。本發(fā)明方法解決了磁共振測(cè)深找水工作中由于強(qiáng)工頻諧波干擾和隨機(jī)噪聲造成的MRS信號(hào)有效提取的難題，獲取的MRS信號(hào)關(guān)鍵特征參數(shù)的擬合誤差較小，同時(shí)本發(fā)明突破了經(jīng)典消噪方法需多通道探測(cè)等其他條件的限制，節(jié)省了大量的人力物力。

附圖說(shuō)明

圖1為本發(fā)明地面核磁共振信號(hào)提取方法流程框圖；

圖2為工頻諧波建模消噪原理及效果圖，圖2(a)是其時(shí)域信號(hào)，MRS信號(hào)淹沒(méi)在噪聲之中；圖2(b)是其功率譜，在諧波的頻點(diǎn)上均存在不同幅度的諧波噪聲；圖2(c)諧波建模后，在測(cè)量數(shù)據(jù)b中減去V_harmonic完成消噪過(guò)程，2(d)諧波建模消噪后，MRS信號(hào)的輪廓呈現(xiàn)衰減趨勢(shì)，在功率譜中只存在MRS信號(hào)的頻點(diǎn)；

圖3為自相關(guān)消噪效果圖；圖3(a)中黑色線是理想的MRS信號(hào)，參數(shù)為e₀＝100nV，f＝2326Hz和灰色線是加入了100nV均方差的高斯白噪聲后的測(cè)量信號(hào)，此時(shí)數(shù)據(jù)的信噪比為0dB，MRS信號(hào)基本淹沒(méi)在噪聲中；圖3(b)是測(cè)量信號(hào)和MRS信號(hào)的頻譜，可見(jiàn)MRS信號(hào)的頻譜峰值是噪聲頻譜均值的10倍。3(c)經(jīng)過(guò)自相關(guān)函數(shù)后，測(cè)量信號(hào)(灰色線)和MRS信號(hào)(黑色線)的結(jié)果；3(d)MRS信號(hào)自相關(guān)后的頻譜峰值與噪聲頻譜；

圖4為自相關(guān)后非線性擬合效果圖，圖4(a)是未經(jīng)過(guò)自相關(guān)函數(shù)的曲線，圖4(b)是經(jīng)過(guò)自相關(guān)函數(shù)的曲線；

圖5為仿真數(shù)據(jù)和數(shù)據(jù)處理結(jié)果圖，圖5(a)含噪MRS信號(hào)；圖5(b)從含噪MRS信號(hào)中減去工頻諧波模型，得到去工頻信號(hào)圖；5(c)自相關(guān)和疊加處理得到的消噪后信號(hào)圖，5(d)MRS信號(hào)的實(shí)部，5(e)MRS信號(hào)的虛部，5(f)MRS信號(hào)實(shí)部的疊加曲線和虛部見(jiàn)；5(g)MRS信號(hào)虛部的疊加曲線；

圖6為實(shí)測(cè)數(shù)據(jù)處理結(jié)果圖；圖6(a-1)為q＝1.4As的MRS信號(hào)；6(a-2)為q＝4.2As的MRS信號(hào)；6(a-3)為q＝8.6As的MRS信號(hào)；

圖6(b-1)為q＝1.4As的信噪比信號(hào)；6(b-2)為q＝4.2As的的信噪比信號(hào)；6(b-3)為q＝8.6As的的信噪比信號(hào)。

具體實(shí)施方式

為了使本發(fā)明的目的、技術(shù)方案及優(yōu)點(diǎn)更加清楚明白，以下結(jié)合實(shí)施例，對(duì)本發(fā)明進(jìn)行進(jìn)一步詳細(xì)說(shuō)明。應(yīng)當(dāng)理解，此處所描述的具體實(shí)施例僅僅用以解釋本發(fā)明，并不用于限定本發(fā)明。

如圖1所示，一種基于工頻諧波建模和自相關(guān)的地面核磁共振信號(hào)提取方法，包括以下步驟：

步驟(1):利用地面核磁共振地下水探測(cè)儀器采集到一組觀測(cè)MRS含噪數(shù)據(jù)；

步驟(2):利用統(tǒng)計(jì)方法(如羅曼諾夫斯基準(zhǔn)則和能量運(yùn)算等算法)判斷是否存在尖峰噪聲，如果存在，去除尖峰噪聲并用插值結(jié)果代替，如果不存在，則測(cè)量數(shù)據(jù)保持不變；

步驟(3):利用諧波建模的方法去除工頻諧波噪聲；

步驟(4):對(duì)步驟(3)的結(jié)果進(jìn)行自相關(guān)和疊加處理減小隨機(jī)噪聲；

步驟(5):對(duì)步驟(4)后的結(jié)果進(jìn)行MRS信號(hào)參數(shù)提??；

如圖2所示，步驟(3)中的工頻建模方法的具體步驟：

步驟3a:設(shè)定工頻諧波基頻f₀的搜索范圍49.9Hz～50.1Hz，f₀以步長(zhǎng)為0.03Hz進(jìn)行粗掃；

步驟3b:工頻諧波噪聲V_harmonic可表示為式(1)的形式

其中，A_n和分別是第n個(gè)諧波的幅度和相位，并有和t是時(shí)間。將表達(dá)式(1)整理成線性矩陣形式Ax＝b，

其中，b＝[V₁,V₂,…,V_P]^T中的V_p(p＝1,2,…P)是t_p時(shí)刻的接收數(shù)據(jù)；x＝[α₁,…,α_N,β₁,…,β_N]^T是諧波系數(shù)，N是諧波個(gè)數(shù)，根據(jù)接收信號(hào)的帶寬(1～3kHz)取N＝100。求解線性方程組(2)得到諧波系數(shù)矩陣x。

步驟3c:分別計(jì)算頻率搜索范圍內(nèi)不同f₀的模型估計(jì)值V_harmonic(f₀)；

步驟3d:計(jì)算V_harmonic(f₀)與測(cè)量數(shù)據(jù)b的差值(誤差)的2范數(shù)：||b-V_harmonic(f₀)||₂；

步驟3e:確定使誤差2范數(shù)取得最小值時(shí)對(duì)應(yīng)的頻點(diǎn)f_min，從而獲得兩側(cè)相鄰的頻點(diǎn)(f_m1,f_m2)；

步驟3f:確定下一次掃描范圍(f_m1,f_m2)，f₀選取步長(zhǎng)為(f_m2-f_m1)/M，其中M＝3～7；

步驟g:判斷是否掃描5次(可以設(shè)置為3-7次，每次搜索的工頻基頻頻率相對(duì)變化量小于0.1％即可)，如果達(dá)到，則根據(jù)誤差2范數(shù)最小時(shí)對(duì)應(yīng)的頻點(diǎn)f_min進(jìn)行工頻諧波建模；否則重復(fù)步驟3c～3g；

步驟h:從步驟(2)的結(jié)果中減去工頻諧波模型V_harmonic(f_min)，從而去除工頻諧波噪聲；

如圖2所示的一組接收信號(hào)數(shù)據(jù)，圖2(a)是其時(shí)域信號(hào)，MRS信號(hào)淹沒(méi)在噪聲之中；圖2(b)是其功率譜，在諧波的頻點(diǎn)上均存在不同幅度的諧波噪聲。諧波建模后，在測(cè)量數(shù)據(jù)b中減去V_harmonic即完成了消噪過(guò)程，如圖2(c)和2(d)所示。諧波建模消噪后，MRS信號(hào)的輪廓呈現(xiàn)衰減趨勢(shì)，在功率譜中只存在MRS信號(hào)的頻點(diǎn)。

本發(fā)明提供的自相關(guān)和疊加處理的具體步驟為：

步驟4a：MRS信號(hào)用V_MRS表示：

其中，e₀表示初始振幅；表示平均弛豫時(shí)間，f表示拉莫爾頻率(與當(dāng)?shù)氐卮艌?chǎng)有關(guān))，表示初始相位。

根據(jù)自相關(guān)原理，將V_MRS代入自相關(guān)函數(shù)的表達(dá)式中：

其中，τ為時(shí)間間隔。由于自相關(guān)函數(shù)是關(guān)于自變量τ的實(shí)偶函數(shù)，所以τ只取其正半軸。表達(dá)式(4)中的積分上下限分別取τ和t_max，t_max是測(cè)量時(shí)間長(zhǎng)度，通過(guò)推導(dǎo)得到V_MRS的自相關(guān)表達(dá)式為：

R(τ)＝RC1(τ)+RC2(τ)+RC3(τ)+RC4(τ)

其中，

由表達(dá)式(5)～(7)可知，MRS信號(hào)自相關(guān)函數(shù)包含四個(gè)分量RC1～RC4，MRS信號(hào)自相關(guān)后的前半部分可以只用RC1代替；后半部分可以用RC1+RC2代替；RC3+RC4的結(jié)果很小，可以忽略。

步驟4b：判斷采集時(shí)間是否大于1s，如果大于1s截取自相關(guān)結(jié)果的前半段用RC1近似；如果小于1s自相關(guān)結(jié)果用RC1+RC2近似。

步驟4c：基于步驟4b的判斷后將多次測(cè)量數(shù)據(jù)b的自相關(guān)結(jié)果R進(jìn)行疊加，即其中N是測(cè)量次數(shù)。N次疊加后信噪比提高倍，進(jìn)一步減小隨機(jī)噪聲。

圖3給出了含噪MRS信號(hào)經(jīng)自相關(guān)和疊加處理的結(jié)果。圖3(a)中黑色線是理想的MRS信號(hào)，參數(shù)為e₀＝100nV，f＝2326Hz和灰色線是加入了100nV均方差的高斯白噪聲后的測(cè)量信號(hào)，此時(shí)數(shù)據(jù)的信噪比為0dB，MRS信號(hào)基本淹沒(méi)在噪聲中。圖3(b)是測(cè)量信號(hào)和MRS信號(hào)的頻譜，可見(jiàn)MRS信號(hào)的頻譜峰值是噪聲頻譜均值的10倍。經(jīng)過(guò)自相關(guān)函數(shù)后，測(cè)量信號(hào)(灰色線)和MRS信號(hào)(黑色線)的結(jié)果見(jiàn)圖3(c)所示。此時(shí)，灰色線與黑色線類似，也呈現(xiàn)明顯的衰減趨勢(shì)，只出現(xiàn)了少量的隨機(jī)噪聲，信噪比為23.75dB，與自相關(guān)之前相比大幅提高。從圖3(d)中的頻譜分析也可以得出，MRS信號(hào)自相關(guān)后的頻譜峰值約為噪聲頻譜均值的100倍，與自相關(guān)之前提高了一個(gè)數(shù)量級(jí)。因此，自相關(guān)函數(shù)能夠有效壓制隨機(jī)噪聲，提高信噪比。

如圖4所示，本發(fā)明提供的參數(shù)提取方法的具體步驟為：

步驟5a:對(duì)疊加的數(shù)據(jù)進(jìn)行希爾伯特變換，希爾伯特變換表達(dá)式為：可將實(shí)信號(hào)R(t)轉(zhuǎn)化為復(fù)信號(hào)其包絡(luò)信號(hào)為再通過(guò)低通濾波，得到兩個(gè)正交分量V_RC＝V_real+i·V_imag；

步驟5b:利用非線性擬合方法對(duì)步驟a得到的兩個(gè)正交分量進(jìn)行參數(shù)提取如式(8)，

獲得MRS信號(hào)的特征參數(shù)A₀和其中再計(jì)算e₀結(jié)果，參數(shù)提取結(jié)果為初始振幅e₀和弛豫時(shí)間

圖3中的數(shù)據(jù)經(jīng)過(guò)Hilbert變換后的兩個(gè)分量見(jiàn)圖4所示，其中圖4(a)是未經(jīng)過(guò)自相關(guān)函數(shù)的曲線，圖4(b)是經(jīng)過(guò)自相關(guān)函數(shù)的曲線。分別利用非線性擬合方法求解得到的兩個(gè)分量的擬合曲線見(jiàn)圖4中黑色實(shí)線和黑色虛線所示。未經(jīng)過(guò)自相關(guān)函數(shù)的數(shù)據(jù)信噪比低，參數(shù)提取結(jié)果分別為e₀＝108.3nV，f＝2326.1Hz；而經(jīng)過(guò)自相關(guān)函數(shù)的數(shù)據(jù)信噪比提高明顯，參數(shù)提取結(jié)果為e₀＝100.5nV，f＝2326.0Hz，誤差分別降低了7.8％、10.0％和0.01％。因此，自相關(guān)函數(shù)能夠提高數(shù)據(jù)的信噪比，從而提高了參數(shù)提取的準(zhǔn)確度。

實(shí)施例1

基于工頻諧波建模和自相關(guān)的地面核磁共振信號(hào)提取方法，本發(fā)明模擬了一組包含16次獨(dú)立采集的測(cè)量數(shù)據(jù)集，參照?qǐng)D1，包括以下步驟：

步驟(1)：利用式(3)構(gòu)造含噪MRS信號(hào)如圖5(a)所示，其參數(shù)為e₀＝100nV，f＝2326Hz和每次采集的諧波噪聲基頻在49.9～50.1Hz隨機(jī)產(chǎn)生，諧波個(gè)數(shù)為100，諧波幅度在200nV內(nèi)隨機(jī)分布；每次采集的隨機(jī)噪聲為200nV的高斯白噪聲；采集數(shù)據(jù)中未加入尖峰噪聲，假設(shè)所有的尖峰噪聲已通過(guò)已有的方法消除。

步驟(2)：利用諧波建模的方法得到工頻諧波模型其中f₀＝50.018125Hz，從含噪MRS信號(hào)中減去工頻諧波模型，得到去工頻信號(hào)，如圖5(b)所示，可見(jiàn)工頻諧波消噪后，數(shù)據(jù)的信噪比大幅提升，平均提高了17.03dB，數(shù)據(jù)幅度減小到±500nV范圍內(nèi)；

步驟(3):對(duì)步驟(2)的結(jié)果進(jìn)行自相關(guān)和疊加處理，減小隨機(jī)噪聲，得到消噪后信號(hào)，如圖5(c)所示，可見(jiàn)自相關(guān)處理后，數(shù)據(jù)呈現(xiàn)明顯的衰減趨勢(shì)，除了信號(hào)初期存在脈沖(隨機(jī)噪聲的自相關(guān)特性)，其他時(shí)刻的數(shù)據(jù)幅度比自相關(guān)之前的數(shù)據(jù)明顯提高而隨機(jī)噪聲的影響較小，因此信噪比有較大程度的提升，為15.64～16.43dB，信噪比平均提高了16.10dB；

步驟(4):對(duì)步驟(3)的結(jié)果進(jìn)行MRS信號(hào)參數(shù)提取，利用希爾伯特變換(參考頻率2325Hz)和低通濾波(截至頻率200Hz)，得到MRS信號(hào)的實(shí)部和虛部見(jiàn)圖5(d)和5(e)所示。建模消噪處理后，由于還存在大量的隨機(jī)噪聲，數(shù)據(jù)的實(shí)部和虛部曲線區(qū)分不明顯，只能觀察到曲線的整體趨勢(shì)稍有起伏(1Hz的震蕩)。相比自相關(guān)處理后，數(shù)據(jù)的實(shí)部和虛部能夠清晰地區(qū)分開(kāi)，兩者都呈現(xiàn)出1Hz的震蕩衰減過(guò)程，再對(duì)16次采集和處理的數(shù)據(jù)進(jìn)行疊加，得到MRS信號(hào)實(shí)部和虛部的疊加曲線見(jiàn)圖5(f)和(g)所示。疊加處理能夠進(jìn)一步降低隨機(jī)噪聲，16次疊加信噪比提高4倍(12.04dB)，利用非線性擬合方法對(duì)兩組數(shù)據(jù)進(jìn)行MRS信號(hào)參數(shù)提取，通過(guò)建模消噪和自相關(guān)處理后的參數(shù)提取結(jié)果為：e₀＝100.45nV，f＝2326.0Hz，相對(duì)誤差分別為：0.45％，-2.14％和0，均控制在±3％以內(nèi)，滿足應(yīng)用要求。

實(shí)施例2

本實(shí)施例以長(zhǎng)春市農(nóng)安縣燒鍋鎮(zhèn)的太平池水庫(kù)旁(該地地磁場(chǎng)為54720nT，對(duì)應(yīng)的Larmor頻率在2330Hz，變化范圍小于2Hz)利用自主研制的地面核磁共振地下水探測(cè)儀器進(jìn)行了野外數(shù)據(jù)采集與數(shù)據(jù)處理實(shí)驗(yàn)。為了從淺到深地對(duì)含水層進(jìn)行探測(cè)，在0.2～8.5As的范圍內(nèi)由小到大按對(duì)數(shù)分布設(shè)置了20組發(fā)射脈沖矩，每組脈沖矩重復(fù)疊加16次，采樣率為50kHz，采集時(shí)間為1s。將第10、16和20個(gè)脈沖矩采集的MRS信號(hào)作為本發(fā)明方法的處理對(duì)象。如圖1所示，基于工頻諧波建模和自相關(guān)的地面核磁共振信號(hào)提取方法，包括以下步驟:

步驟(1):利用地面核磁共振地下水探測(cè)儀器采集到一組觀測(cè)MRS含噪數(shù)據(jù)；

步驟(2):利用統(tǒng)計(jì)方法判斷是否存在尖峰噪聲，如果存在，去除尖峰噪聲并用插值結(jié)果代替，如果不存在，則測(cè)量數(shù)據(jù)保持不變；

步驟(3)：利用諧波建模的方法得到模型其中f₀＝50.018125Hz，從步驟(2)的結(jié)果中減去工頻諧波模型，得到去工頻信號(hào)，如圖6(a-1)為q＝1.4As的MRS信號(hào)；6(a-2)為q＝4.2As的MRS信號(hào)；6(a-3)為q＝8.6As的MRS信號(hào)所示，可見(jiàn)實(shí)測(cè)數(shù)據(jù)中存在大量的諧波噪聲，噪聲水平在2500～3000nV范圍內(nèi)，而MRS信號(hào)幅度在80～375nV，因此信噪比較低。經(jīng)過(guò)工頻諧波消噪處理后，大多數(shù)諧波噪聲顯著降低，諧波基頻在49.975～50.041Hz之間，MRS信號(hào)(2330Hz)的峰值比較明顯，其中q＝1.4As的信號(hào)較大，而q＝4.2As和8.6As的信號(hào)較小，信噪比提高了近20dB；

步驟(4):對(duì)步驟(3)的結(jié)果進(jìn)行自相關(guān)和疊加處理，減小隨機(jī)噪聲，得到消噪后信號(hào)；

步驟(5):對(duì)步驟(4)后的結(jié)果進(jìn)行MRS信號(hào)參數(shù)提取，如圖6(b-1)，6(b-2)，6(b-3))所示，3個(gè)脈沖矩下信噪比增量分別為12.05dB、2.98dB和7.71dB，此時(shí)經(jīng)過(guò)非線性擬合后，參數(shù)提取結(jié)果為：

q＝4.2As時(shí)，e₀擬合結(jié)果為80.7nV，不準(zhǔn)確度為15.9nV，而擬合結(jié)果為235.9ms，不準(zhǔn)確度為38.2ms；q＝1.4As時(shí)，e₀擬合結(jié)果為396.8nV，不準(zhǔn)確度為17.4nV，而擬合結(jié)果為248.6ms，不準(zhǔn)確度為10.8ms；q＝8.6As時(shí)，e₀擬合結(jié)果為175.6nV，不準(zhǔn)確度為13.4nV，而擬合結(jié)果為263.1ms，不準(zhǔn)確度為33.1ms，說(shuō)明SNR的增加，能夠提高非線性擬合結(jié)果的準(zhǔn)確度。

以上所述僅為本發(fā)明的較佳實(shí)施例而已，并不用以限制本發(fā)明，凡在本發(fā)明的精神和原則之內(nèi)所作的任何修改、等同替換和改進(jìn)等，均應(yīng)包含在本發(fā)明的保護(hù)范圍之內(nèi)。