結合TV約束的穿墻雷達擴展目標成像方法

,,2

(1.桂林電子科技大學信息與通信學院,廣西桂林541004;2.廣西無線寬帶通信與信號處理重點實驗室,廣西桂林541004)

0 引言

超寬帶(Ultra Wideband,UWB)穿墻雷達成像(Through-the-Wall Radar Imaging,TWRI)是一種通過電磁波對隱藏目標進行檢測、定位、成像和特征提取的技術,在人質救援、建筑物布局檢測、隱藏目標探測等領域具有廣泛的應用前景[1]。現有研究都是將均勻墻作為背景環境,然而當成像環境是周期性結構的空心墻時,電磁波穿過這種墻體會發生Floquet現象,使用現有的均勻墻環境下的成像算法就會產生多徑虛假像[2]。

針對上述問題,文獻[2]對周期性空心墻后的目標成像問題進行了研究,它強調如果在散射模型中不考慮多徑傳播現象,就會導致成像結果中目標像的扭曲和失真,因此采用基于玻恩近似的線性化散射模型來解決這個問題,但是從最終結果可以看出這種方法的去噪效果并不理想。文獻[3]對周期性空心墻體結構導致的Floquet模式進行了仿真分析,結果證明周期性結構引起的高階空間諧波是導致穿墻成像中雜波和“鬼影”出現的主要原因。在此基礎上,比較了后向投影(Back Projection,BP)成像和傳統的稀疏成像方法。BP成像通過將每條多徑疊加到目標真實位置來消除Floquet現象,但這種方法數據量大、速度慢、主瓣寬且旁瓣多。而傳統稀疏成像雖然精度上有所提高,數據量減小,但是仍然存在速度較慢的問題。不僅如此,上述方法都是針對點目標的成像,考慮到實際情況中更常見的是擴展目標,如果繼續使用點目標模型對擴展目標進行成像,將會出現擴展目標邊緣輪廓不清晰的問題。

全變分法(Total Variation,TV)是一種基于偏微分方程的變分去噪模型,最初應用在圖像稀疏重構中[4]。文獻[5]最先將TV法應用到擴展目標雷達成像中,通過增加一階差分算子D在保留擴展目標邊緣輪廓的同時減少旁瓣。文獻[6]提出一種目標探測和雜波去除方法,它指出通過反復改變天線位置,直接路徑的目標反射會有相同的回波,而間接路徑的回波通常表示為隨機噪聲。因此,該方法可以從各類噪聲中分開目標像點。文獻[5-6]都能夠對擴展目標進行成像,但都是在均勻墻環境下的成像,不能解決空心墻體導致的目標位置偏離問題,同時面臨著需要人工設置參數的問題。

綜上所述,針對穿過周期墻體的成像場景而言,要在保留點目標和擴展目標特性的同時,解決空心墻體導致的目標位置偏離和需要人工設置參數的問題。本文提出一種結合TV約束的自適應擴展目標稀疏成像方法,首先根據電磁波穿過空心墻體的電磁傳播特性將空心墻體等效成三層均勻介質模型,基于貝葉斯推理準則,利用“Integrate-out”方法[7]得到點目標和擴展目標同時存在情況下的最大后驗(Maximum a Posteriori,MAP)估計,然后采用交替循環迭代法對各個未知量進行循環優化直至達到最優值。最后通過仿真和實驗驗證了該算法的有效性。

1 信號模型

如圖1所示,墻體相對介電常數為ε1,厚度為h,其中第一、二、三層的厚度分別為h1,h2,h3,墻體周期為d,空心結構之間的材料厚度為a。L個收發共置的天線構成的天線陣列進行接收和發射信號,天線陣列與墻體平行,發射由J個頻率組成的步進頻信號。當天線發射的電磁波以β角穿過這種特殊的周期結構墻體時,將會在墻體另一側的固定角度上產生與入射方向有關的Floquet模式,其第s階Floquet模式輸出角滿足條件:

式中,k=2πf/c表示自由空間波數。圖1中,天線在每一個頻點發射的電磁波都有很多入射平面波,根據式(1)可以知道,每一個入射波都會產生幾個Floquet透射波,這就會導致計算量的大大增加[3],然而只有直接路徑產生的像是真實目標像,由間接路徑產生的虛假像與真實目標相比旁瓣更大而能量卻較弱,且形狀很不規則,通常將這種間接路徑導致的雜波表示為噪聲[6]。在文獻[8]中將這種復雜的周期性空心墻等效成三層均勻介質模型,如圖2所示。

圖1 電磁傳播模型

圖2 空心墻的等效三層模型

第二層的介電常數表示為ε2=1-g+g·ε1,其中g=a/d表示填充參數。這樣墻體內傳播時延就表示為三層均勻介質內時延的疊加,由此可以得到總的傳播時延:

式中,x p和y p為像素點位的位置,x l和y l為第l個天線的位置,c為電磁波傳播速度,θ≈arctan(|y p-y l|/|x p-x l|)。由此將接收天線接收到的回波信號表示為

式中,η為周圍環境和墻體自身產生的噪聲,σp為第p個像素點的反射系數,τp,l為第l個天線和第p個像素點之間的雙向傳播時延,N x N y為成像區域劃分的像素點數。對式(3)進行矢量化表示,y∈C JL×1是包括所有JL次測量值y(j,l)的列向量:

式(7)中包括了JL的全部測量數據,可以使用較少的數據并從中恢復出σ。用一個MN×JL,(MN?JL)的測量矩陣Φ對y進行采樣[9],得到一個較稀疏的矩陣:

2 成像算法

對于不同天線位置的接收回波,直接路徑產生的目標反射會有相同的回波,而間接路徑產生的多徑回波類似隨機噪聲,由于TV法通常只對噪聲敏感而平滑信號的TV項通常較小,因此本文采用TV法從各類噪聲中分出目標,在去除噪聲的同時保留擴展目標的邊緣信息[6]。然后采用交替循環迭代的方法[10-12],對各個變量進行循環優化,直至求出最優值。這種方法避免了人工設置參數和復雜的噪聲功率求解,簡化了求解過程。

2.1 MAP估計

采用MAP方法,假設噪聲矢量η服從均值為零、協方差矩陣為βI的復高斯分布,則觀測向量的似然函數[13]表示為

由于散射體分布規律的稀疏性,σ通常滿足Laplace分布,并用Jeffrey先驗作為Laplace分布的先驗參數[7],那么σ的先驗信息表示為

式中,第一項表示點目標的散射特性,第二項表示擴展目標的散射特性,α1,α2為Laplace分布的先驗參數,表示對散射體特征的增強程度,D為一階差分算子,[D]i,i=1,[D]i,i+1=-1,1≤i≤N x N y-1,是為了保留擴展目標的邊緣特性。Jeffrey先驗由Fisher信息量I(·)定義[7]:

將式(11)、式(12)和式(13)代入式(10),采用“Integrate-out”方法[7],可得

式中,Γ(·)表示Gamma分布,基于貝葉斯準則,σ的后驗概率表示為

在MAP估計中,通常將p(y)忽略,故σ的M AP估計表示為

為了簡化求解,首先進行取對數操作,即

對上式進行整理,去掉常數項,式(16)可以表示為

2.2 交替循環法

采用傳統的牛頓法[5]來求解J(σ,β)的最優化問題需要的計算量很大,通常目標函數的Hessian矩陣無法保持正定從而使這些方法失效,同時伴隨著正則化參數的選擇困難問題。因此,這里使用交替循環迭代的方法求解反射系數σ和參數β。

首先假定得到了第t次迭代的β(t),求出J(σ,β(t))關于σ的導數并令其為零,可求解出σ(t+1),即

式中,Λ1=diag{1/|σ|i},Λ2=diag{1/|Dσ|i},為了避免出現奇異值,在Λ1,Λ2中加入一個極小值κ,表示為}。式(20)變為

式(25)中,逆運算的求解運算量較大,這里采用文獻[8]的共軛梯度最小二乘法(Conjugate Gradient Least Squares,CGLS)對其進行求解,從而實現快速計算。

接下來,對J(σ(t),β)關于β求導,并令其為零,即

2.3 算法實現

算法流程如表1所示,首先采用BP算法得到初始值σ(0),初始值,最大迭代次數設置為T。循環迭代式(25)和式(27),直到滿足收斂條件。很明顯,σ和β通過式(25)和式(27)的迭代更新,J(σ,β)會隨著t的增加而逐漸減小。文獻[11]中已經證明,交替循環迭代算法是收斂的,收斂條件為‖σ(t+1)-σ(t)‖2/‖σ(t)‖2＜δ,這里的δ為很小的正實數。

表1 算法流程

從算法的流程可以知道,每次迭代過程都需要對σ和β進行更新,其中對σ的求解過程需要進行逆運算,求解式(25)的逆運算計算復雜度為O((N x N y)3)。為了減小計算量,首先采用文獻[12]的“Integrate-out”方法將σ無參數化,在這個過程中,充分利用了目標的先驗信息,提高了成像結果的穩健性。然后采用文獻[11]中的CGLS算法求式(25)的逆運算,CGLS算法在每次迭代只包括兩個矩陣的乘積,因此每次迭代的計算復雜度僅為O(MN×N x N y),而且文獻中提到CGLS算法的迭代次數只與接收信號的信噪比有關,這就在實現快速收斂的同時大大減小了計算量。

3 仿真及實驗結果分析

3.1 仿真結果分析

使用XFDTD電磁仿真軟件進行仿真。非均勻空心墻的介電常數ε1=5,電導率δ=0.01 S/m,厚度為0.2 m,第一、二、三層的厚度分別為0.04,0.12和0.04 m,墻體周期為0.15 m,空心結構之間的材料厚度為0.04 m。由21個收發共置天線組成的天線陣,與墻體平行且距離墻體0.5 m,陣元之間間隔0.05 m,發射頻率2 GHz脈寬1 ns的調制高斯脈沖信號,提取頻段1.5～2.5 GHz的DFT數據,共101個頻點。圖3給出了XFDTD仿真模型,圖3(a)為長0.4 m、寬0.2 m的矩形單目標體,圖3(b)為兩個矩形目標體,各目標與墻體的距離均為1 m。

圖3 XFDTD仿真模型

圖4和圖5分別給出了單目標體和兩個目標體的成像結果,圖像中的虛線方框表示目標位置。從圖中可以看出,圖4(a)和圖5(a)的BP成像方法效果較差,存在大量的雜波,很難識別出目標位置,擴展目標輪廓不能很好保存。圖4(b)和圖5(b)采用傳統的稀疏成像方法,雜波減少,但仍然存在多徑虛假像。而圖4(c)和圖5(c)采用本文的方法成像相比其他兩種效果明顯改善,虛假像得到消除,目標準確且擴展目標輪廓十分清晰,同時本文的優化方法提高了計算速度。

圖4 單目標成像結果

圖5 兩個目標成像結果

為了更直觀地看出各算法性能上的差異,對存在一個目標的情況下,在程序運行時間(Program Running Time,PRT)、目標雜波比(Target-to-Clutter Ratio,TCR)、圖像的熵(Entopy of the Full Image,ENT)和一階差分均值(Mean of First Order Difference,MFOD)四個方面進行比較[14],結果如表2所示。PRT為算法的運算時間,PRT越小,運行時間越短。TCR為目標區域最大像素值和雜波區域平均像素值的比值,TCR值越大,表明目標相對背景雜波就越強,背景雜波抑制越充分,圖像越清晰。ENT為整幅圖的復雜程度,其值越小整幅圖雜波越小、越干凈。MFOD為分布目標的輪廓特性,MFOD值越小,說明檢測區域像素值越接近,區域越平滑。從表中可以看出,BP方法的運算速度最快但背景雜波嚴重,稀疏方法的運算速度很慢,但是圖像相對較清晰,本文方法運算速度介于BP方法和稀疏方法之間,但是其他各方面都有所改善。

表2 PRT、TCR、ENT和MFOD算法比較

3.2 實驗結果分析

使用美國GSSI公司的探地雷達SIR2.0搭建穿墻實驗場景,如圖6(a)所示。使用的是長和寬分別為0.39 m,0.19 m的空心磚,第一、二、三層的厚度分別為0.035,0.12和0.035 m,墻體周期為0.18 m,相對介電常數為5。墻后矩形目標距離墻體1 m左右,目標的長和寬分別為0.7 m和0.2 m。選擇2 GHz的喇叭天線,天線平行于墻體移動21次,合成孔徑長度1 m,在每個位置測量2次,包括有目標和無目標場景。將采集到的數據波形進行平均處理。圖6(b)給出了成像結果,虛線框為目標的真實位置。圖中可以看出使用本文方法能保證在得到目標輪廓信息的同時抑制雜波,證明了本文方法的有效性。

圖6 穿墻實驗結果

4 結束語

本文提出的結合TV約束的擴展目標成像方法,相比較于常規算法解決了空心墻導致的Floquet雜波問題,能夠在空心墻環境下進行高分辨成像,而且該算法能有效對擴展目標進行邊緣成像。同時使用交替循環法進行迭代大大縮減了計算時間。在今后的工作中,將對該方法進行完善使其能應用到更復雜的場景中。

[1]YEKTAKHAH B,SARABANDI K.All-Directions Through-the-Wall Radar Imaging Using a Small Number of Moving Transceivers[J].IEEE Trans on Geoscience and Remote Sensing,2016,54(11):6415-6428.

[2]GENNARELLI G,SOLDOVIERI F.Radar Imaging Through Cinderblock Walls:Achievable Performance by a Model-Corrected Linear Inverse Scattering Approach[J].IEEE Trans on Geoscience and Remote Sensing,2014,52(10):6738-6749.

[3]CHANG P C.Physics-Based Inverse Processing and Multipath Exploitation for Through-Wall Radar Imaging[D].Columbus:The Ohio State University,2011.

[4]ZHAO Yao,LIU Jianguo,ZHANG Bingchen,et al.Adaptive Total Variation Regularization Based SAR Image Despeckling and Despeckling Evaluation Index[J].IEEE Trans on Geoscience and Remote Sensing,2015,53(5):2765-2774.

[5]SUN Bin,JIN Haibo.Feature Enhanced Synthetic Aperture Radar Image Formation[C]∥8th International Conference on Electronic Measurement and Instruments,Xi’an:IEEE,2007:2803-2809.

[6]HANDA H,MANSOUR H,LIU D,et al.Extended Target Localization with Total-Variation Denoising in Through-the-Wall-Imaging[C]∥IEEE 6th International Workshop on Computational Advances in Multi-Sensor Adaptive Processing,Cancun:IEEE,2015:445-448.

[7]OGWORONJO H C,ANDERSON J M M,NGUYEN L.A Parameter-Free MAP Image Reconstruction Algorithm for Impulse-Based UWB Ground Penetrating Radar[C]∥IEEE Global Conference on Signal and Information Processing,Orlando,FL:IEEE,2015:742-746.

[8]BURKHOLDER R J,MARHEFKA R J,VOLAKIS J L.Radar Imaging Through Cinder Block Walls and Other Periodic Structures[C]∥IEEE Antennas and Propagation Society International Symposium,San Diego,CA:IEEE,2008:1-4.

[9]LEIGSNERING M,AHMAD F,AMIN M,et al.Multipath Exploitation in Through-the-Wall Radar Imaging Using Sparse Reconstruction[J].IEEE Trans on Aerospace and Electronic Systems,2014,50(2):920-939.

[10]LIM D,XU L,GIANELLI C,et al.Time-and Frequency-Domain MIMO FLGPR Imaging[C]∥IEEE Radar Conference,Arlington,VA:IEEE,2015:1305-1310.

[11]TAN X,ROBERTS W,LI J,et al.Sparse Learning via Iterative Minimization with Application to MIMO Radar Imaging[J].IEEE Trans on Signal Processing,2011,59(3):1088-1101.

[12]GLENTIS G O,ZHAO K,JAKOBSSON A,et al.Non-Parametric High-Resolution SAR Imaging[J].IEEE Trans on Signal Processing,2013,61(7):1614-1624.

[13]LIU Hongchao,JIU Bo,LIU Hongwei,et al.An Adaptive ISAR Imaging Method Based on Evidence Framework[J].IEEE Geoscience and Remote Sensing Letters,2014,11(6):1031-1035.

[14]晉良念,申文婷,錢玉彬,等.組合字典下超寬帶穿墻雷達自適應稀疏成像方法[J].電子與信息學報,2016,38(5):1047-1054.JIN Liangnian,SHEN Wenting,QIAN Yubin,et al.Adaptive Sparse Imaging Approach for Ultra-Wideband Through-the-Wall Radar in Combined Dictionaries[J].Journal of Electronics&Information Technology,2016,38(5):1047-1054.(in Chinese)