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

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

0 引言

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

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

1 系統(tǒng)模型

1.1 目標(biāo)運(yùn)動模型和量測模型

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

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

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

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

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

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

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

k-1 時刻多目標(biāo)狀態(tài)隨機(jī)有限集Xk-1中的目標(biāo)增廣狀態(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 時刻運(yùn)動狀態(tài)和模型標(biāo)簽，N為所有模型標(biāo)簽組成的離散集。模型轉(zhuǎn)移概率是服從離散馬爾可夫鏈，且存活目標(biāo)狀態(tài)從k-1 時刻轉(zhuǎn)移到k時刻的轉(zhuǎn)移密度為:所以對于增廣狀態(tài)向量(xks,rks)的轉(zhuǎn)移概率可表示為:

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

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

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

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

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

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

2.1 概率假設(shè)密度濾波

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

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

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

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

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

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

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

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

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

式中:

式中:

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

馬爾可夫轉(zhuǎn)移概率密度為:

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

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

得到了更新權(quán)值后，后驗密度可表示為:

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

3 仿真分析

本章給出多機(jī)動目標(biāo)TBD 仿真實驗，仿真中目標(biāo)在3 個模型中相互切換，分別為直線勻速，左轉(zhuǎn)彎運(yùn)動和右轉(zhuǎn)彎運(yùn)動。目標(biāo)的運(yùn)動模型為:

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

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

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

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

將PHD 濾波用于多個弱小目標(biāo)的檢測前跟蹤中有一個前提條件，即圖像中的每一個單元不能同時被多個目標(biāo)影響，也就是說多個目標(biāo)之間不能彼此靠近。假設(shè)共有 3 個目標(biāo)，初始狀態(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為真實目標(biāo)的運(yùn)行軌跡圖。圖2為目標(biāo)位置狀態(tài)的估計圖。由圖2可知，估計的目標(biāo)狀態(tài)都能近似的分布在真實目標(biāo)位置附近。

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

圖1 目標(biāo)真實軌跡Fig.1 Real trajectory of the target

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

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

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

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

4 結(jié)論

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