基于折射波靜校正的未知參數(shù)穿墻成像合成延遲時(shí)算法*

李家強(qiáng)**1a,，劉 然，徐才秀2，陳金立1a,，陳焱博，朱艷萍

(1.南京信息工程大學(xué) a.氣象災(zāi)害預(yù)報(bào)預(yù)警與評(píng)估協(xié)同創(chuàng)新中心;b.電子與信息工程學(xué)院，南京210044；2.江蘇北斗衛(wèi)星導(dǎo)航檢測(cè)中心公司，南京210044)

0 引 言

穿墻成像雷達(dá)(Though-the-Wall Radar Imaging，TWRI)作為一種無(wú)損檢測(cè)技術(shù)，能夠探測(cè)視覺(jué)非透明障礙物后的目標(biāo)，在反恐、救災(zāi)等方面發(fā)揮了重要作用，廣泛應(yīng)用于軍事、民用等領(lǐng)域[1]。

近年來(lái)，在穿墻成像雷達(dá)方面有許多研究致力于對(duì)未知參數(shù)的單層墻體后目標(biāo)進(jìn)行成像。其中文獻(xiàn)[2]討論了相對(duì)介電常數(shù)和墻壁厚度對(duì)目標(biāo)點(diǎn)成像的影響，通過(guò)設(shè)置兩組與墻壁不同距離的天線陣列，對(duì)預(yù)先假設(shè)的幾組墻壁厚度和介電常數(shù)的場(chǎng)景進(jìn)行成像，兩組陣列所成圖像都會(huì)隨著墻體厚度和相對(duì)介電常數(shù)的誤差而產(chǎn)生偏移，其中兩者交叉的點(diǎn)即為真實(shí)目標(biāo)位置。由于該方法需要對(duì)設(shè)置的不同墻體參數(shù)進(jìn)行成像，此過(guò)程計(jì)算復(fù)雜，且墻體參數(shù)的估計(jì)結(jié)果誤差較大。另外，支持向量機(jī)方法(Support Vector Machine,SVM)[3]根據(jù)散射信號(hào)提取出信號(hào)作為輸入、墻體參數(shù)作為輸出得到的函數(shù)關(guān)系，對(duì)墻后目標(biāo)進(jìn)行成像。此方法所估計(jì)出的墻體參數(shù)值易受目標(biāo)位置、大小以及形狀的變化影響。文獻(xiàn)[4]運(yùn)用兩列垂直天線陣列對(duì)穿墻場(chǎng)景進(jìn)行成像，利用中心關(guān)聯(lián)目標(biāo)的位置關(guān)系對(duì)穿墻場(chǎng)景進(jìn)行重聚焦，但該方法易受目標(biāo)的種類和位置的影響。對(duì)于多層墻體成像的研究，文獻(xiàn)[5]主要針對(duì)成像聚焦時(shí)延的計(jì)算問(wèn)題，其每層墻體的參數(shù)都是假設(shè)已知的；文獻(xiàn)[6]通過(guò)構(gòu)造濾波器去除回波域中墻體產(chǎn)生的影響以獲得清晰的成像結(jié)果，該方法每次更新值時(shí)都要進(jìn)行重新計(jì)算并成像，故計(jì)算效率低，且成像效果易受實(shí)際墻體參數(shù)的影響;文獻(xiàn)[7]通過(guò)對(duì)墻體厚度和相對(duì)介電常數(shù)的二維搜索，得到延時(shí)補(bǔ)償值的取值范圍，進(jìn)而對(duì)該取值范圍進(jìn)行遍歷，找到最佳聚焦時(shí)延近似補(bǔ)償由于墻體存在產(chǎn)生的異常傳播時(shí)延。在對(duì)墻體參數(shù)進(jìn)行搜索時(shí)，該方法不能精確地計(jì)算聚焦時(shí)延的取值范圍，導(dǎo)致在遍歷時(shí)迭代次數(shù)過(guò)多。

在實(shí)際應(yīng)用中，墻體一般并非單層均勻介質(zhì)，且墻體參數(shù)未知。由于電磁波在墻體內(nèi)的傳播速度比在空氣中傳播速度慢，導(dǎo)致其回波時(shí)延增長(zhǎng)，成像過(guò)程中目標(biāo)點(diǎn)位置會(huì)產(chǎn)生偏移或出現(xiàn)散焦現(xiàn)象，這給成像工作帶來(lái)巨大挑戰(zhàn)。聚焦時(shí)延的計(jì)算是影響穿墻成像質(zhì)量至關(guān)重要的因素，因此，解決未知參數(shù)多層墻體的時(shí)延計(jì)算問(wèn)題是目前穿墻成像面臨的挑戰(zhàn)之一。為解決未知參數(shù)的多層墻體介質(zhì)，如常見(jiàn)的雙層墻體成像問(wèn)題，本文提出了基于折射波靜校正的合成延遲時(shí)法對(duì)墻后目標(biāo)進(jìn)行精準(zhǔn)地聚焦成像，通過(guò)計(jì)算出每一個(gè)通道內(nèi)的聚焦時(shí)延，進(jìn)而依次對(duì)每個(gè)通道中由于墻體存在產(chǎn)生的異常傳播時(shí)延進(jìn)行一一校正。

1 穿墻成像雷達(dá)場(chǎng)景模型

本文運(yùn)用多發(fā)多收的天線陣列探測(cè)墻后目標(biāo)。如圖1所示，墻體由兩層不同介質(zhì)組成，內(nèi)層墻體為石膏材料，外層墻體為混凝土材料，兩者厚度分別為d1和d2，相對(duì)介電常數(shù)分別為ε1和ε2。天線陣列緊貼墻體依次擺放，發(fā)射信號(hào)為超寬帶窄脈沖信號(hào)[7]，信號(hào)經(jīng)過(guò)墻體的透射和反射至目標(biāo)點(diǎn)處，再由接收天線所接收。

圖1 電磁波發(fā)射接收傳播時(shí)延

本文運(yùn)用基于時(shí)域有限差分法(Finite Difference Time Domain,FDTD)的GprMax2D/3D[8]軟件對(duì)穿墻場(chǎng)景進(jìn)行仿真并獲取回波數(shù)據(jù)。設(shè)發(fā)射信號(hào)為s(t)為負(fù)一階高斯脈沖信號(hào)，其表達(dá)式為

(1)

設(shè)第m(m=1,2,…,M)根發(fā)射天線坐標(biāo)為(xtm,0)，第n(n=1,2,…,N)根接收天線坐標(biāo)為(xrn,0)，目標(biāo)點(diǎn)q坐標(biāo)為(xq,yq)，在不考慮噪聲的影響下，則接收天線接收的回波信號(hào)為aqs(t-τmqn)，其中aq為目標(biāo)的回波幅值，τmqn為電磁波從第m根發(fā)射天線到第n根接收天線的傳播時(shí)延。由圖1可知,

(2)

式中：c為電磁波在空氣中的傳播速度,τtm1與τtm2分別為電磁波在外層和內(nèi)層墻體中發(fā)射天線延遲時(shí)，τrn1與τrn2分別為電磁波在兩層墻體中接收天線的延遲時(shí),l1與l2為電磁波在空氣中的傳播路徑。由后向投影(Back Projection)算法知，成像場(chǎng)景區(qū)域某一像素點(diǎn)p(xp,yp)的像素值為

(3)

式中：τmpn為成像區(qū)域任意點(diǎn)p從第m根發(fā)射天線到第n根接收天線的傳播時(shí)延。

若兩層墻體參數(shù)已知，則可通過(guò)折射定理計(jì)算折射點(diǎn)坐標(biāo)，進(jìn)而計(jì)算出電磁波在墻體內(nèi)的傳播時(shí)延，完成對(duì)墻后目標(biāo)的成像與精準(zhǔn)定位。而在實(shí)際應(yīng)用中，墻體厚度以及相對(duì)介電常數(shù)這些先驗(yàn)信息通常是未知的，這給后續(xù)的成像工作帶來(lái)影響，導(dǎo)致目標(biāo)位置偏移和成像散焦。本文提出基于折射波靜校正的合成延遲時(shí)法，通過(guò)計(jì)算每一個(gè)通道的聚焦時(shí)延，使得目標(biāo)更接近真實(shí)位置，圖像聚焦程度更高。

2 基于折射波靜校正的合成延遲時(shí)算法

2.1 折射波靜校正原理

靜校正技術(shù)屬于地震勘探領(lǐng)域，是一項(xiàng)處理地震資料的技術(shù)。假設(shè)存在兩層不同介質(zhì)，電磁波在下層介質(zhì)中的傳播速度大于在上層介質(zhì)的傳播速度時(shí)，則透射波會(huì)沿水平界面滑動(dòng)，滑行波在上層介質(zhì)中激發(fā)新的波動(dòng)即折射波。基于基本折射方程[9]，折射波靜校正通過(guò)讀取折射波初至?xí)r間[10]計(jì)算介質(zhì)的厚度以及速度信息，從而得到地震波在介質(zhì)中的傳播時(shí)間，完成對(duì)由于地表?xiàng)l件起伏引起的時(shí)差的校正。故應(yīng)用折射波靜校正必須具備兩個(gè)基本條件：一是要存在水平的折射面；二是下層速度要大于上層速度[11]。在本文設(shè)置的穿墻成像場(chǎng)景模型中，由于墻體設(shè)定為平行，故其擁有相對(duì)穩(wěn)定的折射面；另外，由于外層墻體的相對(duì)介電常數(shù)比內(nèi)層墻體相對(duì)介電常數(shù)大，電磁波在介質(zhì)內(nèi)的傳播速度與其相對(duì)介電常數(shù)成反比，故電磁波在外層墻體內(nèi)傳播速度比在內(nèi)層傳播速度小。因此，折射波靜校正理論運(yùn)用于穿墻成像領(lǐng)域能夠計(jì)算電磁波在墻體內(nèi)的傳播時(shí)延，以校正每一通道中由于墻體存在引起的異常傳播時(shí)延誤差。

在本文設(shè)置的穿墻場(chǎng)景中，墻體由兩層不同介質(zhì)組成。當(dāng)電磁波從發(fā)射天線發(fā)出后，一部分透射波會(huì)沿著兩層墻體的界面滑動(dòng)，另一部分透射波在內(nèi)層墻體與空氣的界面上滑動(dòng)，引起墻體振動(dòng)，產(chǎn)生折射波。折射波的折射路徑如圖2所示，在本文模型中折射面有兩個(gè)：內(nèi)層墻體與空氣接觸的界面設(shè)為折射面1，兩層相對(duì)介電常數(shù)不同的墻體接觸的界面設(shè)為折射面2。由基本折射方程知：

圖2 折射波延遲時(shí)路徑

(4)

(5)

2.2 外層墻體電磁波傳播延遲時(shí)

(6)

令

(7)

可得到一系列延遲時(shí)超定方程組，則式(5)可寫為

Gx=b。

(8)

式中：G=[g1,g2,…,gm]T為該超定方程組的系數(shù)矩陣，x為發(fā)射天線以及接收天線延遲時(shí)構(gòu)成的未知數(shù)向量。對(duì)式(8)作差可得

(9)

對(duì)式(9)構(gòu)建最小二乘：

(10)

對(duì)x求導(dǎo)得

(11)

式中：tr[xxTGTG]可寫為tr[xTGTGx]，整理可得

(12)

其中：

(13)

將式(13)代入式(12)中得

(14)

所以有

(15)

則

(16)

則有

x=(GTG)-1b。

(17)

求得的向量x為該超定方程組的解，即電磁波在外層墻體內(nèi)的發(fā)射點(diǎn)傳播延遲時(shí)τm,2和接收點(diǎn)傳播延遲時(shí)τn,2。

2.3 內(nèi)層墻體電磁波傳播時(shí)延

(18)

(19)

對(duì)于相鄰接收天線陣元n與n-1，其截距時(shí)間差[17]可表示為

(20)

則求得的相鄰發(fā)射天線截距時(shí)間差即相鄰發(fā)射天線延遲時(shí)間增量，而相鄰接收天線截距時(shí)間差即為相鄰接收天線延遲時(shí)間增量。

τ=τm,1+τn,1+τm,2+τn,2+τ″ 。

(21)

其中：

(22)

3 實(shí)驗(yàn)與分析

采用GprMax2D/3D建立圖1所示的穿墻仿真場(chǎng)景。天線陣列緊貼外側(cè)墻面擺放，采用3發(fā)13收天線裝置，目標(biāo)為球形理想導(dǎo)體，位于坐標(biāo)(1.1,1.2)處。若不考慮墻體存在，使用后向投影算法對(duì)目標(biāo)場(chǎng)景進(jìn)行直接成像，其成像結(jié)果如圖3所示，可見(jiàn)由于墻體存在引起的異常傳播時(shí)延，成像目標(biāo)位置會(huì)嚴(yán)重偏離真實(shí)坐標(biāo)。

圖3 原始成像效果圖

圖4 接收天線1的回波信號(hào)波形

(a)回波時(shí)刻離散點(diǎn)

如上文所述，首先通過(guò)解式(8)超定方程組求得外層墻體的發(fā)射天線與接收天線延遲時(shí)，再計(jì)算每個(gè)通道內(nèi)的發(fā)射天線與接收天線延遲時(shí)增量來(lái)確定內(nèi)層墻體的延遲時(shí)，最后令其兩者相加則可得電磁波在兩層墻體內(nèi)的傳播延遲時(shí)。圖6給出了本文算法的流程圖。

圖6 合成延遲時(shí)算法流程圖

為驗(yàn)證本文算法在不同墻體兩層介質(zhì)參數(shù)下的有效性，設(shè)置以下三組不同情況墻體參數(shù)值，每組墻體由介質(zhì)1和介質(zhì)2組成，相對(duì)介電常數(shù)和厚度分別為(ε1,d1)和(ε2,d2)，如表1所示。運(yùn)用基于折射波靜校正的合成延遲時(shí)算法對(duì)三種不同墻體參數(shù)情況時(shí)穿墻場(chǎng)景的成像結(jié)果如圖7所示，同時(shí)對(duì)比了運(yùn)用回波域?yàn)V波方法[6]以及自聚焦方法[7]在不同墻體參數(shù)值下的成像情況。與回波域?yàn)V波算法對(duì)比，本文所提算法對(duì)墻后目標(biāo)位置的定位更準(zhǔn)確，且更不易受墻體參數(shù)的影響。自聚焦算法成像遍歷電磁波的傳播時(shí)延的取值范圍，通過(guò)圖像評(píng)價(jià)標(biāo)準(zhǔn)選取成像結(jié)果最佳圖像確定電磁波傳播時(shí)延的值，故當(dāng)墻體參數(shù)改變時(shí)，其成像目標(biāo)位置不易受墻體參數(shù)影響，而本文算法與自聚焦算法對(duì)比，其成像目標(biāo)位置更精準(zhǔn)。

表1 不同墻體參數(shù)取值

圖7 本文算法與兩種算法成像結(jié)果對(duì)比

為進(jìn)一步分析算法的性能，對(duì)本文算法成像結(jié)果的圖像熵進(jìn)行對(duì)比。圖像熵能夠衡量墻后目標(biāo)物體成像的聚焦程度，其值越小則聚焦程度越好。若像素點(diǎn)個(gè)數(shù)為X×Y，p(i,j)表示第i行第j列的像素值，則該圖像的圖像熵ξ可表示為[14]

(23)

根據(jù)式(23) 求得在不同墻體參數(shù)情況下本文所提算法與回波域?yàn)V波算法和自聚焦算法成像結(jié)果的圖像熵值，其對(duì)比如表2所示。

表2 三種算法成像結(jié)果圖像熵對(duì)比

由表2可知，在不同墻體參數(shù)環(huán)境下，本文算法與兩種算法相比，其圖像熵值較小，即本文算法成像結(jié)果在聚焦程度方面比上述兩種算法效果更優(yōu)。

分別對(duì)回波域?yàn)V波算法、自聚焦算法以及本文所提算法進(jìn)行計(jì)算復(fù)雜度對(duì)比分析。對(duì)于回波域?yàn)V波算法，其計(jì)算復(fù)雜度可表示為

C0=mnlgp+mN+3N。

(24)

式中：P為采樣點(diǎn)數(shù)，m、n分別為穿墻成像區(qū)域網(wǎng)格橫向和縱向的劃分個(gè)數(shù)。對(duì)于自聚焦算法，其計(jì)算復(fù)雜度可表示為

C1=LmnP+mN。

(25)

式中：L為自聚焦算法延時(shí)補(bǔ)償值取值范圍的遍歷次數(shù)。而本文所提基于折射波靜校正的合成延遲時(shí)算法，其計(jì)算復(fù)雜度如式(26)所示：

C2=mnlg 2P。

(26)

圖8給出了三種算法在不同采樣點(diǎn)數(shù)情況下計(jì)算復(fù)雜度的曲線圖。由圖8可知，在相同采樣點(diǎn)數(shù)條件下，本文所提算法相較于其他兩種算法，其計(jì)算復(fù)雜度相對(duì)較低。

圖8 三種算法的計(jì)算復(fù)雜度對(duì)比曲線

為驗(yàn)證在噪聲影響下本文所提算法的魯棒性，在設(shè)置同一墻體參數(shù)下運(yùn)用本文算法以及上述兩種算法對(duì)穿墻場(chǎng)景進(jìn)行成像，并給出了不同信噪比條件下該算法的輸入輸出目標(biāo)雜波比(Target Clutter Ratio,TCR)。由圖9可知，相較于自聚焦算法，基于折射波靜校正的合成延遲時(shí)算法在不同信噪比的噪聲情況下能夠輸出更高的目標(biāo)雜波比；相較于回波域?yàn)V波方法，本文所提算法在噪聲信噪比為-3 dB之后能輸出更高的目標(biāo)雜波比。

圖9 噪聲影響下TCR的變化

4 結(jié) 論

針對(duì)未知參數(shù)的多層墻體穿墻雷達(dá)成像目標(biāo)散焦或位置偏移問(wèn)題，本文提出基于折射波合成延遲時(shí)靜校正算法對(duì)電磁波的異常傳播時(shí)延進(jìn)行校正并成像。該算法首先讀取回波信息建立超定方程組計(jì)算外層墻體內(nèi)電磁波的傳播時(shí)延，再通過(guò)合成延遲時(shí)增量估計(jì)內(nèi)層墻體每一通道的電磁波傳播時(shí)延，最后對(duì)每一通道電磁波的異常傳播時(shí)延進(jìn)行一一校正，使得成像目標(biāo)被校正到準(zhǔn)確位置。仿真結(jié)果表明，本文所提算法能夠有效校正墻體存在聚焦時(shí)延的偏差問(wèn)題。與傳統(tǒng)方法相比，本文所提方法能夠更加精確地將目標(biāo)點(diǎn)校正至正確位置，且圖像聚焦效果優(yōu)于傳統(tǒng)算法。理論推導(dǎo)以及程序仿真結(jié)果表明，本文方法在未知墻體參數(shù)穿墻成像的實(shí)際應(yīng)用中具有一定的參考和指導(dǎo)意義。