基于Rao-Blackwellized蒙特卡羅數據關聯的檢測跟蹤聯合優化

陳唯實, 閆軍, 李敬

(中國民航科學技術研究院 機場研究所, 北京 100028)

“低慢小”目標的探測是一個世界性難題，困難主要源于以下3個方面[1-2]：①飛行高度“低”，一般在1 000 m以下，地面雷達觀測時會有大量雜波進入接收機，尤其在建筑密集的城市環境中，探測難度更大。②飛行速度“慢”，一般小于200 km/h，其回波信號處于雜波主瓣區，容易與鳥群等低速雜波混淆。為防止過多虛警，雷達往往設置速度閾值，導致低速的低空小目標被濾除。③雷達散射截面積(RCS)“小”，目標信噪比低，雷達反射信號忽隱忽現，目標不易被發現與識別。隨著無人機技術的日漸成熟，其作為一種典型的“低慢小”目標，嚴重威脅航班安全運行[3-4]。

雷達系統通常包括檢測和跟蹤2個部分。其中，檢測部分將目標從噪聲和雜波環境中提取出來，獲得目標的速度、位置等狀態信息；跟蹤部分對目標下一刻的狀態進行預測，進而形成穩定航跡。在常規雷達系統中，檢測與跟蹤部分通常分別考慮和設計[5-6]。實際上，二者是相互影響的。檢測是跟蹤的前提，良好的檢測是保證跟蹤性能的基礎；同時，跟蹤能進一步完善檢測結果，利用跟蹤獲得的目標動態特征能夠提高檢測能力[7-8]。因此，本文提出一種基于Rao-Blackwellized蒙特卡羅數據關聯的檢測跟蹤聯合優化算法，對雷達目標檢測與跟蹤進行交互處理，將跟蹤結果反饋給檢測部分，提高雷達系統對低空小目標的檢測能力，同時保證下一步的精確跟蹤。

1 檢測跟蹤聯合優化技術

雷達系統在完成目標檢測之后，通常將波門內超過一定閾值的量測數據通過跟蹤算法進行關聯[9]。如將跟蹤結果反饋至檢測中心，進而根據該反饋信息調整波門內的檢測閾值，改善檢測結果，此類算法即可稱為檢測跟蹤聯合優化算法。王云奇[10]結合目標、環境等先驗知識的輔助，量化分析檢測性能與跟蹤性能的關系，實現檢測和跟蹤的耦合。劉紅亮、嚴俊坤等[11-12]提出一種航跡恒虛警的目標檢測跟蹤一體化算法，根據幀虛警概率調整預測波門內的檢測門限，完成目標檢測與跟蹤的聯合優化過程。閆學昭[13]采用DSP等硬件編程方法，搭建了雷達目標檢測與跟蹤交互處理模塊，將跟蹤結果反饋至檢測系統，在提高檢測準確率的前提下實現了精確跟蹤。

綜上所述，這些算法的共同點在于基于一定規則降低波門內的檢測閾值，以犧牲波門內局部虛警率為代價提升局部目標的檢測率[11]。可見，波門的設置至關重要，如果波門面積(體積)過小，則可能將檢測目標排除在外；如果波門面積(體積)過大，則可能引入過多虛警，導致計算機過載。針對此缺點，本文提出一種基于Rao-Blackwellized蒙特卡羅數據關聯的檢測跟蹤聯合優化算法。該算法采用序貫蒙特卡羅(Sequential Monte Carlo,SMC)(又稱為粒子濾波)方法實現多目標跟蹤，由于粒子能夠用于表示多種不同的數據關聯假設，基于SMC方法的多目標跟蹤算法可視為多假設跟蹤(Multiple Hypothesis Tracking, MHT)的推廣[14]。此外，Rao-Blackwellizion方法的應用能提高算法的準確率和效率，采用卡爾曼濾波或擴展卡爾曼濾波等方法估計目標狀態，將SMC方法僅用于數據關聯估計，使得聯合后驗分布由混合高斯表示，比純粹的SMC方法具有更小的方差。同時，本文算法能夠根據粒子的分布范圍確定波門大小，在考慮粒子權重的前提下，利用檢測單元與所有粒子的相對位置對檢測門限進行修正。

2 基于SMC的檢測跟蹤聯合優化算法

本文提出一種基于Rao-Blackwellized蒙特卡羅數據關聯的檢測跟蹤聯合優化算法，包括檢測和跟蹤2個部分，將粒子濾波跟蹤獲取的目標狀態信息反饋至檢測部分，以修正相關檢測單元的檢測門限，實現檢測與跟蹤的聯合優化，算法流程如圖1所示，詳細步驟如下：

圖1 目標檢測跟蹤聯合優化算法流程圖Fig.1 Flowchart of joint optimization algorithm for target detection and tracking

步驟1粒子濾波目標狀態預估。

(1)

計算每個粒子新的非歸一化權重為

(2)

進行權重歸一化處理：

(3)

步驟2檢測單元相似度計算。

基于高斯概率分布函數計算檢測單元uk與每個粒子預估位置的空間相似度：

(4)

式中：Hk和Rk分別為k時刻的量測模型矩陣和量測噪聲矩陣；N(·)表示檢測單元與粒子預估位置的相似度。

結合每個粒子的權重，得到每個檢測單元與所有粒子群預估位置的空間相似度為

(5)

步驟3檢測門限修正。

利用步驟2計算所得的每個檢測單元與所有粒子群預估位置的空間相似度，對每個檢測單元的固定檢測門限θ進行修正[15-16]，修正后的檢測門限表示為

t(uk)=θe-γp(uk)

(6)

式中：γ為門限修正強度參數。

該模型的合理性在于：對于與粒子預估位置相距較遠的檢測單元，其為目標量測的概率較小，則p(uk)值趨近于0，e-γp(uk)值趨近于1，以至于不對閾值進行修正；對于與粒子預估位置相距較近的檢測單元，其為目標量測的概率增大，則p(uk)值增大，且0

步驟4目標檢測。

根據步驟3計算所得的檢測門限，由式(7)判斷檢測單元中是否存在目標：

(7)

式中：rcs(uk)為檢測單元uk處的量測回波強度值，若該值大于或等于檢測門限，則數據確認為目標(E(uk)=1)，由yk表示，反之則無目標(E(uk)=0)。

步驟5數據關聯。

假設馬爾可夫模型為m次，則數據關聯指數ck取決于先前的m個關聯結果ck-m:k-1，假設虛警在量測空間V中均勻分布，卡爾曼濾波的量測相似度計算如下：

(8)

式中：j=1,2，…,T，T為目標數目；KFlh(·)代表卡爾曼濾波器測量相似度估計；Hj,k和Rj,k分別為目標j的量測模型矩陣和量測噪聲矩陣。

對于j=1,2，…,T，有

(9)

數據關聯結果由最優重要性分布決定，其概率分布結果計算如下：

1) 計算非歸一化噪聲關聯概率：

(10)

2) 為每個目標j=1,2，…,T計算非歸一化目標關聯概率：

(11)

3) 歸一化重要性分布：

(12)

步驟6粒子濾波目標狀態修正。

(13)

重新計算粒子權重并估計有效粒子數：

(14)

如有效粒子數過低(如neff

3 實驗結果分析

本節分別針對仿真數據和雷達實測數據，評價本文算法的有效性，評價指標包括目標數Nd、虛警數Nfa、檢測率Pd、虛警率Pfa、工作特征(Receiver Operator Characteristic, ROC)曲線、均方根誤差(Root-Mean-Square Error, RMSE)等。

3.1 仿真數據

(15)

模擬地雜波的瑞利分布概率密度函數為

(16)

式中：b為瑞利系數。

以離散維納過程速度模型建立雜波環境中的單目標運動模型，其中目標的狀態可以寫為

(17)

離散動態可以表示為線性、時不變構造方程：

(18)

式中：qk-1為離散高斯過程白噪聲，矩特征為

(19)

其中：時間步長設定為Δt=0.1；過程噪聲的功率譜密度設定為q=0.1。

為構造雜波測量環境，如果測量值為噪聲，則數據關聯指數ck設定為0，如果為真實目標則設定為1。雜波測量值均勻地分布在空間[-5,5]×[-4,4]中。真實目標的測量模型與附加的高斯噪聲成線性關系。因此，可以將聯合測量相似度表示為

(20)

量測噪聲矩陣為

(21)

數據關聯的先驗值完全獨立，寫為

(22)

說明量測為雜波的概率是pc，量測為目標的概率是1-pc。pc代表了雜波在量測中所占的比例，該值越大，跟蹤的難度也越大。

如圖2所示，經過240步仿真，目標從(-4, -0.2)出發，在0～1.0 s以速度(1, 0)勻速運動，在1.1～3.0 s完成右轉彎，在3.1～3.5 s以速度(0, -1)勻速運動，在3.6～5.5 s完成左轉彎，在5.6～8.0 s以速度(1, 0)勻速運動，在8.1～10.5 s完成左轉彎，在10.6～14.0 s以速度(0, 1)勻速運動，在14.1～16.5 s完成左轉彎，在16.6～19.0 s以速度(-1, 0)勻速運動，在19.1～21.0 s完成左轉彎，在21.1～21.5 s以速度(0, -1)勻速運動，在21.6～23.5 s完成右轉彎，在23.6～24.0 s以速度(-1, 0)勻速運動，直至結束。雜波概率為0.3。針對4類目標RCS模型，按照圖2中的雷達目標運動軌跡，經過1000次蒙特卡羅仿真并取平均值，圖3給出不同雜波概率條件下(pc=0.1, 0.3, 0.5)的ROC曲線和RMSE曲線，對比了采用本文算法對檢測門限進行修正前后的跟蹤結果，粒子數為10。

圖2 雷達目標仿真運動軌跡Fig.2 Simulation of radar target motion trajectory

圖3 雷達目標仿真運動跟蹤效果Fig.3 Simulation results of radar target motion tracking

由圖3中4類模型的ROC曲線可見，本文算法均明顯改善了4類目標模型的檢測效果，其中，模型1的檢測效果稍差，模型2～模型4的檢測效果逐漸改善，這是由目標的散射特性決定的。具體而言，模型1的典型目標為前向觀察的小型噴氣飛機，模型2的典型目標為螺旋槳推進飛機、直升機等，模型3的典型目標為噴氣飛機、大型民用客機，模型4的典型目標為側向觀測的導彈與高速飛行體等。顯然，模型1的檢測難度最大，模型2～模型4的檢測難度逐漸降低。由于目標檢測閾值在修正前后均采用了數據關聯算法剔除雜波，pc值對ROC曲線的影響不大。由圖3中4類模型的RMSE曲線可見，pc值對RMSE曲線的作用明顯，當pc值提高時，RMSE明顯提高，跟蹤效果變差。同時發現，當雜波比例較低時(pc=0.1, 0.3)，量測中的虛警較少，跟蹤精度(RMSE)較高，此時，本文算法利用跟蹤結果對檢測門限進行反饋修正，能夠進一步提高跟蹤精度。但是，當雜波比例較高時(pc=0.5)，雜波和虛警在量測中的比例達到50%，此時由于目標量測較少，導致跟蹤精度(RMSE)降低。由實驗結果可見，如果閾值設置過低(θ<1.5)，檢測門限經過修正后的跟蹤結果反而可能進一步惡化，這可能是由于過低的閾值引入部分虛警導致的。因此，當虛警比例較高時，建議適當提高閾值(θ>1.5)，采用檢測門限反饋修正算法能夠在一定程度上改善跟蹤結果，確保算法收斂。

針對4類目標RCS模型，表1～表4給出了粒子數N為10、50和100的情況下，本文算法給定不同分割閾值時的RMSE值，雜波比例均設定為pc=0.1。通過對以上4組數據的分析，粒子數對檢測結果影響不大，粒子數較多時的檢測結果有時甚至稍遜于粒子數較少的時候。其原因在于本文算法采用Rao-Blackwellization方法將單目標跟蹤與數據關聯分開處理，將SMC(粒子濾波)方法用于數據關聯，采用較少的粒子數，實現雜波與虛警量測中的多目標跟蹤，可視為MHT的推廣。Rao-Blackwellized粒子濾波的理論基礎是某些濾波方程可以閉合的形式計算，其他采用蒙特卡羅采樣，而不是對所有方程都采用采樣方法。Rao-Blackwell的思想實現了較小的估計方差，其可視為用無窮集合去替代有限集合，往往能得到更為準確的結果。

表1 不同粒子數情況下模型1的RMSE值

表2 不同粒子數情況下模型2的RMSE值

表3 不同粒子數情況下模型3的RMSE值

表4 不同粒子數情況下模型4的RMSE值

3.2 實測數據

基于一組S波段非相參雷達數據和一組S波段相參雷達數據，在前期研究成果的基礎上，采用本文算法進行處理，并與現有經典算法進行對比分析，2組數據的采樣周期均為2.5 s。第1組測試數據為S波段非相參雷達采集的圖像序列，圖中目標為一架波音737客機，分辨率1 024×1 024，量程22 km，共120幀。圖4給出了某幀圖像的檢測跟蹤結果。其中，圖4(a)在原始圖像中標定了目標位置，分別采用固定閾值(Fixed Threshold,FT)、單元平均恒虛警(Cell Average-ConstantFalseAlarmRate,CA-CFAR)和本文算法進行目標檢測；圖4(b)為目標跟蹤結果。實驗結果表明，目標檢測閾值設置為θ=70時，FT算法在引入228個虛警的前提下能夠檢測到目標，CA-CFAR算法不僅能檢測到目標，且將虛警數降為25個，基于前期研究成果[18]，結合本文算法，粒子數為10，能夠檢測到小弱目標，且將虛警數減少到2個。圖4(b)為檢測閾值設定為θ=120時，本文算法的目標跟蹤結果，此時120幀圖像的目標檢測累積之和為Nd=119，Nfa=3，所有信息疊加在雷達背景圖像上，3個虛警位置疑為未形成軌跡的低空未知目標。

圖4 S波段非相參雷達目標檢測與跟蹤結果Fig.4 Target detection and tracking results of S-band incoherent radar

為比較各類算法的魯棒性，表5給出了設置不同分割閾值的情況下，各類算法處理120幀雷達圖像獲取的目標數和虛警數總和，包括FT、CA-CFAR、最小選擇恒虛警(Small of-Constant False Alarm Rate,SO-CFAR)和本文算法。形態學處理(Morphological Processing, MP)是一種經典的圖像處理算法，對于剔除雷達圖像中由單個像素組成的雜波信息，效果良好[19]。此處，將MP算法作為所有算法的后處理，將“腐蝕”和“膨脹”處理相結合，既能剔除虛警，又能將屬于同一目標的多個量測區域重新聯通，避免將目標信息誤認為虛警。由表5可見，當檢測到的目標數為114時，FT的虛警數為3 539，CA-CFAR的虛警數為2 469，而SO-CFAR的虛警數為3 541。本文算法在設定不同閾值時，幾乎能檢測到所有目標，虛警數最小值僅為1。總體而言，本文算法最優，CA-CFAR算法優于FT算法，SO-CFAR算法甚至略遜于FT算法，其原因在于SO-CFAR算法通常將局部檢測閾值設置過低，容易導致較高虛警率。

第2組測試數據在某雷達測試外場，由S波段相參雷達采集的輕小型無人機目標數據，該無人機目標沿正北偏西方向逐漸遠離雷達飛行，量程6 km。門限θ為目標RCS值，典型無人機目標的RCS值很小，屬于小弱目標。

表5 S波段非相參雷達目標檢測結果對比

圖5中，采樣周期為360個。表6對比了設定不同θ值時的檢測結果。本文算法將粒子數設定為10，能夠在檢測到339個無人機目標的情況下將虛警數降低到17。相參雷達已經剔除了靜止的背景地物回波，且對虛警具有一定的抑制作用，本文算法的作用主要在于跟蹤目標的同時提高對小弱目標的檢測能力。

圖5 S波段相參雷達無人機目標檢測與跟蹤結果Fig.5 UAV target detection and tracking results of S-band coherent radar

θFT算法本文算法NdNfaNdNfa0．0132721339210．0227321339210．032011733917

4 結 論

本文結合SMC方法計算的粒子范圍確定波門大小，在考慮粒子權重的前提下，利用檢測單元與所有粒子的相對位置對檢測門限進行修正，提高對小弱目標的檢測能力。通過將本文算法與其他經典算法進行對比，得出以下結論：

1) 由仿真結果可見，本文算法對于各種類型雷達目標的檢測結果均有明顯改善，在雜波比例較低時，對跟蹤精度也有一定提高。

2) 本文算法分別對每個目標進行卡爾曼濾波預估、更新和測量相似度估計，比純粹的粒子濾波方法效率更高，以較少的粒子數就能實現對目標的精確跟蹤。

3) 將本文算法與前期研究成果相結合，并應用于S波段非相參雷達采集的圖像序列，檢測結果表明，本文算法能夠改善局部區域內的目標檢測結果，優于其他經典算法。

4) 將本文算法應用于S波段相參雷達數據，能夠提高對低空無人機目標的檢測跟蹤能力。

參考文獻 (References)

[1] 韓崇昭,朱洪艷,段戰勝,等.多源信息融合[M].北京:清華大學出版社,2006:320-365.

HAN C Z,ZHU H Y,DUAN Z S,et al.Multi-source information fusion[M].Beijing:Tsinghua University Press,2006:320-365(in Chinese).

[2] ASLAN M S,SARANL A.Threshold optimization for tracking a nonmaneuvering target[J].IEEE Transactions on Aerospace and Electronic Systems,2011,37(2):2844-2859.

[3] 呂信明.軍用無人機的發展及對策[J].國防科技,2013,34(1):5-7.

LV X M.Military UAV development and countermeasures[J].National Defense Science & Technology,2013,34(1):5-7(in Chinese).

[4] 何煦.對低空小型無人機的對抗方法研究[J].中國新通信,2016,18(15):147-149.

HE X.Research on countermeasures against low-altitude small UAV[J].China New Telecommunications,2016,18(15):147-149(in Chinese).

[5] RICHARDS M A.Fundamentals of radar signal processing[M].New York:McGraw-Hill,2005:1-15.

[6] SKOLNIK M I.Introduction to radar systems[M].New York:McGraw-Hill,2002:5-10.

[7] BAR-SHALOM Y,LI X R,KIRUBARAJAN T.Estimation with application to tracking and navigation:Theory algorithms and software[M].New York:John Wiley & Sons,Inc.,2001:35-42.

[8] LI X R,JILKOV V P.Survey of maneuvering target tracking,Part V:Multiple-model methods[J].IEEE Transactions on Aerospace and Electronic Systems,2005,41(4):1255-1321.

[9] BLACKMAN S S.Multiple hypothesis tracking for multiple target tracking[J].IEEE Aerospace and Electronic Systems Magazine,2004,19(1-2):5-18.

[10] 王云奇.知識輔助檢測跟蹤一體化算法研究[D].成都:電子科技大學,2014:4-12.

WANG Y Q.Research of detection-tracking-integration algorithm with knowledge-aided[D].Chengdu:University of Electronic Science and Technology of China,2014:4-12(in Chinese).

[11] 劉紅亮,周生華,劉宏偉,等.一種航跡恒虛警的目標檢測跟蹤一體化算法[J].電子與信息學報,2016,38(5):1072-1078.

LIU H L,ZHOU S H,LIU H W,et al.An integrated target detection and tracking algorithm with constant track false alarm rate[J].Journal of Electronics & Information Technology,2016,38(5):1072-1078(in Chinese).

[12] 嚴俊坤,劉紅亮,戴奉周,等.一種具有恒虛警性質的檢測跟蹤聯合處理算法[J].電子與信息學報,2014,36(11):2666-2671.

YAN J K,LIU H L,DAI F Z,et al.Joint detection and tracking processing algorithm with constant false alarm rate property[J].Journal of Electronics & Information Technology,2014,36(11):2666-2671(in Chinese).

[13] 閆學昭.雷達目標檢測與跟蹤交互處理[D].大連:大連海事大學,2011:35-40.

YAN X Z.Radar target detection and tracking interactive processing[D].Dalian:Dalian Maritime University,2011:35-40(in Chinese).

[15] MCHUGH J M,KONRAD J,SALIGRAMA V,et al.Foreground-adaptive background subtraction[J].IEEE Signal Processing Letters,2009,16(5):390-393.

[16] STAUFFER C,GRIMSON E.Learning patterns of activity using real-time tracking[J].IEEE Transactions on Pattern Analysis and Machine Intelligence,2000,22(8):747-757.

[17] 黃培康,殷紅成,許小劍.雷達目標特性[M].北京:電子工業出版社,2005:114.

HUANG P K,YIN H C,XU X J.Radar target characteristics[M].Beijing:Publishing House of Electronics Industry,2005:114(in Chinese).

[18] 陳唯實,李敬.基于空域特性的低空空域雷達目標檢測[J].航空學報,2015,36(9):3060-3068.

CHEN W S,LI J.Radar target detection in low-altitude airspace with spatial features[J].Acta Aeronautica et Astronautica Sinica,2015,36(9):3060-3068(in Chinese).

[19] GONZALEZ R C,WOODS R E,EDDINS S L.Digital image processing using MATLAB[M].Englewood Cliffs:Prentice Hall,2009:261-268.