李展輝,楊淼鑫,曹學(xué)峰
(1. 中國(guó)自然資源航空物探遙感中心,北京 100083;2. 河北地質(zhì)大學(xué) 河北省戰(zhàn)略性關(guān)鍵礦產(chǎn)資源重點(diǎn)實(shí)驗(yàn)室,河北 石家莊 050031)
瞬變電磁法對(duì)于地下礦產(chǎn)與水資源勘查以及基礎(chǔ)地質(zhì)調(diào)查具有重要意義。常規(guī)的回線源瞬變電磁法采用感應(yīng)線圈接收發(fā)射電路關(guān)斷后地下介質(zhì)產(chǎn)生的二次感應(yīng)電壓,即接收線圈中磁通量隨時(shí)間的變化率。感應(yīng)電壓場(chǎng)在晚期隨時(shí)間以t-5/2的方式迅速衰減[1],探測(cè)深度普遍限制在1 km以內(nèi)。隨著地表資源的減少,我國(guó)對(duì)深部礦產(chǎn)資源勘查的需求日益迫切,開始大力開展深部礦產(chǎn)資源勘查工作。針對(duì)這一需求,近年來,瞬變電磁的磁感應(yīng)強(qiáng)度B場(chǎng)觀測(cè)越來越受到人們的關(guān)注。在理論上,B場(chǎng)在晚期隨時(shí)間以t-3/2的方式衰減[1],相比感應(yīng)電壓場(chǎng)似乎能在更長(zhǎng)的時(shí)間內(nèi)獲取有效數(shù)據(jù),有利于深部資源的探測(cè)。此外,研究表明,在某些條件下(如電性源瞬變電磁),磁感應(yīng)強(qiáng)度垂直分量Bz場(chǎng)的視電阻率是單一的[2],而感應(yīng)電壓場(chǎng)視電阻率具有多解性,使用Bz場(chǎng)更方便于視電阻率的計(jì)算?;谝陨螧場(chǎng)探測(cè)的優(yōu)勢(shì),一些研究人員研制了磁通門與基于超導(dǎo)量子干涉器(superconducting quantum interference device, SQUID)的磁強(qiáng)計(jì),并開展了實(shí)際測(cè)試與應(yīng)用。澳大利亞的EMIT公司開發(fā)了三軸B場(chǎng)傳感器SMART Fluxgate,已經(jīng)上市售賣并在實(shí)際應(yīng)用中使用;王興春等利用該磁通門探頭在某銅鎳礦區(qū)域開展了探測(cè)效果分析[3];吉林大學(xué)等相關(guān)單位聯(lián)合開發(fā)了基于低溫SQUID的瞬變電磁B場(chǎng)傳感器并開展了試驗(yàn)應(yīng)用,取得了較好的效果[4-5];陳曉東等開發(fā)了高溫SQUID磁強(qiáng)計(jì)并在野外試驗(yàn)應(yīng)用,采集到了精度顯著高于傳統(tǒng)傳感器的B場(chǎng)數(shù)據(jù)[6];Geotech、SkyTEM等航空電磁研制與應(yīng)用公司已經(jīng)分別在各自的部分航空瞬變電磁探測(cè)系統(tǒng)中集成了B場(chǎng)采集系統(tǒng),能夠?qū)崟r(shí)收錄B場(chǎng)[7];Smith與Annan從數(shù)據(jù)本身角度分析了B場(chǎng)相對(duì)感應(yīng)電壓場(chǎng)的優(yōu)勢(shì)[8],并研究了從傳統(tǒng)感應(yīng)線圈中獲取B場(chǎng)數(shù)據(jù)的方法[9];Wolfgram與Thomson在高導(dǎo)地區(qū)使用了航空電磁B場(chǎng)進(jìn)行了分析應(yīng)用[10];崔江偉等研究了接地導(dǎo)線源瞬變電磁感應(yīng)電壓場(chǎng)——B場(chǎng)轉(zhuǎn)換方法,并開展了對(duì)針對(duì)Bz場(chǎng)的視電阻率計(jì)算方法研究,取得了良好的效果[2]。
然而,B場(chǎng)的探測(cè)效果研究大多僅限于此,探測(cè)結(jié)果的橫向分辨率、縱向分辨率、電阻率反演的準(zhǔn)確性等重要參數(shù)少有分析。基于以上現(xiàn)狀,且在目前B場(chǎng)研究與應(yīng)用主要集中在Bz場(chǎng)的基礎(chǔ)上,本文對(duì)Bz場(chǎng)與垂直歸一化感應(yīng)電壓場(chǎng)-dBz/dt的探測(cè)效果進(jìn)行更深入分析。采用回線源瞬變電磁法一維或三維正演算法計(jì)算出一維模型或者三維模型的感應(yīng)電壓場(chǎng)與Bz場(chǎng)數(shù)據(jù),然后針對(duì)合成數(shù)據(jù)開展數(shù)據(jù)對(duì)比、視電阻率對(duì)比、反演結(jié)果對(duì)比等探測(cè)效果分析,對(duì)比感應(yīng)電壓場(chǎng)和Bz場(chǎng)的縱向分辨率、橫向分辨率、電阻率反演結(jié)果的準(zhǔn)確性,為實(shí)際瞬變電磁探測(cè)應(yīng)用提供指導(dǎo)。
本文采用一維或三維正演算法生成合成的感應(yīng)電壓場(chǎng)和Bz場(chǎng)數(shù)據(jù),采用全時(shí)視電阻率計(jì)算方法[2,11]計(jì)算各自的視電阻率,采用一維反演方法[12]進(jìn)行反演對(duì)比分析。
針對(duì)如圖1a所示的層狀模型以及鋪設(shè)于地表之上的水平回線源和接收點(diǎn),可以首先采用水平電偶極子積分的方式計(jì)算出頻率域的磁場(chǎng)響應(yīng),然后采用正弦變換或者余弦變換將頻率域磁場(chǎng)響應(yīng)變換成時(shí)間域的感應(yīng)電壓場(chǎng)或者B場(chǎng)。
圖1 一維層狀模型示意(a)及組成回線源的水平電偶極子AB與接收點(diǎn)C之間的參數(shù)R′與參數(shù)y′示意(b)Fig.1 A schematic figure of the one-dimensional layered model (a) and the delineation of the parameter R′ and y′ between the horizontal electric dipole AB that constitute the loop source and the receiving point C (b)
圖1b中AB為位于地面上方h處的水平電偶極子,電流方向?yàn)锽→A;C為位于地面上方zr處的接收點(diǎn);R′代表接收點(diǎn)C到偶極子AB之間的水平距離;y′代表接收點(diǎn)C到直線AB之間的距離,其正負(fù)符號(hào)由右手原則確定。對(duì)于水平發(fā)射線圈,在接收點(diǎn)C處的頻率域垂直磁場(chǎng)響應(yīng)可表示為[12]
(1)
式中:I為偶極子電流;h為發(fā)射線圈的離地高度;zr為接收點(diǎn)的離地高度;ds為偶極子長(zhǎng)度;λ為漢克爾積分變量;rTE為地表TE模式反射率,由地下地層的電性參數(shù)與厚度參數(shù)組合迭代而成[12]。
求出頻率域磁場(chǎng)后,階躍關(guān)斷源在接收點(diǎn)的歸一化感應(yīng)電壓值Vz(t)或時(shí)域磁感應(yīng)場(chǎng)值bz(t)可以通過正弦變換或余弦變換來實(shí)現(xiàn):
(2)
(3)
受到電子元器件性能限制等因素的影響,實(shí)際發(fā)射波形不可能是理想階躍波形,有一定的開啟與關(guān)斷時(shí)間。為準(zhǔn)確描述接收的感應(yīng)電壓場(chǎng)和Bz場(chǎng),需要計(jì)算全波形的場(chǎng)值。全波形場(chǎng)值可在階躍場(chǎng)值的基礎(chǔ)上通過卷積的方式得到:
(4)
(5)
三維正演采用基于Du Fort-Frankel格式的時(shí)域有限差分方法[13-15]。對(duì)于瞬變電磁法,可以針對(duì)如下Maxwell旋度方程組進(jìn)行二階差分離散:
(6)
式中:E為電場(chǎng);H為磁場(chǎng);μ為磁導(dǎo)率;σ為電導(dǎo)率;γ為確保數(shù)值計(jì)算穩(wěn)定的人為設(shè)置介電常數(shù),隨時(shí)間步長(zhǎng)的變化而變化。
為能夠計(jì)算B場(chǎng),需要將回線源通電時(shí)的靜磁場(chǎng)作為初始磁場(chǎng)加載到正演程序中。為確保磁場(chǎng)散度為零的特征,本文通過磁矢勢(shì)來計(jì)算初始磁場(chǎng)[16],然后加載到正演程序中,同時(shí)計(jì)算出三維模型的B場(chǎng)和感應(yīng)電壓場(chǎng)。三維正演直接計(jì)算出來的是階躍關(guān)斷電流所產(chǎn)生的感應(yīng)電壓場(chǎng)和B場(chǎng),在此基礎(chǔ)上,結(jié)合發(fā)射電流的實(shí)際波形,利用式(4)和式(5)可分別求出全波形的感應(yīng)電壓場(chǎng)和B場(chǎng)。
本文主要對(duì)比感應(yīng)電壓與Bz場(chǎng)的視電導(dǎo)率與反演成像結(jié)果。視電導(dǎo)率的計(jì)算以正演程序?yàn)榛A(chǔ),利用二分法進(jìn)行全時(shí)視電導(dǎo)率計(jì)算[2,11]。本文在視電導(dǎo)率計(jì)算過程中,為準(zhǔn)確獲取感應(yīng)電壓與Bz場(chǎng)對(duì)應(yīng)的視電導(dǎo)率值,二分法搜索的停止閾值設(shè)置為
(7)
式中:dapp為利用視電導(dǎo)率正演所得數(shù)據(jù);dobs為模型的觀測(cè)數(shù)據(jù)。反演則采用Tikhonov正則化約束進(jìn)行反演。反演目標(biāo)函數(shù)可表示為
(8)
式中:d為測(cè)量的數(shù)據(jù);G(m)為正演所得數(shù)據(jù);Wd為數(shù)據(jù)權(quán)重矩陣;m為模型參數(shù);mref為參考模型參數(shù);L為模型約束矩陣;α為正則化參數(shù)。
對(duì)上述目標(biāo)函數(shù)進(jìn)行最小化迭代,可得到關(guān)于第k+1步模型參數(shù)增量Δm的最小二乘方程組:
(9)
式中:J(mk)=?G(mk)/?mk為第k步正演數(shù)據(jù)對(duì)模型參數(shù)的敏感度矩陣。通過求解這一最小二乘方程組,第k+1步模型參數(shù)可表示為mk+1=mk+Δm。直到反演殘差達(dá)到預(yù)定的值,或者不再下降,或者反演次數(shù)達(dá)到預(yù)定的值,反演迭代停止。
針對(duì)本文中的反演問題,由于電導(dǎo)率的變化范圍跨越數(shù)個(gè)量級(jí),感應(yīng)電壓和Bz場(chǎng)的變化范圍也是數(shù)個(gè)量級(jí),為確保反演問題的穩(wěn)定性,一般對(duì)d的對(duì)數(shù)值或者相對(duì)值進(jìn)行擬合,對(duì)電導(dǎo)率的對(duì)數(shù)值進(jìn)行反演。本文在數(shù)據(jù)擬合方面對(duì)d的相對(duì)值進(jìn)行擬合,在反演參數(shù)方面對(duì)電導(dǎo)率的對(duì)數(shù)值進(jìn)行反演,然后轉(zhuǎn)換為電導(dǎo)率。在該條件下,數(shù)據(jù)權(quán)重矩陣Wd定義為
(10)
式中:si為第i個(gè)數(shù)據(jù)di的相對(duì)誤差;Nd為數(shù)據(jù)的數(shù)量。如果某次測(cè)量的絕對(duì)誤差水平為Aerri, 最小相對(duì)誤差水平為Rerri, 那么si可如下確定
(11)
反演均方根(root-mean-square,rms)殘差可定義為
(12)
反演終止條件為rms≤1,或者反演次數(shù)達(dá)到設(shè)定的最大值。由于本文針對(duì)擬合數(shù)據(jù)進(jìn)行反演,rms≤1這一條件總是能夠達(dá)成,因此不再設(shè)置反演次數(shù)上限。
本文在一般地層假設(shè)條件下進(jìn)行反演,即無特定的參考模型,不假定地下電性異常層的數(shù)量、厚度等。因此參考模型mref設(shè)置為空,反演初始模型為電導(dǎo)率為0.001 S/m的均勻半空間模型,模型約束矩陣L為縱向一階導(dǎo)數(shù)約束,正則化參數(shù)α根據(jù)反演情況自適應(yīng)變化。
在均勻半空間條件下,階躍關(guān)斷電流源的Bz場(chǎng)在晚期以t-3/2方式衰減,而感應(yīng)電壓場(chǎng)在晚期以t-5/2方式衰減,因此Bz場(chǎng)在晚期似乎有更大的優(yōu)勢(shì),能夠獲取更晚期信息。然而考慮到實(shí)際觀測(cè)設(shè)備的噪聲水平,Bz場(chǎng)相比感應(yīng)電壓在何種情況下能夠觀測(cè)到更長(zhǎng)時(shí)間需要仔細(xì)分析。目前基于SQUID的磁強(qiáng)計(jì)尚處于試驗(yàn)階段,未進(jìn)行大規(guī)模商業(yè)化應(yīng)用,市場(chǎng)上觀測(cè)B場(chǎng)的主要設(shè)備為澳大利亞EMIT公司的Fluxgate三軸磁通門探頭,探頭精度為3 pT。加拿大Geonics公司的PROTEM系統(tǒng)技術(shù)指標(biāo)達(dá)到了國(guó)際領(lǐng)先,感應(yīng)線圈探頭精度約為0.1 nT/s。因此本文以這兩個(gè)精度閾值為標(biāo)準(zhǔn),考察在不同電導(dǎo)率條件下Bz場(chǎng)和感應(yīng)電壓場(chǎng)分別達(dá)到這兩個(gè)閾值時(shí)的延遲時(shí)間,即精度閾值時(shí)間。
設(shè)有一尺寸為200 m×200 m的單匝發(fā)射線圈放置于一個(gè)均勻半空間模型的地表,在發(fā)射線圈的中心點(diǎn)接收,如圖2所示。均勻半空間模型的電導(dǎo)率值在0.000 2~0.5 S/m之間連續(xù)變化,所采用的發(fā)射電流的波形如圖3所示,為一斜階躍關(guān)斷波形,電流峰值為10 A,關(guān)斷延遲時(shí)間為0.1 ms。時(shí)間零點(diǎn)設(shè)置為電流關(guān)斷至零的時(shí)刻,數(shù)據(jù)采集延遲時(shí)間從此刻開始計(jì)時(shí)。針對(duì)每個(gè)電導(dǎo)率值計(jì)算出感應(yīng)電壓場(chǎng)和Bz場(chǎng),并計(jì)算每個(gè)電導(dǎo)率對(duì)應(yīng)的場(chǎng)值達(dá)到各自探頭精度閾值水平時(shí)的延遲時(shí)間,最后得出感應(yīng)電壓場(chǎng)和Bz場(chǎng)精度閾值時(shí)間隨電導(dǎo)率值變化規(guī)律(圖4)。
圖2 發(fā)射線框尺寸以及接收點(diǎn)的相對(duì)位置Fig.2 The size of the transmitting loop and the relativeposition of the receiving point
圖3 圖2所示發(fā)射線框中的發(fā)射電流波形Fig.3 The waveform of the current in the transmitting loop described in Fig.2
圖4 感應(yīng)電壓場(chǎng)和Bz場(chǎng)精度閾值時(shí)間隨半空間模型電導(dǎo)率的變化Fig.4 The precision threshold delay time of inducedvoltage and Bz varying over the conductivity of the homogeneous half space model
由圖4可見,Bz場(chǎng)的精度閾值時(shí)間與感應(yīng)電壓場(chǎng)的精度閾值時(shí)間呈正比相關(guān),在選取的電導(dǎo)率范圍之內(nèi)都隨著電導(dǎo)率的增大而增大。經(jīng)過進(jìn)一步分析,在電導(dǎo)率值小于0.036 S/m時(shí),Bz場(chǎng)的精度閾值時(shí)間小于感應(yīng)電壓場(chǎng)的精度閾值時(shí)間,而在電導(dǎo)率值大于0.036 S/m時(shí),Bz場(chǎng)的精度閾值時(shí)間大于感應(yīng)電壓場(chǎng)的精度閾值時(shí)間。
圖5給出了電導(dǎo)率為0.036 S/m時(shí)回線中心點(diǎn)的感應(yīng)電壓和Bz衰減曲線,當(dāng)感應(yīng)電壓衰減到0.1 nT/s時(shí),Bz場(chǎng)正好衰減到3 pT。可見對(duì)于本文采用的當(dāng)前主流的探頭精度,以及本文所設(shè)定的常用的TEM測(cè)量方式,只有地下存在電導(dǎo)率高于0.036 S/m(或者電阻率小于27.8 Ω·m)的地質(zhì)體時(shí),Bz場(chǎng)的精度閾值時(shí)間才有可能長(zhǎng)于感應(yīng)電壓場(chǎng)的精度閾值時(shí)間,這一電導(dǎo)率閾值已經(jīng)達(dá)到了塊狀硫化物礦、含水泥質(zhì)黏土、微咸水層等地質(zhì)體電導(dǎo)率范圍。因此,Bz場(chǎng)在探測(cè)良導(dǎo)電異常體方面具有一定的優(yōu)勢(shì),能在更長(zhǎng)的延遲時(shí)間內(nèi)觀測(cè);但是在探測(cè)低導(dǎo)電異常體方面,比如采空區(qū)的探測(cè),效果可能不如感應(yīng)電壓場(chǎng)。
圖5 半空間模型電導(dǎo)率為0.036 S/m時(shí)感應(yīng)電壓和Bz場(chǎng)隨延遲時(shí)間的衰減曲線Fig.5 The attenuation curve of the induced voltageand Bz field of the homogenous model with a conductivity of 0.036 S/m
縱向分辨率對(duì)于地質(zhì)體分層、深度的確認(rèn)至關(guān)重要。本文從視電阻率和一維反演角度詳細(xì)分析感應(yīng)電壓場(chǎng)和Bz場(chǎng)對(duì)地層的縱向分辨能力。
首先針對(duì)高導(dǎo)異常層模型進(jìn)行縱向分辨能力測(cè)試。模型如圖6所示。第一層厚度30 m,電導(dǎo)率為0.01 S/m;第二層為一個(gè)厚度30 m、電導(dǎo)率0.05 S/m的高導(dǎo)異常層;第三層厚度在不同模型中厚度可變,從10 m厚度,以2 m為間隔逐漸變化到208 m,變化趨勢(shì)如圖6b所示,電導(dǎo)率為0.01 S/m; 第四層為一個(gè)厚度為100 m的高導(dǎo)異常層,電導(dǎo)率為0.05 S/m;第五層為一個(gè)電導(dǎo)率為0.01 S/m的半空間。本測(cè)試通過第三層的厚度變化來分析瞬變電磁法感應(yīng)電壓場(chǎng)和Bz場(chǎng)對(duì)第二層和第四層的縱向分辨能力。對(duì)每個(gè)模型采用圖2所示的中心回線裝置進(jìn)行測(cè)量,并將所得數(shù)據(jù)拼接成擬二維數(shù)據(jù)剖面,感應(yīng)電壓場(chǎng)和Bz場(chǎng)的數(shù)據(jù)分別如圖7a和圖7b所示。從數(shù)據(jù)上難以直觀地看出縱向分辨率信息。因此本文對(duì)正演所獲取的數(shù)據(jù)開展了全時(shí)視電導(dǎo)率成像與一維反演成像分析。
圖6 用于縱向分辨能力測(cè)試的高導(dǎo)異常模型示意(a)以及由此生成的100個(gè)參與實(shí)際計(jì)算的模型(b)Fig.6 The model framework containing high conductive layers for longitudinal resolution test (a) and the 100 actural models used in calculation (b)
圖7 針對(duì)圖6所示模型計(jì)算出的感應(yīng)電壓場(chǎng)(a)和Bz場(chǎng)(b)數(shù)據(jù)剖面Fig.7 The induced voltage (a) and Bz data (b) calculated from the models described in Fig.6
以一維正演程序?yàn)榛A(chǔ),開展二分法全時(shí)視電導(dǎo)率計(jì)算分析。感應(yīng)電壓和Bz場(chǎng)全時(shí)視電導(dǎo)率結(jié)果分別如圖8a和圖8b所示。視電導(dǎo)率計(jì)算過程中沒有考慮數(shù)據(jù)的誤差情況,以式(7)作為視電阻率計(jì)算終止條件,獲取高精度的視電阻率值。
圖8 針對(duì)圖7中感應(yīng)電壓場(chǎng)(a)和Bz場(chǎng)(b)計(jì)算得到的全時(shí)視電導(dǎo)率擬二維圖Fig.8 The all-time apparent conductivities calculated from the induced voltage data (a) and Bz data (b) shown in Fig.7
為定量分析第二和第四2個(gè)高導(dǎo)異常層剛好被分辨時(shí)第三層的厚度,本文定義視電導(dǎo)率曲線從單峰值曲線變?yōu)殡p峰值曲線時(shí)對(duì)應(yīng)的模型第三層厚度為視電導(dǎo)率分辨臨界厚度。在此臨界厚度定義下進(jìn)行分析,得出感應(yīng)電壓場(chǎng)視電導(dǎo)率在第25個(gè)模型,即第三層厚度為58 m時(shí),出現(xiàn)了雙峰現(xiàn)象,此時(shí)Bz場(chǎng)視電導(dǎo)率曲線依然是單峰曲線,如圖9a所示;而Bz場(chǎng)直到第35個(gè)模型,即第三層厚度為78 m時(shí),才出現(xiàn)雙峰現(xiàn)象,此時(shí)感應(yīng)電壓場(chǎng)的視電導(dǎo)率具有更明顯的雙峰現(xiàn)象,如圖9b所示??梢?,針對(duì)高導(dǎo)異常層,感應(yīng)電壓場(chǎng)基于視電導(dǎo)率的分辨臨界厚度要小于Bz的分辨臨界厚度,具有更好的縱向分辨能力。同時(shí)也可以看出,Bz場(chǎng)視電導(dǎo)率的峰值時(shí)間一直早于感應(yīng)電壓場(chǎng)視電導(dǎo)率的峰值時(shí)間,Bz場(chǎng)在更早的時(shí)間感知到了高導(dǎo)異常層??梢猿醪脚袛?,Bz場(chǎng)能夠更早地感知到地下高導(dǎo)異常體信息,但是在縱向分辨率上不如感應(yīng)電壓場(chǎng)。
圖9 圖6中模型第三層厚度分別達(dá)到58 m(感應(yīng)電壓視電阻率臨界厚度)(a)和78 m(Bz視電阻率臨界厚度)(b)時(shí)對(duì)應(yīng)的感應(yīng)電壓視電導(dǎo)率和Bz場(chǎng)視電導(dǎo)率曲線Fig.9 The induced voltage all-time apparent conductivity curve and Bz all-time apparent conductivity curve of the modelin Fig.6 when the third layer comes to 58 m (critical resolution thickness of induced voltage apparent conductivity) (a) and 78 m (critical resolution thickness of Bz apparent conductivity) (b), respectively
實(shí)際應(yīng)用中,更重要的是以反演結(jié)果進(jìn)行探測(cè)效果分析,因此進(jìn)一步對(duì)數(shù)據(jù)開展了一維反演分析。在反演過程中,考慮到實(shí)際測(cè)量數(shù)據(jù)的誤差以及儀器探頭自身的精度水平,對(duì)數(shù)據(jù)設(shè)定了一定的數(shù)據(jù)誤差。參考PROTEM系統(tǒng),感應(yīng)電壓場(chǎng)的絕對(duì)誤差水平Aerr為0.1 nT/s,同時(shí)設(shè)定最小相對(duì)誤差水平Rerr為3%;參考EMIT的磁通門探頭,Bz場(chǎng)的絕對(duì)誤差水平為3 pT,同時(shí)設(shè)定最小相對(duì)誤差水平為3%;在此誤差水平下,反演過程中的數(shù)據(jù)權(quán)重矩陣可以用式(10)和式(11)共同計(jì)算出來。分別對(duì)感應(yīng)電壓場(chǎng)與Bz場(chǎng)進(jìn)行反演,結(jié)果分別如圖10a和圖10b所示。與視電導(dǎo)率分析類似,定義反演電導(dǎo)率曲線從單峰值曲線變?yōu)殡p峰值曲線時(shí)對(duì)應(yīng)的模型第三層厚度為反演分辨臨界厚度?;诖伺R界厚度分析,得出感應(yīng)電壓場(chǎng)反演電導(dǎo)率在第6個(gè)模型,即第三層厚度為20 m時(shí),出現(xiàn)了雙峰現(xiàn)象,此時(shí)Bz場(chǎng)反演電導(dǎo)率曲線依然是單峰曲線,如圖11a所示;而Bz場(chǎng)直到第36個(gè)模型,即第三層厚度為80 m時(shí),才出現(xiàn)雙峰現(xiàn)象,此時(shí)感應(yīng)電壓場(chǎng)的反演電導(dǎo)率具有非常明顯的雙峰現(xiàn)象,如圖11b所示。從反演結(jié)果中,可以看到感應(yīng)電壓場(chǎng)在縱向分辨能力上具有更大的優(yōu)勢(shì)。
圖10 針對(duì)圖7中感應(yīng)電壓場(chǎng)(a)和Bz場(chǎng)(b)反演得到的電導(dǎo)率擬二維剖面Fig.10 The conductivities inverted from the induced voltage data (a) and Bz data (b) shown in Fig.7
圖11 圖6中模型第三層厚度分別達(dá)到20 m(感應(yīng)電壓反演臨界厚度)(a)和80 m(Bz反演臨界厚度)(b)時(shí)對(duì)應(yīng)的感應(yīng)電壓反演電導(dǎo)率和Bz場(chǎng)反演電導(dǎo)率曲線Fig.11 The induced voltage inverted conductivity curve and Bz inverted conductivity curve of the model in Fig.6 when the third layer comes to 20 m (critical resolution thickness of induced voltage inversion) (a) and 80 m (critical resolution thickness of Bz inversion) (b), respectively
為進(jìn)一步分析感應(yīng)電壓場(chǎng)和Bz場(chǎng)對(duì)低導(dǎo)異常層的縱向分辨能力,將圖6所示模型中的第二、第四層改為電導(dǎo)率為0.002 S/m的低導(dǎo)異常層,其余模型參數(shù)均不變,模型如圖12所示。首先正演出一維數(shù)據(jù),然后開展視電導(dǎo)率和反演成像分析。感應(yīng)電壓和Bz場(chǎng)數(shù)據(jù)分別如圖13a和圖13b所示,視電導(dǎo)率成像分別如圖14a和圖14b所示。同樣開展視電導(dǎo)率分辨臨界厚度分析,由于回線源瞬變電磁法對(duì)低導(dǎo)異常層的敏感度天然較低,感應(yīng)電壓場(chǎng)直到第32個(gè)模型,即第三層厚度為72 m時(shí),視電導(dǎo)率曲線才出現(xiàn)極小雙谷現(xiàn)象,如圖15a所示;而Bz場(chǎng)則是到第41個(gè)模型,即第三層厚度為90 m時(shí),視電導(dǎo)率曲線才出現(xiàn)極小雙谷現(xiàn)象,如圖15b所示。可見,對(duì)于低導(dǎo)異常層,感應(yīng)電壓場(chǎng)的視電導(dǎo)率分辨能力同樣強(qiáng)于Bz場(chǎng)。不過仍然可以看到Bz場(chǎng)視電導(dǎo)率出現(xiàn)極小波谷的時(shí)間均提前于感應(yīng)電壓場(chǎng)。
圖12 用于縱向能力測(cè)試的低導(dǎo)異常模型示意(a)以及由此生成的100個(gè)參與實(shí)際計(jì)算的模型(b)Fig.12 The model framework containing low conductive layers for longitudinal resolution test (a) and the 100 actural models used in calculation (b)
圖13 針對(duì)圖12所示模型計(jì)算出的感應(yīng)電壓場(chǎng)(a)和Bz場(chǎng)(b)數(shù)據(jù)剖面Fig.13 The induced voltage (a) and Bz data (b) calculated from the models described in Fig.12
圖14 針對(duì)圖13中感應(yīng)電壓場(chǎng)(a)和Bz場(chǎng)(b)計(jì)算得到的全時(shí)視電導(dǎo)率擬二維圖Fig.14 The all-time apparent conductivities calculated from the induced voltage data (a) and Bz data (b) shown in Fig.13
圖15 圖12中模型第三層厚度分別達(dá)到72 m(感應(yīng)電壓視電阻率臨界厚度)(a)和90 m(Bz視電阻率臨界厚度)(b)時(shí)對(duì)應(yīng)的感應(yīng)電壓視電導(dǎo)率和Bz場(chǎng)視電導(dǎo)率曲線Fig.15 The induced voltage all-time apparent conductivity curve and Bz all-time apparent conductivity curve of the model in Fig.12 when the third layer comes to 72 m (critical resolution thickness of induced voltage apparent conductivity) (a) and90m (critical resolution thickness of Bz apparent conductivity) (b), respectively
在反演分析中加入的噪聲與上一個(gè)高導(dǎo)異常模型一致。感應(yīng)電壓和Bz場(chǎng)反演結(jié)果分別如圖16a和圖16b所示。同樣開展反演電導(dǎo)率分辨臨界厚度分析。由于數(shù)據(jù)誤差的加入,感應(yīng)電壓場(chǎng)的分辨臨界模型推后到了第46個(gè)模型,即第三層厚度到100m時(shí),才能夠?qū)蓚€(gè)低導(dǎo)層區(qū)分開來,如圖17a所示;而Bz場(chǎng)受到數(shù)據(jù)噪聲的影響,反演電導(dǎo)率直到第100個(gè)模型也沒有出現(xiàn)雙極小波谷現(xiàn)象,全程未能分辨出兩個(gè)低導(dǎo)異常層,如圖17b所示??梢?,相比感應(yīng)電壓場(chǎng),Bz場(chǎng)對(duì)低導(dǎo)異常層的分辨率同樣更低,并且感應(yīng)電壓場(chǎng)反演得到的電導(dǎo)率異常值更接近于真實(shí)的電導(dǎo)率異常。
圖16 針對(duì)圖13中感應(yīng)電壓場(chǎng)(a)和Bz場(chǎng)(b)反演得到的電導(dǎo)率擬二維剖面Fig.16 The conductivities inverted from the induced voltage data (a) and Bz data (b) shown in Fig.13
圖17 圖12中模型第三層厚度分別達(dá)到100 m(感應(yīng)電壓反演臨界厚度)(a)和208 m(最后一個(gè)模型)(b)時(shí)對(duì)應(yīng)的感應(yīng)電壓反演電導(dǎo)率和Bz場(chǎng)反演電導(dǎo)率曲線Fig.17 The induced voltage inverted conductivity curve and Bz inverted conductivity curve of the model in Fig.12 when the third layer comes to 100 m (critical resolution thickness of induced voltage inversion) (a) and s208 m (the last model) (b), respectively
通過上述高導(dǎo)異常層模型和低導(dǎo)異常層模型的縱向分辨能力分析可知,Bz場(chǎng)相比感應(yīng)電壓場(chǎng)能夠提前感知到地下異常層,但是視電導(dǎo)率和反演電導(dǎo)率縱向分辨率明顯不如感應(yīng)電壓場(chǎng),感應(yīng)電壓場(chǎng)成像對(duì)電導(dǎo)率異常的呈現(xiàn)能力也更強(qiáng)。
橫向分辨率對(duì)于巖體邊界的劃分、礦體范圍的圈定十分重要。以航空TEM測(cè)量為例,目前主流的航空TEM系統(tǒng)水平采樣間隔在5 m左右,旨在利用高測(cè)點(diǎn)密度識(shí)別出目標(biāo)巖體、礦體邊界。為對(duì)比感應(yīng)電壓場(chǎng)與Bz場(chǎng)的橫向分辨率,本文設(shè)計(jì)了如圖18所示的模型。電導(dǎo)率為0.01 S/m的半空間中包含有6個(gè)電導(dǎo)率為0.1 S/m的異常體,編號(hào)從左至右依次為1~6。異常體頂面埋深為10 m,在y方向相當(dāng)于無限延伸。除最右側(cè)的6號(hào)異常體外,其余異常體在x和z(深度)方向的尺寸均為10 m,由于網(wǎng)格劃分的原因,6號(hào)異常體在x方向的尺寸略大于10 m,但不影響分析結(jié)果。異常體之間水平間隔從左至右分別為10 m、20 m、30 m、40 m、55 m。發(fā)射線圈直徑為26 m,匝數(shù)為5匝,發(fā)射電流峰值為175 A,總發(fā)射磁矩為46.5萬Am2,并假定發(fā)射波形為三角波,具體波形如圖19所示,時(shí)間零點(diǎn)為電流關(guān)斷至零的時(shí)刻。接收線圈位于發(fā)射吊艙中心,與發(fā)射線圈共平面。整個(gè)吊艙離地高度為30 m,以5 m采樣間隔沿x方向從-150 m位置飛行到150 m位置。電流關(guān)斷后的數(shù)據(jù)采集延遲時(shí)間范圍為0.01~2.51 ms,時(shí)間門個(gè)數(shù)為30。
圖18 用于橫向分辨能力測(cè)試的三維異常體模型Fig.18 The three-dimensional model for horizontal resolution test
圖19 正演計(jì)算中采用的航空瞬變電磁測(cè)量系統(tǒng)發(fā)射波形(時(shí)間零點(diǎn)為電流關(guān)斷至零的時(shí)刻)Fig.19 The transmitting waveform of the airborne transient electromagnetic survey system used in 3D forward modelling(the zero time is set at when the current descends to zero)
采用三維時(shí)域有限差分進(jìn)行正演計(jì)算,得到感應(yīng)電壓場(chǎng)和Bz場(chǎng)的階躍響應(yīng)數(shù)據(jù),然后通過時(shí)間域卷積得到全波形的感應(yīng)電壓場(chǎng)和Bz場(chǎng)數(shù)據(jù)(圖20)。數(shù)據(jù)上體現(xiàn)出的異常均比較弱,需要進(jìn)一步進(jìn)行電導(dǎo)率成像分析。由于對(duì)于此類三維模型數(shù)據(jù)的一維反演會(huì)產(chǎn)生較明顯的假象[17],而現(xiàn)階段三維反演并不成熟,因此本文在此僅開展全時(shí)視電導(dǎo)率計(jì)算,針對(duì)視電導(dǎo)率進(jìn)行對(duì)比。感應(yīng)電壓數(shù)據(jù)和Bz場(chǎng)數(shù)據(jù)計(jì)算出的全時(shí)視電阻率分別如圖21a和21b所示,圖中標(biāo)記出了4、5、6號(hào)異常體的水平位置。從圖中可以看出,4、5、6號(hào)異常體正上方的異常比較弱,而兩側(cè)的異常比較強(qiáng),形成假異常體。這是由于異常體尺寸相對(duì)發(fā)射線圈較小、分布較淺導(dǎo)致的。根據(jù)瞬變電磁場(chǎng)的煙圈理論,當(dāng)?shù)跖擄w到異常體正上方時(shí),感應(yīng)電磁場(chǎng)向四周擴(kuò)散,正下方的異常體的感應(yīng)強(qiáng)度反而比較弱,相對(duì)的旁側(cè)較近范圍的異常電磁感應(yīng)會(huì)比較強(qiáng)烈。針對(duì)該現(xiàn)象,必須進(jìn)行二維或者三維反演才能將異常體反演到正確的位置上,這一點(diǎn)將在后續(xù)工作中進(jìn)行研究。
圖20 針對(duì)圖18中的模型采用三維正演計(jì)算出的感應(yīng)電壓場(chǎng)(a)以及Bz場(chǎng)(b)剖面Fig.20 The induced voltage (a) and Bz data (b) calculated from the model described in Fig.18 by 3D forward modelling
a—由感應(yīng)電壓數(shù)據(jù)計(jì)算的全時(shí)視電導(dǎo)率擬二維剖面;b—由Bz數(shù)據(jù)計(jì)算的全時(shí)視電導(dǎo)率擬二維剖面;c—t=0.01 ms時(shí)刻感應(yīng)電壓與Bz場(chǎng)視電導(dǎo)率隨水平距離的變化曲線;虛線表示第4、5、6號(hào)異常體在x方向的具體位置;藍(lán)色倒三角形和紅色倒三角形分別為對(duì)應(yīng)視電阻率曲線的峰值a—the all-time apparent conductivities calculated from the induced voltage data;b—the all-time apparent conductivities calculated from the Bz data;c—shows the apparent conductivity curves of the induced voltage and Bz field over x position;the dashed lines delineated the positions of the abnormal bodies 4~6 in x direction;the blue and red inverted triangles labeling the peaks of the corresponding curves圖21 針對(duì)圖20中感應(yīng)電壓數(shù)據(jù)和Bz數(shù)據(jù)計(jì)算得到的全時(shí)視電導(dǎo)率擬二維剖面Fig.21 The all-time apparent conductivities calculated from the induced voltage data and Bz data shown in Fig.20.
從視電導(dǎo)率隨時(shí)間的變化來看,兩種數(shù)據(jù)均在第一個(gè)采樣時(shí)間點(diǎn)獲取到了較強(qiáng)的異常信號(hào),但是Bz場(chǎng)的視電導(dǎo)率在更早的時(shí)間點(diǎn)回到了背景場(chǎng),意味著其更早的探測(cè)到了異常體下方的半空間,再次說明Bz能夠提前感知到地下電性體。從視電導(dǎo)率值上看,感應(yīng)電壓場(chǎng)獲取到了更大的電導(dǎo)率異常值,異常體反映更強(qiáng)烈,而Bz場(chǎng)的視電導(dǎo)率異常值明顯更小。為更準(zhǔn)確地分析橫向分辨能力, 本文以異常響應(yīng)最強(qiáng)的t=0.01 ms(即第一道)時(shí)刻各測(cè)點(diǎn)的視電導(dǎo)率連成的曲線為基礎(chǔ),若兩個(gè)異常體分布在不同的波谷中,則代表這兩個(gè)異常體能被分辨開來。該時(shí)刻的感應(yīng)電壓場(chǎng)和Bz場(chǎng)對(duì)應(yīng)的視電導(dǎo)率隨測(cè)點(diǎn)的變換曲線如圖21c所示。圖中標(biāo)記出了與第4、5、6號(hào)異常體有關(guān)的曲線波峰位置??梢钥闯?,第4、5、6號(hào)異常體恰好落在感應(yīng)電壓場(chǎng)視電導(dǎo)率曲線不同波谷位置,而僅第5、6號(hào)異常體落在了Bz場(chǎng)視電導(dǎo)率的波谷位置,且從波峰到波谷的變化幅值明顯弱于感應(yīng)電壓場(chǎng)。可見,Bz場(chǎng)的橫向分辨率也差于感應(yīng)電壓場(chǎng)。
三維異常體模型的對(duì)比測(cè)試表明,Bz場(chǎng)被再次證明相比感應(yīng)電壓場(chǎng)能夠提前感知到地下異常體,但是橫向分辨率依然差于感應(yīng)電壓場(chǎng),對(duì)異常體電導(dǎo)率的呈現(xiàn)能力亦明顯弱于感應(yīng)電壓場(chǎng)。
假設(shè)由時(shí)域變換到頻域的傅立葉基函數(shù)為e-iωt,那么在頻率域,感應(yīng)電壓場(chǎng)Vz(ω)和Bz(ω)場(chǎng)有如下關(guān)系
Vz(ω)=iωBz(ω) ,
(13)
式中:i為虛數(shù)單位,ω為角頻率。從式(13)中可以看出,Bz場(chǎng)相對(duì)感應(yīng)電壓場(chǎng)具有更豐富的低頻成分,而高頻成分則更少。由于低頻成分相對(duì)高頻成分具有探測(cè)深度大,而分辨率低、體積效應(yīng)更強(qiáng)的特征,這也就決定了Bz相對(duì)感應(yīng)電壓場(chǎng)在相同的時(shí)間點(diǎn)能夠探測(cè)到更深部的信息,或者說能夠提前感知到地下信息,但是分辨率卻較差的特性。本文的一維和三維模型分析再次證明了這一點(diǎn)。尤其是針對(duì)低導(dǎo)異常體,更強(qiáng)的體積效應(yīng)導(dǎo)致原本就比較小的異常信號(hào)很容易被圍巖較強(qiáng)的信號(hào)所淹沒,導(dǎo)致Bz的低導(dǎo)異常體探測(cè)能力非常差。
本文在反演對(duì)比分析中假設(shè)了感應(yīng)電壓場(chǎng)和Bz場(chǎng)具有相同的最小相對(duì)誤差,但是,實(shí)際上Bz場(chǎng)可能具有更強(qiáng)的抗噪能力。因此,為進(jìn)行更準(zhǔn)確的探測(cè)效果分析,需要在實(shí)際測(cè)量中對(duì)數(shù)據(jù)的質(zhì)量進(jìn)行準(zhǔn)確剖析,這也是下一步工作的目標(biāo)之一。
一維模型、三維模型的感應(yīng)電壓場(chǎng)和Bz場(chǎng)數(shù)據(jù)的視電阻率以及反演結(jié)果對(duì)比分析結(jié)果表明:
1) 在本文采用的儀器精度水平以及測(cè)量方式下,只有針對(duì)電導(dǎo)率大于0.036 S/m的高導(dǎo)異常體,Bz場(chǎng)才可能獲取更長(zhǎng)的有效觀測(cè)時(shí)間。
2) 相比感應(yīng)電壓場(chǎng),Bz場(chǎng)能夠提前感知到地下異常體的信息,但是對(duì)電導(dǎo)率的呈現(xiàn)能力弱于感應(yīng)電壓場(chǎng),縱向分辨率、橫向分辨率均較感應(yīng)電壓場(chǎng)差。
3)考慮到Bz場(chǎng)分辨率不足的問題,在實(shí)際測(cè)量過程中,不建議僅測(cè)量Bz場(chǎng),在條件允許的情況下,同時(shí)測(cè)量感應(yīng)電壓和Bz場(chǎng),充分利用感應(yīng)電壓場(chǎng)強(qiáng)異常響應(yīng)和高分辨率的特征,以及Bz場(chǎng)提前感知地下異常體的特征,達(dá)到最優(yōu)的探測(cè)效果。