Joint Target Detection, Tracking and Classification with Forward-Backward PHD Smoothing

TL;DR

This paper introduces a generalized forward-backward PHD smoothing method that incorporates target class features, significantly improving joint detection, tracking, and classification performance in dense clutter environments.

Contribution

It extends the PHD smoothing framework to include mode features like class, enhancing target state estimation and classification capabilities.

Findings

Outperforms existing algorithms in target state estimation.

Reduces OSPA distance by up to 40%.

Improves accuracy in dense clutter environments.

Abstract

Forward-backward Probability Hypothesis Density (PHD) smoothing is an efficient way for target tracking in dense clutter environment. Although the target class has been widely viewed as useful information to enhance the target tracking, there is no existing work in literature which incorporates the feature information into PHD smoothing. In this paper, we generalized the PHD smoothing by extending the general mode, which includes kinematic mode, class mode or their combinations etc., to forward-backward PHD filter. Through a top-down method, the general mode augmented forward-backward PHD smoothing is derived. The evaluation results show that our approach out-performs the state-of-art joint detection, tracking and classification algorithm in target state estimation, number estimation and classification. The reduction of OSPA distance is up to 40%.

Figures

Equations

$$\frac{X_k\text{buildrel}\Delta }{=\left\{x_{k,1},x_{k,2},\cdots x_{k,N_k}\right\}\in \mathcal{F}\left(E_s\right)}$$

$$\frac{Z_k\text{buildrel}\Delta }{=\left\{z_{k,1},z_{k,2},\cdots z_{k,M_k}\right\}\in \mathcal{F}\left(E_o\right)}$$

$${\mathrm{Z}}_{1:k}=\bigcup _{i=1}^k{Z}_i$$

$$D_{k|k-1}(X_k|{\mathrm{Z}}_{1:k-1})=\int f_{k|k-1}(X_k|X)D_{k-1|k-1}(X|{\mathrm{Z}}_{1:k-1})\delta X$$

$$D_{k|k}(X_k|Z_{1:k})=\frac{h_k(Z_k|X_k)D_{k|k-1}(X_k|Z_{1:k-1})}{\int h_k(Z_k|X)D_{k+k-1}(X|Z_{1:k-1})\delta X}$$

$$D_{t|k}(X_t|Z_{1:k})=D_{t|t}(X_t|Z_{1:t})\int f_{t+1|t}(X_{t+1}|X_t)\frac{D_{t+1|k}(X_{t+1}|Z_{1:k})}{D_{t+1|t}(X_{t+1}|Z_{1:t})}\delta X_{t+1}$$

$$D_{k+1|k}(x)={\gamma }_{k+1|k}(x)+\int F_{k+1|k}(x|x^{\prime })\cdot D_{k|k}\left(x^{\prime }\right)d{x}^{\prime }$$

$$D_{k+1|k+1}(x)=D_{k+1|k}(x)\left[1-P_D(x)+P_D(x)\sum _{z\in Z}\frac{L_z(x)}{{\kappa }_k(z)+D_{k+1|k}\left[P_D(x)L_z(x)\right]}\right],$$

$$D_{t|k}(x_t|Z_{1:k})=D_{t|t}(x_t|Z_{1:t})\left[P_s(x_t)\int \frac{D_{t+1|k}(x_{t+1}|Z_{1:k})f_{t+1|t}(x_{t+1}|x_t)}{D_{t+1|t}(x_{t+1}|Z_{1:t})}d{x}_{t+1}+1-P_s(x_t)\right],$$

$$D_{k+1|k}(\ddot{x})={\gamma }_{k+1|k}(\ddot{x})+\int F_{k+1|k}(\ddot{x}|{\ddot{x}}^{\prime })D_{k|k}\left({\ddot{x}}^{\prime }\right)d{\ddot{x}}^{\prime }$$

$$D_{k+1|k+1}(\ddot{x})=D_{k+1|k}(\ddot{x})\left[1-P_D(\ddot{x})+P_D(\ddot{x})\sum _{z\in Z_{k+1}}\frac{L_z(\ddot{x})}{{\kappa }_{k+1}(z)+{\tau }_{k+1}(z)}\right],$$

$$D_{k+1|k}(x,o)={\gamma }_{k+1|k}(x,o)+\sum _{o^{\prime }}\int F_{k+1|k}(x,o|x^{\prime },o^{\prime })D_{k|k}\left(x^{\prime },o^{\prime }\right)d{x}^{\prime }.$$

$$F_{k+1|k}(x,o|x^{\prime },o^{\prime })=P_s(x^{\prime },o^{\prime })f_{k+1|k}(x,o|x^{\prime },o^{\prime })+{\beta }_{k+1|k}(x,o|x^{\prime },o^{\prime })$$

$$D_{k+1|k+1}(x,o)=D_{k+1|k}(x,o)\left[1-P_D(x,o)+P_D(x,o)\sum _{z\in Z_{k+1}}\frac{L_z(x,o)}{{\kappa }_{k+1}(z)+{\tau }_{k+1}(z)}\right]$$

$${\tau }_{k+1}(z)=\sum _o\int P_D(x,o)L_z(x,o)D_{k+1|k}(x,o)dx.$$

$$\begin{array}[]{l}{D_{k+1|k+1}}(x)=\sum\limits_{o}{{D_{k+1|k+1}}(x,o)}=1-\sum\limits_{o}{{P_{D}}(x,o){D_{k+1|k}}(x,o)}+\sum\limits_{z\in{Z_{k+1}}}{\frac{{\sum\limits_{o}{{P_{D}}(x,o){L_{z}}(x,o){D_{k+1|k}}(x,o)}}}{{{\kappa_{k+1}}(z)+{\tau_{k+1}}(z)}}}\end{array}$$

$$D_{t|k}(x_t,o|Z_{1:k})=D_{t|t}(x_t,o|Z_{1:t})\times \left[P_s(x_t,o)\sum _{o^{\prime }}\int \frac{D_{t+1|k}(x_{t+1},o^{\prime }|Z_{1:k})f_{t+1|t}(x_{t+1},o^{\prime }|x_t,o)}{D_{t+1|t}(x_{t+1},o^{\prime }|Z_{1:t})}d{x}_{t+1}+1-P_s(x_t,o)\right]$$

$$D_{t+1|t}(x_{t+1},o^{\prime }|Z_{1:t})=\gamma (x_{t+1},o^{\prime })+\sum _o\int \left[P_s(x_t,o)f_{t+1|t}(x_{t+1},o^{\prime }|x_t,o)+{\beta }_{t+1|t}(x_{t+1},o^{\prime }|x_t,o)\right]D_{t|t}(x_t,o)d{x}_t$$

$$D_{t|k}({\ddot{x}}_t)=D_{t|t}({\ddot{x}}_t)\left[P_s({\ddot{x}}_t)\int \frac{D_{t+1|k}({\ddot{x}}_{t+1})f_{t+1|t}({\ddot{x}}_{t+1}|{\ddot{x}}_t)}{D_{t+1|t}({\ddot{x}}_{t+1})}d{\ddot{x}}_{t+1}+1-P_s({\ddot{x}}_t)\right]$$

$$D_{t+1|t}({\ddot{x}}_{t+1})=\gamma ({\ddot{x}}_{t+1})+\int \left[P_s({\ddot{x}}_t)f_{t+1|t}({\ddot{x}}_{t+1}|{\ddot{x}}_t)+{\beta }_{t+1|t}({\ddot{x}}_{t+1}|{\ddot{x}}_t)\right]D_{t|t}({\ddot{x}}_t)d{\ddot{x}}_t$$

$$D_{k+1|k}(x,c)={\gamma }_{k+1|k}(x,c)+\sum _{c^{\prime }}\int F_{k+1|k}(x,c|x^{\prime },c^{\prime })D_{k|k}\left(x^{\prime },c^{\prime }\right)d{x}^{\prime }$$

$$D_{k+1|k+1}(x,c)=D_{k+1|k}(x,c)\left[1-P_D(x,c)+P_D(x,c)\sum _{\tilde{z}\in {\tilde{Z}}_{k+1}}\frac{L_z(x,c)}{{\kappa }_{k+1}(\tilde{z})+{\tau }_{k+1}(\tilde{z})}\right]$$

$$F_{k+1|k}(x,c|x^{\prime },c^{\prime })=P_s(x^{\prime },c^{\prime })f_{k+1|k}(x,c|x^{\prime },c^{\prime })+{\beta }_{k+1|k}(x,c|x^{\prime },c^{\prime })$$

$${\tau }_{k+1}(\tilde{z})=\sum _c\int P_D(x,c)L_{\tilde{z}}(x,c)D_{k+1|k}(x,c)dx$$

$${f_{k+1|k}}(c|c^{\prime})=\left\{{\begin{array}[]{*{20}{c}}{0,~{}~{}if~{}~{}c\neq c^{\prime}}\\ {1,~{}~{}if~{}~{}c=c^{\prime}}\end{array}}\right.$$

$$\begin{array}[]{l}{D_{k+1|k}}(x,c)={\gamma_{k+1|k}}(x,c)+\sum\limits_{c^{\prime}}{\int{{\beta_{k+1|k}}(x,c|x^{\prime},c^{\prime}){D_{k|k}}\left({x^{\prime},c^{\prime}}\right)dx^{\prime}}}+\int{{P_{s}}(x^{\prime},c){f_{k+1|k}}(x,c|x^{\prime},c)dx^{\prime}}\end{array}$$

$$\begin{array}[]{l}{D_{k+1|k+1}}(x,c)={D_{k+1|k}}(x,c)\times{\kern 1.0pt}\left[{1-{P_{D}}(x,c)+{P_{D}}(x,c)\sum\limits_{\tilde{z}\in{{\tilde{Z}}_{k+1}}}{\frac{{{L_{\tilde{z}}}(x,c)}}{{{\kappa_{k+1}}(\tilde{z})+\sum\limits_{c}{\int{{P_{D}}(x,c){L_{\tilde{z}}}(x,c){D_{k+1|k}}(x,c)}dx}}}}}\right]\end{array}$$

$$D_{k+1|k+1}(x)=\sum _c{D}_{k+1|k+1}(x,c)$$

$$D_{t|k}(x_t,c|Z_{1:k})=D_{t|t}(x_t,c|Z_{1:t})\left[P_s(x_t,c)\sum _{c^{\prime }}\int \frac{D_{t+1|k}(x_{t+1},c^{\prime }|Z_{1:k})f_{t+1|t}(x_{t+1},c^{\prime }|x_t,c)}{D_{t+1|t}(x_{t+1},c^{\prime }|Z_{1:t})}d{x}_{t+1}+1-P_s(x_t,c)\right],$$

$$\begin{array}[]{l}{D_{t+1|t}}({x_{t+1}},c^{\prime}=c|{Z_{1:t}})=\gamma({x_{t+1}},c^{\prime})+\int{{P_{s}}({x_{t}},c){f_{t+1|t}}({x_{t+1}}|{x_{t}},c){D_{t|t}}({x_{t}},c|{Z_{1:t}})d{x_{t}}}\\ {\kern 1.0pt}{\kern 1.0pt}+\sum\limits_{c^{\prime\prime}}{\int{{\beta_{t+1|t}}({x_{t+1}},c^{\prime}=c|{x_{t}},c^{\prime\prime}){D_{t|t}}({x_{t}},c^{\prime\prime}|{Z_{1:t}})d{x_{t}}}}.\end{array}$$

Full text

I Introduction

Multi-target tracking is a challenging problem for surveillance system. It is not only need estimate a set of dynamic target state but also targets numbers and types. Recently, probability hypothesis density (PHD) filter has received much attention [1, 2, 3, 4, 5, 6]. Compared with traditional association-based multi-target tracking approaches like joint probability data association (JPDA) and multiple hypothesis tracking (MHT), PHD filter avoids the data association between the measurement and track due to measurements and false alarms represented as random sets. The PHD filter is capable of dealing with the problem of target birth, spawning, disappearance in dense clutters.

Currently, there are two kinds of PHD implementations: sequential Monte Carlo PHD [3] and Gaussian mixture PHD filter [4]. Compared with Gaussian mixture PHD filter, the sequential Monte Carlo approximation is capable to handle the non-linear non-Gaussian problems. To track maneuvering target, multi-model sequential Monte Carlo PHD filter are introduced in . Smoothing technique uses the observation in the future to improve the current state estimation precision.  [7] derives the PHD smoothing based on finite set statistics method, respectively. It also can be extended to multi-model cases and the multi-model PHD smoothing is obtained. However, Mahler points out that existing works have adopted a bottom-up theoretical approach. That is, they take the PHD filter or the CPHD filter as starting point, and then attempt to generalize it to jump-Markov systems. None have adopted a theoretically top-down approach, which begin with the multi-target Bayes filter as the starting point; generalize it to a multi-target jump-Markov filter; and, then and only then, derive PHD filter equations from this generalized filter. As a result, it is unclear whether any of these proposed jump- Markov PHD filters are fully rigorous from a multi-target-statistics point of view.  [8] derived a JMNS version of the PHD filter on multi-target jump-Markov systems through a top-down method.

The feature or class information of the target is useful to improve the performance of tracking [9]. By incorporating the feature information to particle filter, JPDA and MHT, these methods are more efficient to track closely spaced parallel moving target or crossing moving targets from different classes. Random finite set (RFS) theory also provides an efficient tool to incorporate the feature information. If the feature measurements are considered, RFS-based PHD filter can also be applied to joint detection, tracking and classification [10, 11, 12, 13]. Yang [1] proposed to assign a class-matched PHD-like filter to each type of target, which has a class-dependent kinematic model set to describe the kinematic feature of targets precisely.

PHD smoothing improves the performance of PHD filter. However, to the best of our knowledge, the PHD smoothing with classification information has not been considered yet in literature. Therefore, we try to deal with this issue in this paper and derive the classification-aided PHD smoothing. Specifically, we derived the general feature conditioned forward-backward PHD smoothing through the “top-down” method.

In this paper, we make the following contributions:

• •

The forward-backward PHD-JDTC smoothing is proposed using a ”top-down” approach, which can effectively deal with the multi-target joint detection, tracking and classification problem.

• •

The Sequencial Monte Carlo (SMC) implementation for forward-backward PHD-JDTC smoothing is presented and the stucture of PHD-JDTC smoothing is analyized.

• •

We propose to utilize the signal amplitude of the unknown SNR target for JDTC problem and avoiding the a priori information of the target average SNR. The simulation results show that our approach can make decision on class information effectively and outperform the state-of-art approach — PHD-JDTC filter in target tracking accuracy.

This paper is organized as follows. First, Section II provides a brief review of the forward-backward PHD smoothing based on Random Finite Set theory. It also presents the forward PHD filter with class-dependent kinematic model set. In Section III, PHD filter and smoothing with general feature are derived from a “top-down” approach. In Section  IV, the forward-backward PHD-JDTC smoothing is derived by considering the class information. In addition, its Sequential Monte Carlo implementation is provided. In Section V, simulation case is designed and the evaluation results shows the performance of proposed approach. Finally, conclusion and future work are given in Section VI.

II RFS-based forward-backward PHD smoothing

In multi-target tracking, both the target number and states are random, as well as the number of measurements and the measurements themselves. Therefore the states and measurements could be model by Random Finite Set. This section reviews the forward-backward PHD smoothing, which provides the basis for deriving the forward-backward PHD smoothing with augmented general mode. In Subsection II-A, the Random Finite Set are used to model the multi-target states and observations. The multi-target Bayes filter and smoothing are given. Its first-order approximate, PHD smoothing algorithm is outlined in Subsection II-B.

II-A Multi-target Bayes forward-backward smoothing

Assume that at time $k$, single target state $x_{k}$ belongs to the state space $E_{s}$, i.e. $x_{k}\in E_{s}$, then the multi-target states can be defined as follows:

[TABLE]

Suppose that the single target observation space is $E_{o}$, single target observation at time $k$ is $z_{k}\in E_{o}$ , the multi-states is given by

[TABLE]

where $M_{k}$ and $N_{k}$ are the number of targets and their measurements, respectively, while ${\cal F}({E_{s}})$and ${\cal F}({E_{o}})$ are the finite subsets of $E_{s}$ and $E_{o}$, respectively.

For the measurements from time 1 to $k$, the aggregate of all of them can be expressed as

[TABLE]

In the framework of FISST (Finite Set Statistics), the uncertainty for multi-target states and the corresponding observations are represented by random finite sets.

Forward-backward smoothing consists of forward filtering followed by backward smoothing. In the forward filtering, the posterior density is propagated forward to time k via Bayes recursion. In the backward smoothing step, the smoothed density is propagated backward, from time k¡k, via the backward smoothing recursion. Analogy to the single target Bayes filter predictor, multi-target forward predication is calculated as follows Forward-backward smoothing consists of forward filtering followed by backward smoothing. In the forward filtering, the posterior density is propagated forward to time $k$ via Bayes recursion. In the backward smoothing step, the smoothed density is propagated backward, from time k¡k, via the backward smoothing recursion. Analogy to the single target Bayes filter predictor, multi-target forward predication is calculate as follows,

[TABLE]

where $\int$ and $\delta$ represents random set integral and differential. ${f_{k|k-1}}({X_{k}}|X)$ is the multi-target Makov density.

With the measurements from time $k$, multi-target forward update is given by,

[TABLE]

where ${h_{k}}({Z_{k}}|{X_{k}})$ is the multisource likelihood function.

The smoothed multi-target density is propagated backward, from time $k$ to $k^{\prime}<k$, via the multi-target backward smoothing recursion.

[TABLE]

II-B PHD filter and smoothing

For multi-target Bayes filter, it is intractable to implement in a computational manner. Under the assumption that no target generates more than one measurement and each measurement is generated by no more than a single target, all measurements are conditionally independent of target state, missed detections, and a multi-object Poisson false alarm process, Mahler [10] proposed first-order multi-target moment approximation for multi-target Bayes filter- Probability Hypothesis Density (PHD). Given any region $S$ of single-target state space $X_{0}$, the integral $\int_{S}{{D_{k|k}}(x)dx}$ is the expected number of targets in $S$. In particular, if $S=X_{0}$ is the entire state space then ${N_{k|k}}\buildrel\Delta\over{=}\int{{D_{k|k}}(x)dx}$ is the total expected number of targets in the scene.

Compared with optimal multi-target tracking Bayes recursion, PHD filter is much easier because of its first-order multi-target moment approximation. On the other hand, the computational complexity is small as the integral of PHD filter is performed on the single target space.

PHD forward-backward smoothing could be derived from physical-space approach [14] and standard point process theory [5], respectively.

PHD forward filtering. In prediction step,

[TABLE]

where ${F_{k+1|k}}(x|x^{\prime})={P_{s}}(x^{\prime}){f_{k+1|k}}(x|x^{\prime})+{\beta_{k+1|k}}(x|x^{\prime})$, $\gamma_{k}(x)$ is the intensity of birth targets at time $k$, ${\beta_{k|k-1}}(x|x^{\prime})$ is the intensity of spawned targets. ${P_{s}}(x^{\prime})$ is the survive probability for existing targets.

The state is updatd then

[TABLE]

where $P_{D}(x)$ is the probability of detection for targets and $\kappa_{k}(z)$ is the intensity for clutters.

PHD backward smoothing. Smoothing provides more improved estimation results than filtering due to that it make use of more measurements. There are mainly three kinds of smoothing techniques: fixed-interval smoothing, fixed-point smoothing and fixed-lag smoothing. Fixed lag smoothing estimate the state at time given measurements until for fixed-time lag. Here we consider the fixed lag smoothing for multi-target Bayes backward recursion.

[TABLE]

where

${D_{t+1|t}}({x_{t+1}}|{Z_{1:t}})=\gamma({x_{t+1}})+\int{\left[{{P_{s}}({x_{t}}){f_{t+1|t}}({x_{t+1}}|{x_{t}})+{\beta_{t+1|t}}({x_{t+1}}|{x_{t}})}\right]}{D_{t|t}}({x_{t}})d{x_{t}}$. It should be noted that the backward recursion is initialized with the filtering results at the present time $k$ and stopped at time $k-L$, where $L$ is the time lag of the smoothing algorithm.

III forward-backward PHD filter with general mode

The PHD smoothing algorithm in section II only take into account the kinematic state of targets, which is not able to handle the issue of maneuvering target tracking or joint detection, tracking and classification. To tackle these problems, it is necessary to append the general mode, which might be the kinematic mode, the classification or the combination of them, to the target states. The general mode can be viewed as general jump Markov state, which is just related with the state at last time step while be independent with states of earlier time steps. In this section, the forward-backward PHD smoothing whose states extended with general jump Markov mode are given in details. First the forward extended PHD filter will be outlined in subsection III-A , with the backward extended PHD smoothing recursion derived in subsection III-B.

III-A Forward PHD filter with general mode

Proposition 1

Kinematic mode and class information could be sort as general jump Markov mode.

The single target state space consists of augmented states of the form $\ddot{x}=(x,o)$ , based on which the multi-target state have the form of $\ddot{X}=\left\{{({x_{1}},{o_{1}}),...,({x_{n}},{o_{n}})}\right\}$, where $o$ is the general mode with jump Markov feature. Substitute $x$ with $x_{d}dot$ for the PHD prediction and update equation, then we get the following expressions.

[TABLE]

where ${F_{k+1|k}}(\ddot{x}|\ddot{x}^{\prime})={P_{s}}(\ddot{x}^{\prime}){f_{k+1|k}}(\ddot{x}|\ddot{x}^{\prime})+{\beta_{k+1|k}}(\ddot{x}|\ddot{x}^{\prime})$. The prediction step is

[TABLE]

where ${\tau_{k+1}}(z)=\int{{P_{D}}(\ddot{x}){L_{z}}(\ddot{x}){D_{k+1|k}}(\ddot{x})}d\ddot{x}$.

The corresponding expanded full jump-variable notation will be:

[TABLE]

where

[TABLE]

[TABLE]

[TABLE]

Integrate out mode and the PHD for target we can get multi-target state estimation,

[TABLE]

And the expected target number is ${N_{k+1|k+1}}(x)=\int{{D_{k+1|k+1}}(x)dx}$.

III-B Backward PHD filter with general mode

Proposition 2

General mode extended PHD smoothing can be achieved.

After we get the estimation of target state at time $k$ we can use the do backward smoothing by extending $o$ to $\ddot{x}=(x,o)$

[TABLE]

Where

[TABLE]

Proof: Substitute $x$ with $\ddot{x}$ in Eq. 9, we get

[TABLE]

where

[TABLE]

Then, substitute $\ddot{x}=(x,o)$ to the above equation.

Note that if we use $r$ to replace $o$, we will get the multi-model PHD filter and multi-model PHD smoothing, which is consistant with the work in [14].

IV forward-backward PHD filter with class information

In this section, we will first introduce the PHD-JDTC filter and smoothing and then give the SMC implementation for it.

IV-A PHD-JDTC Filter

The jump Markov PHD filter is mainly designed for the case where the motion pattern of target change. In fact, however, the class information can also be viewed as a jump Markov variable. Specifically, the existing target is not changing with time, which is a special type of jump Markov variable. However, for the birth target, it might be different the parent target. For instance, a missile is launched from aircraft. In [1], the class conditional PHD-JDTC filter is derived based on the addition of point process instensity function. In this section, we will show that PHD-JDTC filter can be easily derived by the ”top-down” approach based on PHD filter with general mode.

Proposition 3

Through the general mode extended method, the joint detection, tracking and classification PHD algorithm will be achieved via a “top-down” way.

Let’s treat the class as a special kind of “mode” which is modeled as Jump-Markov model. If we augment the state $x$ with class information $c$, the class-conditioned PHD filter can be reached by substituting $x$ by $\ddot{x}=(x,c)$ within Equation 9.

[TABLE]

[TABLE]

where

[TABLE]

[TABLE]

Assume that the class of existing target does not change with time, i.e.

[TABLE]

However, for the type of spawning target, it might be a time-variant jump Markov state.

The equation now is in the following form: In the state prediction step,

[TABLE]

Observation is updated as follows,

[TABLE]

The number of target can be calculated by summarize all the classes, i.e.,

[TABLE]

Finally, we get the formulation of PHD-JDTC filter.

IV-B PHD-JDTC smoothing

Proposition 4

The class-conditioned backward PHD filter is,

[TABLE]

where

[TABLE]

Proof: we extend $x$ by adding the class information $c$, i.e. $\ddot{x}=(x,c)$ and substuting to Equation 9. We get

[TABLE]

where

[TABLE]

Assume that the type of existing target does not change with time. Using Equation  26 we can get the PHD-JDTC smoothing.

The final state estimation with smoothing will be

[TABLE]

The number of target after smooth is

[TABLE]

Remark: In the stage of smoothing, different kind of targets could be processed in the class-conditioned PHD smoothings. On the other hand, it can also be seen that there is information interaction in the proposed forward-backward PHD smoothing algorithm. In the forward filtering step, the spawning targets has mutual information exchange in the prediction phase, while different target interact with others by the joint likelihood function in update stage. It also has the information exchange in the backward smoothing due to that the spawning targets might be different with the target they spawned from.

IV-C SMC implementation

Based on the equations above, the forward-backward smoothing can be expressed with an explicit structure, from Fig. 1 we can see the mutual information exchange between targets from different classes.

Next, we will present how to implement PHD-JDTC smoothing with SMC. The filtering step can be divided into three steps: particle prediction, update and resampling.

Assume the state vector of particle with class and model information is

[TABLE]

Step 1: assume at time $k-1$, the PHD is $D_{k-1,c}$. The targets of class $c$ can be represented by $N_{k-1}^{c}$ equal weighted particles, $\{w_{k-1,c}^{s},\ddot{x}_{k-1,c}^{s}\}_{s=1}^{N_{k-1}^{c}}$, the state of targets in different classes will update based on the state transfer model.

For the existing targets, the particles are directly update their state within their own PHD filter, then the predicted class for those particles is

[TABLE]

where $c_{k-1}$ is the class for the original targets. $J_{\beta,c_{k-1}}$ is number of particles for spawn targets. Usually, the recommended distribution is set as same as the transfer probability of spawn targets $p_{\beta,k|k-1}(c|c^{\prime})$. The weight of spawn target is

[TABLE]

For the manuver target, the prediction of model $r_{\beta,k|k-1}^{s}(c_{\beta,k|k-1}^{s})$ will use the same random sampling approach. The predicted particles $\left\{{{x_{c_{\beta,k|k-1}^{s}}}}\right\}_{s=1}^{{J_{{c_{k-1}}}}}$ are sampled in distribution ${\varphi({x_{c_{\beta,k|k-1}^{s}}}|c_{\beta,k|k-1}^{s},{x_{c_{k-1}^{s}}},c_{k-1}^{s})}$ . The weight of particles then becomes

[TABLE]

After the prediction step for all classes of particles, the particles for each class $c$ is represented as $\left\{{\tilde{w}_{k|k-1}^{s},\dddot{x}_{k|k-1,c}^{s}}\right\}_{s=1}^{{{\tilde{N}}_{k|k-1,c}}}$.

Step 2: Particle state upddate

[TABLE]

where

[TABLE]

${g_{k,c,r}}(z_{k,x}^{j}|\dddot{x}_{k|k-1,c}^{s})$ is the state measurement likelyhood function, and ${h_{k,c,r}}(z_{k,f}^{j}|\dddot{x}_{k|k-1,c}^{s})$ is the class measurement likelyhood function.

Step 3:Particle resampling After the state update, the corresponding particles of each category are re-sampled within the category to avoid particle depletion. Finally we get the updated particle collection $\{\bar{w}_{k,c}^{s},\dddot{x}_{k,c}^{s}\}_{s=1}^{N_{k}^{c}}$

Step 4: Particle smoothing. We denote the Probability Hypophisis Density ${D_{t|t}}({x_{t}},{r_{t}},{c_{t}}|{Z_{1:t}})$ as follows,

[TABLE]

For target that belongs to class $c$, we have PHD as

[TABLE]

For steps $t=k-1,...,k-L$

[TABLE]

where

[TABLE]

The output of this step is a new set of particles $\{w_{k-L|k-L}^{(s)},x_{k-L|k-L}^{(s)},r_{k-L|k-L}^{(s)},c_{k-L|k-L}^{(s)}\}_{s=1}^{{L_{k-L}}}$. If we neglect the spawn targets, it is easy to see the PHD-JDTC smoothing is independent with each other across different class of targets.

Step 5: Resampling after smoothing. To avoid the depletion of particles, we resample the particles within each class again. Finally, the ”smoothed” particles belongs to class $c$ of time $k-L$ becomes $\{\bar{w}_{k-L,c}^{s},x_{k-L,c}^{s},r_{k-L,c}^{s}\}_{s=1}^{N_{k-L}^{c}}$

Step 6: State extraction and ordinal estimation. Then, the intra-class particle clustering is performed on different kinds of particles, and the motion state of the multi-target after smoothing is extracted, and all the particle weights are summed together, and the estimated number of target after smoothing is achieved.

Since the smoothing process involves a large number of interactions between particles, the computational complexity becomes higher. In order to achieve real-time requirements, the kd-tree method can be used here. K-d tree method considers the case where a large number of particles are densely distributed in space, and the ”group-to-group” method can greatly increase the computation speed within tolerated error.

V Evaluation

In order to verify the performance of the PHD-JDTC forward-backward smoothing proposed, we compare it with PHD-JDTC filter in target crossover and parallel motion scenarios.

The attribute information of the target in [1] is extracted from the signal-to-noise ratio of the radar measurement signal, which is assumed to be known. In real practice, however, the target true average signal-to-noise ratio is usually unknown or it varies within a certain interval. In [15], the likelihood function of the unknown SNR target is obtained by integrating the interval in the possible signal-to-noise ratio interval of the target. In this paper, we first divide the possible signal-to-noise ratio (SNR) of different categories of targets, and then use the likelihood function of the unknown SNR target as the clutter and the attribute likelihood function of each class target.

Suppose there are two different categories of targets in the surveillance area. The average SNR ratio of the targets is $45dB$ and $25dB$. We first divide the SNR interval into a high SNR interval $[30dB,50dB]$ and a low SNR interval $[10dB,30dB]$. Using the Rayleigh distribution amplitude model, the intensity of the signal is obtained using an envelope detector. The signal amplitude probability density of clutter and different categories of targets are

[TABLE]

[TABLE]

$a$ is the target signal amplitude, $[d_{1},d_{2}]$ defines the possible average signal-to-noise ratio interval of a certain category of targets. The false alarm probability is

[TABLE]

, the target detection probability is

[TABLE]

and $\tau$ is signal detection threshold.

The normalized probability density function over the signal detection threshold $\tau$ can be normalized, and the clutter and attribute likelihood functions of the different classes of targets become:

[TABLE]

[TABLE]

Assume a 2-D scenario with two synchronous sensors, which are located in $[1,0]km$ and $[0,0]km$, respectively. The surveillance area is $[0,2\pi]rad\times[0,15]km$. There are totally four targets which could be classified into two types. Target one (type 1) moves southeastwards with the initial position $[-10,10]km$ at time $k=0s$. At the same time, target two(type 2) and target three (type 1) moves northeastwards in parallel starting from $[-10,-10]km$ and $[-10,-10.5]km$. Target four (type 2) appears in $[1.4,8]km$ at time $120s$ and moves southwards. Target four disappeared at time $480s$ while the rest targets disappeared at $360s$. The trajectories of the targets are shown in Fig. 2.

The target state takes the form of ${X_{k}}=[x_{k},\dot{x}_{k},y_{k},\dot{y}_{k},\omega_{k}]^{T}$. The time evolution of the state $X_{k}$ of target movement is given by ${X_{k}}={F_{k}}{X_{k-1}}+{v_{k}}$. The movement models of the target is described by two models: constant-velocity (CV) model and coordinated turn (CT) model with turn rate unknown. For CV model:

[TABLE]

and ${v_{k,CV}}\sim N(\cdot;0,Q_{k,CV})$, and $dT=6s$ is the sample time step.

[TABLE]

where $l=1{m^{2}}/{s^{3}}$.

For the CT model with turning rate unknown:

[TABLE]

and ${v_{k,turn}}\sim N(\cdot;0,{Q_{k,CT}})$.

[TABLE]

where ${l_{1}}=1{m^{2}}/{s^{3}}$ and ${l_{2}}=0.1ra{d^{2}}/{s^{3}}$. The transition matrix between two motion models is \pi=\!\left[\!\!{\begin{array}[]{*{20}{c}}{0.8}&{0.2}\\ {0.2}&{0.8}\end{array}}\!\!\right].

Assuming that the observation equation for $i$th sensor in polar coordinates is

[TABLE]

The measurements observed by the sensors is shown in Fig. 3. the measurement noise follows Gaussian distribution $w_{k}^{i}\sim N(\cdot;0,R_{w}^{i})$ with the covariance of measurement from two sensors $R_{w}^{1}=R_{w}^{2}=diag([{(300m)^{2}},{\left({{\pi\mathord{\left/{\vphantom{\pi{180}}}\right.\kern-1.2pt}{180}}rad}\right)^{2}}])$.

At time $k$, the intensity of birth target ${\gamma_{k,1}}({x_{k}})$ for type 1 target is

[TABLE]

The intensity of birth target ${\gamma_{k,2}}({x_{k}})$ for type 1 target can be described by

[TABLE]

The intensity of spawn target is

[TABLE]

where ${\lambda_{\beta}}=0.01$,

[TABLE]

[TABLE]

b({c_{k}}|{c_{k-1}})=\left[{\begin{array}[]{*{20}{c}}{0.5}&{0.5}\\ {0.5}&{0.5}\end{array}}\right] The probability of target survival is ${P_{s}}=0.99$. The sequential Monte Carlo approximation implementation have 500 particles for each existing target, while 500 particles are used to initialize each new-born target.

V-A Simulation Results

Fig. 4 shows the estimated trajectories and classes with forward-backward PHD-JDTC smoothing. The proposed algorithm is able to filter out the state of targets, as well as their class.

The OSPA distance for different approaches are provided in Fig. 6. We find that the algorithm’s estimated OSPA distance for different categories of targets is reduced by up to 40% compared with PHD-JDTC. Therefore, the PHD-JDTC forward-backward smoothing proposed outperforms the state-of-art PHD-JDTC filter in tracking accuracy. After smooth, the OSPA distance for target with high SNR ratio is smaller compared with that of low SNR ratio. The reason is that the target with high SNR ratio is more distiguish from clutter and less affected by noise measurements. The target with low signal-to-noise ratio is close to the clutter SNR ratio interval, so the capability to distinguish between target and noise is poor.

The estimation of target number in different categories after smooth is basically the same as that of after smoothing. This is consistent with the theoretical analysis in  [10], smoothing cannot improve the cardinality estimation accuracy. From Fig. 5, we can see that smoothing has hysteresis properties. The disappearance of target occurs beforehand, causing the target number to be inaccurate during target disappears.

Our result also shows that the introduction of category information can effectively deal with the joint target detection, tracking and classification even target has crossover or parallel motions. It is noted that the time complexity of our approach is smaller than the regular PHD smoothing, which does not consider class information.

VI Conclusion and Future Work

In this paper, we derived the generalized PHD smoothing characterized by extending the general mode to PHD filtering and smoothing. Through a “top-down” method, the general mode augmented forward-backward PHD filter is derived. The forward-backward PHD-JDTC filter has explicit structure and its SMC implementation is presented. The evaluation results show the improved performance of our approach in target detection, tracking and classification compared with the state-of-art PHD-JDTC filter. The simulation scenario we considered is a two dimension case. Whereas, it is quite easy to be extended to 3-D by taking into account the pitch angle measuremented by sensors. Furthermore, It is well known that the PHD filter is not able to output the track information, due to the lack of data association. Therefore, it is recommended to incorporate the track label based on our approach to have the track information maintained [16]. Moreover, Gaussian Mixture PHD is the closed form implementation of PHD filter under the assumption of linear and Gaussian noise. The corresponding Gaussian Mixture PHD-JDTC smoothing is worth further study.