Multi-measurement Target Tracking by Using Random Sampling Approach

Multi-measurement Target Tracking by Using Random Sampling Approach Wei-Feng LIU1

Zhong CHAI1

Cheng-Lin WEN1

Abstract: In multi-target tracking field, conventional algorithms supposed that target is a point source and produces at most one measurement. While with the development of modern sensor technology, a target may give multiple measurements. In this paper, we consider that targets have certain geometrical shapes and give multiple measurements and call these targets multi-measurement targets (MmTs). We first build rigid models for the targets in parameter space and then estimate their parameters using the Markov chain sampling approach. Next, we derive the moving state described by target  s centroid with our proposed equivalent measurement. When the number of targets remains unknown, under the Poisson assumption, we use the ratios of Poisson intensities to estimate the number of targets. We also define the probabilistic vectors of type (PVoT) and propose a recursive process for the PVoT. To verify the proposed algorithm, the final experiment proposes three targets, with different shapes and distributions, moving in a 2-dimension plane with constant velocity (CV). The experimental results show that the estimation of target state has an excellent precision and the shape estimation can better and stably reflect the change of target shape. Besides, the target lost rate is around 1.4 % in 500 Monte Carlo (MC) runs. Keywords:

Multi-measurement target (MmT) tracking, information fusion, shape parameter, multiple measurements,

Markov chain sampling

When targets are far from sensor, for example, hundreds of kilometer away, the resolution cell of a radar is much larger than the detected sizes of the targets. In this case, the targets are often assumed to be point sources. The conventional tracking algorithms thus suppose that the detected targets are points and each target produces at most one measurement, for instance, the classical probabilistic data association (PDA)[1] , the joint probabilistic data association (JPDA)[2] , and the multiple hypothesis tracking (MHT)[3] , etc. With the development of modern sensors such as high-resolution radar, far-infrared source detector, and sensor networks, multiple measurements could be received by the sensors. For example, for a fighter some detected sites such as wing, air intake, and front edge of head, are easier to reflect radar wave as shown in Fig. 1 and thus get detected by some advanced radars. Commonly, these kinds of targets are called extended targets. Another instance is when several fighters fly in a formation shown in Fig. 2 and the fighters cannot be discerned. Even though a fighter produces one radar measurement, there are multiple measurements from the formation. It is difficult to confirm Manuscript received October 26, 2011; revised accepted April 28, 2012 Supported by National Natural Science Foundation of China (61175030, 91016020, 61271144, 60934009), Zhejiang Provincial Natural Science Foundation (Y1101218), Hangzhou Dianzi University Scientific Research Foundation (KYS065609051) Recommended by Associate Editor Hai-Bin YU Citation: Wei-Feng Liu, Zhong Chai, Cheng-Lin Wen. Multimeasurement target tracking by using random sampling approach. Acta Automatica Sinica, 2013, 39(2): 164−174 1. Institute of Systems Science and Control Engineering, School of Automation, Hangzhou Dianzi University, Hangzhou 310018, China

the measurement originating from a fighter. These kinds of targets are referred to be group targets. The extended targets and the group targets can be essentially classified into the same problem. Here we called these two types of targets as multi-measurement targets (MmTs). In this paper, we focus on the multi-measurement target tracking. The multi-measurement target (MmT) tracking can be classified into three categories according to various theories. The first is due to the classic Kalman filtering theories. Koch et al. and Feldmann et al. proposed random matrix based approaches which estimate the centroid point while derive the ellipse of the extended target through the random matrix of the second order moment[3−5] . These approaches are easier to understand and have a theoretically solid foundation, but the main work focusses on the single extended target tracking to the best of our knowledge. Specifically, the approaches extend the traditional Kalman filtering of point target to the extended target without clutters and thus neglect the data association. Nevertheless, the data association is a key point for the extended target tracking and requires further research. Richter et al. studied multiple extended targets using the Markov chain Monte Carlo method by combining the models of target birth and death[6] . Baum et al. described the extended target using the random super-surface model and the shape problem for an extended target are changed to the curve problem[7−9] . The second class is based on the random set theory. Mahler proposed the probability hypothesis density (PHD) filter for extended target[10] . Lundquist et al. studied

Wei-Feng LIU et al./ Multi-measurement Target Tracking by Using Random Sampling Approach

the extended targets using the PHD filter. They assumed the measurements coming from the boundary of target[11] . Orguner et al. proposed the Cardinality PHD (CPHD) filter for the extended target[12] . Reference [13] adopted the sequential Monte Carlo (SMC) method of the PHD filter to estimate the states of the group targets. The random set based methods have a splendid paradigm. Their drawbacks lie in the high complexity and abstract concepts. For nonlinear problems, the linearization or nonlinear filtering techniques are generally used and the shape tracking and estimation are usually belonged to the nonlinear problem. The linearization techniques are more fit for the weakly nonlinear issue, while for a strong nonlinear problem such as trigonometric function the linearization technique can be hardly used. The SMC approach is a typical nonlinear method, where the number of particles is proportional to the number of targets. If the number of targets is known, problem may be simpler, otherwise a great number of particles will be produced for the extended PHD filter in dealing with the MmTs (a small number of particles cannot reflect the shape of the MmTs).

Fig. 1

The extended target produces multiple measurements (“+” denotes measurement).

Fig. 2

The group targets have a formation and produce

multiple measurements (“+” denotes measurement).

The third class is the pattern recognition based methods, that is, derivation of the target shape from the measurements of target[13] . These methods need more number of measurements compared with the previous two classes. In most video tracking research, the geometrical shape and MmT are often considered[14−16] , while in some conditions, getting so much measurements may be difficult. This paper fuses the first classic filtering theories and the third pattern recognition approaches. Specifically, when we estimate the

165

state of centroid for the MmT, we adopt the theories of the classic point estimation. While for the shape estimation we combine the methods of pattern recognition. This paper studies the MmT with certain geometrical shape. Its potential applications lie in the fields such as high resolution radar, multiple spectrum, synthetic aperture radar (SAR), infrared tracking, video tracking, and so on. We aim at the following problems: the MmT state (centroid) and parameter estimation, the estimation of number of targets and the computation of the probability of target type. Our main ideas are listed as follows. We describe the measurements of targets by the finite mixture models (FMM) and estimate the parameters of the FMM through the Markov chain Monte Carlo (MCMC) sampling. We import the equivalent measurement (EQM) of an MmT to derive the centroid measurement and estimate the state of the centroid through the EQM. Finally, we estimate the number of MmTs by the ratio of Poisson intensities. The random sampling of the distribution of a random variable is investigated by random trial and further one can derive the statistics of the random variable from the sampling points of the random trial, for example, the familiar Monte Carlo experiment and MCMC. Generally, a proposal distribution is adopted due to the difficultly sampling. The proposal distribution should approximate the original distribution as soon as possible. Nevertheless, the proposal distribution cannot completely reflect the true distribution. A Markov chain approach can sample a great number of particles from a distribution where with high probability. For example, the particle filter (PF) is a typical Monte Carlo method with a single-scan Markov chain[17] . Reference [18] proposed a single-scan and multi-scan MCMC association method to deal with the occlusion of multiple people in video scene. Reference [19] proposed a PF algorithm to estimate the target state by modeling a Markov random field for the target occlusion, where the proposal density is built by the MCMC. We extend the random sampling method of target point source to MmT with certain shape and multiple parameters. Though the parametric vectors have the same form as the state vectors, the nonlinear relationship exists in the parametric vectors. For example, the angle and the side-length in the MmT belong to trigonometric functions. That is, the change of shape parameters are generally nonlinear. Besides, the random sampling is effective for the parametric estimation, data association and state estimation of target centroid. This is another motivation applying the random sampling scheme. This paper is the more complete work of [20]. We define the probabilistic vectors of type (PVoT) for the MmT and give a recursive equation for the PVoT. We also propose the models of two typical MmTs and their estimation including

166

Acta Automatica Sinica, 2013, Vol. 39, No. 2

parameters and centroid states in the paper. The paper is organized as follows. Section 1 gives problem description. Section 2 proposes the tracking algorithm for the MmT with certain geometrical shape. Section 3 shows the proposed random sampling of Markov chain. Section 4 simulates the given algorithm. Last Section 5 concludes this paper.

1

Basic assumptions and problem description

Here we first give the following assumptions for the MmTs: Assumption 1. The MmTs are rigid and have certain geometric configuration. Assumption 2. The multi-measurements follow mixing distributions and the number of the multi-measurements follows Poisson distribution. Assumption 3. The birth places of targets are known but the birth times keep unknown. Assumption 4. The state points of target dynamics and measurement functions are linear Gaussian and the probability of detection is a constant. Remark 1. Assumption 1 shows that the targets are with fixed configuration and thus the number of multimeasurements is Poisson in Assumption 2. It should be noted that the distribution of the number of multimeasurements is different from the distribution of a measurement noise. The former are for multi-measurements while the latter indicates the statistic property for the measurement noise. For example, the former may be uniform while the latter be Gaussian. In Assumption 4, the state points for the MmTs include moving state of centroid point and shape parameters. With all these assumptions, the proposed problem can be formulated as follows: x jk+1 = Ajkx jk + Bkj w jk ,

(1)

y jk+1 = Ckj x jk+1 + v jk+1 , j = 1, · · · , mk ,

(2)

f (zz lk+1 |yy k+1 , x k+1 , πk+1 ) = mk  πkj fkj (zz lk+1 |yy jk+1 , x jk+1 ), j=0 mk 

l = 1, · · · , nk , (3)

πkj = 1, πkj ≥ 0,

j=0

f (|zk+1 |) = Poisson(λk+1 ),

(4)

λk+1 = λ0,k+1 + λ1,k+1 + · · · + λmk+1 ,k+1 ,

(5)

where |zk+1 | is the cardinality (number) of measurement set nk+1 j,c j,c j,c zk+1 = {zz lk+1 }l=1 . The state x jk = (pj,c x,k , p˙ x,k , py,k , p˙ y,k , j T j,c j,c j,c j,c φk ) , here {px,k , p˙x,k , py,k , p˙y,k } is centroidal position and speed for the jth MmT, φjk is shape parameter which will be defined later. λk+1 denotes the total intensity of mea-

surements (including clutter). λ0,k+1 is the intensity for clutter, λ1,k+1 , · · · , λmk+1 ,k+1 are the intensities for individual MmTs. πkj is mixing weight. j is the indicating variable, where j = 0 indicates clutter, y jk+1 is the measurement of centroid. fkj (·|yy jk+1 , x jk+1 ) is the distribution of measurements for an MmT or clutter. w jk and v jk are respectively process and measurement noise, which follow Gaussian distribution with means 0 and covariances Qj and Rj , respectively. The above equations (1) and (2) are respectively dynamic and measurements for the jth MmT. Equation (3) is the multi-measurement function, here we adopt the finite mixture models for its powerful in describing the stochastic distribution. It should be pointed that the centroid measurements cannot be always observed directly. In fact, we derive the measurement set zk+1 and use the set to estimate the centroid measurement y jk+1 .

2

The tracking algorithm for multimeasurement targets

In this section, we first introduce the finite mixture models for multi-measurements and then give the concept of equivalent measurement (EQM). Finally, we propose how to estimate the moving state of a target centroid using the EQM.

2.1

The finite mixture models of multimeasurements

ai } comes from varAssume that measurement set a = {a ious random sources. Each source follows distribution with parameter set φ = {φj }. Then, the measurements can be modeled by the FMMs as follows: ai |φ) = π1 f1 (a ai |φ1 ) + · · · + πm fm (a ai |φm ). f (a

(6)

For a fixed time k, if measurement set a is replaced by zk   = z ik , then the corresponding FMMs have the same appearance like (3). The individual components of the FMMs describe measurements of targets and clutter. For a unknown number of targets, the number of components mk needs to be estimated. Thus, we need to estimate three types of parameters including number of components mk , mixing weights πkj , and distributional parameter φjk . The expectation-maximization (EM) algorithm1 and the Bayesian theorem are often two frameworks. We adopt the latter to derive the parameters. The usual Bayesian esti1 The EM algorithm maximizes the expectation of likelihood function and also has the Bayesian similar form. Nevertheless, this form is the expectation of posterior distribution.

Wei-Feng LIU et al./ Multi-measurement Target Tracking by Using Random Sampling Approach

mation can be formulated as: p(xck , φk |z1:k ) = C −1 f (zk |xck , φk )p(φk |xck , z1:k−1 )p(xck |z1:k−1 ),

(7)

mk ,c j,c j,c j,c j,c x1,c where xck = {x }, x j,c k , · · · , xk k = (px,k , p˙ x,k , py,k , p˙ y,k ), m φ1k , · · · , φ k k }, z1:k = {z1 , · · · , zk }, p(xck , φk |z1:k ) is φk = {φ posterior distribution, C −1 is the normalization, f (·) is likelihood function, p(·) is prior distribution (function). We next focus on these two functions 2.1.1 The likelihood function f (zk |xck , φk ) The multi-measurements from MmTs are often assumed to be independent. This assumption is often reasonable for rigid targets. Therefore, the likelihood function can be derived in the following:

f (zk |xck , φk ) =

nk 

we categorize the MmTs as two types I and II. The type I is straight-side shape shown in Fig. 3, which can be described by three parameters: the state of MmT s centroid xck , the j,q distances of the lines between the centroid and sides {hk j }, and angles of the lines {θkj }. In fact, an angle is enough for the lines due to the rigidity assumption and thus with fixed angles between any two lines. The type II is referred as to be the curved-side shape shown in Fig. 4. We use the Gaussian sum to approximate the curved-side shape. Accordj,q j,q j,q ingly, the parameters of the shape are {πk j , μ k j , Σk j }.

Fig. 3

Some straight line objects whose sides are all straight lines

f (zz ik |xck , φk ) =

i=1 nk mk



167

j xj,c πkj fj (zz ik |x k , φ k ),

(8)

i=1 j=0

where nk is the number of measurements. Here a basic assumption is that a measurement z ik only comes from a random source at most. For any measurement, an indicating information which shows the corresponding relationship of the component and the measurement is lost. The lost in mk formation is described by missing variables e ik = ei,j , k j=0 mk i,j i,j which meet with ek = 0 or 1 and j=0 ek = 1. We supply the variables into observation set zk and derive the complete   observation data z ik , e ik . Therefore, the aforementioned expression (8) is reformulated as follows: f (zk |xck , φk ) =

nk mk  

i,j k

e

πkj

i,j ek

fj

j xj,c (zz ik |x k , φ k ),

(9)

i=1 j=0

Fig. 4

Some curved objects whose sides include at least a curve

1) The parameters for the first type I. As shown in Fig. 5, j,sj T the parameters are defined by h jk = [hj,1 ] and θkj k , · · · , hk j where θk is the angle of straight-side shape. We assume the two parameters follow Gaussian distribution, i.e., θkj ∼ N(ϕjk , Σj,θ k ),

(11)

ρjk , Σj,h N(ρ k ),

(12)

h jk

∼

j j,h where ϕjk , Σj,θ are respectively the prior mean and k , ρ k , Σk covariance.

where the estimation of the missing variable ei,j k is given i,j c by its conditional expectation eˆi,j = E(e |z , x k k , φk ). Let k k i,j i,j i,j c ek be 0 or 1, we have E(ek |zk , xk , φk ) = p(ek |zk , xck , φk ). According to Bayesian formulation: eˆi,j k =

πkj fj (zz ik |ei,j k , φk )

m k

j=0

πkj fj (zz ik |ei,j k , φk )

,

Fig. 5

Parameters for straight line objects

(10)

φjk ). where πkj = p(ei,j k |φ 2.1.2 The prior distribution p(φk |xck , z1:k−1 ) p(xck | z1:k−1 ) Assume that the shape parameters are φ jk = {ζkj,1 , · · · , j,s j,q ζk j }, ζk j ∈ Rqj . It is noted that a target shape can be described by different parameters. For instance, an ellipse shape can be described by its lengths of major axis and minor axis and the center position. Similarly, considering the random property of measurements, we can model the shape using the Gaussian mean μjk and covariance Σjk . For clarity,

2) The second type II. As shown in Fig. 6, assume that j,q mean ξ k j obeys Gaussian distribution: j,qj

ξk j,q

j,qj

μk ∼ N(μ

j,qj

, Σk

),

(13)

j,q

where μ k j and Σk j are the prior mean and variance. j,q j,q Variance Σk j obeys inverse Gamma distribution IG(αk j , j,qj βk ), i.e., j,qj −1

Σk j,qj

where, αk constant.

j,qj

, βk

j,qj

∼ Γ(αk

j,qj

, βk

),

(14)

> 0 are the prior parameters and keep

168

Acta Automatica Sinica, 2013, Vol. 39, No. 2

ˆk = N

mk 

ˆj N k

=

j=1

mk  j=1



 nk  lkj j = ei,j , l k k , λjk i=1

(19)

where [·] means the nearest integer of a float. Fig. 6

Parameters for curved objects

3 With the above likelihood functions and prior distributions, we can derive the posterior distribution and further propose the MCMC sampling scheme.

2.2

The EQM and state estimation for multi-target

In order to estimate the centroid states of the MmTs, the measurements of the MmTs centroid should be derived. Nevertheless, the measurements are usually unavailable. We adopt the conditional expectation of multimeasurements coming from an MmT and call this expectation the EQM[21] . That is = E(yy jk |zk , xk ), yˆj,EQM k

(15)

x1k , is the EQM of the jth MmT, xk = {x where yˆj,EQM k mk · · · , x k }. Directly getting the EQM is difficult, we thus propose the following equation: nk 

yˆj,EQM = k

i eˆi,j k zk

i=1 nk 

i=1

.

(16)

eˆi,j k

Once the EQM is derived, various estimation algorithms such as Kalman filter, extended Kalman filter, particle filter, and unscented filter can be applied. For clarity, we give the following representation: ˆ jk = SE(ˆ x y j,EQM ). k

2.3

Sampling algorithm for the shape parameters

The measurement from a point source reflects the state of the point. Thus, a point source produces a measurement, while the shape parameters of an MmT are not observable. For example, a rectangular target s centroid point {pcx,k , pcy,k }, centroid-side distance, and angle {h1,k , h2,k , θk } cannot be always observed, these parameters can be seen as an implicit state. That is, the obtained measurements do not include these parameters directly, and we need to estimate the parameters from the measurements. For the MmT with rectangular shape in Fig. 7, we need to build the connection between the measurements and the parameter variables. Suppose the measurements in the rectangular MmT are of uniform distribution P t1 , P t2 , P t3 , P t4 }) = f (zk ) = U({P U([pcx,k , pcy,k , hqk , θk ]T ),

(20)

where P t1 ∼ P t4 represent four vertices in a rectangular block. The corresponding parameters are given by [pcx,k , pcy,k , hqk , θk ]T .

(17)

Number of targets estimation

The estimation of number of targets is intractable due to multi-measurements and clutter. Given the distributions of number of targets, we can use them to estimate the number. Under Poisson assumption for the number of targets condition, the below equation is proposed Fig. 7

j

f (Ljk

(λj )t e−λk = t) = k , t!

The parameters for rectangular target

(18)

where Ljk is the number variable of multi-measurement, t is the real number of multi-measurement, λjk is Poisson intensity. The expression shows the Poisson property of the number of multi-measurements originating from the jth MmT. Therefore, we adopt the ratio of estimation of number of multi-measurements lkj and Poisson intensity λjk to evaluate whether an MmT exists (near 1 for an existing target). The total number of targets is given by:

3.1

The parameter jump-Markov process and the detection of target type

Traditional Metropolis-Hastings (MH) sampling scheme adopts the following acceptance-rejection rate to select the candidate sampling point[22] α = min





xk )Q(x xk , x k ) P (x 1,  xk )Q(x xk , x k ) P (x

,

(21)

Wei-Feng LIU et al./ Multi-measurement Target Tracking by Using Random Sampling Approach 

P ik+1|k = PM × P ik ,   diag L1,j , · · · , LM,j k+1 k+1 P ik+1|k , P ik+1 = M,j L1,j k+1 + · · · + Lk+1

xk , x k ) is the proposal density. The candidate where Q(x  sampling point x is selected as: 

x ,

α = 1,

(22)

x ,

α < a, a ∼ U([0, 1]),

(23)

x,

α ≥ a, a ∼ U([0, 1]).

(24)



As for the posterior distribution of the states of targets, the sampling scheme is given by 

α=



xk |z1:k )Q(x xk , x k ) p(x =  xk |z1:k )Q(x xk , x k ) p(x 



where Pr(j1 , j2 ) denotes the probability of the j2 th type (target) jumping to the j1 th type. The matrix PM does not require to be symmetrical. In other words, Pr(j1 , j2 ) is not equal to Pr(j2 , j1 ). Each type of target is described by the following probabilistic vectors of type (PVoT) P ik

=

i [Pk,1 ,···

i , Pk,j ,···

i , Pk,M ]T ,

(26)

i where i is the ith target, j is the jth type of shape. Pk,j is the probability of the ith target belonging to the jth type. When jump matrix PM is identity matrix, the reversible jump MCMC reduces to the MH sampler. Besides, the jump matrix meets the following condition: the sum of each row vector is 1, i.e., M 

Pr(i, j) = 1,

1,j z 1k+1 |x xik+1 ) + · · · + en,j zn xik+1 ). Li,j k+1 |x k+1 = ek p(z k p(z

Pr(i, j) ≥ 0,

(27)

j=1

So we use the max rule to detect the type of each target. Specifically, the maximum in the PVoT is the type of the target. The PVoT can be recursively proceeded by the predicted and update steps. In short, the predicted PVoT is available through the jump matrix and the updated PVoT is given by likelihood sum. The detailed process is as follows: let at step k the PVoT for the ith target be P ik . Then, the predicted PVoT and the updated PVoT are respectively given by the following two equations:

(29)

(30)

P ik+1

Note that the update PVoT should be normalized after each update step and then the next new predicted and update cycle.

3.2

MH sampler is proceeded in the same dimensional space, while the parameters of MmTs are with different dimension. Therefore, we combine Green s reversible jump MCMC sampling scheme[23] . That is, we sample in the same dimensional spaces while import the following jump matrix in different dimensional spaces: ⎡ ⎤ Pr(1, 1) Pr(1, 2) · · · Pr(1, M ) ⎢ ⎥ ⎢ Pr(2, 1) Pr(2, 2) · · · Pr(2, M ) ⎥ ⎢ ⎥ PM = ⎢ ⎥ , (25) .. .. .. .. ⎢ ⎥ . . . . ⎣ ⎦ Pr(M, 1) Pr(M, 2) · · · Pr(M, M )

(28)

where



xk , z1:k−1 )p(x xk |z1:k−1 )Q(x xk , x k ) p(zk |x  . xk , z1:k−1 )p(x xk |z1:k−1 )Q(x xk , x k ) p(zk |x

169

The optimization criterion of parameter estimation

The parameter estimation for MmTs is given under certain optimization criterion. The MH sampler derives the sampling points from the posterior distribution, or maximizing the posterior distribution. In this paper, we propose a new criterion which maximizes the sum of the number of measurements and the likelihood function T  

k  1 i OBN = max (31) Lk + i , Ski > 0, Sk i=1 where Ski is the area of the ith MmT (in fact the likelihood function under the uniform distribution in 2-dimension surface). The area can be gotten by various parameters. The explanation in physically is: to cover as many of points (measurements) using as many smallest area as possible.

3.3

The typical shape I: rectangle modeling and its parameter estimation

Based on the parameters of a rectangle plotted in Fig. 7, the block range of the rectangle can be calculated. Four P t1 , P t2 , P t3 , P t4 } are adopted to represent the vertices {P rectangular MmT. In the block, the measurements follow uniform distribution: 1 1 cS = = . P t1P t2 ||P P t1P t4 | S |P Obviously, we can build the following functions between these four vertices and the shape parameters:   xc + h2 cos(θ) − h1 sin(θ) P t1 = , (32) yc + h2 sin(θ) + h1 cos(θ)   Pt1,x + 2h1 sin(θ) , (33) P t2 = Pt1,y − 2h1 cos(θ)   xc + h1 sin(θ) − h2 cos(θ) , (34) P t3 = yc − h2 sin(θ) − h1 cos(θ)   Pt1,x − 2h2 cos(θ) , (35) P t4 = Pt1,y − 2h2 sin(θ) where P tl = [Ptl,x , Ptl,y ]T , l = 1, 2, 3, 4 are the coordinates of the four vertices. Of course, other parameters can also

170

Acta Automatica Sinica, 2013, Vol. 39, No. 2

be selected and thus the vertices and parameters may have different relationship. Assume that 20 measurements are received from a rectangular MmT, the learning process of the parameters with an initial error (angle and position) in 200 iterative steps are shown in Fig. 8. Despite the nonlinear relationship between the shape parameters and measurements, the dynamics of the shape parameters (or states) against time can be modeled by a linear function according to Assumption 1.

3.4

The typical shape II: ellipse modeling and parameter estimation

Ellipse can be described by its center and major and minor axis, or means and covariance in Gaussian distribution. Fig. 9 shows the estimation process of parameters of four ellipses with 100 measurements in 100 iterative steps. The corresponding parameters are state mean vector ζ k = [pcx,k , p˙ cx,k , pcy,k , p˙ cy,k ] and state covariance Σk and the dynamic process is ζ k+1 = Fζ ζ k + ω k,ζ ,

(38)

Σk+1 = FΣ Σk + ωk,Σ ,

(39)

where Fζ is the linear transfer matrix and FΣ is the rotation matrix. FΣ is identity when no rotation and deformation happen. ω k,ζ is the noise vector, ωk,Σ is noise matrix which is symmetric positive definite.

Fig. 8

Process of parameter estimation for a rectangular target in 200 iterative steps

ξ k+1 = Fξ ξ k + ω k ,  T ξ k = pcx,k , p˙cx,k , pcy,k , p˙cy,k , h1k , h2k , θk , ω k = [ωk,1 , ωk,2 , ωk,3 , ωk,4 , ωk,5 , ωk,6 , ωk,7 ]T , ⎤ ⎡ 1 T 0 0 0 0 0 ⎥ ⎢ ⎢ 0 1 0 0 0 0 0 ⎥ ⎥ ⎢ ⎢ 0 0 1 T 0 0 0 ⎥ ⎥ ⎢ ⎥ ⎢ Fξ = ⎢ 0 0 0 1 0 0 0 ⎥ , ⎥ ⎢ ⎢ 0 0 0 0 1 0 0 ⎥ ⎥ ⎢ ⎥ ⎢ ⎣ 0 0 0 0 0 1 0 ⎦ 0

0

0

0

0

0

(36)

Fig. 9

Process of parameter estimation for elliptical object in 100 iterative steps

(37)

1

where T is the sampling time, ξ k+1 , ξ k are parameter vectors, ω k is noise vector, h1k , h2k are the distances between side and centroid in two directions (see Fig. 7), θk is the horizontal angle. Under the rigid assumption, the shape of a target is fixed. In other words, the parameters of shape can be considered as a conditional distribution under some known values. For instance, the centroid, side length and horizontal angle can be seen as the conditional distribution of the predicted parameters. Here we adopt the conditional Gaussian distribution.

4 4.1

Simulation scenario Scenario configuration

In this section, we track three MmTs using the proposed sampling algorithm to derive the state of the MmTs. Table 1 lists their related prior state information. The three MmTs move in a region with size [50, 150] × [−150, 50] m2 , the probability of detection is PD = 0.9 and the clutter density is γ = 5.0 × 10−4 m−2 (around 20 cluttered measurements). The measurements of the MmT 1 satisfy Gaussian distribution and is with an oval shape. The measurements of Target 2 satisfy uniform distribution and with a rectangle shape, and the measurements of Target 3 uniform and with a circular shape. We assume all the MmTs survive in interval 1 ∼ 100 s.

Wei-Feng LIU et al./ Multi-measurement Target Tracking by Using Random Sampling Approach Table 1

The parameters of multi-measurement targets

Target

1

2

3

Position of birth (m)

(−30, −10)

(60, 20)

(−10, −10)

Poisson intensity λj

40

60

50

Shape

Ellipse

Rectangle

Circle

Parameters

[5, 0.2; 0.2, 8] m2

[5, 5] m

5m

4.2

The estimation and tracking process

Fig. 10 and Fig. 11 illustrate the trajectories in x-y coordinates. We plot all the MmT s measurements and cluttered measurements per 5 steps. From these two figures, we can see that the tracking process is stable. In our 500 Monte Carlo simulations, the tracking lost times is 7 and the lost rate is around 1.4 % (Once the root mean square error (RMSE) in a certain time is larger than threshold 100 in a Monte Carlo run, the lost times plus one). It can be seen from Fig. 12, which shows the parameter estimation per 5 steps, that the estimation approximates the true value listed in Table 1. The estimation of number of targets in Fig. 13 is around 2.8 through averaging 500 MCs. It is smaller than true value because we use ellipse to approximate circular shape, which makes the estimated number of the circular shape smaller and will be investigated below. The estimation fluctuates at the approach point shown in Fig. 14 and Fig. 15. At the approach point, the MmT commonly is easier to be lost. Nevertheless, the PVoT can efficiently detect the shape type of an MmT after they pass the approach point and thus avoid target loss as much as possible. Fig. 14 shows the RMSE of centroid position for the three MmTs. The tracking error is below 1 at most time and keeps stable. The error fluctuates around at time steps 69 s ∼ 70 s of the approach point due to interference. In shape parameter, the ellipse shape in this paper is described by covariance. we use 2-norm of the difference of

Fig. 11

Three tracks in x and y coordinate against time

Fig. 12

Fig. 13

Fig. 10

171

Three multi-measurement targets moving in x-y plane (This figure only shows their centroid positions.)

Estimation for the target shapes

Number of targets against time (500 MC)

two covariances to evaluate the estimated precision of the shape parameter. The upper subgraph of Fig. 15 shows the 2-norm of the estimated error of covariances against time. The 2-norm keeps stable and is below 5. The middle subgraph of Fig. 16 illustrates the estimation of two sidecentroid distance h1k , h2k against time. The tracked target is a square with side length being 5 and the estimated value is

172

Acta Automatica Sinica, 2013, Vol. 39, No. 2

slightly over 5 due to random disturbance and data volume. The bottom subgraph of Fig. 15 gives the radius estimation of Target 3. The true radius is 5 m and the estimated radius is around 4.8 m, which is somewhat smaller than the true value because the random property and finite measurements can hardly represent a circular target. That is, the ellipse shape may be more fitted for the measurements. Under the prior information of the circular shape, we estimate the radius by the area equivalence rule. Specifically, a el√ lipse area is equal to a circle area, we thus have R = ab, where a, b are respectively the minor axis and major axis of the ellipse. This may lead to minor radius estimation. The fluctuations between 69 s ∼ 71 s are similarly due to the target approach.

Fig. 14

P 10 , P 20 , P 30 } and the shape type the initial PVoT is P0 = {P transient matrix PM is as follows: ⎤ ⎤ ⎡ ⎡ 0.6 0.2 0.2 0.7 0.1 0.2 ⎥ ⎥ ⎢ ⎢ P0 = ⎣ 0.2 0.7 0.1 ⎦ , PM = ⎣ 0.2 0.6 0.2 ⎦ 0.1

4.3

0.9

0.2

0.1

0.7

We calculate the probabilities of a group of measurement set belonging to individual target types. Based on maximum rule, it can be seen from Fig. 16 that Targets 1, 2 and 3 respectively belong to Types 1, 2 and 3 with probabilities 0.8, 0.8 and 0.6 in average. The circular target following circle with uniform distribution is similar to rectangle with uniform distribution and ellipse with Gaussian distribution. Therefore, Type 3 has a relatively small value compared to Types 1 and 2. Nevertheless, all the probabilities are dominant, or larger than other types. For example, when a target belongs to Type 1 with probability 0.6, this means the probability of the target belonging to other types (Types 2 and 3) is less than 0.4. The dominant property will help to tell a target type. This can be seen from Fig. 17 that the estimated types of the MmTs are precise and stable.

Estimated RSME for the position (500 MC)

Fig. 16

Fig. 15

0

Estimated error for the shape parameters (500 MC)

Detection of target type

In this simulation, the MmTs include three types: ellipse, rectangle and circle (also seen as ellipse). We respectively label these Types 1, 2, and 3. The initial state (at time step 1) and shape parameters are listed in Table 1. Assume that

Probabilities of target type against time (500MC)

We now illustrate why the probability of type can correctly judge a target type even when two MmTs cross (maybe completely overlap). Assume that an MmT is missed due to the cross. Consider the likelihood update, the probability of the correct type of an MmT will be gradually increased and dominant after the approach point. And thus the detection of a target type can remain a correct type. In this sense, the proposed algorithm can correctly track MmTs after the approach point. To check the computational time, we tested the algorithm in Matlab on a computer with Intel Core 2 CPU, Windows 7 and 2 G Memory. We adopted the iterative process with fixed 30 times in each step. The total CPU time is 4 s in running the whole 100 time steps (or 40 ms per step). So in a data rate of 1 second condition, the real time is

Wei-Feng LIU et al./ Multi-measurement Target Tracking by Using Random Sampling Approach

acceptable even in existing of the 30-iteration procedure. As a comparison, we selected the Gaussian mixture-PHD (GMPHD) algorithm which tracks three point source targets (3 targets, 80 clutters, or the clutter density is 2 × 10−5 m−2 ) and run the GM-PHD algorithm in the same environment, where each step averagely costs above 100 ms.

Fig. 17

Target types against time (500 MC)

173

then estimated the centroid states and shape parameters by multiple measurements, where the centroid state is derived by the EQM and the shape parameters by Markov random sampling scheme. The MmT tracking inevitably involves the detection of the MmT and we proposed the PVoT to detect MmTs and their shape. Besides, we provided recursive equations for the PVoT. We finally validate the proposed algorithm by using three MmTs with different shapes and distributions. The results showed that the given algorithm is better in estimation performance and real time application. When the number of measurements of an MmT is few, for instance, less than 10, the derived measurements can hardly reflect the true shape of the MmTs and thus the estimated shape may have large error. In this case, it is feasible to approximate the shapes of MmTs using ellipses and rectangles. Besides, the Poisson distribution of the number of measurements is based on the rigid assumption. Once this assumption fails or an MmT rotation is perpendicular to the line of sight, the Poisson assumption will be broken and the estimation of the number of targets becomes difficult. Moreover, the size of MmT also changes with the distance between sensor and target. From this point of view, the Poisson assumption of the number of measurements is still strong and additional studies are necessary.

References [1] Bar-Shalom Y. Tracking methods in a multitarget environment. IEEE Transactions on Automatic Control, 1978, 23(4): 618−626 [2] Reid D B. An algorithm for tracking multiple targets. IEEE Transactions on Automatic Control, 1979, 24(6): 843−854 [3] Koch W, van Keuk G. Multiple hypothesis track maintenance with possibly unresolved measurements. IEEE Transactions on Aerospace and Electronic Systems, 1997, 33(3): 883−892

Fig. 18

Flow diagram for the proposed algorithm

Fig. 18 shows the flow diagram for the proposed algorithm. The whole algorithm consists of two steps: the first step is the prediction of parameters, including centroid position and PVoT; the second step is the update for these parameters, which is an iterative process.

5

Conclusion

We proposed an MmT tracking algorithm in this paper. We first classified the MmTs into two types: straight line objects and curved objects. Under the given assumptions, we first modeled the MmTs using random parameters and

[4] Koch J W. Bayesian approach to extended object and cluster tracking using random matrices. IEEE Transactions on Aerospace and Electronic Systems, 2008, 44(3): 1042− 1059 [5] Feldmann M, Fr¨ anken D, Koch W. Tracking of extended objects and group targets using random matrices. IEEE Transactions on Signal Processing, 2011, 59(4): 1409−1420 [6] Richter E, Obst M, Noll M, Wanielik G. Tracking multiple extended objects — a Markov chain Monte Carlo approach. In: Proceedings of the 14th International Conference on Information Fusion. Chicago. Illinois, USA: IEEE, 2011. 314− 321 [7] Baum M, Hanebeck U D. Shape tracking of extended objects and group targets with star-convex RHMs. In: Proceedings of the 14th International Conference on Information Fusion. Chicago, Illinois, USA: IEEE, 2011. 338−345

174

Acta Automatica Sinica, 2013, Vol. 39, No. 2

[8] Baum M, Noack B, Hanebeck U D. Extended object and group tracking with elliptic random hypersurface models. In: Proceedings of the 13th International Conference on Information Fusion. Edinburg, UK: IEEE, 2010. 1−8

[19] Khan Z, Balch T, Dellaert F. MCMC-based particle filtering for tracking a variable number of interacting targets. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2005, 27(11): 1805−1819

[9] Baum M, Hanebeck U D. Random hypersurface models for extended object tracking. In: Proceedings of the 9th IEEE International Symposium on Signal Processing and Information Technology. Ajman, United Arab Emirates: IEEE, 2009. 178−183

[20] Liu Wei-Feng. Research on Multitarget Tracking Algorithm Based on Random Finite Sets and Finite Mixture Models [Ph. D. dissertation], Xi an Jiaotong University, China, 2009. 97−102 (in Chinese)

[10] Mahler R. PHD filters for nonstandard targets I: extended targets. In: Proceedings of the 12th International Conference on Information Fusion. Seattle, WA, USA: ISIF, 2009. 915−921 [11] Lundquist C, Granstr¨ om K, Orguner U. Estimating the shape of targets with a PHD filter. In: Proceedings of the 14th International Conference on Information Fusion. Chicago, Illinois, USA: IEEE, 2011. 49−56 [12] Orguner U, Lundquist C, Granstr¨ om K. Extended target tracking with a cardinalized probability hypothesis density filter. In: Proceedings of the 14th International Conference on Information Fusion. Chicago, Illinois, USA: IEEE, 2011. 65−72 [13] Lian Feng, Han Chong-Zhao, Liu Wei-Feng, Yuan XiangHui. Tracking partly resolvable group targets using SMCPHDF. Acta Automatica Sinica, 2010, 36(5): 731−741 (in Chinese) [14] Rasmussen C, Hager G D. Probabilistic data association methods for tracking complex visual objects. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2001, 23(6): 560−576 [15] Joo S W, Chellpa R. A multiple-hypothesis approach for multiobject visual tracking. IEEE Transactions on Image Processing, 2007, 16(11): 2849−2854 [16] Fleuret F, Berclaz J, Lengagne R, Fua P. Multicamera people tracking with a probabilistic occupancy map. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2008, 30(2): 267−282 [17] Gordon N J, Samlond D J, Smith A F M. Novel approach to nonlinear/non-Gaussian Bayesian state estimation. IEE Proceeding Control Theory and Application, 1993, 140(2): 107−113 [18] Oh S, Russell S, Sastry S. Markov chain Monte Carlo data association for multi-target tracking. IEEE Transactions on Automatic Control, 2009, 54(3): 481−497

[21] Liu W F, Han C Z. Multitarget tracking algorithm based on finite mixture models and equivalent measurement. In: Proceedings of the 11th International Conference on Information Fusion. Cologne, Germany: IEEE, 2008. 1544−1551 [22] Hastings W K. Monte Carlo sampling methods using Markov chains and their applications. Biometrika, 1970, 57(1): 97−109 [23] Green P J. Reversible jump Markov chain Monte Carlo computation and Bayesian model determination. Biometrika, 1995, 82(4): 711−732

Wei-Feng LIU Lecturer at the Institute of Systems Science and Control Engineering, School of Automation, Hangzhou Dianzi University. His research interest covers target tracking, uncertain information processing, and pattern recognition. Corresponding author of this paper. E-mail: [email protected] Zhong CHAI Master student at the Institute of Systems Science and Control Engineering, School of Automation, Hangzhou Dianzi University. He received his bachelor degree from Henan University in 2010. His research interest covers filtering and fusion of measurements with delay. E-mail: [email protected] Cheng-Lin WEN Professor at Hangzhou Dianzi University. His research interest covers multi-source information fusion, technology of system security, and monitoring and fault diagnosis. E-mail: [email protected]