張赟昀 王培杰 陳小斌1，* 詹 艷 韓 冰 王立鳳 趙國澤

1)應(yīng)急管理部國家自然災(zāi)害防治研究院，北京 100085

2)中國地震局地質(zhì)研究所，地震動(dòng)力學(xué)國家重點(diǎn)實(shí)驗(yàn)室，北京 100029

0 引言

極低頻電磁臺(tái)網(wǎng)不僅可監(jiān)測(cè)穩(wěn)定的高信噪比人工源極低頻電磁場信號(hào)，也可監(jiān)測(cè)天然場的電磁信號(hào)(趙國澤等，2010)。然而，天然場電磁信號(hào)包含各種人文噪聲，不同臺(tái)站受到的噪聲干擾情況不同，使得不同臺(tái)站的視電阻率曲線有的頻段數(shù)據(jù)質(zhì)量較好，有的頻段數(shù)據(jù)質(zhì)量較差。數(shù)據(jù)質(zhì)量差的頻段無法為地震短臨前兆提供有用信息，影響了極低頻電磁臺(tái)網(wǎng)數(shù)據(jù)的使用率，因此有必要研究含噪聲干擾情況下的數(shù)據(jù)處理方法。

大地電磁測(cè)深法(Magnetotelluric，MT)(Cagniard，1953)是一種天然源頻率域的測(cè)深方法，探測(cè)深度受限于地下電性結(jié)構(gòu)和觀測(cè)數(shù)據(jù)有效周期等因素，有效周期越長，探測(cè)深度則越大。天然場源信號(hào)弱，存在死頻帶、近場干擾、隨機(jī)噪聲等影響阻抗估計(jì)結(jié)果的因素，在強(qiáng)干擾區(qū)很難獲取十分準(zhǔn)確的電磁響應(yīng)，而穩(wěn)定可靠的阻抗估計(jì)是后續(xù)反演解釋的前提(王家映等，1998)。為得到無偏的大地電磁響應(yīng)，地球物理學(xué)家提出了一系列方法。Sims等(1971)假定天然源大地電磁信號(hào)是平穩(wěn)信號(hào)，磁道在一定區(qū)域內(nèi)為平穩(wěn)分布，且電、磁道噪聲服從高斯分布，可采用最小二乘法進(jìn)行阻抗估計(jì)；Goubau等(1978，1984)和Gamble等(1979a，b)提出了遠(yuǎn)參考道方法(Remote Reference)，通過同時(shí)記錄測(cè)點(diǎn)與參考點(diǎn)的磁場，用互功率譜代替自功率譜以減少相干噪聲的影響；Egbert等(1986)和Chave等(1989，2003)提出穩(wěn)健估計(jì)方法(Robust)，通過降低含噪聲窗口的傅里葉系數(shù)權(quán)重壓制隨機(jī)噪聲。在時(shí)域去噪方面，湯井田等(2008)使用希爾伯特-黃變換(Hilbert-Huang Transformation，HHT)分析大地電磁信號(hào)中的噪聲分布特征，并進(jìn)行干擾壓制，取得了一定的效果。Jin等(2012)應(yīng)用HHT處理海洋電磁數(shù)據(jù)，并與有界影響估計(jì)(Bounded-Influence Remote Reference Processing，BIRRP)(Chaveetal.，2003，2004)的處理結(jié)果進(jìn)行了對(duì)比，發(fā)現(xiàn)HHT得到的結(jié)果更穩(wěn)定。Wang等(2017)提出了STIN(Synthesis of Natural Electric and Magnetic Time-series Using Inter-station Transfer Functions and Time-series from a Neighboring Site)法，通過站間傳輸函數(shù)合成本地道含噪聲時(shí)段的數(shù)據(jù)，避免含噪聲段數(shù)據(jù)對(duì)阻抗估計(jì)的影響。李晉等(2019)結(jié)合變分模態(tài)分解(Variational Mode Decomposition，VMD)和匹配追蹤技術(shù)(Matching Pursuit，MP)壓制音頻大地電磁強(qiáng)干擾。

由于各種新興的數(shù)據(jù)處理方法穩(wěn)定性不足、適用性不強(qiáng)，并沒有得到普遍應(yīng)用，遠(yuǎn)參考道方法和穩(wěn)健估計(jì)方法仍然是當(dāng)前生產(chǎn)中應(yīng)用最廣泛的方法。然而在實(shí)際工作中我們發(fā)現(xiàn)，使用遠(yuǎn)參考道法和穩(wěn)健估計(jì)法相結(jié)合的方式處理在強(qiáng)干擾環(huán)境下采集的數(shù)據(jù)，效果并不理想。其他學(xué)者也針對(duì)這個(gè)問題展開過研究，例如Egbert(1997)利用多臺(tái)站數(shù)據(jù)估算背景噪聲水平，使用低相干噪聲的數(shù)據(jù)估算阻抗，當(dāng)至少存在一個(gè)臺(tái)站不受相干噪聲影響或?yàn)殚g歇相干噪聲源時(shí)，可得到比單臺(tái)站處理更好的阻抗估計(jì)結(jié)果，但是這種方法要求同時(shí)采集多臺(tái)站的數(shù)據(jù)。本文希望找到合適的數(shù)據(jù)處理方案，以提升單臺(tái)站、單遠(yuǎn)參考臺(tái)站數(shù)據(jù)處理結(jié)果的穩(wěn)定性和可靠性。中國東部是經(jīng)濟(jì)發(fā)展繁榮的地區(qū)，人文干擾強(qiáng)烈，本文對(duì)郯廬斷裂帶中南段強(qiáng)干擾區(qū)的大地電磁數(shù)據(jù)進(jìn)行處理，選取了其中部分特征明顯的測(cè)點(diǎn)進(jìn)行分析，并針對(duì)不同的干擾類型提出了一套數(shù)據(jù)處理方法。

1 數(shù)據(jù)處理方法

大地電磁法記錄相互垂直的天然源電磁場分量，在頻率域內(nèi)MT電磁場與阻抗張量滿足如下的線性關(guān)系：

(1)

式(1)中，Ex、Ey、Hx、Hy為4個(gè)水平電磁場分量，Zxx、Zxy、Zyx、Zyy為4個(gè)阻抗張量元素，可簡寫為

E=ZH

(2)

阻抗的特征通常用視電阻率ρ和相位φ來表征。式(3)、(4)中，i，j∈{x，y}：

(3)

(4)

ERH=ZΗRH

(5)

其最小二乘解形式如式(6)，式中，上標(biāo)H為共軛轉(zhuǎn)置符號(hào)：

Z=ERH(HRH)-1

(6)

展開為

(7)

式(7)中，星號(hào)代表共軛，尖括號(hào)代表平均值，即功率譜的平均值。

在大地電磁測(cè)深中，實(shí)測(cè)數(shù)據(jù)通常包括信號(hào)和噪聲：

Em=Es+En

Hm=Hs+Hn

Rm=Rs+Rn

(8)

式(8)中，下標(biāo)m表示實(shí)測(cè)數(shù)據(jù)(measured)，s表示信號(hào)(signal)，n表示噪聲(noise)。實(shí)測(cè)數(shù)據(jù)的張量阻抗最小二乘解為

(9)

由于水平電場變化劇烈，阻抗求解不穩(wěn)定，通常用磁場作為參考道：

(10)

互相關(guān)可突出2個(gè)信號(hào)中相干性強(qiáng)的部分，壓制相干性弱的部分。若2個(gè)信號(hào)中的噪聲不相干，則互相關(guān)可以壓制噪聲；若2個(gè)信號(hào)的噪聲相干性強(qiáng)，互相關(guān)反而會(huì)突出噪聲信號(hào)。因此，最小二乘法阻抗估計(jì)的關(guān)鍵在于選取與Es和Hs相干性強(qiáng)、與En和Hn相干性弱的參考道。當(dāng)觀測(cè)點(diǎn)附近存在強(qiáng)人工源噪聲信號(hào)時(shí)，En和Hn具有強(qiáng)相干性，采用最小二乘法得到的阻抗估計(jì)會(huì)嚴(yán)重失真。Goubau等(1978，1984)和Gamble等(1979a，b)提出遠(yuǎn)參考法，以遠(yuǎn)參考點(diǎn)的磁場作為參考道R參與最小二乘計(jì)算，可有效壓制相干噪聲。

穩(wěn)健估計(jì)(Robust)技術(shù)是一種通過自適應(yīng)加權(quán)來降低飛點(diǎn)對(duì)阻抗估計(jì)影響的方法。Robust將最小二乘法的計(jì)算結(jié)果作為初始阻抗，計(jì)算系統(tǒng)觀測(cè)值與預(yù)測(cè)值之間的誤差，根據(jù)誤差設(shè)置權(quán)函數(shù)對(duì)觀測(cè)值加權(quán)，并利用加權(quán)最小二乘法計(jì)算新的阻抗。重復(fù)上述步驟直至計(jì)算的阻抗不再變化。Larsen等(1996)發(fā)現(xiàn)當(dāng)相干噪聲較大時(shí)，基于Robust技術(shù)的單點(diǎn)最小二乘法阻抗估計(jì)會(huì)給出錯(cuò)誤的結(jié)果，而基于Robust技術(shù)和遠(yuǎn)參考的阻抗估計(jì)結(jié)果則接近真實(shí)阻抗。但Robust技術(shù)的前提是預(yù)測(cè)值接近真實(shí)值，這要求數(shù)據(jù)中干擾較弱的時(shí)間長度遠(yuǎn)大于干擾強(qiáng)的時(shí)間長度。而強(qiáng)人文干擾的地區(qū)顯然不滿足這種條件，Robust處理強(qiáng)干擾持續(xù)時(shí)間長的數(shù)據(jù)往往會(huì)放大干擾信號(hào)而壓制有效信號(hào)。

雖然近場效應(yīng)具有明顯的特征，但對(duì)于近場干擾不太強(qiáng)烈的測(cè)點(diǎn)，其曲線具有與近場類似的特征，視電阻率上翹但角度不到45°，相位下降但不到0°，很難與正常的測(cè)點(diǎn)曲線區(qū)分。若近場影響頻段較小，很容易被誤認(rèn)為是高阻異常。Rhoplus(Parkeretal.，1996)是一種直接反演視電阻率和相位的一維反演方法，該方法考慮視電阻率和相位數(shù)據(jù)的關(guān)聯(lián)性和數(shù)據(jù)的連續(xù)性，可用于大地電磁視電阻率和相位數(shù)據(jù)的一致性檢測(cè)，以判斷曲線是否合理并矯正畸變數(shù)據(jù)(譚捍東等，2004；Sprattetal.，2005；周聰?shù)龋?015)。本文中所有數(shù)據(jù)的處理結(jié)果都已導(dǎo)入大地電磁資料處理與解釋可視化集成系統(tǒng)MT-Pioneer(陳小斌等，2004)中，剔除飛點(diǎn)后與Rhoplus曲線擬合進(jìn)行驗(yàn)證，并以Rhoplus曲線為標(biāo)準(zhǔn)反復(fù)進(jìn)行數(shù)據(jù)編輯。

由于強(qiáng)人文干擾地區(qū)的相干噪聲存在的時(shí)間長，對(duì)全時(shí)段采集的數(shù)據(jù)做Robust處理難以得到準(zhǔn)確的阻抗，因此有必要對(duì)數(shù)據(jù)進(jìn)行分段處理，從中挑選出受干擾較小時(shí)段的數(shù)據(jù)參與疊加。若分段數(shù)較少，很可能每段中都包含強(qiáng)噪聲，不能挑選出不含噪聲時(shí)段的數(shù)據(jù)；若分段數(shù)太多，則低噪聲時(shí)段的數(shù)據(jù)有可能被截?cái)?，也不能獲取較好的阻抗估計(jì)結(jié)果。本文嘗試采用不同的分段數(shù)進(jìn)行處理，分析分段數(shù)對(duì)處理結(jié)果的影響。利用Robust方法處理強(qiáng)干擾數(shù)據(jù)的效果并不穩(wěn)定，本文分析了Robust方法在強(qiáng)干擾區(qū)數(shù)據(jù)處理中的應(yīng)用效果。在強(qiáng)干擾區(qū)數(shù)據(jù)質(zhì)量差的時(shí)段占比很大，需要在進(jìn)行功率譜疊加前予以剔除，本文中所有功率譜都在疊加前經(jīng)過人工挑選。

2 實(shí)測(cè)數(shù)據(jù)采集與處理

2.1 工區(qū)概況及數(shù)據(jù)采集

本文選取中國東部郯廬斷裂帶(合肥—宿遷段)的實(shí)測(cè)數(shù)據(jù)進(jìn)行處理分析。測(cè)區(qū)位于安徽省東北部，該區(qū)大部分為平原，經(jīng)濟(jì)發(fā)展較好。區(qū)內(nèi)村莊城鎮(zhèn)密布，電網(wǎng)縱橫，有風(fēng)力發(fā)電站、太陽能發(fā)電站、電氣化鐵路以及大型工廠，人文干擾強(qiáng)烈。以往采集的數(shù)據(jù)表明，在該區(qū)域很難采集到質(zhì)量好的數(shù)據(jù)。特別是在0.1Hz附近的死頻帶內(nèi)，天然場源信號(hào)弱，近場信號(hào)強(qiáng)，導(dǎo)致在該頻帶內(nèi)難以得到穩(wěn)定的阻抗數(shù)據(jù)。

圖1 大地電磁測(cè)點(diǎn)及遠(yuǎn)參考點(diǎn)分布圖Fig.1 Distribution of MT station and remote reference station.

2.2 遠(yuǎn)參考處理

Shalivahan等(2002)分別使用3個(gè)遠(yuǎn)參考點(diǎn)對(duì)同一測(cè)點(diǎn)進(jìn)行了處理，發(fā)現(xiàn)當(dāng)遠(yuǎn)參考點(diǎn)距離觀測(cè)點(diǎn)較近時(shí)(80km、115km)，遠(yuǎn)參考技術(shù)難以改善死頻帶的處理結(jié)果，當(dāng)遠(yuǎn)參考點(diǎn)距離觀測(cè)點(diǎn)足夠遠(yuǎn)時(shí)(215km)，才能在死頻帶得到穩(wěn)定平滑的阻抗估計(jì)。為改善數(shù)據(jù)質(zhì)量，本文設(shè)置了遠(yuǎn)參考點(diǎn)，采用遠(yuǎn)參考技術(shù)對(duì)工區(qū)內(nèi)的數(shù)據(jù)進(jìn)行處理，并根據(jù)遠(yuǎn)參考技術(shù)的處理效果對(duì)遠(yuǎn)參考點(diǎn)的位置進(jìn)行了調(diào)整。最初的遠(yuǎn)參考點(diǎn)(遠(yuǎn)參一)位于工區(qū)南面安徽省與江西省交界處的山區(qū)內(nèi)(29°55′46″N，117°23′38″E)(圖1)，距離工區(qū)內(nèi)測(cè)點(diǎn)的平均距離約為350km。收取第1輪遠(yuǎn)參考數(shù)據(jù)后，使用該遠(yuǎn)參考點(diǎn)對(duì)工區(qū)內(nèi)的測(cè)點(diǎn)進(jìn)行遠(yuǎn)參考處理，發(fā)現(xiàn)該點(diǎn)自身的數(shù)據(jù)質(zhì)量不佳，存在不連續(xù)的干擾信號(hào)，從圖2 中可以看出該遠(yuǎn)參考點(diǎn)的中低頻段與Rhoplus曲線不能很好地?cái)M合。810測(cè)點(diǎn)參考后數(shù)據(jù)質(zhì)量稍有改善，高頻部分與Rhoplus曲線擬合得更好，中頻近場形態(tài)明顯被壓制，但中低頻部分仍有飛點(diǎn)存在，擬合效果不佳，且與遠(yuǎn)參考點(diǎn)擬合不好的頻帶一致。隨后，我們?cè)谏綎|省臨沂市蒙陰縣(35°31′57″N，118°1′8″E)設(shè)置了第2個(gè)遠(yuǎn)參考點(diǎn)(遠(yuǎn)參二)，該遠(yuǎn)參考點(diǎn)處于工區(qū)北部，距離工區(qū)內(nèi)測(cè)點(diǎn)的平均距離約為260km。遠(yuǎn)參二的低頻部分非常連續(xù)，但0.1～10s處的視電阻率快速上升又急劇下降，雖然曲線平滑，但與Rhoplus曲線擬合不上，相位曲線更是在1～10s處下掉，因此推測(cè)該頻帶可能存在穩(wěn)定干擾。使用該點(diǎn)對(duì)工區(qū)內(nèi)測(cè)點(diǎn)進(jìn)行遠(yuǎn)參考處理，發(fā)現(xiàn)工區(qū)內(nèi)的數(shù)據(jù)有所改善，可能的原因是高頻電磁信號(hào)衰減快，與測(cè)區(qū)干擾相干度低，或遠(yuǎn)參考站的干擾頻帶與測(cè)點(diǎn)的干擾頻帶不一致。但仍有部分測(cè)點(diǎn)在參考后出現(xiàn)異常，如圖2 中305測(cè)點(diǎn)參考后，在0.1～10s范圍內(nèi)視電阻率能完全擬合，但相位擬合不好，推測(cè)參考站的本地噪聲也會(huì)影響遠(yuǎn)參考處理的效果。遠(yuǎn)參三設(shè)置在內(nèi)蒙古鄂托克前旗(107°49′40″N，38°1′27″E)，電磁環(huán)境非常好，幾乎沒有人文噪聲，且距離工區(qū)的平均距離約為1050km，與工區(qū)之間存在相干噪聲的可能性很低，遠(yuǎn)參考處理能很好地改善工區(qū)的數(shù)據(jù)質(zhì)量。圖2 中遠(yuǎn)參三的視電阻率和相位均與Rhoplus曲線擬合得很好，僅在低頻接近10000s處的yx向視電阻率相位未擬合好。509測(cè)點(diǎn)未進(jìn)行遠(yuǎn)參考處理的數(shù)據(jù)在10s處不連續(xù)，處理后數(shù)據(jù)十分連續(xù)，與Rhoplus曲線的擬合情況非常好，且和遠(yuǎn)參考點(diǎn)類似，僅在接近10000s處的yx向視電阻率、相位擬合不佳。

圖2 不同遠(yuǎn)參考站的數(shù)據(jù)處理效果圖Fig.2 Data processing effect of different remote reference stations.紅點(diǎn)、藍(lán)點(diǎn)分別為xy、yx向視電阻率和阻抗相位，橙線、綠線分別為xy、yx向Rhoplus反演曲線，下同

由于 3 個(gè)遠(yuǎn)參考站沒有同時(shí)記錄，僅挑選具有代表性的測(cè)點(diǎn)進(jìn)行比較。不難發(fā)現(xiàn)，所有測(cè)點(diǎn)的數(shù)據(jù)經(jīng)過遠(yuǎn)參考處理后均有改善。一般情況下，只要遠(yuǎn)參考站本身的數(shù)據(jù)干擾較小，經(jīng)過遠(yuǎn)參考處理后的數(shù)據(jù)都會(huì)變好。因此，只要條件允許，都應(yīng)該布設(shè)遠(yuǎn)參考站。下文中的處理過程均是在遠(yuǎn)參考處理的基礎(chǔ)上進(jìn)行的。為簡化敘述，不再對(duì)此作單獨(dú)說明。

2.3 Robust處理

工區(qū)范圍較大，測(cè)點(diǎn)間距大，不同測(cè)點(diǎn)受到干擾的情況不同。為測(cè)試不同干擾情況下Robust技術(shù)的處理效果，將測(cè)點(diǎn)使用Robust技術(shù)處理的結(jié)果與不使用Robust(簡稱非Robust)技術(shù)處理的結(jié)果進(jìn)行對(duì)比。

圖3 709測(cè)點(diǎn)的視電阻率和相位曲線Fig.3 Apparent resistivity and phase curves of MT station 709.

在0.1Hz附近的死頻帶內(nèi)，若沒有近場干擾，信號(hào)整體比較弱，天然場的信號(hào)強(qiáng)度和環(huán)境噪聲的信號(hào)強(qiáng)度處于同一數(shù)量級(jí)。使用非Robust技術(shù)處理時(shí)，對(duì)每段數(shù)據(jù)的譜以相同的權(quán)重進(jìn)行疊加，微弱的天然場信號(hào)會(huì)淹沒在噪聲中。如果使用Robust技術(shù)處理，對(duì)相干性較強(qiáng)的譜以較大的權(quán)重進(jìn)行疊加，能夠突出具有一定相干度的天然場信號(hào)，從而獲得更準(zhǔn)確的阻抗估計(jì)。以709測(cè)點(diǎn)為例(圖3)，從未經(jīng)編輯的曲線中可以看出，709測(cè)點(diǎn)的視電阻率和阻抗相位在高頻部分和低頻部分均展示出基本形態(tài)。對(duì)比Robust處理和非Robust處理的未經(jīng)編輯的曲線可以發(fā)現(xiàn)，經(jīng)Robust處理的曲線形態(tài)更好，0.1Hz附近死頻帶內(nèi)的頻點(diǎn)更接近Rhoplus曲線，而非Robust處理的曲線中飛點(diǎn)現(xiàn)象嚴(yán)重。對(duì)二者的功率譜進(jìn)行精細(xì)挑選后得到最終結(jié)果，可以發(fā)現(xiàn)Robust的處理結(jié)果和Rhoplus曲線能很好地?cái)M合，而非Robust的處理結(jié)果在死頻帶仍存在飛點(diǎn)，不能與Rhoplus曲線擬合。

圖4 209測(cè)點(diǎn)的視電阻率和相位曲線Fig.4 Apparent resistivity and phase curves of MT station 209.

而對(duì)于持續(xù)存在近場干擾源的點(diǎn)，干擾信號(hào)的電場和磁場相干度更高，Robust會(huì)將近場信號(hào)識(shí)別為期望的信號(hào)，處理過程中會(huì)壓制天然場信號(hào)，放大近場干擾信號(hào)，使阻抗估計(jì)偏向于近場信號(hào)的阻抗。工區(qū)內(nèi)靠近大型城市、大型工廠、電氣化鐵路等近場源的測(cè)點(diǎn)曲線都呈現(xiàn)視電阻率45°上升、相位趨于0°的近場特征。以209測(cè)點(diǎn)為例(圖4)，該點(diǎn)SW向1.5km處有村莊，W向1.5km處有工廠，存在潛在的干擾源，使用Robust處理后得到的曲線存在明顯的近場特征，視電阻率曲線在1Hz附近接近45°上升，yx方向的相位曲線接近0°。由于分段進(jìn)行Robust處理后每段的視電阻率和相位均呈現(xiàn)近場形態(tài)，數(shù)據(jù)不具備可編輯的空間，編輯后近場形態(tài)改善不明顯，且中頻仍存在飛點(diǎn)，因此對(duì)該點(diǎn)采用Robust處理只能得到近似于近場的阻抗。采用非Robust方法進(jìn)行處理，即分組功率譜疊加時(shí)不對(duì)每段數(shù)據(jù)賦予權(quán)重，可得到少量無近場特征或近場特征弱的譜，經(jīng)過人工編輯，濾除以近場為主的功率譜，僅保留少量符合大地電磁曲線特征的譜，疊加后可以發(fā)現(xiàn)視電阻率和相位都能與Rhoplus曲線較好地?cái)M合，近場形態(tài)改善明顯。

工區(qū)內(nèi)數(shù)據(jù)處理的難點(diǎn)在于死頻帶和近場干擾這2種情形，同時(shí)也不乏二者同時(shí)存在的測(cè)點(diǎn)，特別是近場頻帶和死頻帶有重合部分的測(cè)點(diǎn)，其數(shù)據(jù)處理十分困難。這種情況下的數(shù)據(jù)中包含近場信號(hào)、有效信號(hào)和隨機(jī)噪聲，其中近場信號(hào)最強(qiáng)，有效信號(hào)和隨機(jī)噪聲較弱。采用Robust處理會(huì)放大相干度高的近場信號(hào)和有效信號(hào)，壓制相干度低的隨機(jī)噪聲，但近場信號(hào)較強(qiáng)，因此阻抗估計(jì)結(jié)果仍以近場阻抗為主；不采用Robust處理則會(huì)讓死頻帶內(nèi)的天然場信號(hào)和隨機(jī)噪聲淹沒在近場信號(hào)中，仍然無法得到準(zhǔn)確的阻抗估計(jì)。針對(duì)這種數(shù)據(jù)，我們把2種方法處理得到的譜放在一起編輯，兼顧二者的優(yōu)點(diǎn)。以測(cè)點(diǎn)307為例(圖5)，分別用2種方法進(jìn)行處理后，均不能得到理想的阻抗估計(jì)，最終編輯得到的曲線中均存在飛點(diǎn)，而將2種處理方法得到的功率譜放在一起進(jìn)行編輯后，曲線在0.1Hz處有明顯改善，且低頻曲線更為平滑。其可能的原因有2點(diǎn)：1)參與挑選和疊加的功率譜增加了一倍，組合得到合理結(jié)果的概率增大，對(duì)結(jié)果的改善存在偶然性；2)不同時(shí)間段的近場信號(hào)、有效信號(hào)和隨機(jī)噪聲強(qiáng)弱存在變化，Robust和非Robust處理得到的譜都參與編輯，近場干擾強(qiáng)時(shí)選擇非Robust處理得到的譜可以降低近場的影響，而在近場干擾弱時(shí)，選擇Robust處理得到的譜有利于從背景噪聲中挑出有效信號(hào)。值得注意的是，307測(cè)點(diǎn)只是個(gè)例，該方法的通用性未經(jīng)考證，僅供參考。

圖5 307測(cè)點(diǎn)穩(wěn)健估計(jì)(Robust)處理、非穩(wěn)健估計(jì)(非Robust)處理及二者疊加處理得到的視電阻率和相位曲線Fig.5 Apparent resistivity and phase curve of MT station 307 obtained by robust processing，non-robust processing and superposition processing.

2.4 分段數(shù)

除采用不同的Robust處理策略外，還可通過增加分段數(shù)來改善數(shù)據(jù)質(zhì)量。死頻帶和近場頻帶都位于中頻帶，中頻帶的數(shù)據(jù)具有足夠長的采樣時(shí)間，增大分段數(shù)會(huì)減小分段處理中每段的時(shí)間，得到更多組功率譜。當(dāng)某時(shí)間段內(nèi)的干擾較小或近場干擾暫停時(shí)，該段數(shù)據(jù)可以得到相對(duì)準(zhǔn)確的阻抗估計(jì)，人工編輯時(shí)將其保留。相對(duì)于每段時(shí)間長、段數(shù)少的處理方式，這種每段時(shí)間短、段數(shù)多的處理方式能保留更多的有效功率譜。以503測(cè)點(diǎn)和504測(cè)點(diǎn)為例(圖6，7)，截取4000min的數(shù)據(jù)，分別采用5段、10段、20段、40段的分段數(shù)，每段對(duì)應(yīng)的時(shí)長分別為800min、400min、200min、100min，經(jīng)過精細(xì)編輯后得到如圖6、7所示的響應(yīng)曲線。不難發(fā)現(xiàn)，隨著分段數(shù)的逐漸增加，響應(yīng)曲線與Rhoplus曲線的擬合越來越好。增大分段數(shù)無疑會(huì)增加人工編輯的工作量，但在電磁環(huán)境復(fù)雜的情況下，為改善數(shù)據(jù)質(zhì)量，增加工作量是值得的。

圖6 503測(cè)點(diǎn)不同分段數(shù)對(duì)應(yīng)的視電阻率和相位曲線Fig.6 Apparent resistivity and phase curves of MT station 503 corresponding to different segment numbers.從左至右對(duì)應(yīng)的每段時(shí)長分別為800min、400min、200min和100min

圖7 504測(cè)點(diǎn)不同分段數(shù)對(duì)應(yīng)的視電阻率和相位曲線Fig.7 Apparent resistivity and phase curves of MT station 504 corresponding to different segment numbers.從左至右對(duì)應(yīng)的每段時(shí)長分別為800min、400min、200min和100min

圖9 305測(cè)點(diǎn)晝夜不同時(shí)長的視電阻率及相位曲線對(duì)比圖Fig.9 Comparison of apparent resistivity and phase curves of MT station 305 with different time and length.

2.5 分時(shí)段處理

人文干擾的強(qiáng)度以天為周期，白天人類活動(dòng)多時(shí)干擾強(qiáng)，夜晚人類活動(dòng)少時(shí)干擾弱。使用305測(cè)點(diǎn)進(jìn)行晝夜數(shù)據(jù)的對(duì)比研究，分別選取10：00—16：00和23：00—5：00作為晝夜數(shù)據(jù)對(duì)比時(shí)段。圖8 對(duì)比了白天和夜晚未編輯與編輯后的數(shù)據(jù)，從圖中可明顯看出夜晚數(shù)據(jù)的質(zhì)量優(yōu)于白天。白天的死頻帶數(shù)據(jù)經(jīng)編輯后仍存在大量飛點(diǎn)，與Rhoplus曲線擬合較差。而夜晚的死頻帶數(shù)據(jù)情形稍好，低頻部分能得到連續(xù)光滑的曲線。為了突出夜晚數(shù)據(jù)質(zhì)量較高的特征，選取了10h白天的數(shù)據(jù)，與5h夜晚的數(shù)據(jù)進(jìn)行對(duì)比，從圖9 中未經(jīng)過人工選譜的原始數(shù)據(jù)可以看出，雖然白天的采集時(shí)間是晚上的2倍，但視電阻率曲線不連續(xù)、相位曲線在10～0.1Hz雜亂無章。經(jīng)人工選譜后白天和夜晚的數(shù)據(jù)均得到改善，但仍是夜晚的數(shù)據(jù)比白天好。

對(duì)305測(cè)點(diǎn)白天(8：00—20：00)和夜晚(20：00—8：00)的數(shù)據(jù)分別截取12h、24h、36h、48h、60h進(jìn)行處理，其中第1段為第1晚的數(shù)據(jù)，后面每段為前面一段加上與之連續(xù)的白天或夜晚12h的數(shù)據(jù)。從圖10 中可以看出，隨著采集時(shí)間的增加，視電阻率曲線的質(zhì)量逐漸變好，60h數(shù)據(jù)的視電阻率曲線與Rhoplus曲線擬合得很好，相位曲線與Rhoplus曲線的擬合尚可。從左往右觀察相位曲線可以發(fā)現(xiàn)，相位曲線的質(zhì)量并不是逐步改善的。從圖10 第2張和第4張小圖中可以看出，加入1個(gè)白天的數(shù)據(jù)后，相位反而變差了。由此可知，由于白天的數(shù)據(jù)干擾較強(qiáng)，對(duì)曲線的整體改善有限，因此在數(shù)據(jù)編輯時(shí)可以有所取舍，若白天數(shù)據(jù)質(zhì)量較差，可在處理時(shí)截去整段數(shù)據(jù)中白天時(shí)段的數(shù)據(jù)。同樣地，建議在進(jìn)行野外數(shù)據(jù)采集時(shí)，應(yīng)盡量多采集夜晚的數(shù)據(jù)，即在下午布設(shè)儀器，根據(jù)周圍的干擾情況確定采集天數(shù)后，在上午結(jié)束采集，這樣可保證夜晚的數(shù)據(jù)必然存在且比白天多。

圖10 305測(cè)點(diǎn)不同時(shí)長數(shù)據(jù)對(duì)應(yīng)的視電阻率及相位曲線圖Fig.10 Apparent resistivity and phase curves of MT station 309 with different duration.

2.6 采集時(shí)長

為了探究采集時(shí)間對(duì)數(shù)據(jù)的影響，我們首先使用合成數(shù)據(jù)研究采集時(shí)長與有效周期的關(guān)系。選取一個(gè)電阻率為 100Ω·m的均勻半空間模型，假設(shè)電場信號(hào)E(ω)穩(wěn)定為1mV/km，根據(jù)H(ω)=E(ω)/Z(ω)=1/Z(ω)，設(shè)計(jì)不同的采樣時(shí)長，用傅里葉逆變換生成電場和磁場的時(shí)間序列，然后對(duì)時(shí)間序列進(jìn)行處理得到視電阻率和阻抗相位。如圖11 所示，由于頻譜泄漏和柵欄效應(yīng)的影響，時(shí)長為8s的數(shù)據(jù)經(jīng)傅里葉變換后其電磁場的振幅譜并不能完全準(zhǔn)確地還原至所能處理的最小頻率0.125Hz處。實(shí)際上，振幅譜在約2Hz處開始偏離理論值，在約1Hz處達(dá)到最低值，相應(yīng)的阻抗在2Hz處開始發(fā)生偏移，1Hz后的阻抗已遠(yuǎn)遠(yuǎn)偏離理論值。采用不同長度的合成數(shù)據(jù)進(jìn)行計(jì)算都能得到一致的結(jié)果。因此我們認(rèn)為，大地電磁數(shù)據(jù)理論上的最大有效周期是數(shù)據(jù)時(shí)長的1/8。

圖11 采樣率為4096Hz、時(shí)長為8s的合成數(shù)據(jù)的處理結(jié)果Fig.11 Synthetic data processing results with a sampling rate of 4096Hz and a length of 8s.a 振幅譜；b 視電阻率；c 相位

為了驗(yàn)證以上的理論推測(cè)，我們選取極低頻電磁臺(tái)網(wǎng)中豐寧臺(tái)的數(shù)據(jù)分析。該臺(tái)建有一個(gè)連續(xù)觀測(cè)的大地電磁觀測(cè)站，可以采集長時(shí)間段的大地電磁數(shù)據(jù)。同時(shí)，該觀測(cè)站的儀器布設(shè)相比較常規(guī)的野外儀器布設(shè)，其電極和磁探頭埋設(shè)更深，數(shù)據(jù)采集更穩(wěn)定。圖12a 中24h的數(shù)據(jù)表明該臺(tái)站的數(shù)據(jù)非常穩(wěn)定，僅通過一個(gè)譜就可以得到連續(xù)的阻抗。截取不同時(shí)長的數(shù)據(jù)進(jìn)行處理，得到有效周期與最短數(shù)據(jù)時(shí)長的對(duì)應(yīng)關(guān)系如圖12b 所示，最短數(shù)據(jù)時(shí)長約為周期的8～10倍。

圖12 豐寧臺(tái)實(shí)測(cè)數(shù)據(jù)的處理結(jié)果Fig.12 Measured data processing result of Fengning station.a 24h數(shù)據(jù)處理得到的視電阻率和相位曲線；b 獲取有效阻抗所需的最短數(shù)據(jù)時(shí)長

從圖10 可見，305測(cè)點(diǎn)的中頻存在干擾，隨著數(shù)據(jù)時(shí)長的增加，數(shù)據(jù)逐漸改善。有效譜的個(gè)數(shù)與頻率成正比，在高頻部分可得到的有效譜個(gè)數(shù)遠(yuǎn)多于低頻部分。隨著觀測(cè)時(shí)長的增加，高頻部分的視電阻率相較于低頻部分改善得更快。因此，對(duì)于中高頻干擾，延長觀測(cè)時(shí)間即可改善數(shù)據(jù)質(zhì)量。

對(duì)比豐寧臺(tái)和305測(cè)點(diǎn)的數(shù)據(jù)可以發(fā)現(xiàn)，不同噪聲條件下的數(shù)據(jù)質(zhì)量相差很大。對(duì)于低噪聲干擾地區(qū)只需記錄10倍有效周期時(shí)長的數(shù)據(jù)，得到一個(gè)有效譜即可；而高噪聲干擾地區(qū)可能需要記錄100倍有效周期甚至更長時(shí)間的數(shù)據(jù)才能得到可靠結(jié)果。因此，在野外觀測(cè)時(shí)，應(yīng)根據(jù)實(shí)際情況合理調(diào)整觀測(cè)時(shí)長，才能在滿足數(shù)據(jù)質(zhì)量要求的前提下以較高的效率完成數(shù)據(jù)采集。

3 結(jié)論

本文結(jié)合實(shí)測(cè)數(shù)據(jù)分析了數(shù)據(jù)處理時(shí)存在的問題，并針對(duì)不同問題提出了數(shù)據(jù)處理技術(shù)方案：

(1)設(shè)置遠(yuǎn)離觀測(cè)點(diǎn)的、低噪聲干擾的遠(yuǎn)參考點(diǎn)是處理強(qiáng)干擾環(huán)境下大地電磁數(shù)據(jù)的必要條件。對(duì)于死頻帶內(nèi)信號(hào)弱的測(cè)點(diǎn)，建議采用Robust技術(shù)進(jìn)行處理；對(duì)于近場干擾強(qiáng)且持續(xù)時(shí)間長的測(cè)點(diǎn)，建議采用非Robust技術(shù)進(jìn)行處理；對(duì)于2種情況都存在的數(shù)據(jù)，可以嘗試同時(shí)使用Robust和非Robust方法得到的功率譜參與挑選。功率譜的人工挑選對(duì)于結(jié)果的影響十分重要，以上3種處理方法得到的功率譜都必須經(jīng)過人工精細(xì)挑選后疊加，才能得到最終的阻抗估計(jì)。

(2)延長采集時(shí)間增加有效譜個(gè)數(shù)，同時(shí)適當(dāng)增加分段數(shù)，減少編輯前的功率譜疊加，并對(duì)大量功率譜進(jìn)行人工挑選可有效改善阻抗估計(jì)結(jié)果。干擾信號(hào)以天為周期，白天干擾強(qiáng)，夜晚干擾弱，在數(shù)據(jù)采集時(shí)應(yīng)盡量保證夜晚數(shù)據(jù)采集時(shí)長，同時(shí)在數(shù)據(jù)處理時(shí)有所取舍，可以舍去部分白天的時(shí)間序列，也可以在功率譜挑選時(shí)去除白天含干擾的功率譜。

(3)理想情況下，最短數(shù)據(jù)采集時(shí)長為最長有效周期的8倍。因此，低噪聲干擾區(qū)記錄10倍有效周期時(shí)長的數(shù)據(jù)就可以獲得比較好的阻抗估計(jì)結(jié)果，而高噪聲干擾地區(qū)可能需要記錄100倍有效周期甚至更長時(shí)間的數(shù)據(jù)才能得到可靠的結(jié)果。在野外觀測(cè)需要根據(jù)實(shí)際噪聲環(huán)境合理設(shè)計(jì)觀測(cè)時(shí)間。

本文的研究區(qū)內(nèi)存在高壓線、發(fā)電站、大型工廠等諸多干擾，但并未單獨(dú)就某種干擾進(jìn)行具體分析。本文所述的Robust、分段數(shù)、分時(shí)段處理等方法適用于干擾持續(xù)存在但部分時(shí)段干擾較弱的情況?？偠灾?，數(shù)據(jù)受干擾的情況不盡相同，數(shù)據(jù)處理的技術(shù)手段也不應(yīng)一成不變。在更好的數(shù)據(jù)處理技術(shù)尚未成熟時(shí)，應(yīng)靈活應(yīng)用已有的技術(shù)對(duì)大地電磁數(shù)據(jù)進(jìn)行處理。