一種應(yīng)用于地下腔體成像的新超寬帶成像算法

姬偉杰，劉平，李杰，韓瀟弘毅

（空軍西安飛行學院，西安 710306）

0 引言

地下目標成像[1-20]是指利用一個二維或三維的成像系統(tǒng)產(chǎn)生地下目標的圖像對地下目標回波的解釋說明，其應(yīng)用范圍非常廣泛，在隧道結(jié)構(gòu)檢測[1-2]、路基檢測及路面沉降探測[3-4]、城市地下管道探測[5-6]、巖溶地質(zhì)勘探[7]、堤壩安全探測[8-9]、排雷作業(yè)[10]等領(lǐng)域都有重要的應(yīng)用。

由于近場成像中的輻射脈沖的寬頻帶特性、背景介質(zhì)以及目標的隨機性和復(fù)雜性，使得定量地描述介質(zhì)分布的電磁參數(shù)存在一定的困難，因此對地下目標成像問題的研究一直是一個熱點問題。在地下目標成像方面，K.Demarest等[11]對多層媒質(zhì)中的散射體建立時域有限差分模型進行了分析；方廣有等[12]在國內(nèi)最早對地下三維目標電磁散射特性進行了研究；A.Sullivan等[13]分別用快速多極子法（Fast Multipole Method，F(xiàn)MM），矩量法（Moment Of Method，MOM）和物理光學（Physical Optics，PO）建立了探測地下啞彈的模型，比較了各種方法的精確度，理論成像結(jié)果和實際系統(tǒng)的成像結(jié)果非常吻合；Y.Bo等[14]采用二維 FDTD研究了探地雷達在實際土壤中的應(yīng)用問題；瑞典學者D.Uduwawala等[15]將FDTD法運用于有耗和色散媒質(zhì)中的探地雷達系統(tǒng)；朱亞平等[16]針對脈沖體制穿墻雷達系統(tǒng)中回波信號結(jié)構(gòu)復(fù)雜、信噪比低、信號檢測困難等問題，提出了一種“小波雙譜”聯(lián)合檢測新算法；楊虎等[17]將多區(qū)域時域偽譜（ Multi-domain Pseudo Spectral Time Domain， MPSTD）算法，以及將特征變量與物理邊界（ Characteristic Variables Physical Boundary ，CV-PB）匹配條件相結(jié)合來模擬探地雷達模型，分別對平坦地表、粗糙地表下不同電磁參數(shù)、不同形狀目標的散射特性進行了分析；魯晶津等[18]采用多重網(wǎng)格法研究了地球電磁三維數(shù)值模擬及應(yīng)用；屈樂樂等[19]研究了壓縮感知理論(Compressive Sensing，CS)在頻率步進探地雷達偏移成像中應(yīng)用。

本文針對地下腔體這一特殊目標，提出基于基本振子線電流輻射模型的改進后向投影算法，應(yīng)用地下孔穴的成像分析?？紤]地下介質(zhì)對成像結(jié)果的影響，分析了地下介質(zhì)的物理模型，構(gòu)建了基于FDTD的半空間全波仿真模型和成像模型，推導了電磁波傳播的時延補償計算公式。分析了目標結(jié)構(gòu)對成像結(jié)果的影響，探討了不同結(jié)構(gòu)地下腔體的成像問題，并對多個腔體存在時情況進行成像仿真，給出多個孔穴的三維像。

1 成像模型

1.1 地下介質(zhì)物理模型

從電磁波的傳播分析及文獻 20的分析可知，地下介質(zhì)的物理性質(zhì)參數(shù)(介電參數(shù)ε rε0、磁導率μ和導電率σ)將決定雷達波的傳播和衰減特征，其中εr為相對介電常數(shù)。一般地下介質(zhì)呈現(xiàn)無磁性或弱磁性，即μ=μ0。因此，大地物質(zhì)的介電性質(zhì)對電磁波傳播、散射的影響最大，眾多學者、專家致力于研究并總結(jié)出了一系列典型介質(zhì)的介電性質(zhì)參數(shù)，詳見表 1。對于野外的土壤，其介電常數(shù)與濕度存在密切關(guān)系，本文采用Wang和Schmugge建立的一種四成份模型[20]來計算土壤的等效介電參數(shù)，即將土壤表示為空氣、固態(tài)土壤、束縛水和自由水四種物質(zhì)的介電混合體，計算公式如下：

表1 部分常見介質(zhì)的介電特征參數(shù)表

其中，mv為土壤中總的體積含水量，mi為土壤的臨界體濕度，p為土壤的積孔率，aε為空氣的介電常數(shù)，rε為巖石的介電常數(shù)，wε為工作頻率下純水的介電常數(shù)，且

1.2 地下腔體的成像模型

探地雷達有多種結(jié)構(gòu)，有單發(fā)射天線單接收天線模式、單發(fā)射天線多接收天線模式和多發(fā)射天線多接收天線模式。單發(fā)單收模式下，可以采用輻射波瓣很窄的天線，探測時同時獲取目標角度和時延信息實現(xiàn)定位，也可以采用輻射波瓣很寬的全向天線，通過移動天線，在不同位置測得目標的回波時延來定位，這種定位方式需要掃描或移動，定位較慢，不利于實時探測。在采用多個接收天線的模式下，可以采用輻射波瓣很寬的全向天線，僅利用各個接收天線及發(fā)射天線的位置和各天線獲得的目標的時延信息實現(xiàn)定位，要實現(xiàn)對目標在二維平面上定位，至少需要兩個接收天線。收發(fā)天線可以是共享的，也可以是分置的，共享時可認為是兩個天線位于同一位置的分置的特殊情況。

本文采用圖1所示的收發(fā)天線共享的多發(fā)多收天線結(jié)構(gòu)，實現(xiàn)地下目標的三維成像。地下目標的電磁散射屬于半空間的電磁散射問題，本文采用FDTD算法進行全波電磁仿真，天線采樣平面位于吸收邊界和總場、散射場邊界之間，直接采樣目標的散射場。

圖1 地下腔體成像模型

2 成像算法

2.1 后向投影算法

如圖1所示，位于采樣平面上的天線發(fā)射時域脈沖信號的表達式為p(t)，系統(tǒng)總采樣點數(shù)為L，第l次天線采樣時所處的位置為O(xl,yl,zl)，l=1,2,3…L。將成像區(qū)域離散，x軸向離散數(shù)目為N，y軸向離散的數(shù)目為M，z軸向離散的數(shù)目為P，則成像區(qū)域中共有N×M×P個點。設(shè)成像區(qū)域中離散點坐標可為O'(xn,ym,zp)，其中，n=1,2,3…N，m=1,2,3...M，p=1,2,3…P。且該離散點上反射系數(shù)為σnmp。則第l個采樣點處接收的回波信號為

其中τl,nmp為天線到采樣點的回波延遲，可表述為

其中，c為真空中的光速，Rl,,nmp為第l個采樣位置到成像區(qū)域離散點的距離。

根據(jù)采樣點目標回波與位置，及成像區(qū)域離散點坐標來求取各離散點散射系數(shù)。后向投影成像算法的處理方法是后向投影疊加，可用公式表示為

式中，u(l,τl,nmp)為第l個采樣點在t時刻的采樣數(shù)據(jù)點。

2.2 成像算法的改進及時延校正

針對后向投影算法的理論基礎(chǔ)基于遠區(qū)場而大多數(shù)無載頻時域脈沖成像來說采樣所得到的數(shù)據(jù)是不滿足這一條件的缺陷，提出一種基于近場數(shù)據(jù)的成像模型，即將采樣數(shù)據(jù)等效為線電流，根據(jù)線電流近場輻射場計算公式求得的輻射場作為新的測量數(shù)據(jù)，再根據(jù)后向投影算法的原理反向投影到成像區(qū)域，獲得成像圖像。

2.3 時延校正

在進行地下腔體成像時，地下介質(zhì)對電磁波的折射會引起傳播路徑的變化，電磁波在地下介質(zhì)中傳播速度也會改變，這導致目標回波到達時間的差異，對這個時間差異進行補償是成像算法的關(guān)鍵。如果時延計算不正確，則成像圖中目標位置會有偏差，還會引起散焦，圖像質(zhì)量下降。假設(shè)地下介質(zhì)的相對介電常數(shù)為1，位于的發(fā)射天線發(fā)射電磁波，經(jīng)墻體折射后照射到目標上的點，然后被目標反射后，再經(jīng)過地層的折射，被采樣天線采樣?？紤]地層的影響后的回波延遲可以修正為

由于由發(fā)射天線到目標的路徑分析及由目標到接收天線的分析是相同的，因而本文只對發(fā)射天線到目標的路徑l1，l2進行分析。由折射定律及相關(guān)物理性質(zhì)得到

其中，n0，n1是空氣和介質(zhì)墻中的折射率。假設(shè)介質(zhì)墻是非磁性的，即μr=1，則

令l1，l2為射線傳播路徑，lm為射線傳播路徑l2在x-z平面上的投影，令

則

式（12）同平方，得

對式（9）進行擴展，得到關(guān)于lm的四次方程如下

當εr≠1時，式（14）簡化為

令

式（15）可以表示為

因而式（20）可以歸納為四階方程的求解，采用費拉里求解方法進行求解，可以求得lm，則射線傳播路徑l1，l2可以表示為

3 仿真實驗分析

為了驗證改進成像算法能夠?qū)Φ叵虑惑w進行有效成像，本文對典型地下目標進行成像驗證，對照射面對成像的影響進行成像分析，對不同結(jié)構(gòu)的地下腔體目標進行了成像。

3.1 時延補償?shù)挠行则炞C

算例1 成像模型及時延補償?shù)姆抡骝炞C。

管道位于地下1.5 m處，半徑為0.5 m，材質(zhì)為金屬，管道埋藏于混凝土中，混凝土的相對介電常數(shù)為 6.0，電導率為0.005S/m，成像區(qū)域位于[-1.7,1.7]m×[-1.7,1.7]m范圍內(nèi)，時間步長為Δt=3.973ps，高斯脈沖參數(shù)τ=25Δt，t0=4τ，成像結(jié)果如圖2所示，圖2(a)的是未進行時延補償?shù)某上窠Y(jié)果，圖2(b)是時延補償后的成像結(jié)果。

圖2 單根地下金屬管道的成像結(jié)果

成像結(jié)果表明：改進后向投影算法能對地下管道進行有效成像，可以較為精確地判別出管道在地層中所處的位置，并能基本的辨別出目標照射面的集合特性；混凝土相對介電常數(shù)同空氣的差異造成目標的成像結(jié)果同真實位置存在偏移，但經(jīng)時延補償后的成像結(jié)果同目標的真實位置相一致，說明時延補償?shù)恼_性。

3.2 目標結(jié)構(gòu)對成像的影響

電磁波照射的目標表面(照射面)主要以三種形式存在，如圖3所示：凸面、凹面、平面。因而對不同照射面進行成像分析能更進一步了解目標成像。

圖3 不同照射面的電磁波反射

算例2 照射面結(jié)構(gòu)對成像的影響。

后向投影算法及其改進算法是基于目標反射信號反向投影到成像區(qū)域，因而采樣信號的波前對成像結(jié)果的影響很大，而影響采樣信號波前最終要的因素是照射面的幾何特性，與陰影面的關(guān)系不大。

圖4 不同照射面結(jié)構(gòu)的目標成像結(jié)果

對于照射面相同、陰影面不同的兩種金屬柱體，成像截面區(qū)域位于[-40,40]mm×[-40,40]mm范圍內(nèi)，時間步長為Δt=0.796ps ，高斯脈沖的特征參數(shù)τ=35Δt，t0=4τ，成像結(jié)果如圖 4所示。成像結(jié)果表明：改進后向投影算法能對照射面為凸面、平面的目標進行準確成像，不能對照射面為凹面的目標進行有效成像，驗證了分析結(jié)果。

3.3 巖、土體中埋場的空洞的成像

算例3 三種不同結(jié)構(gòu)空洞的成像。

空洞的地電模型可以描述為：矩形空洞大小為0.5m×1m，位于水平方向的中心位置，距離混凝土表面1m；拱形空洞大小為1m×1.25m，其中拱高為0.25m，位于水平方向的中心位置，距離混凝土表面 0.75m；橢圓形空洞長軸為 1.0mm，短軸為0.5mm，位于水平方向的中心位置，距離混凝土表面0.75m，?？斩吹慕殡娞匦酝諝庀嗤?，混凝土的相對介電常數(shù)為 8.0，電導率為 0.015S/m，成像算法采用BPDC算法，成像區(qū)域位于[-1.7,1.7]m×[-1.7,1.7]m 范圍內(nèi)，時間步長為Δt=2.575ps，高斯脈沖參數(shù)τ= 25Δt，t0=4τ，成像結(jié)果如圖5所示。

從圖5(a)的成像結(jié)果中可以看出，矩形空洞垂直方向的兩個面都成出像，且成像結(jié)果清晰，這是因為電磁波從地面沿垂直方向傳播至矩形空洞時，空洞內(nèi)是空氣，第一個面回波不明顯，到達矩形空洞第二個面時，分界面發(fā)生突變，反射回波信號較大，而當?shù)诙€面的回波反射信號到達第一個面同樣分界面發(fā)生突變，發(fā)射回波信號同樣比較明顯，且矩形空洞垂直方向的兩個面都為平面，根據(jù)算例2中分析的目標結(jié)構(gòu)對成像結(jié)果的影響，矩形空洞這兩個面都能成像，且成像效果較好；從圖5(b)的成像結(jié)果中可以看出，拱形空洞垂直方向的底面成像結(jié)果清晰，而拱形部分沒有正確成像，這是因為電磁波從地面沿垂直方向傳播至矩形空洞時，空洞內(nèi)是空氣，第一個面回波不明顯，到達矩形空洞第二個面時，分界面發(fā)生突變，反射回波信號較大，而當?shù)诙€面的回波反射信號到達第一個面同樣分界面發(fā)生突變，發(fā)射回波信號同樣比較明顯，但是此時照射面是凹面，根據(jù)算例2中分析的目標結(jié)構(gòu)對成像結(jié)果的影響可知，不能對照射面為凹面的目標進行正確成像，由于散射的能量都匯聚于凹面的焦點處，因而在凹面的焦點處成出一個虛像；從圖5(c)的成像結(jié)果中可以看出，橢圓形空洞垂直方向的兩個面成像都沒有正確成像，這是因為電磁波從地面沿垂直方向傳播至橢圓形空洞時，空洞內(nèi)是 空氣，第一個面回波不明顯，到達橢圓形空洞第二個面時，分界面發(fā)生突變，反射回波信號較大，但此時分界面（照射面）)是凹面，而當?shù)诙€面的回波反射信號到達第一個面同樣分界面發(fā)生突變，發(fā)射回波信號同樣比較明顯，此時照射面分界面（照射面）同樣是凹面，根據(jù)算例2中分析的目標結(jié)構(gòu)對成像結(jié)果的影響可知，不能對照射面為凹面的目標進行正確成像，由于散射的能量都匯聚于凹面的焦點處，因而在凹面的焦點處成出一個虛像。從圖像中能基本判定橢圓形空洞的基本位置，但是難以獲得其它關(guān)于空洞的相關(guān)信息。

圖5 不同結(jié)構(gòu)空洞的成像結(jié)果

算例4 空洞的三維成像。

三維空洞的地電模型可以描述為：模型大小為4×4×4m，正方體空洞大小為0.4×0.4×0.4m，正方體空洞中心位于(0.7m，0.7m，0.7m)，距離混凝土表面1.5m，錐臺空洞底面為圓心位于(-0.15m，0m，-1m)、半徑為0.25m的圓，錐臺空洞上表面為圓心位于(0.1m，0m，-1m)、半徑為 0.05m的圓，錐臺高度為0.25m，距離混凝土表面0.9m，空洞介電特性同空氣相同?；炷恋南鄬殡姵?shù)為8.0，電導率為0.015S/m，成像算法采用改進后向投影算法，時間步長為Δt=2.575ps，高斯脈沖參數(shù)τ= 2 5Δt，t0=4τ，成像結(jié)果如圖6所示。

圖6 三維空洞結(jié)構(gòu)及成像結(jié)果

從圖6的成像結(jié)果中可以看出，改進后向投影算法能夠?qū)Φ叵碌娜S空洞進行較為準確的成像，能準確的判明空洞的位置。通過對不同的地下空洞成像，可以更加全面而細致地了解空洞的準確信息，加深對空洞特征的認識，克服二維雷達探測通常只能探測地質(zhì)體的存在否，但難于對探測對象提供更為準確信息的局限性，提高雷達探測圖譜解釋的可靠性和準確性，為三維探地雷達探測與解釋技術(shù)的開展打下基礎(chǔ)，并可根據(jù)缺陷的分布位置及大小評估襯砌內(nèi)空洞病害的危害程度。

4 結(jié)論

本文提出了一種基于基本振子線電流輻射模型的改進后向投影算法，并應(yīng)用地下孔穴的成像分析?？紤]地下介質(zhì)對成像結(jié)果的影響，分析了地下介質(zhì)的物理模型，構(gòu)建了基于FDTD的半空間全波仿真模型和成像模型，推導了電磁波傳播的時延補償計算公式。分析了目標結(jié)構(gòu)對成像結(jié)果的影響，探討了不同結(jié)構(gòu)地下腔體的成像問題，并能有效地對多個腔體目標進行成像，為堤壩、隧道隱患排查等工程應(yīng)用提供了準確信息，可即時的對工程質(zhì)量進行監(jiān)控，大大減少了工程隱患。

[1]晉良念,歐陽繕,謝躍雷等.基于穩(wěn)健波束形成的超寬帶穿墻成像方法[J].系統(tǒng)工程與電子技術(shù),2011,33(1):208～212.

[2]舒志樂.隧道襯砌內(nèi)空洞探地雷達探測正反演研究[D].重慶大學,2010

[3]覃建波,鄧世坤,李沫.探地雷達在隧道和涵洞工程檢測中的應(yīng)用研究[J].煤田地質(zhì)與勘探,2004,(05):55～58.

[4]葉家瑋,吳鵬,鄭國梁.基于路面雷達的路面結(jié)構(gòu)缺陷檢測方法[J].華南理工大學學報(自然科學版),2004,(09):82～85.

[5]郭秀軍,韓宇,孟慶生等.鐵路路基病害無損檢測車載探地雷達系統(tǒng)研制及應(yīng)用[J].中國鐵道科學,2006,(05):139～144.

[6]趙永輝,謝雄耀,王承.地下管線雷達探測圖像處理及解釋系統(tǒng)[J].同濟大學學報(自然科學版) ,2005,(09):1254～1258.

[7]戰(zhàn)玉寶,張利民,尤春安.探地雷達探測地下管線的研究[J].巖土力學,2004,(S1):133～136.

[8]韋宏鵠,楊順安,劉昌輝.探地雷達在巖土工程應(yīng)用中的進展[J].地質(zhì)科技情報,2005,(S1):133～136.

[9]王國群,何開勝.堤防動物洞穴的探地雷達探測研究[J].巖土力學,2006,(05):838～841.

[10]方廣有,佐藤源之.頻率步進探地雷達及其在地雷探測中的應(yīng)用[J].電子學報,2005,(03):436～439.

[11]Demarest K,Plumb R,Zhubo H.FDTD modeling of scatters in stratified media[J].Antennas and Propagation,IEEE Transactions on,1995,43(10):1164～1168.

[12]方廣有,張忠治,汪文秉.地下三維目標電磁散射特性研究[J].微波學報,1997,(01):8～14.

[13]Sullivan A,Damarla R,Geng N,et al.Ultrawide-band synthetic aperture radar for detection of unexploded ordnance:modeling and measurements[J].Antennas and Propagation,IEEE Transactions on,2000,48(9):1306～1315.

[14]Bo Y,Rappaport C.Response of realistic soil for GPR applications with 2-D FDTD [J].Geoscience and Remote Sensing,IEEE Transactions on,2001,39(6):1198～1205.

[15]Uduwawala D,Gunawardena A.A Fully Three- Dimensional Simulation of a Ground-Penetrating Radar over Lossy and Dispersive Grounds[C].Industrial and Information Systems,First International Conference on.2006.143～146.

[16]朱亞平,沈庭芝,王衛(wèi)江等.穿墻雷達系統(tǒng)中信號檢測的新算法[J].北京理工大學學報,2005,(08):734～738.

[17]楊虎,姜永金,毛鈞杰.平坦/粗糙地表下目標散射特性的MPSTD算法分析[J].微波學報,2007,(06):1～6.

[18]魯晶津.地球電磁三維數(shù)值模擬的多重網(wǎng)格法及應(yīng)用研究[D].中國科學技術(shù)大學,2010.

[19]屈樂樂，方廣有，楊天虹.壓縮感知理論在頻率步進探地雷達偏移成像中的應(yīng)用[J].電子與信息學報,2011,33(1):22～26.

[20]J.R.Wang,T.J.Schmugge.An Empirical Model for the Complex Dielectric Permittivity of Soils as a Function of Water Content [J].IEEE Trans.On Geoscience and Remote Sensing,1980,Vol.18,pp:288-295.