彭凌星， 何現啟， 朱自強

(1.湖南省交通規(guī)劃勘察設計院有限公司， 長沙 410008；2.中南大學地球科學與信息物理學院， 長沙 410083)

地質雷達最初在物探上的應用是用來探測冰川的厚度,后來逐漸拓展為探測地下地質結構的物性特征。1970年后,結合地質雷達和地質鉆孔這一現有通道對地下深處結構的地質信息進行探測逐漸形成了鉆孔地質雷達系統(tǒng)[1-2]。針對地質雷達的研究很多,地質雷達也被應用于多個領域[3-8]。近幾年來,吉林大學的王飛等[9]、中國地質調查局成都地質調查中心的李華等[10-11]、中國科學院武漢沿途所的鐘聲等[12-13]、中南大學的Zhu等[14]也先后對鉆孔雷達進行了研究,從正演模擬到結合實際探測闡述了鉆孔雷達掃描中異常體的圖像特征。而國外的學者也于近幾年來開展了較多對鉆孔雷達的研究,如Grégoire等[15]、James等[16]、Niklas等[17]、Fernández等[18]、Majken等[19]、Dafflon等[20]、Caroline等[21]。這些學者分別采用層析成像、全波場反演、貝葉斯隨機反演等方法對跨孔雷達數據進行了反演成像研究。

其中層析成像法應用最為廣泛,但這些層析成像反演計算中均未考慮到鉆孔雷達天線對反演的影響,例如天線的長度。因此需要將天線的影響加入到數值模擬當中以便提高鉆孔雷達層析成像的精度和準確度?，F根據天線長度對電壓的影響進行公示推導,并對其進行數值模擬驗證。

1 天線收發(fā)信號方程

鉆孔地質雷達發(fā)射和接收天線基本上都為饋電點位于中央部位的柱狀偶極子天線,天線的示意圖如圖1所示。以鉆孔地質雷達的發(fā)射天線為例說明,假定天線的長度坐標為-l～l,Z0為發(fā)射天線的波阻抗,Vg(t)為發(fā)射天線發(fā)出的脈沖電壓。假設天線的收發(fā)距離較遠,且位于天線收發(fā)距離之間的介質為各向同性介質,則接收天線所接到的信號為發(fā)射天線信號遠場的表達式。由于天線收發(fā)距離較遠,且收發(fā)天線位置為互相平行,因此在文中只考慮TE模式(電矢量與傳播方向垂直)的電磁波傳播。

在遠場,可以通過式(1)對發(fā)射的電磁場進行表達[22]:

(1)

式(1)中：r為天線的收發(fā)距離；θ為信號傳播路徑與發(fā)射天線之間的夾角；ω為角頻率；k為波數；ξ為電磁波的波阻抗；I(0,ω)則是發(fā)射天線饋電點處電流(頻率域)。當z=0時；天線有效高度在傳播路徑與發(fā)射天線夾角θ方向上的分量h可以表達為[23]

(2)

圖1 鉆孔雷達發(fā)射天線示意圖Fig.1 Diagram of borehole radar transmitting antenna

而k和ξ可以由相對介電常數ε,相對磁導率μ和相對電導率σ按照式(3)、式(4)計算得到:

(3)

(4)

將式(2)代入式(1)中,可以得到:

(5)

式(5)括號內的表達式就是具有電流為I(z,ω)的無限小的偶極子天線在遠場的電場表達式。式(5)表明鉆孔地質雷達中發(fā)射天線所引發(fā)的遠場電場強度可以通過對無限小的偶極子發(fā)射天線引發(fā)的電場強度進行積分計算得到。

Sato在1991年通過式(6)表述了發(fā)射天線饋電點處的電流(頻率域)與發(fā)射天線電壓脈沖(頻率域)之間的關系[19]:

(6)

式(6)中：Vg(ω)為發(fā)射天線的電壓脈沖(頻率域)；Zin(ω)是發(fā)射天線輸入阻抗。將系數A(z,ω)定義為

(7)

則天線中任意部位的電流強度表達式可以寫為

I(z,ω)=A(z,ω)Vg(ω)

(8)

轉換成時間域內的表達式，則為

I(z,t)=A(z,t)*Vg(t)

(9)

此時,則可以通過計算發(fā)射天線時間域電壓脈沖和系數A(z,ω)的卷積得到發(fā)射天線任何一處的電流時間域表達式。

Sato在1991年提出了遠場的接收電壓頻率域表達式為[19]:

(10)

將式(2)、式(7)代入式(10)中,可以得到電壓沿天線長度的積分表達式:

(11)

電場強度(頻率域)在z方向上的分量:

Ez(z,ω)=Eθ(r,θ,ω)sinθeikzcosθ

(12)

結合式(12)的時間域變換,可以得到式(11)時間域的表達式:

(13)

式(13)表明接收天線的電壓變化可以通過對無限小的接收天線的電壓進行積分計算得到,這與發(fā)射天線引發(fā)的遠場電場強度表達式頗為類似。

根據上述公式的推導結果,得知:發(fā)射天線激發(fā)的電磁波脈沖在遠場所引發(fā)的電場強度分量能夠激發(fā)接收天線產生沿著任意一點往天線兩端傳播的電流。系數A(z,t)則表明發(fā)射天線的電壓脈沖與任一點的電流之間的關系,同時也表明接收天線中任一點的電壓脈沖與發(fā)射天線激發(fā)電流之間的關系。因此,通過Ez(z,t)和A(z,t)的卷積疊加得到發(fā)射天線中電流強度的總和,從而進一步得到接收天線的總電壓。

根據發(fā)射天線激發(fā)電壓與系數A(z,t)的卷積得到發(fā)射天線中任意點電流的分布,從而計算出電場強度；然后根據計算出的接收天線電場強度與系數A(z,t)的卷積得到該點的接收天線電流強度,通過接收天線的長度對天線電場強度進行積分得到接收天線總電流的強度,然后將電流與電阻相乘即可得到該時刻接收天線的電壓。

2 天線收發(fā)信號數值模擬計算

對于偶極子天線,電流與阻抗的關系可以通過式(14)表述:

(14)

假定天線特征阻抗值為Zc,則:

Zin=-iZccot(kl)

(15)

將式(6)和式(15)代入式(14)中,則:

(16)

天線的波數為k,可以通過式(17)進行計算:

(17)

式(17)中：c為真空中電磁波的速度；ε為相對介電常數。

為了對上述公式進行驗證,通過如下例子進行說明。假設天線相對介電常數為9,天線的特征阻抗值為150 Ω,發(fā)射天線的電壓脈沖為高斯脈沖。

高斯脈沖的時間域表達式為

(18)

假定取τ為1 ns,t0為5 ns,發(fā)射天線的電壓時域計算結果如圖2所示。

圖2 激發(fā)電壓的時域波形Fig.2 Time-domain waveform at excitation voltage

根據計算出的電壓,可以得到偶極子天線饋電點處的電流I(0,ω)。假定天線長度為2 m,則式(16)變化為

(19)

計算結果見圖3(取其實數部分成圖)。圖3所示為在饋電點處各個時刻電流的大小。沿天線各處在每個時刻時的電流變化可以根據式(14)進行計算,結果如圖4所示。

圖3 饋電點電流時域波形Fig.3 Time domain waveform of current at the feed point

圖4 天線中任意點電流分布Fig.4 Distribution of antenna current at arbitrary points

由圖4可知,t=t0時,z=0處為電流的最大值,隨后電流強度都在逐漸減小,慢慢變化至0,這就是發(fā)射天線所發(fā)射的信號。由此根據式(1)可以算出遠場的電場強度。

假設發(fā)射天線與接收天線位于同一水平面,即是θ=90°；同時收發(fā)天線之間的距離為6 m,則式(1)可以寫為

(20)

結合式(15)、式(20)可改寫為

(21)

進而對方程求解,得到:

(22)

(23)

計算的E結果如圖5所示。

圖5 觀測點的電場值Fig.5 Value of electric field at observation points

根據上述假設的參數,并將式(7)、式(14)、式(15)代入式(13),可以得到:

(24)

與式(21)類似,可以將式(24)寫為

(25)

(26)

(27)

根據式(27)可以計算出V隨著E的變化,結果如圖6所示。

圖6 接收的電壓波形Fig.6 Received voltage waveform

由圖2～圖6可知,發(fā)射天線激發(fā)的電壓脈沖會產生和電壓一樣波形的電流,從而產生電場,進而被接收天線通過該處觀察到的電場強度轉換成電壓,從而被接收天線接收并顯示。

3 天線數值模擬實例

為了便于計算,模型中的參數設置相對較為簡單。模型的示意圖如圖7(a)、圖7(b)所示。模型的大小為11 m×6 m,天線的頻率為100 MHz,相對介電常數為25,相對磁導率為1,電導率為1 mS/m,天線的長度為1 m。此次數值模擬的天線收發(fā)方式為單發(fā)多收:發(fā)射天線的中心點位于深度5.5 m處,接收天線中心點從0.5 m處每隔0.25 m移動一次直至到達10.5 m處。

發(fā)射天線的激勵脈沖采用微分高斯脈沖,取t0=0.8τ,發(fā)射信號如圖7(c)中所示,接收信號如圖8所示,此時接收信號為一雙曲線狀弧形。

由于接收點過多,數據較密集,從接收的模擬信號中難以看出接收到的信號具體波形。為此,取其中收發(fā)角度最大,也即是接收天線中心點位于深度0.5 m處的信號,如圖9所示。

圖7 天線數值模擬模型Fig.7 Antenna numerical simulation model

圖8 接收的數值模擬信號Fig.8 Received numerical simulation signal

圖9 接收天線位于0.5 m處接收的數值模擬信號Fig.9 The receiving antenna receives the numerical simulation signal at 0.5 m

由圖7(c)、圖9可知,對于天線信號的收發(fā)數值模擬取得了較好的效果。接收信號波形與發(fā)射信號波形基本一致。可以將天線的長度的影響加入到鉆孔地質雷達跨孔數值模擬中。

通過對模型的反演驗證加入天線長度影響后的數值模擬結果。模型中異常體的形狀都較為規(guī)則,為矩形,異常體的磁導率和電導率與背景場相同。背景為介電常數為25的各向均勻介質,異常體的介電常數分別為29(藍色)、21(紅色)；所有介質的電導率均設為1 mS/m；磁導率為真空中的磁導率；天線的總長度為1 m。發(fā)射天線位于圖的左邊,接收天線位于圖的右邊。此次模擬的天線移動方式為多發(fā)多收:發(fā)射天線中心點從0.5 m處每間隔0.25 m移動一次,直到10.5 m處,接收天線也同樣據此方式進行移動。電磁波在天線中傳播的速度為0.12 m/ns。圖10所示為模型的示意圖。

圖10 模型示意圖Fig.10 The model sketch

圖11所示為根據兩種方法得出的合成數據的速度層析成像圖。其中圖11(a)中正演數值模擬采用的是點狀天線；圖11(b)在正演數字模擬中加入了天線的長度。兩圖相比較而言,圖11(b)更為清晰,能夠分辨出模型中所有異常體的位置。采用點狀天線得到的數值模擬結果進行層析成像,可以分辨部分異常體的位置,但并沒有識別出模型中所有異常體。

圖11 速度層析成像結果Fig.11 The tomography results of velocity

4 結論

(1)在考慮到天線長度影響的前提下,根據對鉆孔雷達天線收發(fā)信號方程的推導可以計算出發(fā)射天線和接收天線中電流及電壓的大小。

(2)對于發(fā)射的信號,首先根據激發(fā)的電壓與的卷積得到天線中電流的分布,從而得到計算區(qū)域的電場；對于接收信號接收天線處的電場與卷積得到該無限小的天線電流強度,然后根據天線長度對其進行積分得到天線總電流的大小,與接收天線的負載電阻相乘即可得到該時刻接收到的電壓大小。

(3)根據給定的參數對收發(fā)天線的信號進行理論計算并成圖,計算結果顯示接收信號相比發(fā)射信號并無較大的差異；同時對天線收發(fā)信號進行數值模擬,模擬結果也理論計算結果較為吻合,證明了對天線收發(fā)信號計算方法的準確性。

(4)通過對模型合成數據的層析成像,證明鉆孔雷達數值模擬中加入天線長度的影響后,將會提高鉆孔雷達數值模擬的精度。