李志希,孔令講,賈 勇,趙中興

(電子科技大學(xué)電子工程學(xué)院,四川成都611731)

0 引言

穿墻雷達成像由于其在反恐偵察、抗災(zāi)救援等方面的重要應(yīng)用,近年來逐漸成為一個新興的研究領(lǐng)域[1]。穿墻雷達成像技術(shù)主要包含兩個研究方向,即穿墻人體目標成像[2-3]和建筑布局成像[4-6]。

建筑布局成像因其包含諸多好處近段時間逐步成為一個研究熱點。首先,建筑布局可以作為確定穿墻目標相對位置的重要參照物;其次,建筑布局一旦確定,便能為多徑的提取和目標位置偏移的糾正提供可能性;最后,建筑布局的確立為操作人員進入建筑物內(nèi)部提供便利。目前,世界上存在若干研究機構(gòu)就建筑布局成像進行研究[7-8]。美國陸軍實驗室對一簡單的實際應(yīng)用場景采用兩視角成像融合獲得了一幅粗略的全景圖像。其中,每幅單視角圖像的獲得是在單視角下兩通道圖像簡單相加而來,而各個單通道圖像是利用合成孔徑方式采用后向投影成像算法(Back Projected Algorithm)獲得[9-10]。然而,電磁波信號穿墻傳播引起的信號強度衰減現(xiàn)象帶來多層墻體成像強度的非均勻,這為最終建筑布局的檢測和獲得帶來了較大的困難。

為補償多層墻體成像強度的衰減,本文提出了一種基于霍夫變換門限檢測的補償方法。首先,對預(yù)處理后的實測數(shù)據(jù)采用后向投影成像算法,并用PCA和模糊邏輯圖像增強方法[11]處理,并產(chǎn)生原始圖像。其次,基于固定門限值劃分目標墻體區(qū)域和非目標雜波區(qū)域,門限值的確定采用本文提出的基于霍夫變換直線檢測原理的處理算法,并提出數(shù)據(jù)連通域的概念。再次,獲得各個目標連通域內(nèi)圖像像素強度的最大值并就此得到相應(yīng)的強度補償率。最后,將各個目標連通域內(nèi)圖像的像素統(tǒng)一乘上相應(yīng)的強度補償率,如此便完成了多層墻體成像衰減的補償,能大大提高后續(xù)多視角融合后全景圖像的成像質(zhì)量。

1 成像模型建立

收發(fā)器平行x軸,置于x-y平面。成像區(qū)域可以被理解成M×N數(shù)據(jù)平面(m=1,2,3,…,M和n=1,2,3,…,N),如圖1所示。

圖1 回波模型建立

假設(shè)θ(t)為雷達寬帶傳輸信號,那么成像區(qū)域位于位置(m,n)的像素值可以用信號的傳播時延來表示。在x-z數(shù)據(jù)平面上位于坐標(m,n)目標值ζmn(t)在雷達坐標(p,q)可以表示為

式中,ξ(p,q)為窗函數(shù)(通常是漢明窗),用來抑制旁瓣;τmn(p,q)為時域時延補償。回波信號?(t)是發(fā)射信號的時延乘以傳播過程中的衰減系數(shù)αmn(p,q),即?mn(t)=αmn(p,q)θ(t-τmn(p,q))。令θ(t)=θ(-t),表示與發(fā)射信號匹配的濾波器,那么f mn可表示為

其中的時延補償τmn(p,q)在已知墻體參數(shù)的情況下,可以用固定時延算法、最短距離算法或速度補償算法估計[11]。經(jīng)過上述處理,可以得出區(qū)域的成像:

如此得到的F便是原始輸入圖像矩陣。

2 算法描述

為產(chǎn)生一幅全景布局圖像,現(xiàn)階段一般采用多視角圖像融合的方法。每一幅單視角圖像采用后向投影成像算法獲得,為改善單視角圖像質(zhì)量,抑制雜波,這里采用文獻[11]中提到的基于PCA和模糊邏輯圖像增強方法預(yù)處理圖像。但由信號穿墻傳播衰減引起的多層墻體成像強度的非均勻現(xiàn)象為最終建筑布局的檢測和成像帶來較大影響。對此,為補償多層墻體衰減的非均勻現(xiàn)象,本文提出一種基于霍夫變換門限檢測的補償方法。由于多視角融合算法不是本文討論的重點,不失一般性,這里僅以一幅單視角圖像為例進行討論。假設(shè)雷達在單視角下得到的圖像為{F(x,y)},按幅度歸一化為{I(x,y)}。本文提出的算法包括兩個階段,即檢測階段和補償階段。在檢測階段,算法需要劃分目標墻體區(qū)域和非目標雜波區(qū)域,假設(shè)目標墻體區(qū)域設(shè)為1,非目標雜波區(qū)域設(shè)為0,那么圖像{B(x,y)}可以表示為

此時,問題的關(guān)鍵轉(zhuǎn)變?yōu)槿绾未_立固定門限β,并依此大致劃分目標墻體區(qū)域和非目標雜波區(qū)域。對此本文給出一種基于霍夫變換的門限檢測方法。

2.1 墻體目標檢測階段

在此階段,需要完成固定門限檢測,依此劃分出目標墻體區(qū)域和非目標雜波區(qū)域,因此固定門限的選擇成為了關(guān)鍵。由于多層墻體圖像強度的非均勻,若門限選取過大,則不同墻體無法全部通過門限檢測,出現(xiàn)“漏檢”;若門限選取過小,則部分雜波可能超過檢測門限,出現(xiàn)大量“虛警”。結(jié)合墻體圖像的線性特性,本文提出一種基于霍夫變換的固定門限檢測方法。

為方便說明,這里先介紹霍夫變換的相關(guān)概念。在圖像處理中,霍夫變換是圖像中識別幾何形狀的基本方法之一,有廣泛的應(yīng)用,也有很多改進算法。最基本的霍夫變換是從二值圖像中檢測直線(線段)。

如圖2所示,經(jīng)過霍夫變換,圖像空間中的每個點(x,y)就被映射為一個(r,θ)空間中的正弦曲線。而圖像空間中共線的點所對應(yīng)(r,θ)空間中正弦曲線的交點(r′,θ′)。

圖2 運用霍夫變換直線檢測

于是把在圖像空間中檢測直線的問題轉(zhuǎn)化為在極坐標參數(shù)空間中尋找通過點(r,θ)的最多的正弦曲線數(shù)的問題。

在穿墻實際場景成像中,目標墻體在圖像上近似反映為一系列線段,在圖像上便可用霍夫變換檢測之。但是由于雜波噪聲的存在,在BP成像后的原始圖像中,墻體有可能被噪聲和雜波淹沒,對后續(xù)的墻體檢測帶來困難。針對此情況,這里采用文獻[11]中提到的基于PCA和模糊邏輯圖像增強方法對原始圖像進行預(yù)處理,便于后續(xù)的固定門限檢測。一般情況下,天線合成孔徑路線平行于墻體,且不同墻體之間相互平行,在圖像中反映為平行于X軸的若干條矩形圖案,而矩形可以表征為一系列線段的總和。于是圖像中的目標墻體檢測問題可以轉(zhuǎn)化為直線檢測問題。由于合成孔徑路線平行墻體,且不同墻體相互平行,這在圖像中反映為若干條平行于橫軸的線段。此時傾角θ=90°,于是有

容易知道,理論上,只需對極坐標平面上θ軸坐標為90°附近的點進行反霍夫變換,便能得到不同墻體的真實目標位置。

為了確定合適的門限值,以正確劃分目標墻體區(qū)域和非目標雜波區(qū)域,本文提出基于霍夫變換直線檢測的門限搜索方法。首先,將區(qū)間[βmin,βmax]均勻分成若干等分,βmin代表圖像像素點的像素最小值,βmax代表用Ostu自動門限閾值分割方法獲得的門限值。若以βmin為門限,則所有墻體和雜波能通過門限,產(chǎn)生大量“虛警”;若以βmax為門限,則只有部分墻體能通過門限,產(chǎn)生“漏檢”。因此,最佳門限存在于區(qū)間[βmin,βmax]中。其次,門限值在這些等分值中依次從大到小搜索,直至某一值時,經(jīng)霍夫變換能夠檢測到多層墻體直線,記下此時的門限值。該值即是理論上的最佳固定門限值。在用BP成像和改進的PCA分析方法預(yù)處理得到圖像F后,檢測算法的具體流程如下(這里以兩層墻體成像為例說明):

1)確定最佳門限值的搜索范圍[βmin,βmax]。

2)每搜索一個值βn,便以該值為門限對圖像{F(x,y)}作二值化處理,1為目標區(qū)域,0為非目標區(qū)域,得到圖像{B(x,y)};

3)初始化一塊緩沖區(qū),對應(yīng)于極坐標參數(shù)平面,將其所有數(shù)據(jù)置為0;

4)對于二值化圖像{B(x,y)}每一個值為1的像素點,求出極坐標參數(shù)平面對應(yīng)的曲線,把這些曲線上的所有點的值都加1;

5)對在極坐標參數(shù)平面的傾斜角θ軸上[88°,92°]內(nèi)的所有點作霍夫反變換,得到直角坐標系下的圖像{FH(x,y)};

6)若{FH(x,y)}存在多于一個平行于坐標軸的像素為1的矩形區(qū)域,那么此時的βn即是最佳門限值,此時跳出循環(huán),停止搜索,并記下矩形區(qū)域的位置坐標;否則,令n=n-1,重復(fù)步驟2～6;

7)記下β=βn,β即為最佳門限值。

得到最佳固定門限值β后,便可依此大致劃分目標墻體區(qū)域和非目標雜波區(qū)域,以便后續(xù)強度補償階段的處理。步驟1～7的流程圖如圖3所示。

圖3 固定門限β搜索流程圖

2.2 墻體強度補償階段

本節(jié)介紹的算法是在固定門限β檢測后所進行的墻體強度補償算法。在為原始圖像β給定歸一化判決門限β后,可以得到二值化圖像{B(x,y)}?？梢哉f,根據(jù)二值化圖像{B(x,y)}可以返回到對應(yīng)的原始圖像{F(x,y)},獲得{F(x,y)}上墻體目標像素的位置信息和幅度信息,根據(jù)這些信息進行強度補償。而獲取這些信息,需要引入“數(shù)據(jù)連通域”和“深度優(yōu)先搜索”這兩個概念。

在一個數(shù)據(jù)矩陣里,對于某一個關(guān)注的數(shù)據(jù)a,它的周圍存在8個數(shù)據(jù)a0～7(矩陣邊界數(shù)據(jù)除外),如果這八個數(shù)據(jù)中的某一個a i(i∈0～7)與a的屬性相同,則認為a和a i是連通的。如此類推,這樣的一個以連通的數(shù)據(jù)構(gòu)成的集合,稱為數(shù)據(jù)連通域。在門限檢測后的二值圖像上,通常存在多個這樣的數(shù)據(jù)連通域,這些數(shù)據(jù)連通域記載著目標像素的位置信息,而此時的數(shù)據(jù)連通域包含了全部的墻體目標的位置信息,回溯到原始圖像{F(x,y)}上時,便可只針對數(shù)據(jù)連通域在{F(x,y)}上進行處理。

深度優(yōu)先搜索(DFS)類似于樹的前序遍歷,優(yōu)先往樹狀結(jié)構(gòu)的深度方向搜索這一特性非常適用于圖像數(shù)據(jù)連通域查找這一領(lǐng)域,并且具有邊搜索邊記錄的特點,保證了在搜索的同時,獲得每塊數(shù)據(jù)連通域上對應(yīng)目標的信息。

假設(shè)在對需要處理的單視角圖像上使用DFS得到T塊數(shù)據(jù)連通域,每塊數(shù)據(jù)其幅度最大值為m i,所有數(shù)據(jù)連通域的最大值構(gòu)成最大值集合{m1,m2,m3,…,m T},該集合中最大的元素為M,即,則補償率為

若ratio(i)代表第i個連通域的補償率,那么第i個連通域相應(yīng)原始圖像{F(x,y)}位置坐標為(x,y)的像素點的值在補償后為

該階段的算法步驟如下:

1)對原始圖像{F(x,y)}使用固定門限β,將目標墻體區(qū)域和非目標雜波區(qū)域大致分離,得到二值化圖像{B(x,y)};

2)對{B(x,y)}使用深度優(yōu)先搜索,遍歷整幅圖像后只保留與檢測階段步驟6中得到的矩形區(qū)域位置坐標相同的數(shù)據(jù)連通域,即目標數(shù)據(jù)連通域,并記錄下每塊數(shù)據(jù)連通域的重要參數(shù),如坐標和最大值;

3)利用式(8)對目標數(shù)據(jù)連通域上需要進行幅度補償?shù)南袼剡M行處理。

其中步驟2中深度優(yōu)先搜索的處理流程如圖4所示。

整個流程可以由圖5的框圖表示。

3 實測數(shù)據(jù)結(jié)果驗證

雷達系統(tǒng)發(fā)射信號類型為步進頻,實驗采用“走-停-走”的合成孔徑方式進行場景成像,孔徑合成路徑平行于墻體,實驗場景設(shè)置如圖6所示,這里采用兩層墻體成像為例。

圖4 深度優(yōu)先搜索處理流程

圖5 本文提出的算法流程圖

圖6中黑點表示合成孔徑路線天線的不同位置方位,發(fā)射信號步進頻起始頻點fstart=1 GHz,終止頻點fstop=2 GHz,頻階Δf=2 M Hz,中心頻率fc=1 GHz,單頻持續(xù)時間T0=0.5μs,采用BP成像算法,并對成像結(jié)果采用基于PCA和模糊邏輯圖像增強方法[1]處理圖像,得到原始圖像{F(x,y)},如圖7(a)所示。

由圖7(a)可以看出,由于第一層墻體wall1對信號的衰減,造成第二層墻體wall2的成像強度遠低于第一層墻體,對后續(xù)的建筑布局重構(gòu)帶來困難。圖7(b)是采用本文提出的基于霍夫變換門限檢測的多層墻體成像強度衰減補償方法的處理結(jié)果?？梢钥闯?采用本文提出的算法,不僅完成了對第二層墻體wall2的檢測,也較好地完成了對其的強度補償。

圖6 實驗場景設(shè)置

圖7 處理前后對比圖像

為量化體現(xiàn)本文算法的優(yōu)越性,這里以均方誤差MSE和峰值信噪比PSNR為例說明。均方誤差越小,峰值信噪比越高說明圖像與參考圖像的差異越小,處理效果越好。首先給出背景抽取圖像以作標準參考圖像,如圖8所示。

圖8 背景抽取圖像

記背景抽取圖像矩陣為F bs,則均方誤差MSE和峰值信噪比PSNR可以表示為

表1 圖像處理前后的均方誤差和峰值信噪比

由表1可以看出,運用本算法后圖像的均方誤差降低、峰值信噪比提高,說明圖像運用本算法后更接近真實場景圖像,定量上體現(xiàn)了該算法的優(yōu)越性。

4 結(jié)束語

本文提出一種基于霍夫變換門限檢測的多層墻體成像強度補償算法,用于穿墻建筑布局成像中,該算法在雜波抑制、墻體補償、提高信噪比等方面有很大優(yōu)勢,對改善圖像質(zhì)量方面有顯著效果。該算法分為兩個階段,分別是墻體目標檢測和墻體強度補償。在墻體目標檢測階段,關(guān)鍵是確定合適的固定門限值,以區(qū)分目標墻體和非目標雜波。對此本文采用基于霍夫變換直線檢測的算法搜索最佳門限值。在墻體強度補償階段,根據(jù)目標連通域的參數(shù)特性,獲得相應(yīng)連通域的補償率,并對墻體強度進行補償。實測數(shù)據(jù)處理結(jié)果驗證了本文提出的算法的有效性和可行性。

[1]Yoon Y S,Amin M G.Through-the-Wall Radar Imaging Using Compressive Sensing Along Temporal Frequency Domain[C]∥2010 IEEE International Conference on Acoustics Speech and Signal Processing(ICASSP),Dallas,TX:[s.n.],2010:2806-2809.

[2]Smith G E,Mobasseri B G.Robust Through-the-Wall Radar Image Classification Using a Target-Model Alignment Procedure[J].IEEE Trans on Image Processing,2012,21(2):754-767.

[3]Debes C,Zoubir A M,Amin M G.Enhanced Detection Using Target Polarization Signatures in Throughthe-Wall Radar Imaging[J].IEEE Trans on Geoscience and Remote Sensing,2012,50(5):1968-1979.

[4]Lagunas E,Amin M G,Ahmad F,et al.Sparsity-Based Radar Imaging of Building Structures[C]∥2012 Proceedings of the 20th European Signal Processing Conference(EUSIPCO),Bucharest:[s.n.],2012:864-868.

[5]Lagunas E,Amin M G,Ahmad F,et al.Determining Building Interior Structures Using Compressive Sensing[J].Journal of Electronic Imaging,2013,22(2):021003.[6]Aftanas M,Drutarovsky M.Imaging of the Building Contours with Through the Wall UWB Radar System[J].Radio Engineering,2009,18(3):258-264.

[7]Nguyen L,Ressler M,Sichina J.Sensing Through the Wall Imaging Using the Army Research Lab Ultra-Wideband Synchronous Impulse Reconstruction(UWB SIRE)Radar[C]∥Proceedings of SPIE,Orlando,FL:[s.n.],2008:1-4.

[8]Ao Zhuojun,Kong Lingjiang,Jia Yong,et al.Multi-Channel Target Cooperative Detection in Through-the-Wall-Radar Imaging Based on Random Sparse Array[C]∥CIE International Conference on Radar,Chengdu,China:[s.n.],2011:176-178.

[9]Cui G,Kong L,Yang J.A Back-Projection Algorithm to Stepped-Frequency Synthetic Aperture Through-the-Wall Radar Imaging[C]∥1st Asian and Pacific Conference on Synthetic Aperture Radar,Huangshan:[s.n.],2007:123-126.

[10]趙彧,黃春琳,粟毅,等.超寬帶穿墻探測雷達的反向投影成像算法[J].雷達科學(xué)與技術(shù),2007,5(1):49-54.ZHAO Yu,HUANG Chun-lin,SU Yi,et al.BP Imaging Algorithm for UWB-TWDR[J].Radar Science and Technology,2007,5(1):49-54.(in Chinese)

[11]Riaz M M,Ghafoor A.Principle Component Analysis and Fuzzy Logic Based Through Wall Image Enhancement[C]∥Progress In Electromagnetics Research,[s.l.]:[s.n.],2012:461-478.