侯利明,連峰,王偉

(1.西安交通大學(xué)智能網(wǎng)絡(luò)與網(wǎng)絡(luò)安全教育部重點(diǎn)實(shí)驗(yàn)室,710049,西安;2.陜西省組合與智能導(dǎo)航重點(diǎn)實(shí)驗(yàn)室,710068,西安;3.中國電子科技集團(tuán)公司第二十研究所導(dǎo)航重點(diǎn)實(shí)驗(yàn)室,710068,西安)

近20年來,隨機(jī)有限集(RFS)理論[1-2]被廣泛地應(yīng)用于多目標(biāo)跟蹤領(lǐng)域。與傳統(tǒng)的多目標(biāo)跟蹤算法[3-5]相比,基于RFS理論的多目標(biāo)跟蹤算法不需要明確量測(cè)與目標(biāo)之間復(fù)雜繁瑣的對(duì)應(yīng)關(guān)系,尤其適用于多目標(biāo)跟蹤問題。目前,基于RFS的多目標(biāo)跟蹤濾波器主要有概率假設(shè)密度(PHD)濾波器[6]、勢(shì)概率假設(shè)密度(CPHD)濾波器[7]、多目標(biāo)多伯努利(MeMBer)濾波器[1]和δ-廣義標(biāo)記多伯努利(δ-GLMB)濾波器[8-9]。特別地,δ-GLMB濾波器不僅在跟蹤精度方面優(yōu)于其他RFS濾波器,且能自動(dòng)生成目標(biāo)航跡,是首個(gè)理論可證的具有精確閉式解的多目標(biāo)貝葉斯跟蹤濾波器。由于δ-GLMB濾波器的遞推過程是嚴(yán)格封閉的,因此它對(duì)目標(biāo)個(gè)數(shù)及其狀態(tài)的估計(jì)精度高于其他同類濾波器,但同樣面臨計(jì)算量大的問題。該濾波器的計(jì)算復(fù)雜度理論上與目標(biāo)的個(gè)數(shù)呈線性關(guān)系,與量測(cè)的個(gè)數(shù)呈3次方關(guān)系[10],但在其實(shí)現(xiàn)過程中仍需要花費(fèi)較大的計(jì)算代價(jià)用于假設(shè)排序問題。文獻(xiàn)[8-9]給出了該濾波器在線性高斯模型下的高斯混合(GM)實(shí)現(xiàn)和非線性模型下的序貫蒙特卡洛(SMC)實(shí)現(xiàn)方法,其中,由于SMC-δ-GLMB濾波算法受粒子數(shù)和粒子退化等問題的影響,盡管可以通過增加粒子數(shù)以獲得較高的跟蹤精度,但大量的粒子同樣意味著需要消耗更大的計(jì)算量,難以實(shí)現(xiàn)δ-GLMB濾波器對(duì)目標(biāo)的快速跟蹤。在實(shí)際的多目標(biāo)跟蹤應(yīng)用中,非線性模型則更為普遍,如雷達(dá)和聲吶等量測(cè)均為非線性[11],因此研究非線性模型下δ-GLMB濾波器的多目標(biāo)跟蹤問題具有非常重要的意義。

GM-δ-GLMB濾波器在線性高斯模型下具有閉式解,并且可以將GM-δ-GLMB濾波器與擴(kuò)展卡爾曼(EK)[11]和無味卡爾曼(UK)[12]相結(jié)合用于解決非線性模型下的多目標(biāo)跟蹤問題。由于EK算法只是簡(jiǎn)單地將非線性系統(tǒng)模型進(jìn)行局部線性化處理,存在嚴(yán)重的模型描述誤差,當(dāng)模型的非線性程度較高時(shí),該濾波算法可能會(huì)導(dǎo)致被跟蹤目標(biāo)丟失,且在對(duì)非線性函數(shù)進(jìn)行Taylor級(jí)數(shù)展開時(shí),需要求解Jacobians矩陣和Hessians矩陣(二階EKF),而這2種矩陣的計(jì)算容易導(dǎo)致數(shù)值不穩(wěn)定,且計(jì)算量較大,更有甚者,它們根本就不存在[13]。UK算法則是通過選擇一組確定的采樣點(diǎn)來近似非線性函數(shù)的均值和方差,精度可達(dá)到二階,但其采樣點(diǎn)及權(quán)系數(shù)一般由狀態(tài)維數(shù)和憑慣例或經(jīng)驗(yàn)而設(shè)定的參數(shù)共同來決定,很難達(dá)到較高精度的跟蹤效果。關(guān)于EK和UK濾波算法的上述限制條件,直接阻礙了δ-GLMB濾波器向非線性多目標(biāo)跟蹤應(yīng)用領(lǐng)域的推廣。

本文將積分卡爾曼(QK)[14-15]非線性濾波算法與GM-δ-GLMB濾波器相結(jié)合,提出了具有更加嚴(yán)謹(jǐn)?shù)臄?shù)學(xué)理論基礎(chǔ)且適用性更強(qiáng)的QK-GM-δ-GLMB濾波器。為了驗(yàn)證該算法的有效性,通過在不同的雜波密度和檢測(cè)概率條件下的仿真實(shí)驗(yàn),將所提算法的濾波精度和實(shí)時(shí)性能分別與EK-GM-δ-GLMB、UK-GM-δ-GLMB和SMC-δ-GLMB 3種濾波器做了較詳細(xì)的對(duì)比。實(shí)驗(yàn)結(jié)果表明:QK-GM-δ-GLMB濾波器的跟蹤精度是4種濾波器中最優(yōu)的,盡管其時(shí)間消耗有所增加,但與另外兩種δ-GLMB濾波器的GM實(shí)現(xiàn)算法均保持在同一量級(jí),這在許多實(shí)際應(yīng)用中是可以接受的。

1 QK-δ-GLMB濾波器的GM實(shí)現(xiàn)

假設(shè)每個(gè)目標(biāo)均服從如下的非線性運(yùn)動(dòng)和量測(cè)模型,即

xk=f(xk-1)+wk-1

(1)

zk=h(xk)+vk

(2)

式中:k為觀測(cè)時(shí)刻;xk為k時(shí)刻目標(biāo)的狀態(tài)向量;zk為k時(shí)刻傳感器的觀測(cè)向量;f(·)和h(·)分別為已知的非線性狀態(tài)轉(zhuǎn)移函數(shù)和非線性觀測(cè)函數(shù);wk-1和vk分別表示相互獨(dú)立的過程噪聲和量測(cè)噪聲,并假設(shè)它們服從均值為零、協(xié)方差矩陣分別為Qk-1和Rk的高斯分布。

積分卡爾曼濾波算法的主要思想是基于Gauss -Hermite數(shù)值積分規(guī)則求取一組帶權(quán)重的積分點(diǎn),利用這些積分點(diǎn)求取密度函數(shù)的均值和協(xié)方差估計(jì)[14-15]。假設(shè)向量x～N(x;0,Inx),則nx維的積分規(guī)則可由下式表示

(3)

通過構(gòu)造一個(gè)具有零對(duì)角元素的對(duì)稱三對(duì)角矩陣J來求取其特征值,其中,εl表示矩陣J的第l個(gè)特征值,且該矩陣其他元素可由下式求得

(4)

根據(jù)上述所求取的一組帶權(quán)重積分點(diǎn),可以推導(dǎo)出基于QK的δ-GLMB濾波器的高斯混合實(shí)現(xiàn)過程,具體預(yù)測(cè)步驟如下。假設(shè)k時(shí)刻的多目標(biāo)后驗(yàn)概率密度具有如下的δ-GLMB RFS形式[9]

(5)

式中:ξ為直到k時(shí)刻的所有關(guān)聯(lián)映射對(duì)應(yīng)于標(biāo)簽集合I的關(guān)聯(lián)歷史,ξ=(θ1,…,θk),θk表示k時(shí)刻航跡標(biāo)簽到傳感器量測(cè)的關(guān)聯(lián)映射;X表示帶標(biāo)簽的多伯努利狀態(tài)集合;數(shù)對(duì)(I,ξ)表示當(dāng)前時(shí)刻所有航跡的標(biāo)簽集合I∈F(L)與直到當(dāng)前時(shí)刻的標(biāo)簽-量測(cè)關(guān)聯(lián)歷史ξ∈Ξ之間的一種對(duì)應(yīng)關(guān)系;w(I,ξ)表示該假設(shè)為真的概率;F(L)表示由標(biāo)簽空間L的子集所構(gòu)成的集合;離散空間Ξ表示航跡標(biāo)簽到量測(cè)的關(guān)聯(lián)過程。直到k時(shí)刻的所有關(guān)聯(lián)歷史ξ條件下標(biāo)簽為l對(duì)應(yīng)目標(biāo)航跡的廣義標(biāo)記多伯努利分量的空間分布密度函數(shù)由J(ξ)(l)個(gè)高斯項(xiàng)組成

(6)

由δ-GLMB分布的共軛先驗(yàn)性[8]可知,k+1時(shí)刻的多目標(biāo)預(yù)測(cè)分布密度也具有δ-GLMB形式,即

(7)

(8)

式中:pB(x,l)表示新生目標(biāo)的狀態(tài)空間密度;HL(l)與HB(l)分別表示存活目標(biāo)和新生目標(biāo)的示性函數(shù)。標(biāo)簽為l∈L的存活航跡的空間分布密度具有如下的高斯混合形式

(9)

(10)

(11)

(12)

為簡(jiǎn)單起見,假設(shè)目標(biāo)航跡存活概率獨(dú)立于標(biāo)簽l和狀態(tài)x,即pS(x,l)=pS,航跡消亡概率為qS(x,l)=1-pS(x,l)=1-pS。

設(shè)新生目標(biāo)概率密度可表示為如下的LMB RFS形式

(13)

(14)

(15)

(16)

式中:I+?L+=L∪B表示預(yù)測(cè)目標(biāo)的標(biāo)簽空間,而新生目標(biāo)航跡的假設(shè)wB(J)則與具體的應(yīng)用場(chǎng)景相關(guān)。

用求得的積分點(diǎn)bl和相應(yīng)權(quán)重wl來近似式(9)中第j個(gè)高斯項(xiàng)的均值和協(xié)方差,具體實(shí)現(xiàn)步驟如下。

步驟1對(duì)式(6)中組成標(biāo)簽為l的航跡中的第j個(gè)高斯項(xiàng)的協(xié)方差矩陣進(jìn)行Cholesky分解

(17)

步驟2計(jì)算組成標(biāo)簽l的第j個(gè)高斯項(xiàng)的積分點(diǎn)

(18)

步驟3非線性狀態(tài)方程傳遞積分點(diǎn)

(19)

步驟4計(jì)算預(yù)測(cè)狀態(tài)均值

(20)

步驟5計(jì)算狀態(tài)預(yù)測(cè)協(xié)方差矩陣

(21)

(22)

則更新后的后驗(yàn)密度仍然為δ-GLMB形式

δI+(L(X))[p(ξ,θ)(·|Z)]X}

(23)

式中Θ表示航跡標(biāo)簽到量測(cè)關(guān)聯(lián)θ的映射空間,即L+→{0,1,…,|Z|},0表示漏檢。此時(shí)的后驗(yàn)標(biāo)簽集合與預(yù)測(cè)航跡標(biāo)簽集合相等,仍為I+。更新假設(shè)的概率為

(24)

令k+1時(shí)刻標(biāo)簽l的航跡關(guān)聯(lián)映射為θ(l),若θ(l)>0,表示該航跡被傳感器收到的量測(cè)集合zθ(l)更新;若θ(l)=0,則表示沒有量測(cè)與之關(guān)聯(lián),該航跡漏檢。總的量測(cè)更新廣義似然函數(shù)可表示為

(25)

假設(shè)檢測(cè)概率獨(dú)立于標(biāo)簽l和狀態(tài)x,即pD(x,l)=pD,漏檢概率qD(x,l)=1-pD(x,l)=1-pD。

對(duì)于更新航跡,標(biāo)簽為l的航跡關(guān)聯(lián)映射滿足θ(l)>0,式(25)的廣義似然函數(shù)可具體表示為

(26)

將式(22)與式(26)相乘并進(jìn)行歸一化,即可得航跡l更新后的δ-GLMB空間分布密度

(27)

(28)

歸一化常數(shù)

(29)

(30)

式中:q(ξ,j)(z(θ(l)),l)表示目標(biāo)(x,l)生成量測(cè)z(θ(l))的似然概率;κ(z(θ(l)))表示量測(cè)過程中的雜波密度函數(shù)。

對(duì)于漏檢航跡,因沒有量測(cè)與之關(guān)聯(lián),故θ(l)=0,此時(shí)的廣義似然函數(shù)可重新表示為

(31)

(32)

將式(22)與式(31)相乘并除以式(32),可得標(biāo)簽l的航跡漏檢時(shí)的空間分布密度,此時(shí)的假設(shè)更新概率、狀態(tài)均值和狀態(tài)協(xié)方差矩陣分別對(duì)應(yīng)于該航跡的預(yù)測(cè)空間密度式(22)中的各參數(shù),即

(33)

(34)

(35)

步驟1對(duì)式(22)中組成航跡l的第j個(gè)高斯項(xiàng)的協(xié)方差矩陣進(jìn)行Cholesky分解

(36)

步驟2計(jì)算標(biāo)簽l的第j個(gè)高斯項(xiàng)的積分點(diǎn)

(37)

步驟3非線性量測(cè)方程傳遞積分點(diǎn)

(38)

步驟4計(jì)算預(yù)測(cè)量測(cè)均值

(39)

步驟5計(jì)算量測(cè)預(yù)測(cè)協(xié)方差矩陣

S(ξ,j)=

(40)

步驟6計(jì)算互協(xié)方差矩陣

(41)

步驟7計(jì)算卡爾曼增益

K(ξ,j)=G(ξ,j)[S(ξ,j)]-1

(42)

步驟8均值更新

(43)

步驟9協(xié)方差矩陣更新

(44)

2 仿真驗(yàn)證

2.1 仿真場(chǎng)景及參數(shù)設(shè)定

在二維平面內(nèi),設(shè)觀測(cè)區(qū)域范圍為[0°,180°]×[0 m,2 000 m],觀測(cè)時(shí)長為50 s,傳感器采樣間隔T=1 s,最多可能同時(shí)有8個(gè)目標(biāo)隨機(jī)出現(xiàn)在該區(qū)域的任何位置,且目標(biāo)出現(xiàn)和消失的時(shí)刻都是不確定的。目標(biāo)真實(shí)運(yùn)動(dòng)軌跡如圖1所示。

圖1 極坐標(biāo)下多目標(biāo)的真實(shí)運(yùn)動(dòng)軌跡

目標(biāo)的運(yùn)動(dòng)模型為CT模型[4],且觀測(cè)數(shù)據(jù)存在角度和距離的觀測(cè)噪聲,k時(shí)刻目標(biāo)的狀態(tài)向量xk=[px,k,vx,k,py,k,vy,k,ωk]T由位置分量(px,k,py,k)、速度分量(vx,k,vy,k)和轉(zhuǎn)彎速率ωk組成。CT模型的狀態(tài)轉(zhuǎn)移方程可表示為

xk=F(ωk-1)xk-1+Gkwk-1

(45)

F(ωk-1)=

(46)

(47)

含方位角和徑向距離噪聲的非線性量測(cè)方程為

(48)

目標(biāo)的存活概率pS,k=0.99,檢測(cè)概率pD,k=0.98。雜波模型服從強(qiáng)度為κk=λcVu(z)的泊松RFS,其中λc表示平均每幀雜波的個(gè)數(shù),V為傳感器監(jiān)控區(qū)域的面積,u(·)為該區(qū)域內(nèi)的均勻概率分布。

采用剪切和合并技術(shù)限制高斯項(xiàng)數(shù)目的增長,設(shè)航跡的截?cái)喔怕书撝礣′=10-5,高斯項(xiàng)權(quán)重的閾值T″=10-5,合并閾值U=4,航跡最大值Mmax=300,每條航跡對(duì)應(yīng)高斯項(xiàng)個(gè)數(shù)的最大值Jmax=100,其具體剪切和合并過程與文獻(xiàn)[16]所述一致。

實(shí)驗(yàn)通過100次蒙特卡洛(MC)仿真驗(yàn)證了所提算法的有效性,采用最優(yōu)子模式分配(OSPA)距離DOSPA[17]來評(píng)價(jià)幾種多目標(biāo)跟蹤算法的性能。本實(shí)驗(yàn)中,將OSPA表達(dá)式的誤差調(diào)節(jié)參數(shù)c和階次參數(shù)p分別設(shè)置為100、2。與已有的EK-GM-δ-GLMB、UK-GM-δ-GLMB和SMC-δ-GLMB濾波算法關(guān)于平均OSPA距離和平均單步運(yùn)行時(shí)間在不同檢測(cè)概率和不同的雜波強(qiáng)度條件下做了比較,其中,在SMC實(shí)現(xiàn)過程中,每個(gè)假設(shè)航跡被分配的最大粒子數(shù)為10 000,最小粒子數(shù)為2 000。

仿真實(shí)驗(yàn)運(yùn)行環(huán)境軟件和硬件參數(shù)為Matlab R2017b,Windows 7 64 bit,Intel Core i5-4570 CPU 3.20 GHz,RAM 4.00 GB。

2.2 實(shí)驗(yàn)結(jié)果分析

圖2給出了4種算法在檢測(cè)概率pD為0.98、雜波密度λc為3.2×10-3(rad·m)-1的實(shí)驗(yàn)場(chǎng)景中經(jīng)100次MC仿真實(shí)驗(yàn)后DOSPA隨時(shí)間變化曲線的對(duì)比。從圖2可以看出,4種算法的DOSPA具有相同的變化趨勢(shì),說明在存在目標(biāo)消失和新生目標(biāo)出現(xiàn)的場(chǎng)景中,它們均可以準(zhǔn)確跟蹤目標(biāo),但在跟蹤精度上有所差別。

圖2 100次MC仿真實(shí)驗(yàn)后DOSPA的對(duì)比

由圖2可知,各濾波器跟蹤精度依QK-GM-δ-GLMB、EK-GM-δ-GLMB、UK-GM-δ-GLMB、SMC-δ-GLMB順序依次降低,具體原因分析如下。

(1)在δ-GLMB濾波器的更新過程中涉及到量測(cè)信息與航跡標(biāo)簽之間的關(guān)聯(lián)問題,而這種關(guān)聯(lián)問題存在較大的不確定性,隨著量測(cè)信息的不斷增加,量測(cè)-標(biāo)簽之間的關(guān)聯(lián)假設(shè)呈指數(shù)增長,造成大量粒子被用來逼近這些小權(quán)重的關(guān)聯(lián)假設(shè),使得大量的計(jì)算資源被用于粒子的重采樣過程,從而導(dǎo)致跟蹤精度與計(jì)算時(shí)間效率都較低。SMC-δ-GLMB算法的跟蹤精度在很大程度上取決于粒子數(shù)的大小,隨著粒子數(shù)的增加,SMC-δ-GLMB算法的實(shí)現(xiàn)精度也會(huì)隨之升高,但同時(shí)意味著需要更長的計(jì)算時(shí)間,不易實(shí)現(xiàn)快速跟蹤。

(2)UK-GM-δ-GLMB算法的濾波參數(shù)值一般是憑經(jīng)驗(yàn)和慣例而得到的,針對(duì)復(fù)雜和高狀態(tài)維數(shù)的多目標(biāo)跟蹤問題往往難以湊效。

(3)在本次實(shí)驗(yàn)中,EK-GM-δ-GLMB算法的跟蹤精度僅次于QK-GM-δ-GLMB算法,這是因?yàn)獒槍?duì)弱非線性問題,EK算法是首要考慮的方法,當(dāng)非線性程度較高時(shí),很容易出現(xiàn)濾波發(fā)散問題,從而導(dǎo)致被跟蹤目標(biāo)的丟失。

(4)本文所提的QK-GM-δ-GLMB濾波器具有較高的跟蹤精度,因?yàn)槠渌袨V波參數(shù)均是基于Gauss-Hermite積分規(guī)則得到,所以得到的參數(shù)更加準(zhǔn)確,從而可實(shí)現(xiàn)較高精度的目標(biāo)跟蹤。

為進(jìn)一步詳細(xì)比較4種算法的目標(biāo)跟蹤精度和實(shí)時(shí)性,在表1和表2中分別列出了不同場(chǎng)景中各濾波器經(jīng)100次MC仿真實(shí)驗(yàn)后的DOSPA和單步運(yùn)行時(shí)間的比較結(jié)果。

表1給出了當(dāng)雜波密度λc為0.8×10-3(rad·m)-1時(shí),在檢測(cè)概率pD為0.80、0.90、0.98的條件下4種濾波器的DOSPA和單步運(yùn)行時(shí)間的實(shí)驗(yàn)數(shù)據(jù)。由表1可以看出,在相同的雜波密度場(chǎng)景下,各濾波器的DOSPA和單步運(yùn)行時(shí)間都隨檢測(cè)概率的增加而減小,其中,在高檢測(cè)概率場(chǎng)景中,QK-GM-δ-GLMB算法的跟蹤精度最高,而在較低檢測(cè)概率場(chǎng)景中QK-GM-δ-GLMB與EK-GM-δ-GLMB算法同時(shí)具有最高跟蹤精度,且實(shí)時(shí)性優(yōu)于SMC-δ-GLMB算法。

表2給出了在檢測(cè)概率pD為0.98的條件下,當(dāng)雜波密度λc分別為0.8、1.6、3.2×10-3(rad·m)-1時(shí)4種濾波器的DOSPA和平均單步運(yùn)行時(shí)間的實(shí)驗(yàn)數(shù)據(jù)。由表2可以看出,在同一檢測(cè)概率條件下,各濾波器的DOSPA和單步運(yùn)行時(shí)間都隨著雜波密度的增強(qiáng)而增加,其中,在不同的雜波強(qiáng)度下,QK-GM-δ-GLMB算法具有最高的跟蹤精度,同樣其實(shí)時(shí)性表現(xiàn)處于UK-GM-δ-GLMB算法和SMC-δ-GLMB算法之間。當(dāng)雜波密度較高時(shí),SMC-δ-GLMB算法的計(jì)算代價(jià)是巨大的。

表1 不同檢測(cè)概率條件下的DOSPA和平均單步運(yùn)行時(shí)間

表2 不同雜波密度下的DOSPA和平均單步運(yùn)行時(shí)間

綜上所述,本文所提的QK-GM-δ-GLMB濾波算法在不同的場(chǎng)景中都能夠以較高精度對(duì)目標(biāo)進(jìn)行準(zhǔn)確跟蹤,而其對(duì)應(yīng)的時(shí)間消耗在實(shí)際應(yīng)用中處于可接受的范圍內(nèi),其適用場(chǎng)景也得到進(jìn)一步擴(kuò)展。

3 結(jié) 論

在非線性高斯多目標(biāo)運(yùn)動(dòng)場(chǎng)景下,給出了QK-GM-δ-GLMB濾波算法的遞推過程,QK-GM-δ-GLMB濾波算法具有良好的濾波性能,是一種適應(yīng)性強(qiáng)、估計(jì)精度高的跟蹤算法,其時(shí)間開銷在實(shí)際的多目標(biāo)跟蹤應(yīng)用中是完全可接受的,所以QK-GM-δ-GLMB濾波算法具有一定的工程應(yīng)用價(jià)值。

由于QK-GM-δ-GLMB濾波算法在非線性場(chǎng)景中具有較高的跟蹤精度,下一步將致力于將其應(yīng)用于機(jī)動(dòng)目標(biāo)、擴(kuò)展目標(biāo)和群目標(biāo)等較復(fù)雜的多目標(biāo)跟蹤問題。