基于代價參考粒子濾波器組的多目標檢測前跟蹤算法

盧 錦 馬令坤 呂春玲 章為川 Sun Chang-Ming

低信噪比情況下的多目標跟蹤是雷達、聲吶、紅外探測、圖像處理等領域面臨的難題之一,其目的是從觀測數據中估計時變的多目標狀態,包括目標數量、位置等信息.傳統先檢測后跟蹤方法難以有效處理低信噪比目標[1-2].為提高低信噪比目標的檢測和跟蹤性能,學者們提出檢測前跟蹤(Trackbefore-detect,TBD)方法.TBD 方法從未經門限檢測的原始觀測數據中同時檢測和跟蹤低信噪比目標[3-8].實現多目標TBD 的算法主要包括粒子濾波[9-11]、概率假設密度(Probability hypothesis density,PHD)[12-13]、多目標伯努利濾波器[14]等.

文獻[9]針對第2 個目標自第1 個目標產生的雙目標跟蹤與檢測問題,提出基于粒子濾波的多目標TBD 方法,該方法是文獻[6]中單目標粒子濾波的多目標檢測前跟蹤算法(Particle filter based track-before-detect,PF-TBD)方法的拓展.但隨著目標數量的增加,該計算計算量急劇增加,效率較低.文獻[10]提出一種多模型多目標TBD 方法.該方法需假設最大目標數量,首先,估計所有可能的目標存在組合或模式及其對應的聯合后驗概率密度和目標模型的概率;然后,整合不同模式下目標對應的狀態向量,獲得所有目標的狀態估計.但隨著目標數量的增加,該方法會出現算法結構復雜、計算量大問題.文獻[11]提出一種基于粒子濾波的多目標TBD 算法,該方法假設最大的目標數量,用同樣數量的樣本表示每個目標的狀態和存在情況,再聯合估計所有目標的聯合后驗概率密度,據此估計各個時刻的目標狀態和數量.隨著目標數量的增加,該方法中的粒子數將成倍增加,因此也存在計算量大的問題.

PHD 和多目標伯努利濾波器是基于隨機有限集(Random finite set,RFS)的多目標貝葉斯濾波器[12-14]的近似.基于RFS 的多目標貝葉斯濾波器將多目標狀態和觀測視為RFS,估計此RFS 的后驗概率密度函數.文獻[12]提出基于PHD 的單傳感器TBD 方法(Probability hypothesis density based track-before-detect,PHD-TBD),實現紅外圖像多目標跟蹤.文獻[14]針對圖像觀測,提出基于多伯努利濾波的TBD 算法.

然而,無論是基于粒子濾波的TBD 方法,還是基于RFS 的TBD 方法,都將多目標視為整體,估計多個目標的聯合狀態的后驗概率密度函數.隨著目標數量的增加,必然會出現算法結構復雜、計算量大問題.針對上述問題,在單傳感器條件下,本文提出一種基于代價參考粒子濾波器組的多目標TBD方法(Cost-reference particle filter bank bas-ed multi-target track-before-detect,CRPFB-MTBD).該方法是文獻[15]基于代價參考粒子濾波器組的單目標TBD 方法的拓展.CRPFB-MTBD 首先將多目標檢測與跟蹤問題轉化為多個單目標檢測與跟蹤問題,序貫地估計各個目標的狀態序列;然后,通過比較各個目標航跡之間的距離來刪減多余目標,從而最終確定目標數量和航跡.該方法無需假設目標數量,算法結構與目標數量無關,可并行執行,運行時間極短.仿真結果表明,與基于粒子濾波和RFS 的TBD 算法相比,本文方法可有效估計時變的多目標數量和狀態,且運行時間極短.

1 狀態空間模型

1.1 目標模型

假設k時刻有Nk個目標以特定的強度水平在x-y平面運動,離散運動模型為:

式中,vl,k表示k時刻第l個目標的過程噪聲服從零均值協方差為Q的高斯分布,且在幀間、分辨單元間,各個目標間相互獨立;k表示觀測時刻,k=1,···,K;xl,k表示k時刻第l個目標的狀態向量:

式中,[·]T表示向量的轉置.(xl,k,yl,k) 和分別表示k時刻第l個目標的位置和速度,l=1,···,Nk,k時刻可能同時存在Nk個目標.F表示轉移矩陣,定義如下:

式中,△T表示采樣時間,?表示克羅內克積.Q定義如下:

式中,q1表示目標運動的噪聲水平.由于本文方法無需目標強度信息,故在目標模型中并未包含對目標強度信息的估計,但傳感器模型中包含目標強度.為方便起見,定義第l個目標在k時刻的強度為Il,k.

1.2 傳感器模型

紅外傳感器提供某一監測區域的二維圖像序列,每幅圖像包含N×M個分辨單元(像素).一個分辨單元對應一個△x×△y的矩形觀測區域,第 (i,j) 個分辨單元對應的觀測區為(i△x×j△y),i=1,···,N,j=1,···,M.

以△T為間隔記錄測量圖像,k時刻,第(i,j)個分辨單元的測量強度為:

式中,Sk(l)∈{0, 1} 表示k時刻第l個目標出現或消失狀態;表示k時刻第l個目標對于第 (i,j) 個分辨單元的強度水平的貢獻,k時刻第(i,j)個分辨單元的強度水平是所有目標的強度水平之和;是第 (i,j) 個分辨單元的背景噪聲,假設其在幀間、像素間相互獨立且服從零均值、方差為σ2的高斯分布.本文采用傳感器擴散模型,近似為式(6),表示強度為Il,k、位置為(xl,k,yl,k) 的點目標對于第 (i,j) 個分辨單元的強度貢獻:

式中,Σ 是已知參數,表示傳感器的模糊系數.

信噪比(Signal-to-noise ratio,SNR)定義為:

綜上,從1 到K時刻的觀測序列為: 1≤i ≤N,1≤j ≤M,1≤k ≤K},多目標跟蹤的目標是從觀測序列ZK中估計目標的數量,并估計每個目標的狀態序列,先驗信息為 0≤xl,k ≤N△x,0≤yl,k ≤M△y.

2 基于代價參考粒子濾波器組的多目標檢測前跟蹤算法

針對線性或分段線性運動模型,文獻[15]提出CRPFB,實現了單目標的快速檢測和跟蹤.本節基于CRPFB,提出未知數量的多目標TBD 方法.

2.1 代價參考粒子濾波器組

基于式(1)和式(5)的狀態空間模型,CRPFB基本結構如圖1 所示.CRPFB 包含Ms個并行的代價參考粒子濾波器(Cost-reference particle filter,CRPF).每個CRPF 采用不同的先驗信息,獲得不同估計結果比較各個CRPF的累積代價將累積代價最大(或最小,與累積代價的定義有關)的CRPF 的估計結果作為CRPFB 的估計結果.CRPFB 中并行CRPF 數量Ms是一個較難確定的參數,類似于粒子濾波中的樣本數量(粒子數).因此,CRPFB 性能會隨Ms的增加而增加,但不會一直增加.本文通過對比不同Ms對算法性能的影響,來選擇合適的Ms數量.

圖1 CRPFB 的基本結構Fig.1 Basic structure of CRPFB

步驟1.先驗信息計算.圖1 中藍色背景框部分表示第ms個CRPF 的先驗信息,其計算原理及過程如下.

基于式(1)中的線性運動模型(系統噪聲為零均值),第l個目標在k時刻的狀態為其在x方向的運動速度的期望為常數:

進一步地,可從各個時刻第l個目標在x方向的位置來估計:

由式(9)可得,當xl,1=xms,xms ～U(0,N△x),U(·,·) 表示均勻分布,考慮到 0≤xl,K ≤N△x,可得:

同理,當yl,1=yms,yms ～U(0,M△y),考慮到0≤yl,K ≤M△y,第l個目標在y方向的速度均值的范圍為:

基于式(10)和式(11)的假設和結論,將速度的均值近似為速度,則第l個目標的先驗信息如表1所示.

表1 第 l 個目標的先驗信息Table 1 Apriori information for the l -th target

圖2 對比了原始先驗信息(0≤xl,k ≤N△x,0≤yl,k ≤M△y)與表1 先驗信息.在圖2 中,矩形框為觀測區域,對應第l個目標(包括每個目標)原始先驗信息,陰影部分面積為表1 先驗信息.由圖2 可以看出,基于上述結論,當初始時刻目標位置確定時,目標的先驗信息精度可大幅提高.

圖2 原始先驗信息與表1 先驗信息的對比Fig.2 Comparison of original prior information and the prior information in table 1

步驟2.代價參考粒子濾波器.假設Ms個不同的目標初始位置,根據表1 可獲得Ms種精確的先驗信息.文獻[15]指出,當Ms足夠大時,CRPF 的樣本數可取1.將觀測序列ZK和Ms種先驗信息輸入Ms個 CRPF 中,可獲得Ms種估計結果.各個CRPF在運行中相互獨立,只在濾波結束后比較累積代價,因此CRPFB 具有完全的并行結構,其運行時間僅由單個CRPF 運行時間和累積代價的比較過程決定.

圖1 中橙色背景框部分為第ms個CRPF 的濾波過程,僅有1 個目標,樣本數設置為1 時,CRPF的估計過程如下.

步驟3.CRPFB 的狀態估計.圖1 中的綠色背景框部分即為CRPFB 的狀態估計過程,其具體過程如下.

Ms個CRPF 共獲得Ms個估計結果{,···,}.比較Ms個估計結果的累積代價,將累積代價最大的CRPF 的狀態序列作為CRPFB 的估計結果,對應的累積代價記為Ccum:

與一般粒子濾波算法相比,CRPFB 中的CRPF 僅有一個樣本,故無需重采樣過程,可并行運行大量具有不同初始狀態的CRPF,實現對目標狀態序列的估計.

2.2 基于CRPFB-MTBD 算法

基于上述的CRPFB,本文提出CRPFB-MTBD 算法.該算法包括3 個連續過程,如圖3 所示.1)獲得所有可能目標的狀態序列.當全部觀測數據輸入后,序貫地應用CRPFB,獲得所有可能目標的狀態序列;2)關聯過程(估計目標數量).比較所有可能目標航跡間的歐氏距離,刪除距離相近航跡中累積代價較小者,從而確定觀測時段內出現的總體目標數量;3)判斷各個目標出現的具體時刻.根據各狀態序列累積代價的幅度變化特點,確定目標出現的具體時刻.

圖3 CRPFB-MTBD 算法基本框架Fig.3 Basic structure of CRPFB-MTBD

1)獲得所有可能狀態序列

圖4 是基于CRPFB-MTBD 算法的步驟1,估計所有可能的目標航跡.在該過程中,首先輸入觀測數據ZK,執行第1 次CRPFB.若CRPFB 的累積代價大于門限V1,則從ZK中減去CRPFB 的估計觀測,再次輸入CRPFB.此過程可執行若干次,直到CRPFB 的累積代價小于等于V1.當虛警概率Pfa給定,V1可通過對大量無目標的觀測累積代價排序獲得.圖4 中,表示第l次CRPFB 的估計結果,l個目標的估計觀測之和,Cl,cum表示第l個CRPF的估計狀態序列的累積代價.若上述過程共經過l次CRPFB 濾波,則該步驟輸出l個估計結果{,},即為所有可能目標的狀態序列.表示第1 到第

圖4 估計可能的目標狀態序列Fig.4 Estimation of all possible targets' state sequences

2)關聯過程(估計目標數量和相應的狀態序列)

圖5 是基于CRPFB-MTBD 算法的第2 個步驟,估計觀測時段內總的目標數量.將1)中獲得的所有可能目標的狀態序列 {,,···,} 排列組合,計算各估計結果之間的歐氏距離.d()即為第li個估計結果與第π(li)個估計結果的歐氏距離,如式(14)所示.當2 個的狀態序列的歐氏距離小于等于給定門限V2時,刪除累積代價較小的的狀態序列.門限V2會隨應用場景和分辨率等的不同而變化,因此難以確定.本文通過仿真實驗,對比不同V2時算法的性能,來確定門限V2值.剩余航跡的個數即為觀測時段內出現的目標數量.圖5 中,即為lest個目標的可能的狀態序列:

圖5 估計目標數量Fig.5 Estimation of target number

3)判斷各個目標出現的具體時刻

基于CRPFB-MTBD 算法假設各個目標在各個時刻都存在,而在實際中,各個目標出現的時刻可能不同.因此,本步驟判斷各個目標存在的具體時刻,圖6 是判斷各個目標出現和消失的具體時刻.

圖6 判斷各個目標存在的具體時刻Fig.6 Determination of the specific moments when each target existing

首先,輸入2)估計到的第l個目標的估計狀態和累積代價Cl,其中Cl={Cl,1,Cl,2,···,Cl,K},l=1, 2,···,lest,Cl,k為第l個目標在第k時刻的代價.

然后,比較Cl,cum與 min{Cl,cum}、Cl,cum與max{Cl,cum},以確定第l個目標出現和消失的時刻.其中表示第l個目標在第k時刻的累積代價.min{Cl,cum} 和 max{Cl,cum} 分別表示Cl,cum中的最小值和最大值.當 |Cl,cum-min{Cl,cum}|小于等于V3時,沒有目標回波累積,即k1,k2,···,kd時刻目標不存在;當 |Cl,cum-max{Cl,cum}|小于等于V3時,目標回波也沒有累積,即kg+1,···,K時刻目標未出現.

最后,輸出第l個目標的存在時刻kd+1,···,kg和相應的狀態序列.

經過上述3 個連續過程,即可基于CRPFB 估計未知數量的目標和狀態序列.其中,過程1)和3)可并行執行,能夠極大縮短計算時間.

3 仿真實驗及分析

仿真實驗包括2 個部分: 1)比較本文基于CRPFB-MTBD 算法與基于傳統粒子濾波的多目標檢測前跟蹤算法(Particle filter based multi-target track-before-detect,PF-MTBD)[11]、PHD-TBD[12]和基于伯努利濾波的TBD 方法(Bernoulli based track-before-detect,Bernoulli-TBD)[14]的性能;2)對影響CRPFB-MTBD 性能的因素進行分析.

3.1 CRPFB-MTBD 與現有方法的性能對比

在仿真實驗1)中,通過檢測和跟蹤變化數量的目標來比較CRPFB-MTBD、PF-MTBD、PHDTBD 和Bernoulli-TBD 在目標數量估計、目標狀態序列估計、運行速度等方面性能.其中,目標數量估計和目標狀態序列估計通過最優次模式分配準則(Optimal subpattern assignment metric,OSPA)[16-17]來評判.在多目標跟蹤中,OSPA 評估觀測時段內估計的目標狀態序列集合與真實目標狀態序列集合之間的差距,包括OSPA 位置誤差、OSPA勢誤差和兩者的組合OSPA 總體誤差.具體計算公式和Matlab 代碼可參見文獻[16-17].此外,通過相同平臺,對算法的運行速度與各算法的平均單次運行時間進行比較.

仿真環境設置如下: 假設觀測區域內存在3 個目標,其初始狀態分別為x1,1=[17, 0, 13,-0.13]T,x2,1=[4, 0.2, 4, 0.2]T,x3,1=[3, 0.2, 17, 0]T,q1=0.001.目標強度由式(7)的SNR 決定,q2=0.010.整個檢測和跟蹤過程持續時間為60 s,采樣時間△T=1 s,像素分辨單元△x=△y=1 m,監測區域大小為 3 0 m×30 m 的序列圖像,傳感器模糊系數Σ=0.7.假設背景噪聲為方差σ2=1 的高斯噪聲,目標的運動方程如式(1)所示.目標1 在第5 幀進入監視區域,在第60 幀消失;目標2 在第15 幀進入監視區域,在第55 幀消失;目標3 在第10 幀進入監視區域,在第40 幀消失.圖7 是當SNR=10 dB時,第2 幀(無目標)、第14 幀(1 個目標)、第35 幀(3 個目標)和第55 幀(2 個目標)的一次觀測圖.由圖7 可以看出,當SNR=10 dB 時,很難直接從圖像上判斷目標的數量和位置.

各算法參數和條件設置如下: 1) PF-MTBD.設最大目標數量為5,每個目標的樣本數為10 000;2) PHD-TBD.每幀中搜索新生目標的粒子數為2 500,每個期望目標對應的樣本數為2 000,其他參數設置與文獻[12]相同;3) Bernoulli-TBD.設最大目標數量為100,每個伯努利濾波器中生存樣本數量為5 000,新生樣本數量為2 000,存在概率門限為0.5;4) CRPFB-MTBD.CRPF 數量Ms=3 0 000,每個CRPF只采用1個粒子,當Pfa=10-3時,門限V1=101.097 9,V2=6.000 0,V3=8.884 0.在信噪比分別為8 dB、6 dB 時,分別進行100 次蒙特卡洛仿真,對比4 種方法的OSPA 結果.OSPA 參數為c=5,p=2.

圖8、圖9 分別表示當SNR 為6 dB、8 dB,3個目標時,4 種方法的OSPA.圖8(a)和圖9(a)表示OSPA 總體誤差隨時間變化情況,圖8(b)和圖9(b)為OSPA 位置誤差隨時間變化情況,圖8(c)和圖9(c)為OSPA 勢誤差隨時間變化情況.由圖8 和圖9 可以看出,CRPFB-MTBD 的目標數量和目標狀態估計性能均優于PF-MTBD、PHD-TBD 和Bernoulli-TBD.由圖9 可以看出,CRPFB-MTBD 方法對于位置較近的目標有一定的區分能力.圖10 和圖11為當SNR=8 dB、3 個目標時,CRPFB-MTBD 方法的目標狀態序列估計結果和目標數量估計結果.

圖8 當SNR=6 dB,3 個目標時,4 種方法的OSPAFig.8 Comparison of OSPAs resulted from 4 algorithms when SNR=6 dB and 3 targets

圖9 當SNR=8 dB,3 個目標時,4 種方法的OSPAFig.9 Comparison of OSPAs resulted from 4 algorithms when SNR=8 dB and 3 targets

圖10 當SNR=8 dB,3 個目標時,CRPFB-MTBD 的目標狀態估計結果Fig.10 State estimation of CRPFB-MTBD when SNR=8 dB and 3 targets

圖11 當SNR=8 dB,3 個目標時,CRPFB-MTBD 的目標數量估計結果Fig.11 Estimation of target number provided by CRPFB-MTBD when SNR=8 dB and 3 targets

表2 為仿真實驗1)的4 種算法平均單次運行時間比較.其中,CRPFB-MTBD 進行了10 次CRPFB 循環(在本次實驗中,將CRPFB 的序貫運行次數設置為10 次;在實際運用中,CRPFB 的序貫運行次數不確定).由表2 可以看出,CRPFB-MTBD 的運行時間遠小于其他3 種算法.原因是PFMTBD 和PHD-TBD 將所有目標視為整體,目標數量越多,所需樣本數越大,運算量越大,且難以并行執行,故耗時長.雖然在Bernoulli-TBD 中,各個伯努利濾波器可并行執行,但每個濾波器中的樣本數量較多,故Bernoulli-TBD 的運行時間也較長.CRPFB-MTBD 算法即使每次CRPFB 中包含的CRPF 數量達到30 000 個,但30 000 個CRPF 可并行執行,每個CRPF 只需要1 個樣本,故耗時短.

表2 4 種算法的平均單次運行時間 (s)Table 2 Average single running time of 4 algorithms (s)

3.2 CRPFB-MTBD 性能影響因素分析

仿真實驗2)分析目標的數量、并行CRPF 數量和門限對CRPFB-MTBD 算法性能的影響.

仿真環境和參數設置如下: 目標數量增加至5 個.5 個目標的初始狀態分別為x1,1=[7, 0, 23, 0.1]T,x2,1=[30, 0.15, 4, 0.2]T,x3,1=[13, 0.2, 17, 0]T,x4,1=[40,-0.1, 40,-0.1]T,x5,1=[20, 0.1, 10, 0.2]T.CRPFB-MTBD 參數設置與仿真實驗1)相同.

1)目標數量對CRPFB-MTBD 算法性能的影響.圖12 為當SNR=6 dB,1、3、5 個目標時,100次CRPFB-MTBD 仿真的平均OSPA.由圖12 可以看出,隨著目標數量的增加,CRPFB-MTBD 的性能有所下降.

圖12 當SNR=6 dB 時,目標數量對CRPFB-MTBD性能的影響Fig.12 Impact of target number on the performance of CRPFB-MTBD when SNR=6 dB

2)并行CRPF 數量對CRPFB-MTBD 算法性能的影響.圖13 和圖14 分別比較3 個和5 個目標情況下,當SNR=6 dB,CRPF 的數量分別為1 000、2 000、3 000、5 000、10 000、20 000、30 000 時,100次仿真實驗的平均OSPA.由圖13 和圖14 可以看出,隨著CRPF 的數量的增加,CRPFB-MTBD 性能逐步提升.當CRPF 的數量超過5 000 時,性能逐漸穩定;隨CRPF 數量的增加,仍有小幅提高.但針對某個具體問題,CRPF 數量如何設置,仍有待進一步研究,此問題類似粒子濾波中粒子數量的設置.

圖13 當SNR=6 dB,3 個目標時,CRPF 數量對CRPFB-MTBD 性能的影響Fig.13 Impact of the number of CRPFs on the performance of CRPFB-MTBD when SNR=6 dB and 3 targets

圖14 當SNR=6 dB,5 個目標時,CRPF 數量對CRPFB-MTBD 性能的影響Fig.14 Impact of the number of CRPFs on the performance of CRPFB-MTBD when SNR=6 dB and 5 targets

3) 門限對CRPFB-MTBD 算法性能的影響.在本次仿真實驗中,CRPFB-MTBD 算法的參數設置與仿真實驗1)相同.如第2.2 節所述,CRPFBMTBD 包含V1、V2、V3三個檢測門限.其中V1和V3可根據觀測給定虛警估計,V2是經驗值.門限V2至關重要,一方面用于確定觀測時段內出現的目標數量,另一方面用于確定目標的可能狀態序列.在仿真實驗1)中,V2設置為6.圖15、圖16 比較了當3 個和5 個目標情況下,V2為4, 6, 8, 10 時,100 次仿真的平均OSPA.由圖15 和圖16 可以看出,針對本文設置的仿真情況,當V2=6 時,算法的勢估計性能(即目標數量估計性能) 最好;當V2=10時,算法的位置估計誤差(即目標狀態序列估計性能)最好.綜合考慮勢誤差和位置誤差,故設置V2=6.

圖15 當SNR=6 dB,3 個目標時,門限 V2 對CRPFB-MTBD 性能的影響Fig.15 Impact of threshold V2 on the performance of CRPFB-MTBD when SNR=6 dB and 3 targets

圖16 當SNR=6 dB,5 個目標時,門限 V2 對CRPFB-MTBD 性能的影響Fig.16 Impact of threshold V2 on the performance of CRPFB-MTBD when SNR=6 dB and 5 targets

仿真實驗結果表明,與現有經典方法相比較,本文提出的CRPFB-MTBD 算法的估計性能、運算速度均有明顯提高.實驗結果表明,隨著目標數量的增加和信噪比的下降,CRPFB-MTBD 算法性能會有所下降.此外,并行CRPF 的數量和門限V2的選擇,也會對CRPFB-MTBD 的性能有所影響.

4 結束語

針對圖像序列中多目標檢測問題,基于CRPFB提出一種新穎的低信噪比多目標估計方法.與現有方法相比較,該方法將多目標估計問題轉化為多個單目標檢測和估計問題,打破了現有的將多目標視為整體的思路.與現有幾種經典算法進行對比實驗,結果表明: 與性能最佳的PHD-TBD 相比,CRPFBMTBD 的目標數量估計和位置估計精度提高了約50%,如圖8 和圖9 所示;與運行速度最快的Bernoulli-TBD 相比,CR-PFB-MTBD 運行速度提高了約600 倍,如圖2 所示.由不同影響因素的仿真實驗結果可以看出: 1)隨著目標數量的增加,CRPFBMTBD 的估計性能會下降;2)并行CRPF 數量對CRPFB-MTBD 性能影響很大,只有在CRPF 的數量足夠大時,CRPFB-MTBD性能才會趨于穩定;3)門限V2的選擇直接影響CR-PFB-MTBD 性能.上述圖8～ 圖16 和表2 仿真實驗結果表明,針對在圖像序列中檢測和跟蹤較低信噪比的多目標問題,本文提出的CRPFB-MTBD 算法以極短的運行時間獲得了檢測性能、估計性能的可觀改善.但CRPFB-MTBD 算法采用仿真分析方法確定并行CRPF 數量和門限V2等參數,這不利于該算法的廣泛應用.在后續研究中,將進一步研究這些參數與應用場景的關系,進而確定參數的選取原則.