文獻(xiàn)標(biāo)識(shí)碼: A
DOI:10.16157/j.issn.0258-7998.2015.11.038
中文引用格式: 晉良念,錢玉彬,申文婷,等. 基于改進(jìn)OMP的超寬帶穿墻雷達(dá)稀疏成像方法[J].電子技術(shù)應(yīng)用,2015,41(11):135-139.
英文引用格式: Jin Liangnian,Qian Yubin,Shen Wenting,et al. Sparse imaging for ultra-wideband through-the-wall radar based on modified OMP algorithm[J].Application of Electronic Technique,2015,41(11):135-139.
0 引言
穿墻雷達(dá)成像(Through the Wall Radar Imageing,TWRI)使用電磁波穿透墻體和內(nèi)部封閉的建筑物進(jìn)行監(jiān)測(cè),以確定建筑結(jié)構(gòu)布局,區(qū)分建筑物內(nèi)部的活動(dòng)情況,檢測(cè)、識(shí)別和跟蹤運(yùn)動(dòng)目標(biāo)。TWRI技術(shù)在搜救、墻后目標(biāo)檢測(cè)和城區(qū)環(huán)境的監(jiān)視等方面?zhèn)涫荜P(guān)注[1-4]。TWRI系統(tǒng)通常利用超寬帶信號(hào)和大孔徑天線陣列提供高的距離向和方位向分辨率,能夠在較短的采樣時(shí)間內(nèi)利用較小的存儲(chǔ)空間獲得高質(zhì)量圖像,與之對(duì)應(yīng)的后向投影一類算法需要大的空間和時(shí)間(或頻率)采樣數(shù)據(jù)才能獲得高質(zhì)量成像。
近年來,壓縮感知理論的興起,國(guó)內(nèi)外學(xué)者將該理論很好地應(yīng)用到TWRI中,提出了一些高分辨稀疏成像方法[5-6]。它們首先將成像區(qū)域劃分為有限個(gè)網(wǎng)格,當(dāng)目標(biāo)恰好落在網(wǎng)格點(diǎn)上,稀疏成像效果很好。一旦目標(biāo)偏離了網(wǎng)格點(diǎn),基矩陣不匹配,從而引起成像模糊,甚至散焦。為此,文獻(xiàn)[5]通過建立off-gird稀疏表示模型,提出了一種基于L21范數(shù)懲罰項(xiàng)的FISTA算法,該算法雖能較準(zhǔn)確地恢復(fù)目標(biāo)散射系數(shù),但是算法的性能受人工參數(shù)設(shè)置的影響,在實(shí)際應(yīng)用中很難正確地選擇。文獻(xiàn)[6]針對(duì)穿墻雷達(dá)偏離網(wǎng)格目標(biāo)成像中提出了一種貝葉斯推理方法,克服了FISTA算法性能受人工參數(shù)設(shè)置的影響。但是,該方法在成像區(qū)域劃分網(wǎng)格數(shù)過大時(shí)處理速度會(huì)很慢,所以需要尋求新的方法提高處理速度。
本文提出一種改進(jìn)的正交匹配追蹤方法(GBTOMP算法)實(shí)現(xiàn)偏離網(wǎng)格目標(biāo)超寬帶穿墻雷達(dá)稀疏成像。該方法以正交匹配追蹤(Orthogonal Matching Pursuit,OMP)算法為基礎(chǔ),先通過最速上升梯度優(yōu)化方法搜索到目標(biāo)的真實(shí)位置,逐漸減小網(wǎng)格偏移量,修正模型中的感知矩陣,再由修正的感知矩陣通過最小二乘方法計(jì)算目標(biāo)散射系數(shù),但在此計(jì)算過程中,會(huì)引入冗余下標(biāo)。本文利用貝葉斯假設(shè)檢驗(yàn)去冗余,這樣既可以減小重構(gòu)過程的運(yùn)算量,還可以增強(qiáng)重構(gòu)的精度和抗噪性能。仿真和實(shí)驗(yàn)結(jié)果驗(yàn)證了本文方法的有效性和可行性。
1 稀疏表示模型
一般的穿墻傳播模型可以等效為如圖1所示的兩層傳播模型[7]。假設(shè)均勻介質(zhì)墻體的厚度為d,相對(duì)介電常數(shù)為?著。將成像區(qū)域在方位向和距離向上離散成Nx Ny個(gè)網(wǎng)格點(diǎn),對(duì)應(yīng)的像素值用?滓k,l,k=1,…,Nx,l=1,…,Ny表示,當(dāng)網(wǎng)格點(diǎn)有目標(biāo)時(shí)k,l不為0,反之k,l等于0。為了方便表示,將k,l按照順序重新排列成Nx Ny維的目標(biāo)散射系數(shù)列向量,令所有元素的網(wǎng)格位置用矩陣表示,其中的元素用p=[xp,yp]T表示,對(duì)應(yīng)的像素值用p表示,則第n個(gè)天線回波信號(hào)的離散數(shù)據(jù)模型為:
式中,yn=[yn(t1),…,yn(tM)]T,An為M×Nx Ny的矩陣,第m行的元素為[s(tn,m-n,1),…,s(tn,m)],這里n,p為第n個(gè)天線到第p個(gè)網(wǎng)格點(diǎn)之間的傳播時(shí)延。
根據(jù)圖1所示的傳播模型,n,p可以近似表示為:
式中,xn、yn分別表示天線的坐標(biāo)位置,n,p為傳播的入射角,當(dāng)成像區(qū)域距離天線較遠(yuǎn)時(shí)可以近似為:
將N個(gè)天線到成像區(qū)域中所有像素點(diǎn)的傳播時(shí)延按照式(2)和式(3)來計(jì)算,并把所有天線的接收數(shù)據(jù)進(jìn)行堆疊,由此所構(gòu)造的離散數(shù)據(jù)模型為:
然而,實(shí)際情況中目標(biāo)很可能偏離預(yù)設(shè)網(wǎng)格,如圖1所示,這種偏移會(huì)影響A′的精確計(jì)算。定義成像空間所有網(wǎng)格的偏移量矩陣為,A′關(guān)于的函數(shù)為A′。另外,對(duì)于穿墻成像場(chǎng)景來說滿足稀疏特性,因此可以利用較少的測(cè)量值恢復(fù)。令表示從y中隨機(jī)選擇選擇Q1個(gè)天線單元和每個(gè)天線單元隨機(jī)選擇Q2個(gè)數(shù)據(jù)得到的Q1Q2個(gè)數(shù)據(jù)矢量,并定義這種隨機(jī)選擇形式構(gòu)成的測(cè)量矩陣用表示[8],由此得到大小為Q1Q2×1的測(cè)量數(shù)據(jù)(或稀疏成像)模型(為了書寫方便將A為:
根據(jù)式(5)的稀疏成像模型可以看出,感知矩陣A有關(guān),因此求解式(4)反演需要事先準(zhǔn)確估計(jì)?專參數(shù)。當(dāng)且僅當(dāng)A滿足約束等距性質(zhì)時(shí),可以通過求解如下優(yōu)化問題將高維信號(hào)從低維的測(cè)量值中準(zhǔn)確重構(gòu)出來。
由于最小0范數(shù)問題是一個(gè)非確定性多項(xiàng)式難問題,目前的重構(gòu)算法都是基于0范數(shù)的變換或逼近處理,其中貪婪算法以其運(yùn)算速度快,實(shí)現(xiàn)簡(jiǎn)單的優(yōu)勢(shì)被用于穿墻雷達(dá)成像中。
2 GBTOMP的稀疏成像方法
2.1 GBTOMP算法
本文提出的GBTOMP算法是在正交匹配追蹤算法[9]的迭代過程中采用梯度優(yōu)化最速上升算法校正網(wǎng)格偏移量,然后采用貝葉斯假設(shè)檢驗(yàn)[10]剔除冗余下標(biāo)。該算法整體上可以分為兩大步:第一,根據(jù)殘差與感知矩陣的內(nèi)積最大值找到最匹配原子位置(預(yù)設(shè)網(wǎng)格位置)。由于目標(biāo)可能會(huì)偏離預(yù)設(shè)的網(wǎng)格位置,因此在該位置需要進(jìn)行梯度優(yōu)化的最速上升,以搜索到目標(biāo)的真實(shí)位置,從而以目標(biāo)真實(shí)位置建立新的感知矩陣,并由該矩陣?yán)米钚《朔ㄇ蟮蒙⑸湎禂?shù)向量的近似值,同時(shí)通過交替迭代更新下標(biāo)集求解來得到目標(biāo)散射系數(shù)。第二,OMP算法在有噪聲和沒有噪聲的情況下都會(huì)存在冗余下標(biāo)問題。為解決此問題,GBTOMP算法最后采用貝葉斯假設(shè)檢驗(yàn)理論給出判決準(zhǔn)則對(duì)第一步輸出的預(yù)選下標(biāo)集中原子進(jìn)行篩選,得到等于或者大于信號(hào)真實(shí)下標(biāo)集估計(jì)后,再用最小二乘算法進(jìn)行重構(gòu),以改善算法的性能,增強(qiáng)其抗噪能力。GBTOMP算法流程如下。
初始化:構(gòu)造測(cè)量向量,網(wǎng)格偏移矩陣,感知矩陣A。
(1)采用OMP算法找到匹配原子對(duì)應(yīng)的空間網(wǎng)格位置,利用最速上升梯度優(yōu)化方法對(duì)偏離網(wǎng)格進(jìn)行估計(jì),再利用最小二乘法對(duì)目標(biāo)散射系數(shù)進(jìn)行粗估計(jì),輸出預(yù)選下標(biāo)集P和粗重構(gòu)目標(biāo)散射系數(shù)X。
(2)從預(yù)選下標(biāo)集P出發(fā),通過貝葉斯假設(shè)檢驗(yàn)準(zhǔn)則設(shè)置一個(gè)適當(dāng)門限Thj對(duì)P中的元素進(jìn)行篩選,將大于門限的元素都保留在最終的下標(biāo)集F中。
輸出:根據(jù)最終下標(biāo)集F,利用最小二乘法得到目標(biāo)散射系數(shù)最終的估計(jì)值。
2.2 網(wǎng)格偏移量的估計(jì)
在GBTOMP算法的迭代過程中需要估計(jì)網(wǎng)格偏移量,其代價(jià)函數(shù)的構(gòu)造是關(guān)鍵。OMP算法是以測(cè)量數(shù)據(jù)和感知矩陣A構(gòu)造代價(jià)函數(shù)的。為了能夠使得估計(jì)準(zhǔn)確,并使處理過程更加穩(wěn)定,本文使用OMP算法相關(guān)函數(shù)的平方開方項(xiàng)來表示[11]:
首先初始化迭代更新次數(shù)k=1,殘差向量r。計(jì)算感知矩陣A所有列和殘差向量r的內(nèi)積值,找到內(nèi)積值最大時(shí)所對(duì)應(yīng)的列下標(biāo)pk,接著進(jìn)入搜索目標(biāo)真實(shí)位置過程,滿足搜索條件|Ji-Ji-1|<a結(jié)束,這里的i代表搜索過程中的迭代次數(shù),a代表很小的一個(gè)正數(shù)。在此過程中,從?仔中取出第pk列賦給?仔0,以該點(diǎn)坐標(biāo)(xp,yp)計(jì)算出與隨機(jī)選擇Q1個(gè)天線組成新集合的第q個(gè)索引號(hào)天線位置之間的雙程傳播時(shí)延為:
為求解代價(jià)函數(shù)J關(guān)于所選擇坐標(biāo)(xp,yp)的梯度,需要事先求解?子q,p關(guān)于(xp,yp)的偏導(dǎo),分兩種情況進(jìn)行:
(1)當(dāng)xq=xp時(shí),即q,p=0,q,p簡(jiǎn)化為:
(2)當(dāng)xq≠xp時(shí),即q,p≠0,q,p關(guān)于xp和yp的偏導(dǎo)數(shù)為:
假設(shè)輻射脈沖波形為高斯脈沖s(t),則第q個(gè)天線和所選目標(biāo)坐標(biāo)(xp,yp)對(duì)應(yīng)的感知矩陣Aq、B是大小為Q2×1的矩陣,它們的第q2行分別為:
式中是當(dāng)前(xp,yp)對(duì)應(yīng)的感知矩陣,僅僅是A的一個(gè)子集,這樣處理是為了提高運(yùn)行速度,減小存儲(chǔ)空間。B1x和B1y表示對(duì)矩陣分別關(guān)于xp和yp的偏導(dǎo)數(shù),它們的具體表達(dá)式與s(t)有關(guān)。以梯度搜索方法更新所選空間的坐標(biāo)為:
式中是搜索步長(zhǎng),本文中以網(wǎng)格尺寸為參考的很小的正數(shù)。由此得到更新后的重新計(jì)算代價(jià)函數(shù)Ji,搜索迭代次數(shù)i=i+1,一旦代價(jià)函數(shù)滿足條件|Ji-Ji-1|<a就終止搜索過程。終止后,將所選的網(wǎng)格空間位置?仔1所對(duì)應(yīng)的感知矩陣保存在A2矩陣中,同時(shí)將所選的感知矩陣列下標(biāo)pk添加到預(yù)選下標(biāo)集P中,由A2通過最小二乘法計(jì)算出第k次粗估計(jì)的散射系數(shù)xpre,并更新殘差向量。
繼續(xù)搜索下標(biāo)并校正網(wǎng)格偏移量,如此循環(huán)迭代,k=k+1,直到殘差向量r滿足條件,粗估計(jì)結(jié)束。將k次循環(huán)求得的散射系數(shù)xpre對(duì)應(yīng)其下標(biāo)集P賦給Nx Ny維列向量X。
X{P}=xpre(18)
2.3 貝葉斯假設(shè)檢驗(yàn)去冗余
為剔除OMP算法輸出的預(yù)選集P中的冗余下標(biāo),根據(jù)貝葉斯假設(shè)檢驗(yàn)?zāi)P蚚10]構(gòu)造似然函數(shù),其似然比檢驗(yàn)公式為:
根據(jù)判決準(zhǔn)則可知,當(dāng)(zj-mj)2>Thj時(shí),下標(biāo)j對(duì)應(yīng)的散射系數(shù)?滓j的值為非零,更新此時(shí)的j為最終下標(biāo),反之,下標(biāo)j對(duì)應(yīng)的散射系數(shù)?滓j為零,則從預(yù)選下標(biāo)集中剔除此時(shí)的j值。對(duì)預(yù)選下標(biāo)集中元素進(jìn)行篩選,最終剩下的真實(shí)下標(biāo)集為F,接下來利用最小二乘法進(jìn)行目標(biāo)散射系數(shù)重構(gòu)完成整個(gè)GBTOMP算法。
3 仿真與實(shí)驗(yàn)結(jié)果分析
3.1 仿真結(jié)果分析
利用GPRMAX產(chǎn)生仿真數(shù)據(jù)[12],發(fā)射信號(hào)采用高斯脈沖,脈沖中心頻率為1 GHz。墻體厚度為0.20 m,相對(duì)介電常數(shù)為6,電導(dǎo)率為0.03 S/m。
場(chǎng)景1:天線陣列由16個(gè)收發(fā)共置單元組成,隨機(jī)選擇其中的10個(gè)天線單元,陣元間距為0.1 m。成像區(qū)域是[-1 m,1 m]×[0.3 m,1.3 m],距離向和方位向的網(wǎng)格尺寸均為0.05 m。假設(shè)成像區(qū)域中有三個(gè)點(diǎn)目標(biāo),其坐標(biāo)分別為(-0.414 m,0.514 m),(0.025 m,0.725 m),(0.325 m,0.965 m)。圖2為OMP方法和本文GBTOMP方法成像結(jié)果。從圖2(a)可以看出,OMP方法成像出現(xiàn)嚴(yán)重散焦,同時(shí)帶來虛假像,從圖2(b)可以看到,估計(jì)的目標(biāo)位置和預(yù)設(shè)目標(biāo)位置完全吻合,同時(shí)目標(biāo)像聚焦程度高,很明顯本文方法性能優(yōu)于傳統(tǒng)OMP方法。
場(chǎng)景2:天線陣列由45個(gè)收發(fā)共置單元組成,隨機(jī)選擇其中的10個(gè)天線單元,陣元間距為0.04 m,與墻體之間的距離為0.2 m。圖3給出了2維GPRMAX模型,圖中左邊的長(zhǎng)方形目標(biāo)中心位于(0.7 m,3.25 m),長(zhǎng)為0.4 m,寬為0.3 m,右邊的圓形目標(biāo)體的中心位置位于(1.25 m,3.05 m),半徑為0.2 m。GPRMAX網(wǎng)格單元尺寸為0.005 m,采樣時(shí)窗為30 ns。假設(shè)回波信號(hào)由有目標(biāo)時(shí)減去無目標(biāo)時(shí)的信號(hào)進(jìn)行模擬。
圖4給出了兩個(gè)擴(kuò)展目標(biāo)的成像結(jié)果,成像區(qū)域的方位向和距離向的網(wǎng)格尺寸設(shè)為0.03 m。圖中的矩形和圓形虛線框代表兩目標(biāo)的真實(shí)位置。從圖4可以看出,傳統(tǒng)OMP方法的成像效果較差,雖然能夠大致識(shí)別目標(biāo)位置,但是引入了很多的虛假目標(biāo)像,而本文方法GBTOMP成像效果明顯改善,目標(biāo)像清晰。
3.2 實(shí)驗(yàn)結(jié)果分析
使用美國(guó)GSSI公司的探地雷達(dá)SIR-20搭建穿墻實(shí)驗(yàn)場(chǎng)景,如圖5所示。實(shí)驗(yàn)墻體為實(shí)心磚墻,墻體厚度0.20 m,相對(duì)介電常數(shù)6.4。選用1 GHz的喇叭天線,架高1.2 m,貼著墻體。水平移動(dòng)喇叭天線共掃描21個(gè)測(cè)點(diǎn),數(shù)據(jù)處理時(shí)隨機(jī)選擇10個(gè)測(cè)點(diǎn),測(cè)點(diǎn)間距為0.05 m,合成孔徑長(zhǎng)度為1 m,要求在每個(gè)天線孔徑測(cè)試點(diǎn)處測(cè)量2次,包括有人和無人的場(chǎng)景。人體目標(biāo)高度為1.80 m,體型寬度約0.40 m,站在墻后1 m處。SIR-20系統(tǒng)數(shù)據(jù)采集的參數(shù)設(shè)置:每道采樣1 024點(diǎn),每秒60道,時(shí)間窗15 ns。
首先將SIR-20系統(tǒng)采集的數(shù)據(jù)波形做平均處理,然后經(jīng)過去噪和雜波相消、自動(dòng)增益控制等信號(hào)處理得到較好的目標(biāo)回波,所有成像都是以天線陣列中心為坐標(biāo)原點(diǎn)。圖6給出實(shí)驗(yàn)數(shù)據(jù)成像結(jié)果比較,圖中的矩形虛線框代表目標(biāo)的真實(shí)位置。從圖中可以看出,傳統(tǒng)OMP方法目標(biāo)成像散焦現(xiàn)象較嚴(yán)重,目標(biāo)區(qū)域不能確定,且周圍還有較多的虛假像。而本文GBTOMP方法成像則不同,雜波相對(duì)較少,且目標(biāo)位置清晰可見,由此看來通過校正目標(biāo)偏離網(wǎng)格能夠很好解決傳統(tǒng)稀疏成像所出現(xiàn)的散焦現(xiàn)象。
4 結(jié)論
本文提出的GBTOMP稀疏成像方法,交替迭代得到目標(biāo)真實(shí)空間位置和目標(biāo)散射系數(shù),目標(biāo)位置的校正保證了基矩陣的匹配,實(shí)現(xiàn)目標(biāo)的高質(zhì)量稀疏成像。仿真和實(shí)測(cè)數(shù)據(jù)處理結(jié)果表明:與OMP成像相比,本文GBTOMP方法成像效果在目標(biāo)位置準(zhǔn)確度還是目標(biāo)聚焦程度有明顯的提高。但是,由于以點(diǎn)稀疏模型對(duì)擴(kuò)展目標(biāo)進(jìn)行成像,存在圖像不夠平滑,也不連續(xù)等問題,下一步工作將對(duì)擴(kuò)展目標(biāo)的稀疏成像方法展開研究。
參考文獻(xiàn)
[1] TIVIVE F H C,BOUZERDOUM A,AMIN M G.A subspaceprojection approach for wall clutter mitigation in through--wall radar imaging[J].IEEE Transactions on Geoscienceand Remote Sensing,2015,53(4):2108-2122.
[2] LI G,BURKHOLDER R J.Hybrid matching pursuit for distributed through-wall radar imaging[J].IEEE Transactionson Antennas and Propagation,2015,63(4):1701-1711.
[3] 蔣留兵,韋洪浪,許騰飛,等.EEMD生命探測(cè)雷達(dá)人體數(shù)量識(shí)別技術(shù)[J].電子技術(shù)應(yīng)用,2014,40(5):122-125.
[4] LIU J,CUI G,JIA Y,et al.Sidewall detection using mul-tipath in through-wall radar moving target tracking[J].IEEEGeoscience and Remote Sensing Letters,2015,12(6): 1372-1376.
[5] XIA S,LIU F.Off-grid compressive sensing through-the-wall radar imaging[C].SPIE Defense+Security.International Society for Optics and Photonics,2014:90771F-90771F-8.
[6] 晉良念,錢玉彬,劉慶華,等.超寬帶穿墻雷達(dá)偏離網(wǎng)格目標(biāo)稀疏成像方法[J].儀器儀表學(xué)報(bào),2015,36(4).
[7] JIN T,CHEN B,ZHOU Z.Imaging-domain estimation of wall parameters for autofocusing of through-the-wall SAR imagery[J].IEEE Transactions on Geoscience and Remote Sensing,2013,51(3):1836-1843.
[8] LAGUNAS E,AMIN M G,AHMAD F,et al.Joint wall mitigation and compressive sensing for indoor image reconstruction[J].IEEE Transactions on Geoscience and Remote Sensing,2013,51(2):891-906.
[9] 李少東,裴文炯,楊軍,等.貝葉斯模型下的OMP重構(gòu)算法及應(yīng)用[J].系統(tǒng)工程與電子技術(shù),2015,37(2).
[10] 甘偉,許錄平,蘇哲,等.基于貝葉斯假設(shè)檢驗(yàn)的壓縮感知重構(gòu)[J].電子與信息學(xué)報(bào),2011,33(11):2640-2646.
[11] GURBUZ A C,TEKE O,ARIKAN O.Sparse ground-penetrating radar imaging method for off-the-grid target problem[J].Journal of Electronic Imaging,2013,22(2):021007-021007.
[12] XIE J L,XU J L.Ground penetrating radar-based experimental simulation and signal interpretation on roadway roof separation detection[J].Arabian Journal of Geosciences,2014:1-8.