燕 歡， 李 晉,2， 李春陽(yáng)

(1.湖南師范大學(xué) 物理與信息科學(xué)學(xué)院，長(zhǎng)沙 410081; 2.中南大學(xué) 有色金屬成礦預(yù)測(cè)與地質(zhì)環(huán)境監(jiān)測(cè)教育部重點(diǎn)實(shí)驗(yàn)室， 地球科學(xué)與信息物理學(xué)院，長(zhǎng)沙 410083)

子空間增強(qiáng)在大地電磁信噪分離中的應(yīng)用研究

燕 歡1， 李 晉1,2， 李春陽(yáng)1

(1.湖南師范大學(xué) 物理與信息科學(xué)學(xué)院，長(zhǎng)沙 410081; 2.中南大學(xué) 有色金屬成礦預(yù)測(cè)與地質(zhì)環(huán)境監(jiān)測(cè)教育部重點(diǎn)實(shí)驗(yàn)室， 地球科學(xué)與信息物理學(xué)院，長(zhǎng)沙 410083)

為了提高強(qiáng)干擾地區(qū)大地電磁測(cè)深數(shù)據(jù)質(zhì)量、保留大地電磁信號(hào)低頻段的有用信息，壓制大尺度強(qiáng)干擾已刻不容緩。針對(duì)礦集區(qū)典型的強(qiáng)干擾類(lèi)型，引入子空間增強(qiáng)方法對(duì)大地電磁信號(hào)與強(qiáng)干擾進(jìn)行信噪分離研究。通過(guò)計(jì)算機(jī)模擬常見(jiàn)的大尺度方波和充放電三角波干擾，討論了子空間增強(qiáng)算法中幀長(zhǎng)與特征值判別閾值的最優(yōu)選取范圍，給出了算法流程，并對(duì)礦集區(qū)實(shí)測(cè)大地電磁數(shù)據(jù)進(jìn)行了信噪分離處理。實(shí)驗(yàn)結(jié)果表明：該方法是切實(shí)可行的，能有效剔除時(shí)間域序列中的大尺度強(qiáng)干擾，近源效應(yīng)得到了有效壓制；經(jīng)處理后，重構(gòu)的大地電磁信號(hào)中包含了更為豐富的細(xì)節(jié)成分，視電阻率曲線(xiàn)更為光滑、連續(xù)，低頻段的大地電磁數(shù)據(jù)質(zhì)量得到了明顯改善。

大地電磁； 信噪分離； 子空間增強(qiáng)； 強(qiáng)干擾

0 引言

大地電磁測(cè)深(Magnetotelluric，簡(jiǎn)稱(chēng)MT)法是由Tikhonov和Cagniard[1-2]提出的以天然交變電磁場(chǎng)為場(chǎng)源的電磁勘探方法。由于天然電磁場(chǎng)信號(hào)弱、頻帶寬，因此極易受到各種噪聲的干擾，是典型的非線(xiàn)性、非平穩(wěn)信號(hào)[3-5]。隨著人類(lèi)文明的不斷進(jìn)步和經(jīng)濟(jì)的高速發(fā)展，人文電磁干擾日益嚴(yán)重，尤其在礦集區(qū)，電話(huà)線(xiàn)路、無(wú)線(xiàn)通訊基站、采礦用的大功率直流電機(jī)車(chē)等引起強(qiáng)能量的電磁干擾，導(dǎo)致原本微弱的大地電磁信號(hào)幾乎完全被湮沒(méi)在這些“大尺度”干擾中。諸多現(xiàn)代信號(hào)處理方法(如Hilbert-Huang變換[6]、廣義S變換[7]、方差比維納濾波[8]、數(shù)學(xué)形態(tài)濾波[9]、同步時(shí)間序列[10]等),被應(yīng)用到大地電磁噪聲壓制，均在一定程度上改善了大地電磁測(cè)深數(shù)據(jù)質(zhì)量，推動(dòng)了大地電磁測(cè)深法在各種噪聲中的發(fā)展。

子空間增強(qiáng)是一種非線(xiàn)性、非平穩(wěn)信號(hào)分析方法，上世紀(jì)90年代該方法已被應(yīng)用于語(yǔ)音增強(qiáng)領(lǐng)域。由于語(yǔ)音信號(hào)的噪聲方差可以估計(jì)，根據(jù)子空間增強(qiáng)算法原理，帶噪語(yǔ)音信號(hào)的向量空間可分解為信號(hào)子空間和噪聲子空間，通過(guò)剔除噪聲子空間，再利用信號(hào)子空間即可恢復(fù)較為純凈的語(yǔ)音信號(hào)。子空間增強(qiáng)法現(xiàn)已被逐步推廣到故障診斷、地震資料處理等領(lǐng)域[11-12]。由于大地電磁信號(hào)是典型的非線(xiàn)性、非平穩(wěn)信號(hào)，為此借助語(yǔ)音增強(qiáng)、地震資料處理等領(lǐng)域廣泛應(yīng)用的子空間增強(qiáng)法，提出分離礦集區(qū)微弱大地電磁信號(hào)和典型強(qiáng)干擾的研究思路，分析了關(guān)鍵參數(shù)的最優(yōu)選取方案，給出了具體的算法處理流程。

1 子空間增強(qiáng)基本原理

子空間分解的降噪方法早期由Dendrinos[13]提出，該算法利用帶噪信號(hào)組成具有Toeplitz結(jié)構(gòu)的矩陣，隨后對(duì)該矩陣進(jìn)行奇異值分解；根據(jù)最小二乘估計(jì)，忽略較小的奇異值后可將降維得到的數(shù)據(jù)根據(jù)Toeplitz矩陣結(jié)構(gòu)恢復(fù)出原信號(hào)。然而，子空間增強(qiáng)法具有一定的局限性，只能處理白噪聲、不能處理有色噪聲。在此基礎(chǔ)上，Ephraim等[14]做了進(jìn)一步完善，利用信號(hào)協(xié)方差矩陣的特征值分解(Eigen value decomposition，EVD)，將帶噪信號(hào)的向量空間分解為噪聲子空間和信號(hào)子空間兩個(gè)相互正交的子空間再進(jìn)行語(yǔ)音增強(qiáng)處理。隨后，Jensen等[15]將基于SVD(Singular value decomposition，SVD)的子空間算法擴(kuò)展到有色噪聲領(lǐng)域。

子空間增強(qiáng)具有控制信號(hào)失真和噪聲殘留的特征，其基本原理是將受干擾的帶噪信號(hào)映射到兩個(gè)子空間①信號(hào)子空間；②噪聲子空間，通過(guò)數(shù)學(xué)方法將噪聲子空間全部置零，濾除信號(hào)子空間中所包含的噪聲干擾，再利用信號(hào)子空間恢復(fù)出較為純凈的原始信號(hào)[16]。

以一維離散信號(hào)為例，假定噪聲是加性的，且與純凈信號(hào)相互獨(dú)立，則帶噪信號(hào)可表示為：

y=x+n

(1)

式中：y、x和n分別表示帶噪信號(hào)、純凈信號(hào)和加性噪聲的向量表示,對(duì)其進(jìn)行分幀，可得到相應(yīng)的K維向量，其帶噪信號(hào)、純凈信號(hào)和噪聲的協(xié)方差矩陣Ry、Rx和Rn關(guān)系如下：

Ry=E[yyT]=E[xxT]+E[nnT]=

Rx+Rn

(2)

式中:E[·]表示均值運(yùn)算。

假設(shè)純凈信號(hào)估計(jì)為：

(3)

式中:H表示K×K維的線(xiàn)性估計(jì)器矩陣，則該估計(jì)器的真實(shí)值與估計(jì)值的誤差ε為：

εx+εn

(4)

式中：εx和εn分別表示信號(hào)的失真向量和殘余噪聲向量。子空間增強(qiáng)即通過(guò)線(xiàn)性估計(jì)投影到信號(hào)子空間，從而得到原始信號(hào)的較好估計(jì)。

tr(HRxHT-HRx-RxHT+Rx)

(5)

tr(HRnHT)

(6)

通過(guò)求解以下時(shí)域約束條件方程，可獲得優(yōu)化的線(xiàn)性估計(jì)器：

(7)

式中:σ2為正常數(shù)。在可接受的噪聲殘差水平σ2下，該估計(jì)器能最小化信號(hào)失真。

上述方程的解為：

Hopt=Rx(Rx+μRn)-1

(8)

式中:μ為拉格朗日乘子。

由Rx的特征值分解為：

Rx=UΛxUT

(9)

式中：U表示Rx的歸一化特征向量矩陣；Λx表示由Rx的特征值組成的對(duì)角矩陣：

(10)

Hopt=UΛx(Λx+μUTRnU)-1UT

(11)

當(dāng)噪聲n為方差σn的加性白噪聲時(shí)，Rn可表示為：

Rn=σnI

(12)

當(dāng)噪聲n為有色噪聲時(shí)，往往用對(duì)角矩陣Λn來(lái)近似UTRnU：

(13)

為此，Hopt可優(yōu)化為式(14)。

Hopt=UΛx(Λx+μΛn)-1UT

(14)

2 仿真信號(hào)分析

2.1 子空間增強(qiáng)關(guān)鍵參數(shù)分析

觀測(cè)礦集區(qū)海量大地電磁數(shù)據(jù)的時(shí)間域波形可知，大尺度方波和充放電三角波干擾通常出現(xiàn)在各種采樣數(shù)據(jù)中，是最為典型且嚴(yán)重的噪聲干擾類(lèi)型[17]。為了驗(yàn)證文中所提方法的可行性，計(jì)算機(jī)模擬大尺度方波和充放電三角波進(jìn)行仿真實(shí)驗(yàn)。大尺度方波通過(guò)選取Matlab自帶的block信號(hào)疊加4dB的白噪聲來(lái)模擬(圖1)。大尺度充放電三角波通過(guò)選取幅值為5%、峰谷距為10、間距為100的充放電三角波疊加隨機(jī)噪聲來(lái)模擬(圖2)。

圖1 模擬大尺度方波干擾Fig.1 Simulated large-scale square interferences

圖2 模擬大尺度充放電三角波干擾Fig.2 Simulated large-scale charging and discharging triangular interferences

由于子空間增強(qiáng)算法中分幀幀長(zhǎng)和協(xié)方差矩陣特征值的判別閾值對(duì)濾波精度尤為關(guān)鍵，為此引入波形相似度對(duì)該關(guān)鍵參數(shù)進(jìn)行性能評(píng)價(jià)。相似度的定義為式(15)。

(15)

式中:f(n)、g(n)分別表示兩個(gè)離散序列。NCC的值在-1～1之間，-1表示變換前后波形反向，“0”表示正交，“1”表示完全相同；NCC越接近“1”，則表示波形之間越相似。

從圖3可知，大尺度方波的幀長(zhǎng)取值在0～8范圍內(nèi)相似度逐漸上升；當(dāng)幀長(zhǎng)取值在8～18時(shí)相似度趨于平穩(wěn)，且達(dá)到最大值0.9；幀長(zhǎng)取值在18以后，相似度逐漸下降。大尺度充放電三角波的幀長(zhǎng)取值在0～10范圍時(shí)，相似度逐漸上升；在10～18時(shí)，相似度接近于0.8；幀長(zhǎng)為18以后，相似度逐漸下降。分析圖4可知，大尺度方波的判別閾值在10～100時(shí)，相似度在0.75～0.8之間；大尺度充放電三角波具有同樣的變化趨勢(shì)，其相似度在0.9～0.95之間。由于模擬的大尺度強(qiáng)干擾相比實(shí)測(cè)大地電磁信號(hào)要單一得多，而判別閾值本身即為噪聲殘留與信號(hào)失真的折中，且每一幀噪聲的特征值變化沒(méi)有實(shí)測(cè)信號(hào)復(fù)雜。綜上分析，當(dāng)幀長(zhǎng)選取10～18、判別閾值選取40～60時(shí)，濾波效果最優(yōu)。

圖3 不同幀長(zhǎng)濾波性能對(duì)比Fig.3 Filtering performance comparison of different frames

圖4 不同判別閾值濾波性能對(duì)比Fig.4 Filtering performance comparison of different criterion thresholds

2.2 算法流程

根據(jù)子空間增強(qiáng)基本原理，將該算法應(yīng)用到模擬大尺度信噪分離，其具體步驟如下：

1)將模擬帶噪信號(hào)分段，每段為一幀，幀長(zhǎng)在最優(yōu)范圍內(nèi)選取。

2)構(gòu)造協(xié)方差矩陣，對(duì)其進(jìn)行特征值分解。

3)假設(shè)大尺度噪聲與微弱信號(hào)互不相關(guān)，設(shè)置最優(yōu)判別閾值，分離出噪聲子空間和信號(hào)子空間。

4)將噪聲子空間的特征值置零，重構(gòu)原始信號(hào)。

圖5所示為模擬大尺度干擾經(jīng)子空間增強(qiáng)后的去噪效果圖。文中幀長(zhǎng)為12，判別閾值為40。

從圖5可知，疊加在大尺度干擾上的信號(hào)極其微弱，經(jīng)文中所提方法處理后，提取的大尺度干擾輪廓較為光滑、連續(xù)，重構(gòu)信號(hào)中較好地保留了更為豐富的細(xì)節(jié)成分。

圖5 模擬大尺度干擾去噪效果圖Fig.5 Denoising effect of simulated large-scale interferences(a)模擬大尺度干擾(方波);(b)子空間分離出的大尺度干擾(方波);(c)重構(gòu)信號(hào)(方波);(d)模擬大尺度干擾(充放電三角波);(e)子空間分離出的大尺度干擾(充放電三角波);(f)重構(gòu)信號(hào)(充放電三角波)

3 實(shí)測(cè)數(shù)據(jù)處理

3.1 時(shí)間域去噪效果

為了驗(yàn)證子空間增強(qiáng)法的實(shí)用性，對(duì)礦集區(qū)受強(qiáng)干擾嚴(yán)重的大地電磁數(shù)據(jù)進(jìn)行噪聲壓制。由于礦集區(qū)大地電磁數(shù)據(jù)量龐大，文中給出某測(cè)點(diǎn)在同一時(shí)刻電道(EX、EY)和磁道(HY、HX)中四道時(shí)間域片段的去噪效果，如圖6所示。

從圖6可知，由于該測(cè)點(diǎn)所處地區(qū)噪聲干擾源眾多，導(dǎo)致時(shí)間域波形中出現(xiàn)大量能量強(qiáng)、頻帶寬的大尺度類(lèi)方波、類(lèi)充放電三角波等噪聲干擾；波形的整體形態(tài)呈不連續(xù)狀態(tài)、跳躍性強(qiáng)，且電道EX/EY和磁道HY/HX出現(xiàn)干擾的時(shí)刻具有一定的相關(guān)性，微弱的大地電磁有用信號(hào)幾乎被完全湮沒(méi)。經(jīng)子空間增強(qiáng)法處理后，大尺度干擾的輪廓特征幾乎被完整地提取出來(lái)，分離出的大地電磁微弱信號(hào)在剔除基線(xiàn)漂移的同時(shí)，保留了有用信號(hào)更多的細(xì)節(jié)信息。由于實(shí)測(cè)大地電磁信號(hào)所含噪聲種類(lèi)眾多，即使是同種形態(tài)的噪聲其幅值、間距也是多種多樣的。在如此復(fù)雜的干擾環(huán)境下，子空間增強(qiáng)算法仍能較清晰地提取出疊加在原始大地電磁信號(hào)上的大尺輪廓，說(shuō)明該方法在時(shí)間域?qū)Υ蟪叨雀蓴_的分離是可行的。

3.2 視電阻率曲線(xiàn)

選用加拿大鳳凰公司生產(chǎn)的V5-2000大地電磁測(cè)深儀采集的數(shù)據(jù)進(jìn)行分析，該儀器采樣頻率通常為高頻(2 560 Hz)、中頻(320 Hz)和低頻(24 Hz)，數(shù)據(jù)的存儲(chǔ)格式為T(mén)SH、TSL，其中TSH文件記錄了2 560 Hz和320 Hz采樣率的中高頻數(shù)據(jù)，TSL文件記錄了24 Hz采樣率的低頻數(shù)據(jù)[18]。由于高頻和中頻是交替采集，每5 min采集1次高頻或中頻，其中只有1 s的高頻數(shù)據(jù)和連續(xù)8 s的中頻數(shù)據(jù)，所含信息量較少；低頻數(shù)據(jù)為全時(shí)間段采集，能記錄豐富的大地電磁數(shù)據(jù)，且更能反映測(cè)點(diǎn)深部的電性結(jié)構(gòu)信息，為此我們重點(diǎn)對(duì)低頻信號(hào)(TSL文件)進(jìn)行去噪處理。

圖7所示為廬縱礦集區(qū)B40993測(cè)點(diǎn)的原始時(shí)間域波形(EX、EY、HX、HY)和視電阻率曲線(xiàn)。

從圖7(a)可知，原始數(shù)據(jù)電道EX、EY中出現(xiàn)大量大尺度類(lèi)方波干擾，磁道HX、HY中出現(xiàn)大量大尺度類(lèi)充放電三角波干擾。分析圖7(b)可知，原始數(shù)據(jù)的視電阻率曲線(xiàn)整體形態(tài)連續(xù)性差。在大于40 Hz時(shí)，視電阻率曲線(xiàn)比較平穩(wěn)，且變化趨勢(shì)一致；在0.05 Hz～50 Hz，視電阻率曲線(xiàn)呈45°左右漸近線(xiàn)快速上升，0.05 Hz左右視電阻率值達(dá)到最大值接近于100 000 Ω·m，表現(xiàn)為典型的近源效應(yīng)；在0.005 Hz～0.05 Hz，視電阻率曲線(xiàn)突然下降至100 Ω·m左右，誤差棒增大，并出現(xiàn)不同程度的突跳與畸變。分析原始數(shù)據(jù)的時(shí)間域波形和視電阻率曲線(xiàn)可知，該測(cè)點(diǎn)受到了各種復(fù)雜干擾源的影響，其數(shù)據(jù)已不能客觀反映真實(shí)的地電信息。雖然功率譜篩選具有較強(qiáng)的適用性，但對(duì)礦集區(qū)高強(qiáng)度近源干擾也無(wú)能為力。

圖6 實(shí)測(cè)大地電磁數(shù)據(jù)去噪效果Fig.6 Denoising effect of measured magnetotelluric data(a)原始實(shí)測(cè)信號(hào)EX;(b)子空間分離出的大尺度干擾EX;(c)重構(gòu)大地電磁信號(hào)EX;(d)原始實(shí)測(cè)信號(hào)EY;(e)子空間分離出的大尺度干擾EY;(f)重構(gòu)大地電磁信號(hào)EY;(g)原始實(shí)測(cè)信號(hào)HY;(h)子空間分離出的大尺度干擾HY;(i)重構(gòu)大地電磁信號(hào)HY;(j)原始實(shí)測(cè)信號(hào)HX;(k)子空間分離出的大尺度干擾HX;(l)重構(gòu)大地電磁信號(hào)HX

圖8為子空間增強(qiáng)法得到的時(shí)間域波形和視電阻率曲線(xiàn)。其中，圖8(a)為經(jīng)子空間增強(qiáng)法處理后，重新還原至V5 System 2000觀測(cè)的時(shí)間域波形。

從圖8(a)可知，充斥在時(shí)間域中的大尺度類(lèi)方波和類(lèi)充放電三角波已被基本剔除、基線(xiàn)漂移得到了有效壓制，處理后的時(shí)間域波形基本連續(xù)、無(wú)明顯跳變。分析圖8(b)可知，經(jīng)子空間增強(qiáng)法處理后，0.05 Hz～50 Hz呈近45°上升的近源趨勢(shì)完全消失，曲線(xiàn)變得光滑、連續(xù)、同步；在0.05 Hz～5 Hz，視電阻率最大值下降了近3個(gè)數(shù)量級(jí)，近15個(gè)頻點(diǎn)的數(shù)據(jù)質(zhì)量得到了明顯改善。分析對(duì)比圖7可知，時(shí)、頻域的處理更加真實(shí)、合理、置信，其結(jié)果已逐步接近于測(cè)點(diǎn)本身所固有的地下介質(zhì)電性結(jié)構(gòu)。

圖9所示為在子空間增強(qiáng)法基礎(chǔ)上做簡(jiǎn)單的功率譜篩選得到的視電阻率曲線(xiàn)。

分析圖9可知，經(jīng)簡(jiǎn)單功率譜篩選后，視電阻率曲線(xiàn)在0.000 5 Hz～0.05 Hz的突跳頻點(diǎn)得到了較好恢復(fù)，曲線(xiàn)下掉的趨勢(shì)抬升顯著，且誤差棒減??；視電阻率曲線(xiàn)的整體形態(tài)更為光滑、連續(xù)，視電阻率曲線(xiàn)相對(duì)穩(wěn)定。實(shí)驗(yàn)結(jié)果表明，經(jīng)子空間增強(qiáng)法處理后，受近源干擾嚴(yán)重的測(cè)點(diǎn)其低頻段的整體數(shù)據(jù)質(zhì)量得到了明顯改善。

圖7 測(cè)點(diǎn)B40993原始大地電磁數(shù)據(jù)Fig.7 Original MT data at site B40993(a)時(shí)間域波形;(b)視電阻率曲線(xiàn)

圖8 子空間增強(qiáng)方法Fig.8 The subspace enhancement method(a)時(shí)間域波形;(b)視電阻率曲線(xiàn)

圖9 經(jīng)簡(jiǎn)單功率譜篩選的視電阻率曲線(xiàn)Fig.9 Apparent resistivity curve through simple power spectrum filtering

4 結(jié)論

1)將子空間增強(qiáng)算法應(yīng)用于礦集區(qū)大地電磁信噪分離，通過(guò)模擬典型的強(qiáng)干擾類(lèi)型，分析了算法中影響濾波精度的幀長(zhǎng)及特征值判別閾值的最優(yōu)選值范圍，給出了算法流程。

2)通過(guò)模擬信號(hào)仿真和實(shí)測(cè)數(shù)據(jù)處理，驗(yàn)證了方法的可行性。受近源干擾嚴(yán)重的測(cè)點(diǎn)經(jīng)子空間增強(qiáng)法處理后，時(shí)間域的大尺度干擾和基線(xiàn)漂移得到了有效壓制，視電阻率曲線(xiàn)的整體形態(tài)更為光滑、連續(xù)，低頻段的大地電磁數(shù)據(jù)質(zhì)量得到了明顯改善。子空間增強(qiáng)法對(duì)提升礦集區(qū)大地電磁測(cè)深深部勘探能力，具有重要的實(shí)際應(yīng)用價(jià)值。

[1] TIKHONOV A N. On determining electrical characteristics of the deep layers of the Earth’s crust[J]. Dok1. Akad. Nauk. SSSR, 1950, 73(2): 295-297.

[2] CAGNIARD L.Basic theory of the magnetotelluric method of geophysical prospecting[J].Geophysics,1953,18(3):605-635.

[3] 朱威, 范翠松, 姚大為, 等. 礦集區(qū)大地電磁噪聲場(chǎng)源分析及噪聲特點(diǎn)[J]. 物探與化探, 2011, 35(5):658-662. ZHU W, FAN C S , YAO D W,et al. Noise source analysis and noise characteristic study of MT in an ore concentration area[J]. Geophysical and Geochemical Exploration, 2011, 35(5): 658-662. (In Chinese)

[4] 王書(shū)明, 王家映. 大地電磁信號(hào)統(tǒng)計(jì)特征分析[J]. 地震學(xué)報(bào), 2004, 26(6): 669-674. WANG S M, WANG J Y. Analysis on statistic characteristics of magnetotelluric signal[J]. Acta Seismologica Sinca, 2004, 26(6): 669-674. (In Chinese)

[5] 王書(shū)明, 王家映. 大地電磁信號(hào)非最小相位性的討論[J]. 地球物理學(xué)進(jìn)展, 2004, 19(2): 216-221. WANG S M, WANG J Y. Discussion on the non-minimum phase of magnetotelluric signals[J]. Progress in Geophysics, 2004, 19(2): 216-221. (In Chinese)

[6] 湯井田, 化希瑞, 曹哲民, 等. Hilbert-Huang變換與大地電磁噪聲壓制[J]. 地球物理學(xué)報(bào), 2008, 51(2):603-610. TANG J T, HUA X R, CAO Z M, et al. Hilbert-Huang transformation and noise suppression of magnetotelluric sounding data[J]. Chinese Journal Geophysics, 2008, 51(2):603-610. (In Chinese)

[7] 景建恩, 魏文博, 陳海燕, 等. 基于廣義S變換的大地電磁測(cè)深數(shù)據(jù)處理[J]. 地球物理學(xué)報(bào), 2012, 55 (12): 4015-4022. JING J E, WEI W B, CHEN H Y, et al. Magnetotelluric sounding data processing based on generalized S transformation[J]. Chinese Journal Geophysics，2012, 55(12): 4015-4022. (In Chinese)

[8] KAPPLE K N. A data variance technique for automated despiking of magnetotelluric data with a remote reference[J]. Geophysical Prospecting, 2012, 60(1): 179-191.

[9] 湯井田, 李晉, 肖曉, 等. 數(shù)學(xué)形態(tài)濾波與大地電磁噪聲壓制[J]. 地球物理學(xué)報(bào), 2012, 55(5): 1784-1793. TANG J T, LI J, XIAO X, et al. Mathematical morphology filtering and noise suppression of magnetotelluric sounding data[J]. Chinese Journal Geophysics, 2012, 55(5): 1784-1793. (In Chinese)

[10]王輝, 魏文博, 金勝, 等. 基于同步大地電磁時(shí)間序列依賴(lài)關(guān)系的噪聲處理[J]. 地球物理學(xué)報(bào), 2014, 57(2): 531-545. WANG H, WEI W B, JIN S, et al. Removal of magnetotelluric noise based on synchronous time series relationship[J]. Chinese Journal Geophysics, 2014, 57(2): 531-545. (In Chinese)

[11]唐貴基，龐彬，劉尚坤. 基于奇異差分譜和平穩(wěn)子空間分析的滾動(dòng)軸承故障診斷[J]. 振動(dòng)與沖擊, 2015, 34 (11): 83-87. TANG G J, PANG B, LIU S K. Fault diagnosis of rolling bearings based on difference spectrum of singular value and stationary subspace analysis[J]. Journal of Vibration and Shock, 2015, 34 (11): 83-87. (In Chinese)

[12]陸文凱, 丁文龍, 張善文, 等. 基于信號(hào)子空間分解的三維地震資料高分辨率處理方法[J]. 地球物理學(xué)報(bào), 2009, 52 (8): 2174-2181. LU W K, DING W L, ZHANG S W, et al. A high resolution processing technique for 3_D seismic data based on signal subspace decomposition[J]. Chinese Journal Geophysics, 2009, 52(8): 2174-2181. (In Chinese)

[13]DENDRINOS M, BAKAMIDIS S, CARAYANNI G. Speech enhancement from noise: A regenerative approach[J]. Speech Communication, 1991, 10(2): 45-57.

[14]EPHRAIM Y, VAN TREES H L. A signal subspace approach for speech enhancement[J]. IEEE Transactions on Speech and Audio Processing, 1995, 3(4): 251-266.

[15]JENSEN S H, HANSEN P C, HANSEN S D, et al. Reduction of broad-band noise in speech by truncated QSVD[J]. IEEE Transactions on Speech and Audio Processing, 1995, 3(6): 439-448.

[16]李晉, 湯井田, 王玲, 等. 基于信號(hào)子空間增強(qiáng)和端點(diǎn)檢測(cè)的大地電磁噪聲壓制[J]. 物理學(xué)報(bào), 2014,63(1): 019101. LI J, TANG J T, WANG L, et al. Noise suppression for magnetotelluric sounding data based on signal subspace enhancement and endpoint detection[J]. Acta Physica Sinica, 2014, 63(1): 019101. (In Chinese)

[17]湯井田, 劉子杰, 劉峰屹, 等. 音頻大地電磁法強(qiáng)干擾壓制試驗(yàn)研究[J]. 地球物理學(xué)報(bào), 2015, 58(12): 4636-4647. TANG J T, LIU Z J, LIU F Y, et al. The denoising of the audio magnetotelluric data set with strong interferences[J]. Chinese Journal Geophysics, 2015, 58(12): 4636-4647. (In Chinese)

[18]代小威. 基于V5-2000格式MT時(shí)間序列處理與功率譜估計(jì)及軟件開(kāi)發(fā)[D]. 武漢：長(zhǎng)江大學(xué), 2015. DAI X W. Processing and power spectrum estimation of MT time serizes in V5-2000 format and software development[D]. Wuhan: Yangtze University, 2015. (In Chinese)

Application research of subspace enhancement for magnetotelluric signal to noise separation

YAN Huan1, LI Jin1,2, LI Chunyang1

(1.Institute of Physics and Information Science, Hunan Normal University, Changsha 410081, China; 2.School of Geosciences and Info-Physics, Key Laboratory of Metallogenic Prediction of Non-Ferrous Metals and Geological Environment Monitor, Ministry of Education, Central South University, Changsha 410083, China)

In order to improve the quality of magnetotelluric sounding data and retain the useful information of low frequency band for magnetotelluric, suppressing the large-scale strong interference is imperative. Aimed at typical strong interference types of ore concentration area, we introduced the subspace enhancement method to study signal to noise separation of magnetotelluric data and strong interference. Through simulated large-scale square wave and charge and discharge triangle wave by computer, we discussed the optimal range of frame and criterion threshold of eigenvalue in subspace enhancement, and the concrete steps of the algorithm are presented. Then, the measured magnetotelluric data of ore concentration area were processed by signal to noise separation. Experiment results show that the method is feasible, which can effectively eliminate the large-scale strong noise interference in the time domain and the near source interference of ore concentration area has been effectively suppressed. At the same time, the apparent resistivity curve is more smoothly and continuously. Moreover, the quality of magnetotelluric sounding data is improved significantly in the low frequency band.

magnetotelluric sounding data; signal to noise separation; subspace enhancement; strong interference

2016-06-14 改回日期：2016-12-15

國(guó)家自然科學(xué)基金資助項(xiàng)目(41404111)；湖南省自然科學(xué)基金資助項(xiàng)目(2015JJ3088)；中國(guó)博士后科學(xué)基金資助項(xiàng)目(2015M570687) ;湖南省研究行科研創(chuàng)新項(xiàng)目(CX2017B224)

燕歡(1991-)，女，碩士，主要從事礦集區(qū)大地電磁強(qiáng)干擾壓制研究，E-mail：aviva163@126.com。

李晉(1981-)，男，博士后，副教授，碩士生導(dǎo)師，主要從事信號(hào)處理及電磁勘探研究， E-mail：geologylj@163.com。

1001-1749(2017)03-0346-08

P 631.2

A

10.3969/j.issn.1001-1749.2017.03.08