海上低空突防群目標(biāo)跟蹤的IMM-Bayesian實現(xiàn)

王婷婷，繆禮鋒，程 然

（中國航空工業(yè)集團(tuán)公司雷華電子技術(shù)研究所，江蘇無錫214063）

在現(xiàn)代對海戰(zhàn)爭中，群目標(biāo)跟蹤問題受到了國內(nèi)外學(xué)者的廣泛關(guān)注，尤其海上低空突防編隊群目標(biāo)具有檢測率低、機(jī)動性強(qiáng)、戰(zhàn)術(shù)性合并分裂等特征，在實戰(zhàn)中靈活性更強(qiáng)、威脅性更大、對抗難度也更高。由于海上低空突防編隊內(nèi)的個體目標(biāo)較為密集，海雜波干擾較為嚴(yán)重，受傳感器分辨率的限制，完全分辨編隊內(nèi)個體目標(biāo)較為困難[1]。在這種情況下，試圖對群內(nèi)所有個體分別進(jìn)行跟蹤變得不現(xiàn)實，這時可以將群目標(biāo)看作一個整體進(jìn)行跟蹤[2]。在運(yùn)動過程中，基于特定的突防戰(zhàn)術(shù)或目的，海上低空突防編隊群目標(biāo)隨時會發(fā)生轉(zhuǎn)彎等整體機(jī)動，還會出現(xiàn)分裂、合并等編隊群目標(biāo)特有的機(jī)動模式，這種情況下，群內(nèi)目標(biāo)回波的相對位置發(fā)生縮放、旋轉(zhuǎn)等變換[3-5]，導(dǎo)致群目標(biāo)跟蹤更為復(fù)雜。

基于群整體跟蹤的思想，目前最經(jīng)典的群跟蹤的方法主要有兩種：一類是Mahler在隨機(jī)集框架下提出的概率密度假設(shè)（PHD）算法[6-8]，遞推的傳遞目標(biāo)狀態(tài)的分布信息，不需要數(shù)據(jù)關(guān)聯(lián)，即從目標(biāo)狀態(tài)分布中同時提取目標(biāo)的個數(shù)和狀態(tài)估計，但沒有考慮群目標(biāo)形狀特征，也無法直接生成目標(biāo)航跡，并且算法復(fù)雜，不利于實際工程應(yīng)用。還有一類是KOCH等人提出的一種群目標(biāo)Bayesian遞推算法[9]，可以在跟蹤群目標(biāo)的質(zhì)心狀態(tài)同時，對群目標(biāo)的形狀進(jìn)行估計。但是它目前僅用于跟蹤無雜波條件下的單個群目標(biāo)。隨后，國內(nèi)外學(xué)者又對這種算法進(jìn)行了改進(jìn)[10-16]，例如，F(xiàn)eldmann對群目標(biāo)跟蹤Bayesian算法進(jìn)行了改進(jìn)，考慮了量測不確定性對目標(biāo)形狀估計的影響。

本文在經(jīng)典Bayesian算法的基礎(chǔ)上，針對海上低空突防編隊群目標(biāo)的3種典型機(jī)動模式（拐彎、分裂、合并），將交互式多模型（Interacting Multiple Model，IMM）算法與群目標(biāo)Bayesian遞推算法相結(jié)合，并利用隨機(jī)矩陣描述群目標(biāo)的形狀擴(kuò)展?fàn)顟B(tài)，采用IMMBayesian算法對機(jī)動群目標(biāo)進(jìn)行有效跟蹤的同時，用橢圓擬合群目標(biāo)的形狀。

1 群目標(biāo)模型

將海上低空突防群目標(biāo)看作一個整體，假設(shè)在tk時刻群質(zhì)心的運(yùn)動狀態(tài)向量為，包括了位置、速度、加速度等空間維數(shù)。表示tk時刻目標(biāo)的位置表示tk時刻目標(biāo)的速度表示tk時刻目標(biāo)的加速度。Bayesian遞推算法用表示群質(zhì)心的運(yùn)動狀態(tài)向量，并增加了一個正定隨機(jī)矩陣Xk來描述該群的形狀特征，該方法將目標(biāo)的形狀用一個橢圓來近似，橢圓方程可以用二次型方程表示為：

式（1）中：xc為橢圓中心；正定矩陣Xk的特征值為橢圓各軸的長度，其所對應(yīng)的特征向量為其橢圓各軸的方向向量。

該方法的性質(zhì)完備地描述了空間或平面內(nèi)一個橢圓的全部信息。

用狀態(tài)(xk,Xk)描述tk時刻群目標(biāo)的狀態(tài)，依據(jù)經(jīng)典Kalman濾波理論，對群目標(biāo)質(zhì)心建立如下動力學(xué)方程：

式（2）中：Fk為狀態(tài)轉(zhuǎn)移矩陣；wk為過程噪聲，是零均值的高斯白噪聲，即wk～N(0,Δk|k-1)。

式中：?為kronecker積，采用Van Keuk模型[17]；

其中，θ為機(jī)動飛行相關(guān)的時間常數(shù)。

群目標(biāo)量測模型為：

式（4）中：Hk為量測矩陣；vkj為量測噪聲，是零均值的高斯白噪聲，vk～N(0 ,λXk+Rk)，其中，λ為標(biāo)量，用于描述Xk對量測的影響程度，Rk為量測噪聲協(xié)方差；nk為tk時刻某群目標(biāo)得到的量測個數(shù)。

2 群目標(biāo)跟蹤的Bayesian方法

2.1 群目標(biāo)跟蹤預(yù)測部分

利用先前的信息，對被跟蹤目標(biāo)的下一時刻的信息進(jìn)行估計，即

式中，Pk-1|k-1為tk-1時刻狀態(tài)估計誤差協(xié)方差。

那么，一步預(yù)測：

所以，目標(biāo)的運(yùn)動狀態(tài)預(yù)測為：

式（10）中存在一個未被限制的參量vk-1|k-1，由逆Wishart分布相關(guān)知識可知，此參量與Xk-1有關(guān)，它在每一時刻的取值必須在實際跟蹤中進(jìn)行合理的設(shè)定或遞歸處理。通常設(shè)定vk-1|k-1和Xk-1|k-1與上一次預(yù)測后的值相關(guān)。p(Xk|Zk-1)同樣滿足逆Wishart分布，并有形式：

考慮一個關(guān)于vk-1|k-1和vk|k-1的衰減模型，這種衰減也與采樣區(qū)間T相關(guān)，即

存在未被限定的參量τ作為衰減因子常量，需要在實際跟蹤中進(jìn)行合理的設(shè)定。當(dāng)τ=∞時，代表群目標(biāo)的形狀不發(fā)生變化。

假設(shè)預(yù)測步前后Xk-1和Xk的期望值相等。由逆Wishart分布的知識，有：

式（13）中，d為矩陣X的維數(shù)。

那么，形狀擴(kuò)展?fàn)顟B(tài)的預(yù)測為：

2.2 群目標(biāo)跟蹤更新部分

對被跟蹤目標(biāo)的相關(guān)信息的更新，即

引入似然函數(shù)：

式（16）中，nk為當(dāng)前時刻的群目標(biāo)的量測個數(shù)，其與當(dāng)前的動力學(xué)狀態(tài)和形狀狀態(tài)無關(guān)。

p(nk|xk,Xk)可以看作一個常數(shù)，并且

利用預(yù)測部分得到的結(jié)果，由Bayesian先驗概率公式，可得：

所以：

目標(biāo)的運(yùn)動狀態(tài)更新為：

式（20）、（21）中：

式（24）中：

所以，

目標(biāo)形狀擴(kuò)展的更新為：

至此，對于被跟蹤的群目標(biāo)的一個可用的跟蹤算法已被建立起來。

3 IMM-Bayesian算法

在廣義Bayesian算法的基礎(chǔ)上，采用具有Markov切換系數(shù)的IMM算法，推導(dǎo)了IMM-Bayesian算法對海上低空突防編隊群目標(biāo)進(jìn)行跟蹤。在tk時刻，利用IMM進(jìn)行目標(biāo)狀態(tài)估計的計算時，考慮每個模型濾波對于目標(biāo)的擴(kuò)展部分，器都有可能成為當(dāng)前有效地系統(tǒng)模型濾波器，每個濾波器的初始條件都是基于前一時刻各條件模型濾波結(jié)果的合成。IMM-Bayesian算法對海上低空突防編隊群目標(biāo)跟蹤的具體遞推步驟如下：

1）模型條件初始化。模型條件重初始化是在假定第j個模型在當(dāng)前時刻有效的條件下，與其匹配的濾波器的輸入由上一時刻各濾波器的估計混合而成。

對于j=1,2,…,r，重初始化的運(yùn)動狀態(tài)、形狀擴(kuò)展?fàn)顟B(tài)與協(xié)方差陣按混合估計分別為：

2）模型條件濾波。模型條件濾波是在給定重初始化的狀態(tài)和協(xié)方差陣的前提下，在獲得新的量測zk之后，進(jìn)行狀態(tài)估計更新。

對于j=1,2,…,r，狀態(tài)預(yù)測為：

對于j=1,2,…,r，分別計算量測預(yù)測殘差與其協(xié)方差陣：

在Gaussian假設(shè)下，似然函數(shù)可以計算如下：

濾波增益陣、狀態(tài)估計更新、形狀擴(kuò)展?fàn)顟B(tài)更新和狀態(tài)估計更新誤差協(xié)方差陣如下：

式（42）中：

3）模型概率更新。對j=1,2,…,r，計算模型概率：

4）輸出交互。給出運(yùn)動狀態(tài)估計、形狀擴(kuò)展?fàn)顟B(tài)估計和誤差協(xié)方差陣，分別為：

4 轉(zhuǎn)彎模式下群目標(biāo)跟蹤實現(xiàn)

假設(shè)海上低空突防群目標(biāo)的質(zhì)心按圖1軌跡勻速運(yùn)動，速度為 300 m/s，在采樣區(qū)間k=13～15、k=25～29、k=39～44進(jìn)行機(jī)動轉(zhuǎn)彎，分別轉(zhuǎn)過 45°、90°和90°，其他時間段均為300 m/s的勻速運(yùn)動。量測噪聲服從分布，其中α=1/4。將表征機(jī)動群目標(biāo)形狀特征的正定對稱矩陣初始化為X=diag([ 2 000,2000])，Rk=diag([ 25,25])。在每次掃描過程中，量測數(shù)據(jù)的個數(shù)服從均值為λ=12的泊松分布。群目標(biāo)形狀更新時間延遲參數(shù)τ=8T。

圖1 群目標(biāo)運(yùn)動軌跡Fig.1 Trajectory of group target

對于給定的場景，采用3個模型。第1個模型為CV模型；第2個模型為CT模型，ω1=-0.01πrad/s；第3個模型為CT模型，ω2=-0.0125πrad/s。則群目標(biāo)跟蹤結(jié)果、模型匹配結(jié)果以及跟蹤誤差如圖2～4所示。

圖3 模型概率Fig.3 Model probability

圖4 群目標(biāo)跟蹤誤差Fig.4 Tracking error of group target

通過仿真可看出，采用IMM-Bayesian方法對在海上低空突防編隊群目標(biāo)轉(zhuǎn)彎的情況下可以達(dá)到較好的跟蹤效果，圖3給出了不同采樣時刻的模型匹配結(jié)果，可以看出選用的3個模型可以較好地完成模型匹配，模型匹配與真實運(yùn)動過程模型基本吻合。從圖4可以看出，分別采用多模型進(jìn)行群目標(biāo)跟蹤平均誤差較小，跟蹤精度較高，但在目標(biāo)機(jī)動位置誤差會增大。

5 群目標(biāo)的分裂與合并實現(xiàn)

針對海戰(zhàn)場探測區(qū)域中出現(xiàn)的低空突防群目標(biāo)進(jìn)行分裂或合并的情況，根據(jù)最近鄰的思想，將得到的量測進(jìn)行分類，再根據(jù)群目標(biāo)的跟蹤算法對其實施有效跟蹤。在第k個掃描周期，采用分類器將得到的所有量測進(jìn)行分類，分類個數(shù)記作nClass(k)。與前一時刻的nClass(k-1)比較。若nClass(k)=nClass(k-1)，則說明檢測到的群目標(biāo)沒有進(jìn)行分裂或合并，先采用關(guān)聯(lián)算法進(jìn)行k與k-1時刻各分群的關(guān)聯(lián)，關(guān)聯(lián)之后分別進(jìn)行跟蹤；若nClass(k)≠nClass(k-1)，則說明檢測到的群目標(biāo)進(jìn)行了分裂（滿足nClass(k)＞nClass(k-1)）或合并（滿足nClass(k)＜nClass(k-1)），則需對k時刻的各個分群分別進(jìn)行初始化。

群目標(biāo)分離場景設(shè)置：假設(shè)海上低空突防編隊群目標(biāo)質(zhì)心初始位置[0,0]，開始以300 m/s的速度沿45°方向運(yùn)動。t=39T時，群目標(biāo)中部分目標(biāo)出現(xiàn)分離，形成分群1和分群2。分群1始終保持勻速直線運(yùn)動，分群2以(1 . 873×2π/360) rad/s的角速度勻速向右偏轉(zhuǎn)，當(dāng)t＞64T后，分離出來的目標(biāo)2機(jī)動停止，保持勻速直線運(yùn)動。

圖5 群目標(biāo)分離Fig.5 Splitting trajectory of group target

圖6 群目標(biāo)分離過程Fig.6 Splitting process of group target

仿真可以看出，IMM-Bayesian算法可在海上低空突防編隊群目標(biāo)分裂時進(jìn)行穩(wěn)定跟蹤，隨機(jī)矩陣可對群目標(biāo)形狀進(jìn)行擬合，圖6清晰地展現(xiàn)了群目標(biāo)分裂過程中目標(biāo)的跟蹤情況，當(dāng)出現(xiàn)2個橢圓時，說明群目標(biāo)已完成分裂，之后當(dāng)成2個群目標(biāo)進(jìn)行穩(wěn)定跟蹤。

群目標(biāo)合并場景設(shè)置：假設(shè)檢測到海上有2個分群目標(biāo)，分群1初始位置為[0,3488]，分群2初始位置為[0,0]，開始均以300 m/s的速度沿x軸方向運(yùn)動。t=39T時，分群1和分群2開始機(jī)動，在39T＜t＜54T時間段，分群1以(1 .5×2π/360) rad/s的角速度勻速向右偏轉(zhuǎn)；分群2以(- 1.5×2π/360) rad/s的角速度勻速向左偏轉(zhuǎn)。在 54T＜t＜69T時間段，分群1以(-1.5×2π/360) rad/s的角速度勻速向左偏轉(zhuǎn)；分群2以(1 . 5×2π/360) rad/s的角速度勻速向右偏轉(zhuǎn)。之后，分群1和分群2機(jī)動停止，保持勻速直線運(yùn)動。

圖7 群目標(biāo)合并軌跡Fig.7 Merging trajectory of group target

圖8 群目標(biāo)合并過程Fig.8 Merging process of group target

仿真可以看出，IMM-Bayesian算法可在海上低空突防編隊群目標(biāo)合并時進(jìn)行穩(wěn)定跟蹤，隨機(jī)矩陣可實施估計群目標(biāo)形狀，圖8展現(xiàn)了群目標(biāo)合并過程中目標(biāo)的跟蹤情況，當(dāng)2個橢圓合并成1個橢圓時，說明群目標(biāo)已完成合并，之后當(dāng)成1個群目標(biāo)進(jìn)行穩(wěn)定跟蹤。

6 結(jié)論

針對拐彎、分裂、合并等典型海上低空突防編隊群目標(biāo)機(jī)動模式，將IMM算法與Bayesian遞推算法結(jié)合，完成對海上低空突防群目標(biāo)整體的跟蹤，并利用隨機(jī)矩陣描述群目標(biāo)形狀特征，采用橢圓擬合群目標(biāo)形狀。其中針對群目標(biāo)出現(xiàn)的分裂合并現(xiàn)象，采用基于最近鄰算法思想的跟蹤算法對群目標(biāo)進(jìn)行穩(wěn)定跟蹤。本算法只針對單個傳感器對海上群目標(biāo)的跟蹤，又由于用橢圓僅是粗略地對群目標(biāo)形狀進(jìn)行估計，需進(jìn)一步研究使其能更加精確的估計群目標(biāo)的形狀。

參考文獻(xiàn)：

[1]連峰，王婷婷，韓崇昭.多個不可分辨目標(biāo)群的聯(lián)合檢測與估計誤差界[J].西安交通大學(xué)學(xué)報，2015，49（11）：89-95.LIAN FENG，WANG TINGTING，HAN CHONGZHAO.Error bounds of joint detection and estimation for multiple unresolved target-groups[J].Journal of Xi’an Jiaotong University，2015，49（11）：89-95.（in Chinese）

[2]耿文東.基于群目標(biāo)幾何中心的群起始算法研究[J].系統(tǒng)工程與電子技術(shù)，2008，30（2）：269-272.GENG WENDONG.Study of group-initialization method based on group-target center of geometry[J].Systems Engineering and Eletronics，2008，30（2）：269-272.（in Chinese）

[3]王海鵬，唐田田，王子玲，等.面向群目標(biāo)典型機(jī)動的多傳感器精細(xì)跟蹤方法[J].船舶電子工程，2016，36（8）：76-80.WANG HAIPENG，TANG TIANTIAN，WANG ZILING，et al.Multi-sensor refined tracking algorithm of typical group targets maneuvering models[J].Ship Electronic Engineering，2016，36（8）：76-80.（in Chinese）

[4]耿文東，王元欽，董正宏.群目標(biāo)跟蹤[M].北京：國防工業(yè)出版社，2014：100-122.GENG WENDONG，WANG YUANQIN，DONG ZHENGHONG.Group-target tracking[M].Beijing：National Defense Industry Press，2014：100-122.（in Chinese）

[5]耿文東.基于PDA的群目標(biāo)合并與分離方法研究[J].無線電工程，2007，37（2）：24-30.GENG WENDONG.A study on group target merging and splitting method based on PDA[J].Radio Engineering，2007，37（2）：24-30.（in Chinese）

[6]MAHLER R.Multi-target Bayes filtering via first-order multi-target moments[J].IEEE Transactions on Aerospace and Electronic Systems，2003，39（4）：1152-1178.

[7]MAHLER R.PHD filters of higher order in target number[J].IEEE Transactions on Aerospace and Electronic Systems，2007，43（4）：1523-1543.

[8]張慧，徐暉，王雪瑩，等.一種基于橢圓隨機(jī)超曲面模型的群目標(biāo)高斯混合PHD濾波器[J].光學(xué)學(xué)報，2013，33（9）：1-10.ZHANG HUI，XU HUI，WANG XUEYING，et al.A Gaussion mixture PHD filter for group targets tracking based on ellipse random hypersurface models[J].Acta Optica Sinica，2013，33（9）：1-10.（in Chinese）

[9]KOCH W.Bayesian approach to extended object and cluster tracking using random matrices[J].IEEE Transactions on Aerospace and Electronic Systems，2008，44（3）：1042-1059.

[10]FELDMANN M，F(xiàn)R?NKEN D.Tracking of extended objects and group targets using random matrices—a new approach[C]//Proceedings of the 11thInternational Conference on Information Fusion.Cologne：IEEE，2008：1-8.

[11]FELDMANN M，F(xiàn)R?NKEN D，KOCH W.Tracking of extended objects and group targets using random matrices[J].IEEE Transactions on Signal Processing，2011，59（4）：1409-1420.

[12]FELDMANN M，F(xiàn)RANKEN D.Advances on tracking of extended objects and group targets using random matrices[C]//Proceedings of the 12thInternational Conference on Information Fusion.Seattle，Washington，USA：IEEE，2009：1029-1036.

[13]韓玉蘭，朱洪艷，韓崇昭.采用隨機(jī)矩陣的多擴(kuò)展目標(biāo)濾波器[J].西安交通大學(xué)學(xué)報，2015，49（7）：98-104.HAN YULAN，ZHU HONGYAN，HAN CHONGZHAO.A multi-target filter based on random matrix[J].Journal of Xi’an Jiaotong University，2015，49（7）：98-104.（in Chinese）

[14]WIENEKE M，KOCH W.Probabilistic tracking of multiple extended targets using random matrices[C]//Proceedings of SPIE.2010，7698：1-12.

[15]BAUM M，NOACK B，HANEBECK U D.Mixture random hypersurface models for tracking multiple extended objects[C]//Proceedings of the IEEE Conference on Decision and Control.Piscataway，NJ，USA：IEEE，2011：3166-3171.

[16]王婷婷，何科峰，程然.具有形狀信息的多個群目標(biāo)跟蹤算法[J].雷達(dá)科學(xué)與技術(shù)，2017，15（5）：531-536.WANG TINGTING，HE KEFENG，CHENG RAN.A tracking algorithm for multiple group-targets with shape information[J].Radar Science and Technology，2017，15（5）：531-536.（in Chinese）

[17]VAN KEUK G，BLACKMAN S S.On phased-array radar tracking and parameter control[J].IEEE Transactions on Aerospace&Electronic Systems Aes，1993，29（1）：186-194.

[18]HAFF L R.Identities for the inverse Wishart distribution with computational results in linear and quadratic discrimination[J].The Indian Journal of Statistics：Series B，1982，44（3）：245-258.