基于深度迭代網絡的穿墻雷達成像方法

王玉皞 張玥 周輝林 劉且根 蔡琦

（南昌大學信息工程學院,南昌 330031）

引 言

穿墻雷達成像(through-the-wall radar imaging，TWRI)是一種非接觸和非破壞性電磁感知技術，由于電磁波具有穿透非金屬建筑材料的能力，可實現墻后或封閉環境中隱蔽目標的探測、定位、跟蹤以及成像，廣泛應用于災難救援、安防和城區巷戰等民用和軍事領域[1-3].在TWRI 中，由于墻體反射而產生的強雜波往往會淹沒來自墻后或建筑物內靜止目標的反射信號，影響成像精度和后續目標檢測性能.目前為止已經有許多墻體雜波抑制方法被提出，如空間濾波[4]、子空間分解[5]和稀疏信號表示[6]等.但這些方法大部分采用兩階段策略，分別解決墻體雜波抑制和圖像重建問題，其圖像的重構質量完全取決于墻體雜波濾除的效果.

文獻[7]提出了一種迭代優化方法，將雷達數據分解為對應雜波分量的低秩矩陣和對應目標分量的稀疏矩陣，然后采用迭代優化算法抑制墻體雜波并同時得到目標圖像.與傳統的兩階段TWRI 方法相比，該方法精簡了TWRI 過程，并且更高效地抑制墻體雜波以及重構目標圖像.進一步地，文獻[8]在聯合低秩和稀疏模型上引入了全變分正則化器，將原始全變分先驗推廣為聯合微分稀疏的全變分正則化，以保證目標區域的連續性并有效抑制圖像鬼影.然而，這些基于迭代的成像方法在模型建立和優化求解過程中對于各種正則化參數、閾值的選取沒有具體公式指導，都是采用試錯法選擇最佳參數，且該參數與場景密切相關，存在計算代價高、對初值敏感及重構速度慢等問題.

近年來，深度神經網絡由于具有從數據中學習特征的能力，在信號處理領域及圖像視覺相關領域都取得了極大的突破.在雷達領域，深度學習也受到了廣泛的關注.針對單通道步進頻率連續波穿墻雷達，文獻[9]引入了ResNet 深度卷積神經網絡，解決封閉建筑空間中多個移動人體目標計數的問題.文獻[10]從數據驅動的角度展開了相應的基于深度網絡的探地雷達前向數值模擬研究，然而只是將深度網絡作為給定輸入及輸出條件下數據驅動的函數逼近器，并沒有研究更具有物理意義的求解框架.文獻[11]提出一種去噪自編碼器的深度網絡，在不需要已知墻體參數和建筑物布局條件下，實現穿墻雷達的圖像重建.文獻[12]提出了一種基于全卷積網絡的自適應目標檢測方法，通過提取多尺度特征來實現TWRI 中弱散射目標成像.生成式對抗網絡(generative adversarial nets，GAN)作為無監督學習的新興方法，在TWRI 領域也體現了一定優勢，通過生成器和鑒別器的相互博弈能夠有效地抑制TWRI 中的多徑鬼影[13-14].這些網絡都為純數據驅動，在雜波抑制和成像方法上相較于傳統重建方法具有一定優勢，但也存在依賴大數據集、訓練時間長、泛化能力差等問題.文獻[15]則提出了一種具有低秩投影的稀疏自編碼器，通過解決增廣拉格朗日乘數優化問題來確定網絡的權重，能夠有效抑制墻體雜波并恢復目標信號.

有別于上述純數據驅動和物理模型約束的深度網絡，本文提出一種物理模型驅動的深度迭代網絡，實現穿墻雷達室內目標成像.該網絡結合數據驅動與領域知識，將墻體雜波抑制以及目標圖像重建問題表示為聯合低秩與稀疏分解驅動的正則化優化問題，然后基于輪換迭代的求解策略將該問題分解為兩個子問題，分別求解低秩分量及稀疏分量的子優化問題，并推導其更新公式，將傳統迭代優化算法的每步優化求解過程映射到深度迭代網絡結構中的每個階段，形成融合物理模型驅動的可學習深度迭代網絡求解框架.該方法無需手動調節算法參數和預先選擇正則化函數，與傳統迭代方法相比較，本文提出的深度學習方法在推理階段，重構速度快的同時，能高效抑制墻體雜波并且精確重構墻體后方的目標圖像.

1 聯合低秩與稀疏分解的優化問題

本文采用沖擊脈沖雷達體制進行目標探測，雷達系統的收發一體天線陣列平行于墻體放置，收發一體天線陣元自左向右均勻步進前進，TWRI 系統探測示意圖如圖1 所示.圖1 中x軸方向表示方位向，y軸方向表示距離向，步進方向上的天線陣元移動位置數為M，在每個天線陣元位置上雷達系統收發的頻率個數為N，則系統在第m(m=1,2,···,M)個天線陣元位置和第n(n=1,2,···,N)個頻點測量到的雷達回波信號數據yn,m可以表示為

圖1 TWRI 系統探測示意圖Fig.1 TWRI system detection diagram

式中：σw為墻體的反射系數；fn為第n個工作頻點；τw為第m個天線陣元位置和墻體之間的傳輸時延；P為探測場景中墻后的目標個數；σp為第p個目標的反射系數；τm,p為第m個天線陣元位置和第p個目標之間的傳輸時延.當天線沿水平方向均勻步進時，天線與墻體的距離保持不變，墻體回波僅與收發天線的相對位置有關，與其絕對位置無關[16].

為了進行圖像重建，將待成像區域沿距離向和方位向劃分成Q=I×J個均勻空間網格，其中I和J分別代表距離向和方位向的離散網格數，每個離散網格數都代表一個成像像素點.則待重建圖像經過按列堆疊轉換為Q×1維 場景反射系數向量s，其第q(q=1,2,···,Q)個 元素sq表示第p個目標占用的第q個像素點的反射系數 σp.

在每個天線陣元測量位置得到的N個工作頻點的信號和vn,m可以分別堆疊成N×1維向量和vm. 其中，第m個天線陣元測量位置上的目標信號可以表示為

式中，Ψm為N×Q維矩陣，其第n行第q列的元素可以表示為

τm,q表示第m個天線陣元位置和第q個成像網格之間的傳輸時延.

圖2 給出了電磁波在第m個天線位置和第q個成像網格之間的傳播路徑.Tx 表示發射機，Rx 表示接收機.lm,q,air2,t，lm,q,wall,t，lm,q,air1,t分別表示電磁波在發射天線和第q個成像網格之間經過各介質的距離；而lm,q,air2,r，lm,q,wall,r，lm,q,air1,r則分別表示電磁波在接收天線和第q個成像網格之間經過各介質的距離.因此，第m個天線陣元位置和第q個成像網格之間的傳播時延 τm,q表示為[17]

圖2 電磁波傳播路徑Fig.2 Electromagnetic wave propagation path

式中：c為電磁波在自由空間的傳播速度；v=為電磁波在墻體內的傳播速度，εw為墻體的相對介電常數.

將M個天線陣元位置的目標信號進行堆疊可得到總的目標回波信號yt為

再將M個天線陣元位置的N×1維 向量ym、和vm堆疊成MN×1維向量，可得出式(1)的向量形式為

在穿墻雷達探測場景中墻體位置相對固定，每個天線位置上的墻雜波相關性較強，因此yw可認為是一個包含墻體回波信號的低秩分量.

基于壓縮感知(compressive sensing，CS)理論[18]，通過測量矩陣 Φ來表示壓縮測量向量B和整個信號向量y之間的關系：

式中：Φ為K×MN維的測量矩陣，由從MN×MN的單位矩陣中隨機選取K行組成.將式(7)代入式(8)可得：

式中：Ω=ΦΨ；n=Φv.因此，根據式(9)及各向量的特性，雷達探測目標圖像重建問題可以通過建模一個低秩聯合稀疏的最小化問題來求解：

式中：‖·‖為核范數，表示yw的奇異值和；‖·‖1,2為混合的l1,2范數，其值等于s每一行的l2范數和；λ1和λ2分別對應低秩分量yw和稀疏分量s的正則化參數.

2 基于可迭代的深度網絡方法

2.1 輪換迭代求解策略

對式(10)的問題模型采用輪換迭代求解策略，將其轉化為兩個優化子問題，分別為低秩分量的最小化問題求解和稀疏分量的最小化問題求解：

對于子問題(11)，可采用奇異值閾值 (singular value thresholding,SVT) 法求解低秩分量的近端映射[19]：

對于子問題(12)，可采用軟閾值法求解稀疏分量的近端映射：

在yw和s的傳統迭代求解過程中，迭代收斂的條件取決于選取合適的正則化參數 λ1和λ2，而 λ1和λ2通常需要不斷試錯調整后手動設定，且奇異值閾值和軟閾值操作的閾值也需提前給定，這將導致迭代計算過程繁瑣且難以得到最佳取值.因此，為了克服上述不足，本文采用深度學習的框架來展開優化該迭代算法，上述參數在網絡訓練中將作為可學習參數不斷更新以達到最優值.

2.2 可學習迭代展開網絡結構

在深度學習的框架下，展開一個迭代算法能夠加速算法的收斂，而迭代算法提供了一種自然的遞歸架構，可解決特定的問題.在深度學習中可將迭代算法的每一次迭代看作深度學習網絡中的一層，然后將這些層連接起來，通過訓練這種網絡架構來提高迭代的收斂性，顯著減少迭代次數[21].因此，本文基于深度網絡架構將2.1 節中的迭代算法展開，提出相應的遞歸神經網絡模型，第k次迭代對應前饋網絡的第k層.在神經網絡框架中，使用卷積層來代替矩陣乘法.因此，根據式(14)和(16)來形成一個遞歸神經網絡，將依賴 Φ和 Ω的六個矩陣均替換為適當大小的二維卷積核···,由此可以得出多層前饋網絡中第k+1層的方程：

式中，*表示一個卷積算子.在網絡訓練過程中，這些卷積核是可學習的.利用卷積層空間不變性的優勢，可顯著減少學習參數的規模，從而加速訓練過程.

圖3 給出了展開算法的網絡架構圖.在網絡訓練過程中，將獲取到的雷達回波數據經過測量矩陣 Φ壓縮采樣后，得到CS 測量向量B作為網絡的輸入，同時初始化=0、s0=0后進行可學習迭代網絡重建，網絡中的每一層對應一次迭代求解，經過第K層后得到最終的低秩分量輸出和稀疏分量輸出sK，而稀疏輸出sK即為雷達探測目標重建圖像.

圖3 展開算法的網絡架構Fig.3 Network architecture of the untoued algorithm

圖4 為單層網絡的架構圖，每層網絡都包含卷積層和激活層，卷積層中使用六個卷積核(,···,)與卷積層的輸入進行交替卷積運算，卷積層輸出xy對應式(17)中的對應式(18)中的在卷積層后連接激活層，激活層中進行奇異值和閾值操作，獨立學習每層的閾值系數，分別得出對應低秩分量和稀疏分量的迭代解.

圖4 單層網絡架構圖Fig.4 Architecture of single-layer network

在網絡訓練過程中，第k層奇異值閾值的實際閾值由給出，軟閾值操作的實際閾值由·αs·max(sk)給 出，其中：σ(x)=1/(1+e-x)為 sigmoid 函數；αy和αs為給定的固定標量，αy=0.4，αs=1.8；而每層的和在網絡訓練中是可學習的.將目標真值標簽圖像G作為網絡的標簽數據，G可通過模擬被探測目標的分布得出.整個訓練網絡采用反向傳播的監督方式，損失函數定義為網絡預測的圖像sK與標簽圖像G之間的均方誤差之和：

式中：T為數據集的數量；sK(B,β)為具有可學習參數的網絡稀疏輸出，K為網絡層數.

本網絡基于PyTorch 神經網絡框架實現，Pytorch通過自動微分功能執行反向傳播.整個網絡使用自適應矩估計(adaptive moment estimation,ADAM)優化器進行訓練，學習率設置為0.001，訓練數據集為電磁波傳播仿真軟件gprMax 仿真后得到的200 組不同場景下的雷達回波數據.整個網絡共有12 層，前3 層使用的卷積核大小為5×5×1、步長為(1，1，1)、填充為(2，2，0)和一個偏置項，而后9 層使用的卷積核大小為3×3×1、步長為(1，1，1)、填充為(1，1，0)和一個偏置項.表1 為展開迭代深度網絡算法流程圖.

表1 展開迭代深度網絡算法流程圖Tab.1 Flowchart of the expand iterative deep network algorithm

3 實驗結果分析

為驗證本文方法的有效性及準確性，采用基于時域有限差分法的軟件gprMax 進行雷達探測場景的仿真.在本文實驗中，gprMax 的激勵源為中心頻率1 GHz 的Ricker 波，收發一體天線陣元在掃描過程中移動間隔為0.02 m，測線位置上陣元數M=43，每個測量位置的頻點數N=64.成像區域大小設置為1 m×1 m，并將成像區域劃分為51×51 個成像網格.雷達回波數據中加入信噪比(signal to noise ration,SNR)為20 dB 的復加性高斯白噪聲，再通過測量矩陣 Φ來進行降采樣，從64×43=2 752 個頻域樣本數據中隨機選擇1 200 個數據進行成像.圖5 探測場景中探測目標P=1，墻體設置為相對介電常數 εw=6的單質混凝土，厚度為0.1 m，墻后目標為半徑0.13 m 圓柱形物體，其中心位置位于墻體后方0.54 m 處.圖6探測場景中探測目標P=2，墻體設置為相對介電常數 εw=6的單質混凝土，厚度為0.1 m，墻后目標分別是兩個半徑為0.12 m 和0.1 m 的圓柱形物體，其中心位置分別位于墻體后方的0.54 m 及0.3 m 處，兩個目標物體的中心位置相距0.56 m.兩圖中白色實線標注的為目標所在區域.

圖5 單目標仿真及成像結果對比Fig.5 Single target simulation and imaging results comparison

圖6 雙目標仿真及成像結果對比Fig.6 Double target simulation and imaging results comparison

為研究本文方法的可靠性和準確性，采用經典的后向投影(back projection，BP)算法、CS 算法、低秩聯合稀疏迭代重建的方法和本文方法進行成像比較，結果如圖5(c)～(f)及圖6(c)～(f)所示.從兩種場景的成像結果都可以看出：BP 成像的結果存在較多的虛像，成像精度不夠高；CS 成像和低秩聯合稀疏迭代重建的結果能夠較為準確地實現墻后目標的定位，但是成像的輪廓與真值標簽圖仍存在一定差距；而本文算法不僅能夠準確實現目標的中心定位，同時也增強了目標的輪廓特性，更加接近目標重建真值圖，明顯優于前述成像結果.

為定量比較以上4 種算法的成像重建性能，采用目標雜波比(target-to-clutter ratio,TCR)、峰值SNR(peak SNR，PSNR)和程序運行時間(program running time,PRT)作為定量評價標準.TCR 用來評估圖像的聚焦度，TCR 值越高，圖像的聚焦度越高；PSNR 用來評估圖像的失真度，PSNR 值越高，圖像的失真越少；PRT 用來評估算法運行時間，PRT 值越小，PRT 越短.

TCR(單位dB)定義為

式中：Nt和Nc分別表示成像結果中的目標區域和雜波區域所對應的像素個數；At和Ac分別表示成像結果中的目標區域和雜波區域；Iq表示第q個像素點的像素值.

PSNR(單位dB)定義為

式中：G(i,j)表 示真值標簽圖像；H(i,j)表示重建目標圖像.

表2 及表3 分別給出了單目標場景和雙目標場景下4 種成像算法成像重建性能對比結果.反復進行100 次實驗后求平均值，得出4 種算法的TCR值、PSNR 值和PRT 值.可以看出，低秩聯合稀疏成像的TCR 值較高，但其PSNR 值較低，PRT 值最高，說明該方法圖像存在失真且耗時長.本文方法的TCR 值、PSNR 值和PRT 值都為最佳，說明本文所提方法的成像質量最高，且更能滿足實時性要求.

表2 單目標重建性能對比Tab.2 Comparison of single target reconstruction performance

表3 雙目標重建性能對比Tab.3 Comparison of double target reconstruction performance

圖7 分別給出了雷達回波數據中加入SNR 為0～60 dB 的復加性高斯白噪聲下重建目標圖像的TCR 值及PSNR 值.可以看出，不論是單目標場景還是雙目標場景，本文所提出的方法在低SNR 的情況下，重構圖像依然具有較高的TCR 值和PSNR 值，表明能夠有效地抑制墻體雜波，重構出高精度的目標圖像，且在SNR 大于20 dB 時，本文方法的圖像重建性能接近穩定.

圖7 不同噪聲水平下的TCR 和PSNRFig.7 TCR and PSNR at different noise levels

4 結 論

本文提出了一種基于深度迭代網絡的TWRI 方法，該方法將成像問題轉化為聯合低秩與稀疏分解驅動的優化問題，并將求解算法展開到深度網絡結構中，替代原有的迭代求解方法形成了物理模型驅動的深度網絡求解框架.實驗結果分析表明，該方法能夠實現比其他方法更好的性能，兼顧了圖像重建精度和計算復雜度.在后續工作中，將繼續探究該框架下的雷達成像算法，致力于將該方法更加完善地應用于更真實、更復雜的成像場景中.