薛秋條,寧巧嬌,吳孫勇,2,蔡如華,伍雯雯

(1.桂林電子科技大學 數(shù)學與計算科學學院，廣西 桂林 541004;2.廣西密碼學與信息安全重點實驗室，廣西 桂林 541004)

0 引言

多機動目標跟蹤是目標通過在不同的模型之間的轉換，利用模型轉移概率矩陣來實現(xiàn)的一種方法，因為基于單一模型的跟蹤是不能準確預測目標的運動軌跡的[1]。跳躍馬爾可夫系統(tǒng)(jump Markov system，JMS)模型是解決機動目標跟蹤的有效方法[2-4]。JMS 通過使用多個描述目標運行軌跡的模型并行實現(xiàn)機動目標跟蹤。對于多目標的跟蹤，Mahler 利用隨機有限集理論(random finite sets, RFS)，提出了有限集統(tǒng)計學理論(finite set statistics, FISST)[5]，隨后提出了概率假設密度(probability hypothesis density，PHD)濾波[6]。PHD濾波是用一階統(tǒng)計矩近似后驗概率密度函數(shù)，該算法將復雜的多目標狀態(tài)空間運算轉化為單目標狀態(tài)空間運算，有效避免了復雜數(shù)據(jù)聯(lián)合問題。國內(nèi)外學者將PHD、cardinalized probability hypothesis density(CPHD)和Bernoulli 濾波使用點量測與JMS 模型結合，實現(xiàn)對多機動目標的跟蹤[7-10]。

檢測前跟蹤(track before detect，TBD)是解決低信噪比條件下弱小目標檢測與跟蹤的有效方法[11]。TBD 不用事先設定檢測門限對目標有無進行檢測，而是直接利用傳感器量測數(shù)據(jù)進行目標跟蹤，可以最大限度地保留目標的信息。近年來，國內(nèi)外學者將TBD技術與隨機有限集濾波進行結合，有效地實現(xiàn)對多個弱小目標的檢測與跟蹤問題[12-14]。對于機動弱小目標的檢測與跟蹤問題，國內(nèi)外學者采用多模的隨機有限集濾波與TBD 技術進行結合，也能有效地實現(xiàn)[15-16]。本文針對低信噪比條件下機動弱小目標的檢測與跟蹤問題，提出基于JMS-SMC-PHD 濾波的檢測前跟蹤算法。該算法以SMC-PHD 濾波為基礎，在目標跟蹤的預測階段利用跳躍馬爾可夫系統(tǒng)模型實現(xiàn)對目標狀態(tài)的轉移，在權值更新階段直接使用紅外傳感器量測數(shù)據(jù)實現(xiàn)對所有粒子的權值進行更新。最后通過仿真實驗證明，所提算法可以有效地實現(xiàn)多個弱小機動目標的檢測與跟蹤。

1 系統(tǒng)模型

1.1 目標運動模型和量測模型

針對多目標的跟蹤問題，假設在k(k=1,2,…)時刻，在觀測區(qū)域內(nèi)存在Nk個目標，并用一個隨機有限集表示為:

式中:E(χ)為目標狀態(tài)的所有有限子集。目標的狀態(tài)集包含目標的存活、新生和衍生部分，k時刻目標s(k=1,2,…,Nk)的運動方程可表示為:

紅外傳感器以一定的幀頻對監(jiān)視區(qū)域進行監(jiān)測，得到一幅關于監(jiān)測區(qū)域的二維圖像，圖像包括nXm的感應像元每一個分辨單元對應一個矩形區(qū)域ΔxXΔx，且每一分辨單元(i,j)的中心為(iΔx,jΔx)，i=1,…,n，j=1,…,m。則分辨單元(i,j)的觀測強度zk(i,j)可以表示為:

式中:ωk(i,j)為零均值高斯白噪聲，即并且各個單元之間是相互獨立的。C表示目標擴散影響區(qū)域其中為離目標狀態(tài)xks最近的整數(shù)單元。hk(i,j)(xks)為目標xks對分辨單元(i,j)的貢獻強度，一般采用點擴散函數(shù)的形式:

式中:∑為傳感器引入的模糊參數(shù)，一般是已知的。點擴散函數(shù)產(chǎn)生的強度大小以目標xks影響最大單元向周圍單元減小，所以目標的影響區(qū)域是有限的，在最靠近目標的單元取得最大值。從而k時刻的量測可以表示為:以及從1～k時刻的所有量測數(shù)據(jù)可以表示為:Z1:k={zl:l=1,…,k}。

1.2 跳躍馬爾可夫系統(tǒng)模型

跳躍馬爾可夫系統(tǒng)通過一組參數(shù)化的狀態(tài)空間模型去描述，各個狀態(tài)空間模型之間轉移是服從馬爾可夫鏈[5]。JMS 提供了一種不再是通過單一模型的方法去模擬機動目標，是解決機動目標跟蹤的有效方法。

k-1 時刻多目標狀態(tài)隨機有限集Xk-1中的目標增廣狀態(tài)(xk-1s,rk-1s)∈E(χ)=RnXN(s=1,…,Nk-1)要么以概率為Ps,k(xk-1s,rk-1s)存活到下一個時刻，要么以概率為1-Ps,k(xk-1s,rk-1s)死亡，其中xk-1s∈Rn和rk-1s∈N分別表示k-1 時刻運動狀態(tài)和模型標簽，N為所有模型標簽組成的離散集。模型轉移概率是服從離散馬爾可夫鏈，且存活目標狀態(tài)從k-1 時刻轉移到k時刻的轉移密度為:所以對于增廣狀態(tài)向量(xks,rks)的轉移概率可表示為:

紅外傳感器量測的似然函數(shù)可表示為:

式中:H1表示存在Nk個目標，H0表示不存在目標。pN(zk(i,j))表示分辨單元(i,j)處背景噪聲的概率密度函數(shù)，則表示目標加噪聲的概率密度函數(shù)。對于高斯背景噪聲條件下，上式可以表示為:

由于目標影響區(qū)域是有限僅能對其所在單元和鄰近區(qū)域有顯著影響，因此似然函數(shù)可進一步描述為:

2 基于JMS-SMC-PHD 濾波的TBD 算法

根據(jù)貝葉斯理論，濾波的實現(xiàn)可以分為預測階段和更新階段，預測是實現(xiàn)目標狀態(tài)的轉移，更新則是在預測的基礎上利用當前量測來做出修正。根據(jù)目標的狀態(tài)轉移服從一階馬爾科夫過程，所以預測和更新過程可以分別描述為:

式中:fk|k-1(Xk|Xk-1)是狀態(tài)轉移概率密度;gk(Zk|Xk)為似然函數(shù)。最后再利用后驗概率密度去估計目標的狀態(tài)。

2.1 概率假設密度濾波

PHD 濾波是貝葉斯濾波的一種智能近似，通過傳遞后驗概率密度的一階統(tǒng)計矩來代替后驗概率密度函數(shù)，可以有效避免數(shù)據(jù)聯(lián)合問題。

k時刻的預測密度和后驗密度分別用vk|k-1和vk表示，則PHD 濾波的預測可描述為:

式中:fk|k-1(?|x′)為已知目標狀態(tài)x′轉移到k時刻的轉移概率密度;pS,k|k-1(x′)為目標x′存活的概率;βk|k-1(?|x′)為目標x′在k時刻衍生RFS 的強度;γ(?)為新生RFS 的強度。

PHD的更新利用k時刻多目標的量測信息和預測的概率密度，所以后驗密度可描述為:

式中:pD,k(x)為k時刻的檢測概率;Zk為k時刻多目標的量測;κk(?)為雜波RFS 的密度。

2.2 JMS-SMC-PHD 濾波實現(xiàn)TBD

PHD 的序貫蒙特卡羅實現(xiàn)也稱PHD 的粒子濾波實現(xiàn)，是用一組帶權值的粒子去近似估計目標的后驗概率密度。跳躍馬爾可夫系統(tǒng)下的 PHDvk-1(xk-1,rk-1|Z1:k-1)用一組帶權值ωk-1(s)的粒子表示為:

式中:δ(?)為狄拉克函數(shù);rk-1為模型標簽;Lk-1是存活粒子數(shù)。然后通過預測和更新階段實現(xiàn)濾波:

1)預測:包含了新生粒子的產(chǎn)生和持續(xù)存活粒子的傳遞，新生粒子是由k-1 時刻的量測集zk-1={zk(i,j),i=1,…,n,j=1,…m}通過降序排列后，強度較大的前M個分辨單元所產(chǎn)生的Jk個新生的粒子。由于在強度較大的位置是真實目標所在位置的概率比較大，所以這樣產(chǎn)生的粒子更加接近真實目標。假設qk( ?|xk-1,Z1:k)和pk(?|Z1:k)分別為存活目標和新生目標的重要密度函數(shù)。則預測PHD 可近似為:

式中:

式中:

其中，衍生密度函數(shù)為:

馬爾可夫轉移概率密度為:

2)對上述得到的粒子根據(jù)量測進行粒子權值的更新，對于在量測zk中，僅有那些鄰近的分辨單元受到較大影響，因此似然比可計算為:

假定由狀態(tài)xk引起的有效量測為:從而對于狀態(tài)的權值更新為:

得到了更新權值后，后驗密度可表示為:

5)提取目標的狀態(tài):根據(jù)估計目標數(shù)zk利用k_means 方法對重采樣后的粒子進行聚類劃分，從而提取目標的狀態(tài)。

3 仿真分析

本章給出多機動目標TBD 仿真實驗，仿真中目標在3 個模型中相互切換，分別為直線勻速，左轉彎運動和右轉彎運動。目標的運動模型為:

vk為零均值高斯白噪聲，其協(xié)方差為:

式中:q1=2，q2=1 分別表示目標運動和強度的過程噪聲大小。模型r=1 時，為勻速直線運動模型，模型旋轉角度ω=0 rad/s;模型r=2 時，模型旋轉角度ω=10 rad/s;模型r=3 時，模型旋轉角度ω=-10 rad/s。馬爾可夫轉移概率矩陣和新生RFS 模型概率密度分別為:

傳感器以間隔T=1 s 接收30 幀圖像，其他參數(shù)值為:Δx=Δy=1，n=m=20，σ=1，∑=0.7。目標的存活概率為Ps=0.98。采用平面中M=20 個最大強度單元作為新生目標的潛在位置，速度均服從U[-0.4,0.4]的均勻分布，強度值服從U[I-1,I+1]的均勻分布，不考慮衍生情況即bk=[0;0;0]。采用最優(yōu)子模式分配度量(optimal subpattern assignment,OSPA)[17]來檢驗目標的跟蹤精度，OSPA 距是描述集合間的誤差的一種有效方法，其最大的優(yōu)勢是同時考慮了兩個集合間的勢誤差(檢測精度)和位置誤差(跟蹤精度)。OSPA 距的參數(shù)為:階數(shù)p=1，截斷誤差c=10。紅外傳感器強度量測的信噪比為:

式中:σ為量測圖中噪聲強度協(xié)方差。

將PHD 濾波用于多個弱小目標的檢測前跟蹤中有一個前提條件，即圖像中的每一個單元不能同時被多個目標影響，也就是說多個目標之間不能彼此靠近。假設共有 3 個目標，初始狀態(tài)為x1=[5,0.3,1,0.4,I]T，從1～25 s 時刻存活。x2=[19,-0.4,9,0.5,I]T，從3～30 s 時刻存活。x3=[5,0.3,18,-0.4,I]T，從6～27 s 時刻存活，其中I=15。圖1為真實目標的運行軌跡圖。圖2為目標位置狀態(tài)的估計圖。由圖2可知，估計的目標狀態(tài)都能近似的分布在真實目標位置附近。

使用Matlab R2013a 運行100 次蒙特卡羅實驗。圖3為目標狀態(tài)估計的平均OSPA 誤差。本文方法(JMS-PHD-TBD)與目標模型為單一模型的概率假設密度濾波的檢測前跟蹤方法(SING-PHD-TBD)作比較。由圖可知，由于在初始時刻，目標的初始狀態(tài)選取為空集的情況，所以算法在開始還不能估計出目標的狀態(tài)，而在目標出現(xiàn)的時候，估計都出現(xiàn)了一定的延時。但JMS-PHD-TBD 的估計誤差總體上要優(yōu)于SING-PHD-TBD。而單模型估計機動目標也能估計出目標的狀態(tài)是由于新生粒子時采用的是高強度位置產(chǎn)生粒子機制，這樣能保證真實目標在所采樣的粒子中的概率較高。并且由于弱小目標速度比較小，模型的轉換影響會稍微較小，但這也是不容忽視的。

圖1 目標真實軌跡Fig.1 Real trajectory of the target

圖2 目標位置狀態(tài)估計Fig.2 Target location status estimation

圖3 平均OSPA 誤差估計Fig.3 Mean OSPA error estimation

圖4為100 次蒙特卡羅實驗得到的平均目標數(shù)目的估計圖。由圖可知，在新目標出現(xiàn)時，目標數(shù)目的估計出現(xiàn)了延時，這也是導致圖3中估計的OSPA 誤差出現(xiàn)較大峰值的原因。JMS-PHD-TBD 方法對與目標的數(shù)目估計更加準確。出現(xiàn)目標跟蹤丟失的情況比較少，更加能準確地估計目標的狀態(tài)。

圖4 平均目標數(shù)目估計Fig.4 Average target cardinality estimation

4 結論

本文對跳躍馬爾可夫系統(tǒng)模型下的概率假設密度濾波的檢測前跟蹤方法進行深入的分析和研究。在序貫蒙特卡羅概率假設密度濾波的基礎上，融合跳躍馬爾可夫系統(tǒng)，用于解決多機動弱小目標的檢測與跟蹤。仿真實驗表明，通過跳躍馬爾可夫系統(tǒng)對模型進行狀態(tài)轉移并直接對單幀量測數(shù)據(jù)更新權值的處理，可以有效地實現(xiàn)多機動弱小目標的檢測與跟蹤。