基于GPR數值模擬正演的地下介質分析方法

宋正淳，邱哲理，張雪冰，孫 超，胡 昕，覃育祥

吉林建筑大學 測繪與勘查工程學院,長春 130118

0 引言

探地雷達[1]探測地下介質多為0.5 m～5 m區間內,具有無損探測、探測速度快且探測過程連續、分辨率高、操作方便靈活、探測費用低等優點.因其優點故目前廣泛用于考古、礦產勘查、巖土工程勘探、工程質檢、建筑結構檢測以及軍事勘查等眾多領域中.探地雷達其工作原理是利用天線激發高頻脈沖電磁波傳至地下有耗介質中,地下各層有耗介質各不相同且介質具有不均勻性、各向異性、強衰減性，故接收的電磁波信號波動也不同.利用其電磁信號的差異性,可將一道道回波信號接收的振幅與時間、頻率等參數來分析和推斷被探測土體的地下構造并數據處理繪制成圖像.而本文研究的gprMax[2]軟件則是將電磁波信號模擬仿真,仿真信號波動數據處理組成B-scan圖[3],屆時可通過圖形研究地下不同介質的位置及形態.其方法是基于時域有限差分(Finite difference time domain,縮寫為FDTD)[4]的麥克斯韋[5]方程求解電磁信號強度.FDTD方法是將仿真區域拆分成一個個空間網格,用數值有限差分方法逼近麥克斯韋方程,并將時間離散化,再根據已知的初始邊界條件和介質,通過數據計算處理,將一段時間內的反射信號數據進行整合,即可將仿真模型的二維或三維的數據可視化.介電常數存在差異的兩個界面會反射電磁波,觀察生成的B-scan圖的電磁波現象,通過兩組實驗對比觀察其極化[6]方向和相位[7]變化的規律.

1 探地雷達數據模擬

1.1 探地雷達數據模擬基本原理

從宏觀上來說,所有的電磁現象均可由麥克斯韋方程解釋.以下均是一階偏微分方程,用來表示基本電磁場量與其源的關系:

×E=-?B?t

(1)

(2)

(3)

(4)

式中，E是電場強度,V/m;H是磁場強度,V/m;B是磁通密度,Wb/m2;D是電位移,C/m;J是電流密度,A/m2;q是電荷密度,C/m3;t是時間,s.除上述式子外,其媒質本構關系式是:

D=εE，B=μH，J=σE

(5)

式中，ε是相對介電常數;μ是磁導率,H/m;σ是電導率,S/m.基于麥克斯韋方程組下,根據模型的位置及時間的條件下求解上述方程式,得到探地雷達模擬數據.將模擬的電磁信號數據分析處理即可得到A-scan,B-scan,C-scan圖像,B-scan圖像即是本文主要分析研究的重點,B-scan圖像可將模擬地質的地下模型剖切面用二維表示出來.

1.2 電磁波在界面的反射系數與相位變化

電磁波被探地雷達的天線激發信號射入介質界面上會發生反射和透射現象,其原理類似于光入射到水面和玻璃面.在GPR中利用反射器的傳統技術可以將速度走時分析與CMP數據采集方法相結合,但確定雷達波的速度則需用反射波的傳播時間來確定.其方法存在缺陷,電磁波進入介面后的回波信號深度有局限性,不能確定最低反射面以下介質的速度.電磁波從天線中釋放出其電磁信號可能會產生極化,極化類型由天線的設計確定.天線的設計可以是線性的,圓形的或橢圓形的.絕大多數的探地雷達發出的波是線性的極化,線性的極化又可分為垂直極化和平行極化.極化類型決定電磁波如何反射,電磁波在界面的反射系數可定義為反射波的振幅(或能量)Er與入射波的振幅Ei之比.需要注意的是,垂直極化和平行極化的反射系數有所不同,實際探地雷達數據采集多以垂直極化為主,垂直極化的反射系數可表示為R⊥[8]:

R⊥=ErEi=ε1cosθ1-ε2-ε1sin2θ1ε1cosθ1+ε2-ε1sin2θ1

(6)

式中,ε1,ε2分別對應第一和第二層介質的相對介電常數;θ1是入射角度.

電磁波入射到不同介質時,一部分會發生反射現象,另一部分電磁波則會發生透射現象透過分界面.電磁波在界面的透射系數可定義為Et是透射波的振幅(或能量)與Ei是入射波的振幅之比.

T⊥=EtEi=2ε1cosθ1ε1cosθ1+ε2-ε1sin2θ1

(7)

式中,T⊥是垂直極化的透射系數.

理想完美電導體(Perfect Electric Conductor,PEC)不發生透射現象,是因為PEC導體的相對介電常數無窮大,當第二層介質為PEC時,上式中分母無窮大,分子是除PEC材料以外的其他非導體任意材料,透射系數T⊥趨近于零.

反射系數和透射系數之間滿足公式:

1+R⊥=T⊥

(8)

電磁波的頻率由發射的波源決定,探地雷達中多選擇Richer型波源,頻率不變.電磁波從天線釋放出信號波,信號波在同種介質中傳播波長不發生變化.當電磁波從介質進入其他介質中時波長會減小,振幅也會因阻尼效應減小;當電磁波靠近PEC時,電磁波不會進入PEC中而是會將電磁波反射.共享相同頻率的多個采樣道電磁波相互作用產生相位變化,相位變化是由兩個或多個采樣道電磁波波形的振幅波峰和波谷位置關系確定,通過距離、時間、度數測量相位.介面中某一點的波阻抗公式為:

η=EH=(με)0.5

(9)

式中,η是電場與磁場的復振幅之比定義為波阻抗,Ω.

2 模型建模與結果分析

因本文主要討論的是探地雷達的高頻脈沖電磁波入射到地下的電磁波相位及極化方向,而設置兩層介質界面后第二層界面的入射波是來自于第一層界面的透射波,其信號波動不太明顯,故擬構建兩組濕沙土樣為背景介質的模型(如圖1,圖2所示),并分別設置2組不同尺度、位置和深度的金屬和空洞異常體.

圖1 小尺度圓柱異常體正演模型,圓柱直徑均為75 mmFig.1 The forward model of small-scale cylindrical anomalies at the cylinder diameter of 75 mm

圖2 大尺度圓柱異常體正演模型,圓柱直徑均為300 mmFig.2 The forward model of large-scale cylindrical anomalies at the cylinder diameter of 300 mm

2組模型均由表層0.03 m的空氣和0.6 m深的濕沙構成;異常體均為圓柱形異常體,左側2個異常體采用PEC金屬電導體,右側2個異常體均為空洞.模型一中異常體直徑均為75 mm,深度依次為400 mm,350 mm,400 mm,350 mm;模型二中異常體直徑均為300 mm,深度依次為350 mm,250 mm,350 mm,250 mm.在gprMax正演模擬中,選擇激發源均為主頻為900 MHz的Ricker子波,激發源與接受源間距均為50 mm.

不同深度與直徑的兩組管線模型的gprMax正演結果如圖3,圖4所示.

圖3 小尺度圓柱異常體正演剖面,圓柱直徑均為75 mmFig.3 The forward section of small-scale cylindrical anomaly at the diameter of 75 mm

圖4 大尺寸圓柱異常體正演剖面,圓柱直徑均為300 mmFig.4 The forward section of large-scale cylindrical anomaly at the diameter of 300 mm

電磁波傳播到存在介電常數差異的界面時會在相應位置會呈現波形,紅色與藍色分別代表正向和負向振幅,圖中顏色越深代表其電場強度越強.正演模型的關鍵過程如下(由于兩組原理類似故僅對圖4過程進行描述),其中左側PEC異常體上方的電磁波傳播大致可描述為2個過程.

過程1 以圖4為例.探地雷達從天線發射的高頻電磁脈沖(波形約可用Ricker子波表示),其波形為一組由正向-負向-正向-負向振幅組成的子波(圖中由紅-深藍-深紅-淺藍表示,見2 ns處直達波),通過較薄的空氣層傳入濕沙介質中.電磁波在濕沙介質中傳播至PEC金屬異常體表面,在剖面約時間深度6 ns處可見明顯的雙曲線(參見紅色箭頭),其顏色依次是淺藍-深紅-深藍-淺紅表示.通過對比直達波可知極性發生了反轉.

過程2 電磁波經PEC異常體時反射到濕沙界面和空氣界面之間,一部分電磁波透射到空氣介層,另一部分又反射到PEC異常體表面處(多次波).因為一部分電磁波透射到了空氣介質中,圖4中黑色箭頭處電磁波產生衰減,但其振幅依次為微紅-淺藍-淺紅-微藍,隨后再次傳播的多次波因振幅衰減在圖中已無法識別.圖2中左二的藍色圓柱完美電導體由于位置距離空氣界面層較遠,所以圖中可識別的電磁波只有一次.上述兩過程中的反射系數R⊥,透射系數T⊥以及反射與透射電場強度之比Er∶Et參見表1.

對于右側空洞異常體,電磁波傳播可大致描述為3個過程：

過程1 電磁波由空氣界面透過到濕沙介面之后運動到左測空心圓柱異常體頂部,有一部分透射到了空洞空氣介質管內,另一部分則反射到了濕沙層頂部,回波被接收天線接收,見紅色箭頭處雙曲線同相軸,按接收時間順序其電場強度由淺紅-藍-紅-淺藍表示.

過程2 之前一部分電磁波透射到空洞內部后,電磁波運動到空洞底部后一部分電磁波經反射,回波被接收天線接收,見藍色箭頭位置,其顏色變化依次是微紅-淺藍-淺紅-微藍.

過程3 當電磁波傳播到濕沙與表層空氣界面時,一部分經反射又傳播到空洞頂部,部分電磁波經反射后又回到地表被接收天線接收,在圖4黑色箭頭處產生電磁波,電場強度顏色變化依次是微紅-微藍-微紅-微藍.左二與左四情況相同均因埋深較深且沒有進行衰減補償處理,所以后續多次波不明顯.

通過圖片對比發現電磁波在濕沙界與PEC界面顏色不發生改變,正負極振幅不發生改變.電磁波由濕沙界面到圓柱空心管或由空心管到濕沙界面時,顏色發生改變,正負極振幅發生相位變化.其原因是完美電導體沒有介電常數,不發生透射,只有反射,極化方向不發生改變.空心管則是由高介電常數向低介電常數透射電磁波,極化方向發生變化,磁場強度會由衰減.反射系數與透射系數的模型正演關系見表1,表2.

表1 金屬異常體反射系數與透射系數的模型正演Table 1 Forward model of reflection coefficient and transmission coefficient of metal abnormal bodies

表2 空心異常體反射系數與透射系數的模型正演Table 2 Forward model of reflection coefficient and transmission coefficient of hollow abnormal bodies

3 結語

本文結合gprMax進行探地雷達正演模擬,探討電磁波在地下介質中的運動路徑和相位變化規律,通過gprMax對不同材質、不同深度、不同大小的地下異常體進行正演模擬并對比,就此算出關鍵反射界面的反射系數和投射系數.通過已有實驗數據的介電常數進行仿真模擬地下金屬異常和空洞,對比驗證了其結果與經典電磁學中反射系數、投射系數以及相位變化的計算結果一致,證實了結合gprMax正演模擬的地下目標體解釋的可行性.