An optimal two-stage algorithm for highly maneuvering targets tracking

ARTICLE IN PRESS Signal Processing 89 (2009) 532–547

Contents lists available at ScienceDirect

Signal Processing journal homepage: www.elsevier.com/locate/sigpro

An optimal two-stage algorithm for highly maneuvering targets tracking Ali Karsaz, Hamid Khaloozadeh  Department of Electrical Engineering, K.N. Toosi University of Technology, Dr. Shariati Street, Teharn, Iran

a r t i c l e i n f o

abstract

Article history: Received 6 February 2008 Received in revised form 3 October 2008 Accepted 6 October 2008 Available online 29 October 2008

A new general two-stage algorithm which extends the target dynamic model to higher order time-derivatives of acceleration is proposed. The conventional input estimation techniques assume constant acceleration level and there are uncovered a generalized acceleration modeling. In contrast, the augmented algorithms, which are based on the jerk modeling, are computationally expensive. In this paper, an innovative scheme is developed to overcome these drawbacks by using a new generalized dynamic modeling of acceleration and an optimal two-stage modified Kalman filter. The proposed scheme is based on the basic Fisher and Bayesian uncertainty models. The optimality of the proposed two-stage modified Kalman filter is proved. Comparisons with two well known approaches using Monte Carlo simulation show that the proposed scheme has a computational advantage over the augmented algorithms and also more significant improvement than the input estimation techniques. & 2008 Elsevier B.V. All rights reserved.

Keywords: Maneuvering target tracking Fisher and Bayesian uncertainty models Jerk model Two-stage Kalman filter

1. Introduction In the airborne vehicles, complex and highly maneuvering targets tracking (MTT) with more accuracy is desirable [1]. Many approaches for tracking a maneuvering target are covered with varying degrees of sophistication in several books on target tracking [1–3]. The motion profile of an arbitrary highly maneuvering target would have a number of non-zero acceleration time-derivatives, which all of them should be included in target motion modeling. The well known ‘‘augmented approach’’ or ‘‘full state’’ methods solve the MTT problems, considering the input parameters as a part of an augmented state to be estimated [4,5]. For example, a successful model of target motion based on the augmented approach is called jerk model and includes position time-derivatives up to the third order [5]. However, ‘‘reduced state’’ methods [6] do not augmented the state, and usually yield a better

 Corresponding author.

E-mail address: [email protected] (H. Khaloozadeh). 0165-1684/$ - see front matter & 2008 Elsevier B.V. All rights reserved. doi:10.1016/j.sigpro.2008.10.016

performance. The full state estimators suffer from two major drawbacks which stem from assuming constant input [7] and are computationally expensive. Decisionbased methods are the most natural approaches to the MTT problem [7]. The key step in the decision based approaches is that the state estimation process is directed by decision regarding target maneuver [7]. Decision-based techniques are divided into three classes which are known as equivalent-noise methods [8,9], input detection and estimation (IDE) algorithms [10–12], and switching model approaches [13]. The equivalent-noise methods assume the filter correction can take the simple form of increasing the Kalman filter (KF) gain, equivalent to placing more weights on the measurements. The IDE algorithm consists of three steps; input estimation, maneuver detection and state estimate correction. The IDE algorithm was first developed by Chan et al., in [11]. They used a simplified batch least square algorithm to estimate the acceleration magnitude when a maneuver is detected. It should be noted that the equivalent-noise approach does not rely on direct estimation of the unknown input, whereas the IDE approach explicitly

ARTICLE IN PRESS A. Karsaz, H. Khaloozadeh / Signal Processing 89 (2009) 532–547

estimates the unknown input. Although IDE approach is attractive—because it intends to relax restrictive assumptions about maneuver dynamics modeling—it suffers from inaccessibility to a priori knowledge of maneuvering target dynamics which is necessary for proper acceleration estimation [7]. In [10], Wang et. al., presented a constant acceleration (CA) method, which the predicted states for the maneuvering target are related to the corresponding states without maneuvering via appropriate filters. Therefore, when target moves with nonconstant acceleration, the performance of the estimation is degraded. In [12], the unknown input was defined as a sum of elementary time functions. Although this input modeling is more general than the constant-input modeling in the original IDE algorithm, the performance will degrade if any input dynamics exists. The modified input estimation (MIE) algorithm proposed in [14] and the enhanced input estimation (EIE) algorithm proposed in [15] are two modifications of the original IDE algorithm by relaxing some assumptions. Switching model as another approach to decision based methods, switches between maneuvering and non-maneuvering models. For example the variable state dimension (VSD) algorithms, which use a constant velocity (CV) and a CA models for target dynamics to achieve good performance, are most general and flexible than other methods [16]. Model-based adaptive filtering techniques have evolved to the interacting multiple model (IMM) algorithm [17], in which the change of the plant is modeled as a Markovian parameter having a transition probability. IMM trackers hypothesize two or more maneuvering modes and assume that the mode changes are modeled by a hidden Markov process. The measurements are filtered through each mode to produce a set of state estimates conditioned on the hypothesized maneuvering mode. Then, the outputs are combined as a weighted sum where the weights are proportional to the model likelihoods [18]. The most basic IMM [19] has one low acceleration mode and one high acceleration mode and assumes that the Markov chain transition probabilities are stationary and known. More complicated IMMs can have non-stationary transition probabilities, autocorrelated maneuvers, and adaptive mode sets [18]. They perform well when the modes accurately represent the true accelerations and are relatively robust to small modeling errors. Their computational complexity increases linearly with the number of maneuvering modes [20]. Recently, neural networks and fuzzy logic schemes have been used for maneuvering targets with intelligent adaptation rules [21–23]. Friedland [24] introduced a method of separating estimation of the unknown input from the dynamic variables and Blair used this method in the MTT problem [25]. The basic idea was to decouple the augmented Kalman filter (AKF) into two-stage filters in order to reduce computation and memory requirements [26–29]. Recently, Hsieh and Chen [26,27] derived a two-stage KF for a more general case in which the bias is expressed by the first-order auto-regressive model. The two-stage

533

filtering method suggested for MTT problem in [25] suffers from two major drawbacks. These drawbacks arise from constant acceleration and the input observability assumptions. In this paper, we present a new modified KF algorithm based on the Fisher uncertainty modeling for highly MTT. A partitioned dynamic model for target acceleration vector is proposed based on the augmented jerk model of [5]. It is shown that the maneuver tracking algorithm proposed in [10] and [25] are special case of our proposed method. The motivation of our proposed method is generating a two-stage structure to obtain the optimal performance when the input term is not observable through the measurements. The rest of the paper is organized as follows. The models of uncertainty and some background materials on maneuvering target modeling are briefly reviewed in Section 2. The new scheme formulations are derived in Section 3. The simulation results of the proposed method, the CA and the augmented approaches are presented and compared in Section 4. 2. Models of uncertainty and MTT formulation The two basic uncertainty models, which have been considered in this paper, are the Bayesian and Fisher models. These models are specific cases of the state space structure-white process [30]. A standard Fisher model for static case (off-line model) is defined as Z ¼ HX þ V

(1)

where the state vector X is a completely unknown vector variable, Z is the observation vector, V is stochastic vector variable (observation or measurement noise) with the covariance matrix R and H is measurement matrix. Minimum variance estimation of state vector X for a given observation vector Z can be obtained from the following equation: X^ ¼ ðHT R1 HÞ1 HT R1 Z

(2)

The Bayesian models are one of the most important and common used models of uncertainty. The complete definition of the Bayesian, discrete time model for linear dynamic systems is summarized as follows: Xðn þ 1Þ ¼ FðnÞXðnÞ þ GðnÞwðnÞ ZðnÞ ¼ HðnÞXðnÞ þ vðnÞ

(3)

where X(U) is the system state, X(0) is the initial condition, Z(U) is the observation vector, w(U) is the white system uncertainty and v(U) is the white observation uncertainty with the following covariance matrices: Efvðn1 ÞvT ðn2 Þg ¼ Rðn1 Þdðn1  n2 Þ, Efwðn1 ÞwT ðn2 Þg ¼ Q ðn1 Þdðn1  n2 Þ, EfXð0Þ X T ð0Þg ¼ c and EfXð0Þg ¼ 0; EfwðnÞg ¼ EfvðnÞg ¼ 0; for all n, EfvðÞwðÞT g ¼ 0, where d(U) denotes the Kronecker delta function and the operator E{U} indicates the expected value. The matrices F(n), G(n) and H(n) in (3) are assumed to be known function of step time n. In the Bayesian case, stochastic probabilistic models are used for all the uncertainties. Thus, X(0), v(U) and w(U) are modeled as zero mean uncertainty random processes. The desired

ARTICLE IN PRESS A. Karsaz, H. Khaloozadeh / Signal Processing 89 (2009) 532–547

latter observation Zðn þ 1Þ. The solution of the filtering problem, after some manipulation, leads to the standard Kalman filter (SKF). The SKF is an efficient and unbiased ^ þ ˜ þ 1Þ ¼ Zðn þ 1Þ  Zðn filter, such that the sequence Zðn ^ 1jnÞ ¼ Zðn þ 1Þ  Hðn þ 1ÞXðn þ 1jnÞ which is the residue of observation is a stochastic zero mean white process. i.e., ˜ ˜ 1 ÞZðn ˜ 2 ÞT g ¼ Rz dðn1  n2 Þ. EfZðnÞg ¼ 0, EfZðn The system model of a target moving with CV in a straight line is different from that of a target moving with acceleration or maneuver. It is assumed that the target moves in a plane, which is the two-dimensional case, such as a helicopter with fixed elevation. The state equation for the non-maneuvering model is given by (3) where XðnÞ ¼ ½xðnÞ

vx ðnÞ

yðnÞ

vy ðnÞT

(4)

where x(n), y(n) and vx(n), vy(n) denote the position and the velocity of the target along the x and y directions, respectively. The expression for G(n), F(n) and H(n) as functions of the update time T (T is the time interval between two consecutive measurements) are 2 3 1 T 0 0 2 3T 6 7 T 2 =2 T 0 0 60 1 0 07 6 7 6 7 FðnÞ ¼ 6 7; GðnÞ ¼ 4 5 ; 2 60 0 1 T 7 0 0 T =2 T 4 5 0 0 0 1 " # 1 0 0 0 (5) HðnÞ ¼ 0 0 1 0 Also the measurement equation is given by (3), where v(U) is the radar measurement noise. The maneuvering model deals with the acceleration as an additive term: Xðn þ 1Þ ¼ FðnÞXðnÞ þ CðnÞUðnÞ þ GðnÞwðnÞ

(6)

T

where U(n) ¼ [ux(n)uy(n)] and C(n) ¼ G(n)  U(n) is the acceleration of target which is modeled as an unknown variable. In contrast of non-maneuvering target which the ˜ þ 1Þ is zero, for maneuvering target case this mean of Zðn sequence is no longer zero mean, using the SKF [7]. Many approaches for target maneuvering detection, including maneuver onset and termination, have been investigated in several literatures [10,31–34]. The maneuvering detection algorithm which is implemented in this paper is the same as [33]. 3. The proposed maneuvering target tracking scheme The conventional IDE approaches are un-model based maneuvering algorithms. This is the main drawback of these methods [7]. In contrast, other conventional augmented techniques include higher order time-derivative up to the third order, which overcome this problem, are computationally expensive. In this paper the proposed approach intends to relax restrictive assumptions concerning the maneuver dynamics modeling and using a new optimal two-stage algorithm. For this purpose a generalized modeling is developed such that the additive term has a new dynamic equation. The additive term has a

distinct dynamic equation, including the jerk dynamics in [5], and the early work of Singer [35]. This generalized acceleration model (GAM) can be extended to the higher order time-derivative in terms of acceleration.

3.1. Generalized acceleration dynamic model based on Gauss–Markov target motion modeling It is assumed that the target jerk, j(t) is a zero-mean first-order stationary Markov process with auto-correlation Rj ðtÞ ¼ EfjðtÞjðt þ tÞg ¼ s2j expðajtjÞ or equivalently, power spectrum Sj ðwÞ ¼ 2as2j =ðw2 þ a2 Þ. Where sj is the variance of the target jerk and a is the reciprocal of the jerk time constant t ¼ 1/a. Thus it depends on how long the maneuvering lasts [5,35]. It is assumed that the process j(t), is the state process of a linear time-invariant system j(t) ¼ aj(t)+wj(t). Where wj(t) is zero-mean white noise with constant power spectral density Swj ¼ 2as2j . The discrete-time state space representation for twodimensional case of jerk model based on augmented assumption as suggested in [5] is d X Aug ðtÞ ¼ GAug ðtÞX Aug ðtÞ þ GAug ðtÞwðtÞ dt

(7)

where 2 2 3 x 0 6 6 7 60 6 x_ 7 6 6 7 6 6 7 60 6 x€ 7 6 6 7 6 6 7 60 6x7 6 6 7 X Aug ðtÞ ¼ 6 7; GAug ðtÞ ¼ 6 60 6y7 6 6 7 6 6 7 60 6 y_ 7 6 6 7 6 6 7 60 6 y€ 7 4 4 5 y 0 2 3T 0 6 7 607 6 7 6 7 607 6 7 6 7 617 6 7 GAug ðtÞ ¼ 6 7 607 6 7 6 7 607 6 7 6 7 607 4 5 ___

form of the filtering solution in the dynamic system modeling is a recursive relationship equation which ^ þ 1jn þ 1Þ in the terms of XðnjnÞ ^ expressing Xðn and the

___

534

1

0

0

0

0

0

0

1

0

0

0

0

0

0

1

0

0

0

0

0

a

0

0

0

0

0

0

0

1

0

0

0

0

0

0

1

0

0

0

0

0

0

0

0

0

0

0

0

3

0

7 0 7 7 7 0 7 7 7 0 7 7 7, 0 7 7 7 0 7 7 7 1 7 5 a

(8)

1 The corresponding discrete-time system (8) is obtained as X Aug ðn þ 1Þ ¼ F Aug ðnÞX Aug ðnÞ þ W Aug ðnÞ

(9)

where, X Aug ðnÞ ¼ ½xðnÞ vx ðnÞ ux ðnÞ jx ðnÞ yðnÞ vy ðnÞT and W Aug ðnÞ ¼

Z

tnþ1

GAug ðt nþ1 ; lÞGAug ðlÞwðlÞ dl

(10)

tn

For highly maneuvering target cases, the product aT is small. Thus, targets have sustained jerk level at jerk time

ARTICLE IN PRESS A. Karsaz, H. Khaloozadeh / Signal Processing 89 (2009) 532–547

period and FAug(n) 2 1 6 60 6 60 6 6 60 lim F Aug ðnÞ ¼ 6 60 aT!0 6 6 60 6 6 40 0

where

goes to T

T 2 =2

T 3 =6 2

0

0

0

0 0

0 0

0 0

1 0

T 1

T =2 T

0

0

1

0

0

0

0

0

0

1

T

T 2 =2

0

0

0

0

1

T

0

0

0

0

0

1

0

0

0

0

0

0

0

3

7 0 7 7 0 7 7 7 0 7 7 T 3 =6 7 7 7 T 2 =2 7 7 7 T 5

The process noise covariance matrix EfW Aug ðnÞW TAug ðnÞg has been obtained: " # L44 044 , lim Q Aug ðnÞ ¼ 2as2j aT!0 044 L44 3 2 7 T =252 T 6 =72 T 5 =30 T 4 =24 6 7 6 T 6 =72 T 5 =20 T 4 =8 T 3 =6 7 6 7 7 L44 ¼ 6 6 5 7 T 4 =8 T 3 =3 T 2 =2 7 6 T =30 4 5 T 3 =6 T 2 =2 T T 4 =24

X m ðnÞ ¼ ½xðnÞ vx ðnÞ yðnÞ vy ðnÞT

(18)

U Gam ðnÞ ¼ ½ux ðnÞ jx ðnÞ uy ðnÞ jy ðnÞT

(19)

AGam ðnÞ ¼ eGGam ðtnþ1 tn Þ ¼ eGGam T

(20)

WU Gam ðnÞ ¼

Q Aug ðnÞ ¼

(12)

Obviously, this augmented jerk model has a higher dimension than the Singer acceleration model and the method suggested in [10] and [25] (for more details see [5]). If aT ! 1 and sampling rate T is fixed, the discretetime model of (7) behaves like as a CA model as it was obtained in [5]. In this paper, the highly maneuvering targets which usually have the long maneuvering time constant (or small aT) are discussed. The continuous time state space representation for acceleration dynamic equation of the GAM system is considered as below: d U Gam ðtÞ ¼ GGam ðtÞU Gam ðtÞ þ GGam ðtÞwGam ðtÞ dt

___

where the parameters of generalized acceleration UGam(t), GGam(t) and GGam(t) using (8) are 2 0 1 0 6 6 0 a 0 h iT 6 U Gam ðtÞ ¼ x€ x y€ y ; GGam ðtÞ ¼ 6 60 0 0 4 0 0 0 2 3 0 6 7 617 6 7 GGam ðtÞ ¼ 6 7 607 4 5

(13)

0

3

7 0 7 7 7, 1 7 5 a

(14)

1

(15) U Gam ðn þ 1Þ ¼ AGam ðnÞU Gam ðnÞ þ W U Gam ðnÞ

(16)

The measurement equation for the GAM system is (17)

GGam ðtnþ1 ; tÞGGam ðtÞwGam ðtÞ dt

0

0 2

0

2

1 6 60 AGam ðnÞ ¼ 6 60 4 0

a2

0

1 6 60 6 þ  ¼ 6 60 4

(21)

p

0

q

0

0

1

0

0

p

0

q

0

0

1

0

0

0

0

3

0

0

0

a3

7 07 7 7 p7 5 q

3

7 07 7 p7 5 q

(22)

where T2 T3 ð1  eaT Þ þ a2   ¼ ; 2 3! a 2 3 T T  a3 þ    ¼ eaT q ¼ 1  aT þ a2 2 3! p¼T a

(23)

WGamU(n) is a zero-mean sequence with the following covariance: U U T EfW U Gam ðnÞW Gam ðnÞ g ¼ Q Gam ðn1 Þdðn1  n2 Þ Z tnþ1 Z tnþ1 ¼ GGam ðtnþ1 ; tÞGGam ðtÞEfwGam ðtÞwTGam ðlÞg tn

X m ðn þ 1Þ ¼ F Gam ðnÞX m ðnÞ þ C Gam ðnÞU Gam ðnÞ þ W XGam ðnÞ

t nþ1

tn

Where the superscript ‘m’ denotes the maneuvering target state. Despite the assumption in Refs. [25,26], and [29], the acceleration vector UGam(n) is not observable through the measurement process in (17). FGam(n) is equal to F(n) as defined in (5). Through a series expansion of the (20) the state transition matrix is obtained: 3 2 3 2 1 0 0 0 0 1 0 0 7 6 7 6 6 0 1 0 0 7 6 0 a 0 0 7 7 6 7 6 AGam ðnÞ ¼ 6 7T 7þ6 60 0 1 07 60 0 0 1 7 5 4 5 4 0 0 0 a 0 0 0 1 2 3 2 3 0 a 0 0 0 0 0 a2 6 7 6 7 3 6 0 a2 0 0 7 T 2 6 0 a3 0 0 7 6 7 6 7T þ6 þ6 7 7 2 7 3! 6 0 0 0 a 7 2 60 0 0 a 5 4 5 4

model

The corresponding discrete-time reduced state of GAM system using (9) and (13) is obtained as

Zðn þ 1Þ ¼ Hðn þ 1ÞX m ðn þ 1Þ þ Vðn þ 1Þ

Z

1 (11)

___

535

tn

 GTGam ðlÞGTGam ðt nþ1 ; lÞ dt dl

(24)

The generalized matrix AGam(n) and the generalized covariance matrix Q U Gam ðnÞ are obtained for highly maneuvering targets when aT ! 0 as below: 2 3 1 T 0 0 60 1 0 07 6 7 lim AGam ðnÞ ¼ 6 (25) 7 40 0 1 T 5 aT!0 0

0

0

1

ARTICLE IN PRESS 536

A. Karsaz, H. Khaloozadeh / Signal Processing 89 (2009) 532–547

2

T 3 =3 6 2 6 T =2 6 lim Q U ðnÞ ¼ 2asj 6 6 0 aT!0 Gam 4 0

T 2 =2

0

T

0

0

T 3 =3

0

T 2 =2

0

3

7 0 7 7 7 T 2 =2 7 5 T

for each sample time as (26)

¼ F Gam ð0ÞX m ð0Þ þ C Gam ð0Þ½A1 Gam ð0ÞU Gam ð1Þ U X  A1 Gam ð0ÞW Gam ð0Þ þ W Gam ð0Þ

In the proposed method, positions and velocities are updated based on the relations (15) and (16) if considering the generalized matrix CGam(n) as 2

T 2 =2 6 6 T 6 C Gam ðnÞ ¼ 6 6 0 4 0

T 3 =6

0

T 2 =2

0

0

T 2 =2

0

T

0

3

7 0 7 7 7 T 3 =6 7 5 T 2 =2

X m ð1Þ ¼ F Gam ð0ÞX m ð0Þ þ C Gam ð0ÞU Gam ð0Þ þ W XGam ð0Þ

(27)

¼ Xð1Þ þ C Gam ð0ÞA1 Gam ð0ÞU Gam ð1Þ U  C Gam ð0ÞA1 Gam ð0ÞW Gam ð0Þ

X m ð2Þ ¼ F Gam ð1ÞX m ð1Þ þ C Gam ð1ÞU Gam ð1Þ þ W XGam ð1Þ ¼ F Gam ð1ÞXð1Þ þ F Gam ð1ÞC Gam ð0Þ½A1 Gam ð0ÞU Gam ð1Þ U X  A1 Gam ð0ÞW Gam ð0Þ þ C Gam ð1ÞU Gam ð1Þ þ W Gam ð1Þ

Eq. (9) is used for updating the parameters in the augmented method as suggested in [5]. The quantities W XGam ðnÞ and V(n) are zero-mean sequences with the following covariance: EfVðn1 ÞVðn2 ÞT g ¼ Rðn1 Þdðn1  n2 Þ

(32)

Since AGam(n) is a non-singular matrix defined in (25), UGam(0) in (32) is replaced by (16). It is assumed that Xm(0) ¼ X(0), so Xð1Þ ¼ F Gam ð0ÞXð0Þ þ W XGam ð0Þ ¼ F Gam ð0Þ X m ð0Þ þ W XGam ð0Þ.

¼ Xð2Þ þ ½F Gam ð1ÞC Gam ð0ÞA1 Gam ð0Þ þ C Gam ð1ÞU Gam ð1Þ U  F Gam ð1ÞC Gam ð0ÞA1 Gam W Gam ð0Þ

¼ Xð2Þ þ ½F Gam ð1ÞC Gam ð0ÞA1 Gam ð0Þ 1 U þ C Gam ð1Þ½A1 Gam ð1ÞU Gam ð2Þ  AGam ð1ÞW Gam ð1Þ

(28)

U  F Gam ð1ÞC Gam ð0ÞA1 Gam W Gam ð0Þ 1 ¼ Xð2Þ þ ½F Gam ð1ÞC Gam ð0ÞA1 Gam ð0Þ þ C Gam ð1ÞAGam ð1ÞU Gam ð2Þ

XU T EfW XGam ðn1 ÞW U Gam ðn2 Þ g ¼ Q Gam ðn1 Þdðn1  n2 Þ 3 2 7 T =252 T 6 =24 0 0 7 6 6 T 6 =72 T 5 =20 0 0 7 6 7 X 6 7 Q Gam ðnÞ ¼ 2asj 6 7 0 0 T 7 =252 T 6 =72 7 6 4 5 0 0 T 6 =72 T 5 =20

1 U  ½F Gam ð1ÞC Gam ð0ÞA1 Gam ð0Þ þ C Gam ð1ÞAGam ð1ÞW Gam ð1Þ U  F Gam ð1ÞC Gam ð0ÞA1 Gam ð0ÞW Gam ð0Þ

(33)

For arbitrary sample time n, we have (29)

m

X ðn þ 1Þ ¼ Xðn þ 1Þ þ MGam ðn þ 1ÞU Gam ðn þ 1Þ n X PðiÞW U  Gam ðiÞ

(34)

i¼0

The cross covariance matrix between W XGam ðnÞ and WU Gam ðnÞ considering (12) is XU T EfW XGam ðn1 ÞW U Gam ðn2 Þ g ¼ Q Gam ðn1 Þdðn1  n2 Þ 3 2 5 4 T =30 T =24 0 0 6 7 6 T 4 =8 T 3 =6 0 0 7 6 7 XU 6 7 Q Gam ðnÞ ¼ 2asj 6 7 5 4 0 T =30 T =24 7 6 0 4 5 0 0 T 4 =8 T 3 =6

where MGam and P are defined as follows: M Gam ðn þ 1Þ ¼ ½F Gam ðnÞMGam ðnÞ þ C Gam ðnÞA1 Gam ðnÞ, n ¼ 2; 3; :::. M Gam ð1Þ ¼ C Gam ð0ÞA1 Gam ð0Þ

(30)

PðiÞ ¼ ½Pði  1Þ þ ½

ni1 Y

(35)

F Gam ðn  jÞC Gam ðiÞA1 Gam ðiÞ;

j¼0

i ¼ 1; 2; ::; n

The acceleration vector can update independently by reduced order Eq. (16) in the proposed method. The maneuvering model of conventional methods such as [10,12] and [25] is a special case of the new input model in (16), when jx(n) ¼ jy(n) ¼ 0.

Pð0Þ ¼ ½

n 1 Y

F Gam ðn  jÞC Gam ð0ÞA1 Gam ð0Þ

(36)

j¼0

Using zero-mean property of W U Gam ðnÞ and one-step prediction logic from dynamic equation (34) the predicted state is obtained as follows: m ^ þ 1jnÞ þ M Gam ðn þ 1ÞU ^ Gam ðn þ 1jnÞ X^ ðn þ 1jnÞ ¼ Xðn

(37)

3.2. Optimal partitioned state estimation algorithm The non-maneuvering model can be obtained by replacing UGam(n) ¼ 0 in (15) as below: Xðn þ 1Þ ¼ F Gam ðnÞXðnÞ þ W XGam ðnÞ

(31)

The state of the target without maneuvering is denoted by X(n). The state of the maneuvering target model (15) using non-maneuvering state in Eq. (31) can be calculated

After determining the presence of the maneuvering, one can switch from non-maneuvering algorithm to maneuvering algorithm. The residues of non-maneuver and maneuver target observations are defined, respectively, ^ þ 1jnÞ ˜ þ 1jnÞ ¼ Zðn þ 1Þ  Hðn þ 1ÞXðn Zðn

(38)

m m Z˜ ðn þ 1jnÞ ¼ Zðn þ 1Þ  Hðn þ 1ÞX^ ðn þ 1jnÞ

(39)

ARTICLE IN PRESS A. Karsaz, H. Khaloozadeh / Signal Processing 89 (2009) 532–547

In the following, we derive a recursive algorithm to estimate UGam(n) in order to minimize the covariance matrix of acceleration vector estimation, is denoted by P U Gam ðnÞ. Using (37), (38) and (39), we represent UGam(n+1) in the equation, which relates the acceleration state to the measurement residue as below: ^ þ 1jnÞ ˜ þ 1jnÞ ¼ Zðn þ 1Þ  Hðn þ 1ÞXðn Zðn m

¼ Zðn þ 1Þ  Hðn þ 1ÞðX^ ðn þ 1jnÞ  M Gam ðn þ 1ÞU Gam ðn þ 1ÞÞ

and measurement covariance matrix <ðn þ 1Þ defined as below: <ðn þ 1Þ 82 9 3 m < = Z˜ ðn þ 1jnÞ m 5½Z˜ ðn þ 1jnÞ  U˜ Gam ðn þ 1jnÞ ¼E 4 : U˜ Gam ðn þ 1jnÞ ; P P 2 3 ðn þ 1jnÞ  ðn þ 1jnÞ Zm Z m U Gam 6 7 6 7 ¼6 P P 7 T 4 ðn þ 1jnÞ ðn þ 1jnÞ 5 Z m U Gam

m ¼ Zðn þ 1Þ  Hðn þ 1ÞX^ ðn þ 1jnÞ

m

(40) ˜ þ 1jnÞ ¼ Hðn þ 1ÞM Gam ðn þ 1ÞU Gam ðn þ 1Þ Zðn (41)

It is clear from the non-maneuvering algorithm which is ˜ þ 1jnÞ exists a SKF, that the non-maneuvering residue Zðn and easily obtainable while the maneuvering residue m Z˜ ðn þ 1jnÞ is actually not available. Note that (41) is in the standard form Z(n) ¼ H(n)X(n)+v(n), and can be viewed as an observation model of UGam(n+1) and m measurement noise Z˜ ðn þ 1jnÞ which is a zero mean P white random process with covariance matrix Z m ðn þ 1jnÞ [7]. Suppose that U˜ Gam ðn þ 1jnÞ denotes the acceleration state residue: ^ Gam ðn þ 1jnÞ U˜ Gam ðn þ 1jnÞ ¼ U Gam ðn þ 1Þ  U

P ðn þ 1jnÞ Zm 6 6 P <ðn þ 1Þ ¼ 4 ðn þ 1jnÞT 

˜ þ 1jnÞ Zðn ^ Gam ðn þ 1jnÞ U

"

# ¼

Hðn þ 1ÞM Gam ðn þ 1Þ 2

þ4

I m Z˜ ðn þ 1jnÞ

U˜ Gam ðn þ 1jnÞ

# U Gam ðn þ 1Þ

3 5

ðn þ 1jnÞ

Z m U Gam

P

ðn þ 1jnÞ

3 7 7 5

(45)

U Gam

Using the one-step prediction logic from dynamic Eq. ^ Gam ðn þ 1jnÞ is calculated from U ^ Gam ðnjnÞ: (16), U ^ Gam ðn þ 1jnÞ ¼ AGam ðnÞU ^ Gam ðnjnÞ U

(46)

In the following, we derive a recursive algorithm to P estimate U Gam ðn þ 1Þ to minimum U Gam ðn þ 1jn þ 1Þ using the standard Fisher estimation described by (2): ^ Gam ðn þ 1jn þ 1Þ ¼ U

½M Gam ðn þ 1ÞT Hðn þ 1ÞT I P ðn þ 1jnÞ Zm 6 6 6 P 4 ðn þ 1jnÞT 2

" 

P



Z m U Gam

(42)

(1) Use the one-step prediction logic to change the ^ Gam ðn þ 1jn þ 1Þ from U ^ Gam ðn þ problem to calculate U ˜ þ 1jnÞ. 1jnÞ and Zðn ^ Gam ðn þ (2) Solve a Fisher estimation problem where U ˜ þ 1jnÞ are two available variables con1jnÞ and Zðn sidered on an unknown vector U Gam ðn þ 1Þ such as below:

P



Z m U Gam

The desired form of the filtering solution is a difference ^ Gam ðn þ equation (discussed in Section 2) expressing U ^ Gam ðnjnÞ on Zðn ˜ þ 1Þ as state by (41). 1jn þ 1Þ in terms of U From Eqs. (41) and (42) and by some manipulation as below, we can illustrate a standard Fisher model algorithm for maneuvering case. The following two major steps for accelerate estimation are assumed:

"

U Gam

2

þ Hðn þ 1ÞM Gam ðn þ 1ÞU Gam ðn þ 1Þ ¼ Z˜ ðn þ 1jnÞ þ Hðn þ 1ÞM Gam ðn þ 1ÞU Gam ðn þ 1Þ

m þ Z˜ ðn þ 1jnÞ

537

Hðn þ 1ÞM Gam ðn þ 1Þ

ðn þ 1jnÞ

Z m U Gam

P

ðn þ 1jnÞ

U Gam

31 7 7 7 5

#!1

I

 ½M Gam ðn þ 1ÞT Hðn þ 1ÞT I P P 2 3 ðn þ 1jnÞ  ðn þ 1jnÞ 1 Zm Z m U Gam 6 7 6 7 6 P P 7 4 ðn þ 1jnÞT ðn þ 1jnÞ 5 Z m U Gam

" 

˜ þ 1jnÞ Zðn

U Gam

#

^ Gam ðn þ 1jnÞ U

X m m ðn þ 1jnÞ ¼ EfZ˜ ðn þ 1jnÞZ˜ ðn þ 1jnÞT g

(47)

(48)

Zm m Where the covariance matrix of Z˜ ðn þ 1jnÞ, which is the estimated covariance matrix of the maneuver target P observation, denoted by Zm ðn þ 1jnÞ must be calculated. From (39) X X ðn þ 1jnÞ ¼ Hðn þ 1Þ ðn þ 1jnÞHðn þ 1ÞT þ Rðn þ 1Þ Zm

(43)

where from (42) the acceleration error covariance matrix is defined as below: X ðn þ 1jnÞ ¼ EfU˜ Gam ðn þ 1jnÞU˜ Gam ðn þ 1jnÞT g (44) U Gam

Therefore, (43) is a standard Fisher model (1) with m measurement noise vector ½Z˜ ðn þ 1jnÞ  U˜ Gam ðn þ 1jnÞT

(49) P m m where ðnjnÞ ¼ E½ðX m ðnÞ  X^ ðnjnÞÞðX m ðnÞ  X^ ðnjnÞÞT  is the state estimation covariance matrix. It is seen that the acceleration state residue U˜ Gam ðn þ 1jnÞ is time-correlated m with the maneuvering observation residue Z˜ ðn þ 1jnÞ, P denoted by Zm U Gam ðn þ 1jnÞ. X m ðn þ 1jnÞ ¼ EfZ˜ ðn þ 1jnÞU˜ Gam ðn þ 1jnÞT g (50) Z m U Gam

538

Table 1 The proposed method (PM) and the augmented method (AM) equations. Dynamic equations State vectors Matrices

vx ðnÞ ux ðnÞ jx ðnÞ yðnÞ vy ðnÞ uy ðnÞ jy ðnÞT 2 3 3 2 3 2 3 T 2 =2 T 3 =6 0 0 T 0 0 1 T 0 0 1 0 T 6 7 2 6 T 7 60 1 0 07 60 07 0 0 7 T =2 1 0 07 6 7 6 6 7 7 7 ; AGam ðnÞ ¼ 6 7 ; C Gam ðnÞ ¼ 6 7 ; HðnÞ ¼ 6 7 6 0 40 0 1 T 5 40 15 0 1 T5 0 T 2 =2 T 3 =6 7 4 5 2 0 0 0 1 44 0 0 0 1 44 0 0 24 0 0 T T =2 44 3 2 3 2 AM 0 0 1 T T 2 =2 T 3 =6 0 0 1 0 T 7 6 60 07 60 1 0 0 7 T T 2 =2 0 0 7 6 7 6 7 6 7 60 0 60 07 1 T 0 0 0 0 7 6 7 6 7 6 7 6 7 60 0 0 0 0 1 0 0 0 0 7 7 6 F Aug ðnÞ ¼ 6 ; HAug ðnÞ ¼ 6 7 7 60 0 2 3 60 17 0 0 1 T T =2 T =6 7 6 7 6 7 6 60 07 2 7 60 0 7 6 0 0 0 1 T T =2 7 6 7 6 7 6 4 0 05 40 0 0 0 0 0 1 T 5 0 0 28 0 0 0 0 0 0 0 1 88 3 2 3 2 2 PM T 7 =252 T 6 =72 T 3 =3 T 2 =2 0 0 T 5 =30 7 6 6 0 22 7 6 2 6 4 7 6 T =72 T 5 =20 7 6 T =2 6 T =8 T 0 0 7 6 7 6 XU 7; Q U ðnÞ ¼ 2asj 6 Q XGam ðn1 Þ ¼ 2asj 6 7; Q Gam ðnÞ ¼ 2asj 6 6 Gam 6 6 0 6 0 T 7 =252 T 6 =72 7 0 T 3 =3 T 2 =2 7 7 6 5 4 4 5 4 022 T T 6 =72 T 5 =20 0 0 T 2 =2 0 3 2 7 6 5 4 AM T =252 T =72 T =30 T =24 7 6 6 " # 6 T =72 T 5 =20 T 4 =8 L44 044 T 3 =6 7 7 6 T 2 ; L44 ¼ 6 5 EfW Aug ðnÞW Aug ðnÞg ¼ Q Aug ðnÞ ¼ 2asj 7 4 3 2 044 L44 6 T =30 T =8 T =3 T =2 7 5 4 4 3 2 T =6 T =2 T T =24 AM X Aug ðnÞ ¼ ½xðnÞ 2 1 6 60 F Gam ðnÞ ¼ 6 40

PM

T 4 =24

0

3

T =6

0 5

0

T =30

0

T 4 =8

0

3

7 7 7 7; T =24 7 5 T 3 =6 0

4

EfV ðn1 ÞV ðn2 ÞT g ¼ Rðn1 Þdðn1  n2 Þ

^ þ 1jn þ 1Þ ¼ F Gam ðnÞXðnjnÞ ^ ^ PM Xðn þ Kðn þ 1Þ½Zðn þ 1Þ  Hðn þ 1ÞF Gam ðnÞXðnjnÞ P Kðn þ 1Þ ¼ ðn þ 1jn þ 1ÞHT ðn þ 1ÞR1 ðn þ 1Þ P P P P P ðn þ 1jn þ 1ÞÞ ¼ ðn þ 1jnÞ  ðn þ 1jnÞHT ðn þ 1Þ½Rðn þ 1Þ þ Hðn þ 1Þ ðn þ 1jnÞHT ðn þ 1Þ1 Hðn þ 1Þ ðn þ 1jnÞ P P ðn þ 1jnÞ ¼ F Gam ðnÞ ðnjnÞF TGam ðnÞ þ Q XGam ðnÞ ^ Gam ðnjnÞ þ K Gam ðn þ 1Þ½Zðn ^ Gam ðnjnÞ ^ Gam ðn þ 1jn þ 1Þ ¼ AGam ðnÞU ˜ þ 1jnÞ  M Gam ðn þ 1ÞT Hðn þ 1ÞT AGam ðnÞU U P P P K Gam ðn þ 1Þ ¼ 2 U Gam ðn þ 1jnÞM Gam ðn þ 1ÞT Hðn þ 1ÞT ½3Hðn þ 1ÞM Gam ðn þ 1Þ UGam ðn þ 1jnÞM Gam ðn þ 1ÞT Hðn þ 1ÞT þ Zm ðn þ 1jnÞ1 P P P P P P T T T T 1 Hðn þ 1ÞM Gam ðn þ 1Þ UGam ðn þ 1jnÞ U Gam ðn þ 1jn þ 1Þ ¼ U Gam ðn þ 1jnÞ  4 U Gam ðn þ 1jnÞM Gam ðn þ 1Þ Hðn þ 1Þ  ½3Hðn þ 1ÞM Gam ðn þ 1Þ U Gam ðn þ 1jnÞM Gam ðn þ 1Þ Hðn þ 1Þ þ Z m ðn þ 1jnÞ P P m T U ^ ^ ^ U Gam ðn þ 1jnÞ ¼ AGam ðnÞ U Gam ðnjnÞAGam ðnÞ þ Q Gam ðnÞ, X ðn þ 1jnÞ ¼ Xðn þ 1jnÞ þ M Gam ðn þ 1ÞU Gam ðn þ 1jnÞ AM X^ Aug ðn þ 1jn þ 1Þ ¼ F Aug ðnÞX^ Aug ðnjnÞ þ K Aug ðn þ 1Þ½Zðn þ 1Þ  HAug ðn þ 1ÞF Aug ðnÞX^ Aug ðnjnÞ P K Aug ðn þ 1Þ ¼ Aug ðn þ 1jn þ 1ÞHTAug ðn þ 1ÞR1 ðn þ 1Þ P P P P P T T 1 HAug ðn þ 1Þ Aug ðn þ 1jnÞ Aug ðn þ 1jn þ 1Þ ¼ Aug ðn þ 1jnÞ  Aug ðn þ 1jnÞHAug ðn þ 1Þ½Rðn þ 1Þ þ HAug ðn þ 1Þ Aug ðn þ 1jnÞH Aug ðn þ 1Þ P P P T T ðn þ 1jnÞ ¼ F Aug ðnÞ ðnjnÞF ðnÞ þ Q Aug ðnÞ, ðnjnÞ ¼ E½ðX Aug ðnÞ  X^ Aug ðnjnÞÞðX Aug ðnÞ  X^ Aug ðnjnÞÞ  Aug

Aug

Aug

Aug

ARTICLE IN PRESS

MTT filters

AM X Aug ðn þ 1Þ ¼ F Aug ðnÞX Aug ðnÞ þ W Aug ðnÞ, ZðnÞ ¼ HAug ðnÞX Aug ðnÞ þ vðnÞ PM X m ðnÞ ¼ ½xðnÞ vx ðnÞ yðnÞ vy ðnÞT , U Gam ðnÞ ¼ ½ux ðnÞ jx ðnÞ uy ðnÞ jy ðnÞT

A. Karsaz, H. Khaloozadeh / Signal Processing 89 (2009) 532–547

Covariance matrices

m PM X m ðn þ 1Þ ¼ F Gam ðnÞX m ðnÞ þ C Gam ðnÞU Gam ðnÞ þ W XGam ðnÞ, U Gam ðn þ 1Þ ¼ AGam ðnÞU Gam ðnÞ þ W U Gam ðnÞZðn þ 1Þ ¼ Hðn þ 1ÞX ðn þ 1Þ þ Vðn þ 1Þ

ARTICLE IN PRESS A. Karsaz, H. Khaloozadeh / Signal Processing 89 (2009) 532–547

Since (50) is not generally equal to zero, the CA approach proposed in [10], which has been considered P is a sub-optimal algorithm. Z m U Gam ðn þ 1jnÞ ¼ 0, P The extra computation of Z m U Gam ðn þ 1jnÞ is the main reason for the computational complexity in the augmented state methods. Therefore, if this term can be eliminated, one can reduce the computational complexity. In the following, we derive an expression to relate P ðn þ 1jnÞ to the acceleration error covariance Z m U GamP matrix UGam ðn þ 1jnÞ. Using the Eqs. (17), (34), (37) and (39), gives

539

¼ Hðn þ 1ÞX˜ ðn þ 1jnÞ þ Hðn þ 1ÞM Gam ðn þ 1ÞU˜ Gam ðn þ 1jnÞ n X PðiÞW U  Hðn þ 1Þ Gam ðiÞ m Z˜ ðn þ 1jnÞ ¼ Hðn þ 1ÞX˜ ðn þ 1jnÞ

þ Hðn þ 1ÞM Gam ðn þ 1ÞU˜ Gam ðn þ 1jnÞ n X PðiÞW U  Hðn þ 1Þ Gam ðiÞ

(52)

i¼0

It should be noted that the term Hðn þ 1ÞX˜ ðn þ 1jnÞ

m m Z˜ ðn þ 1jnÞ ¼ Hðn þ 1ÞX˜ ðn þ 1jnÞ

˜ þ 1Þ in (38). Using the Eq. (52), in (52) is not equal to Zðn P the cross covariance matrix Z m U Gam ðn þ 1jnÞ can be

m ¼ Hðn þ 1Þ½X m ðn þ 1Þ  X^ ðn þ 1jnÞ

¼ Hðn þ 1Þ½Xðn þ 1Þ þ MGam ðn þ 1ÞU Gam ðn þ 1Þ n X ^ PðiÞW U  Gam ðiÞ  Hðn þ 1Þ½Xðn þ 1jnÞ

calculated as X m ðn þ 1jnÞ ¼ Ef½Z˜ ðn þ 1jnÞU˜ Gam ðn þ 1jnÞT g

i¼0

Z m U Gam

^ Gam ðn þ 1jnÞ þ M Gam ðn þ 1ÞU ^ þ 1jnÞ ¼ Hðn þ 1Þ½Xðn þ 1Þ  Xðn

¼ Ef½Hðn þ 1ÞX˜ ðn þ 1jnÞ þ Hðn þ 1ÞM Gam ðn þ 1ÞU˜ Gam ðn þ 1jnÞ n X T ˜ PðiÞW U  Hðn þ 1Þ Gam ðiÞU Gam ðn þ 1jnÞ g

þ Hðn þ 1ÞM Gam ðn þ 1Þ½U Gam ðn þ 1Þ n X ^ Gam ðn þ 1jnÞ  Hðn þ 1Þ PðiÞW U U Gam ðiÞ

i¼0

(53)

i¼0

2

(51)

i¼0

x 104

4

x 104 Actual position Proposed method

3

CA method Augmented method

-2

y (m)

x (m)

0

-4

2

Actual position Proposed method

-6

1

CA method Augmented method

-8

0

20

40 60 Time (sec)

80

0

100

140

20

40 60 Time (sec)

80

100

80

100

80 Proposed method

120

Proposed method

CA method

CA method

60

Augmented method

100

RMSE of y (m)

RMSE of x (m)

0

80 60 40

Augmented method

40

20

20 0

0

20

40 60 Time (sec)

80

100

0 0

20

40 60 Time (sec)

Fig. 1. The actual value and the estimation of x, y positions and RMS errors estimations by the proposed, the CA and the augmented methods.

ARTICLE IN PRESS 540

A. Karsaz, H. Khaloozadeh / Signal Processing 89 (2009) 532–547

P Since EfX˜ ðn þ 1jnÞU˜ Gam ðn þ 1jnÞT g and Ef½ ni¼1 PðiÞW U Gam ðiÞ U˜ Gam ðn þ 1jnÞT g are zero, then, X X ðn þ 1jnÞ ¼ Hðn þ 1ÞM Gam ðn þ 1Þ ðn þ 1jnÞ (54) m



X

ðn þ 1jnÞM Gam ðn þ 1ÞT Hðn þ 1ÞT

U Gam

þ

X ðn þ 1jnÞ1

(56)

Zm

U Gam

Z U Gam

m However, ½Z˜ ðn þ 1jnÞ  U˜ Gam ðn þ 1jnÞT is still a zeromean white process, thus (47) is an optimal solution

X

ðn þ 1jn þ 1Þ ¼

U Gam

4

ðn þ 1jnÞM Gam ðn þ 1ÞT Hðn þ 1ÞT

 ½3Hðn þ 1ÞMGam ðn þ 1Þ X  ðn þ 1jnÞMGam ðn þ 1ÞT Hðn þ 1ÞT

modified as below (see Appendix A):

U Gam

^ þ 1jn þ 1Þ ¼ U ^ Gam ðn þ 1jnÞ Uðn

X þ ðn þ 1jnÞ1 Hðn þ 1ÞM Gam ðn þ 1Þ

^˜ þ 1Þ ˜ þ 1Þ  Zðn þ K Gam ðn þ 1Þ½Zðn ^ Gam ðn þ 1jnÞ þ K Gam ðn þ 1Þ½Zðn ˜ þ 1Þ ¼U

m



^ Gam ðn þ 1jnÞ  Hðn þ 1ÞM Gam ðn þ 1ÞU

Z X

ðn þ 1jnÞ

(57)

U Gam

(55) where the Kalman gain of maneuver estimation is denoted

X

ðn þ 1jnÞ ¼ AGam ðnÞ

U Gam

X

ðnjnÞATGam ðnÞ þ Q U Gam ðnÞ

(58)

U Gam

X

by KGam(n+1):

ðn þ 1jnÞ ¼ Hðn þ 1ÞM Gam ðn þ 1Þ

Z m U Gam

X

ðn þ 1jnÞ

(59)

U Gam

ðn þ 1jnÞM Gam ðn þ 1ÞT Hðn þ 1ÞT Table 1 shows the summarized equations for both the proposed and the augmented methodologies.

U Gam

 ½3Hðn þ 1ÞM Gam ðn þ 1Þ

0

1400 Actual velocity

1200

-500

Proposed method CA method

vy (m/sec)

1000 vx (m/sec)

X U Gam

Kalman algorithm for acceleration state estimation is

K Gam ðn þ 1Þ ¼ 2

ðn þ 1jnÞ

U Gam

based on the Fisher model (see Appendix A). The new

X

X

-1000 -1500

Actual velocity

600 400

Proposed method

-2000

Augmented method

800

CA method

200

Augmented method

-2500

0 0

20

40 60 Time (sec)

80

100

0

100

40 60 Time (sec)

80

100

80

100

50 Proposed method

Proposed method

CA method

80

RMSE of vy (m/sec)

RMSE of vx (m/sec)

20

Augmented method

60 40 20 0

CA method

40

Augmented method

30 20 10 0

0

20

40 60 Time (sec)

80

100

0

20

40 60 Time (sec)

Fig. 2. The actual value and the estimation of vx(n), vy(n) and RMS errors of x and y velocities estimations by the proposed, the CA and the augmented methods.

ARTICLE IN PRESS A. Karsaz, H. Khaloozadeh / Signal Processing 89 (2009) 532–547

Remark 1. Table 1 shows, the proposed two-stage scheme involves the computation effort by two partitioned 4  4 subsystem matrices. In contrast, the augmented method deals with a 8  8 system matrix that is more computationally expensive than the proposed method. The proposed scheme also eliminates the extra computation P of the cross covariance matrix Z m UGam ðn þ 1jnÞ. Remark 2. The acceleration filter is turned on and off by a maneuvering detector. In some situations, which maneuvering is not declared, the proposed scheme uses only a 4  4 subsystem, but the augmented method still uses the 8  8 system.

541

state vector includes the acceleration and jerk parameters. Therefore, the state vector for two-dimensional case would have the eight components. The sampling time is T ¼ 0.1 s and target maneuvering is applied at ninth second (90th sample). The initial conditions are selected similar for the CA method and the augmented method as well as the proposed method. We consider the target initial conditions for state and acceleration vectors as below: X m ð0Þ ¼ ½2165 m

 80 m=s

U Gam ð0Þ ¼ ½0 g 0 g=s 0 g 0 g=s

1250 m

25 m=sT

T

X Aug ð0Þ ¼ ½2165 m  80 m=s 0 g 0 g=s 1250 m 25 m=s 0 g 0 g=s

4. Simulation results To evaluate the proposed algorithm, an example of target which turns, in two-dimensional space is simulated such as a ship or an aircraft with constant elevation. In this simulation example, the performance of the new algorithm for the maneuvering target tracking as the generalized acceleration model (GAM) has been compared with the work suggested in [10] as a CA method and the work of [5] as an augmented method. The maneuvering target in this simulation example has a medium jerk level value. As mentioned before in the augmented method the

The target begins to maneuver as U Gam ð9Þ ¼ ½0 g  0:07 g=s 0 g 0:04 g=sT for 9 sptp90 s. The system noise covariance matrices are selected for the CA, the augmented and the proposed methods, respectively, as mention in (12) and (26), (29) and (30), respectively. The measurement pstandard deviations of x and y ffiffiffiffiffiffi positions are: sx ¼ 10 10 m; sy ¼ 20 m. Thus, the noise covariance matrix is   1000 0 RðnÞ ¼ 0 400

4

0

Actual Acceleration Proposed method

3

-1 uy (m/sec2)

ux (m/sec2)

CA method

-2 -3 -4

Actual Acceleration

Augmented method

2 1

Proposed method

-5

0

CA method

-6

Augmented method

-1 20

40 60 Time (sec)

80

100

0

Proposed method

3

2

2.5

40 60 Time (sec)

80

100

Proposed method CA method

Augmented method

RMSE uy (m/sec2)

RMSE of ux (m/sec2)

CA method

20

2 1.5 1 0.5 0

Augmented method

1.5 1 0.5 0

0

20

40

60

Time (sec)

80

20

40 60 Time (sec)

80

Fig. 3. The actual value and the estimation of acceleration in x and y directions and corresponding RMS errors by the proposed method compared with the CA and the augmented methods.

ARTICLE IN PRESS 542

A. Karsaz, H. Khaloozadeh / Signal Processing 89 (2009) 532–547

for all methods. The jerk time period is t ¼ 81 s so a ¼ 0.01223 s1. The root mean square error (RMSE) index is used for results evaluation. Fig. 1 shows the actual values and the estimations of x, y and RMS errors of x and y positions by the proposed, the CA and the augmented methods.

Fig. 2 shows the actual values and the estimations of vx(n), vy(n) and the RMS errors of x and y velocities by the proposed method compared with the CA and the augmented methods. The actual values and the estimations of the acceleration in the x and y directions and their corresponding RMS errors can be seen in Fig. 3.

Jerk (m/sec3)

0.1 0.08 0.06

Actual Jerk

0.04

Proposed method Augmented method

0.02 0 0

10

20

30

40 50 time (sec)

60

70

80

90

RMSE of jerk (m/sec3)

0.06 Proposed method

0.05

Augmented method

0.04 0.03 0.02 0.01 0 0

10

20

30

40 50 time (sec)

60

70

80

90

Fig. 4. The actual value and the estimation of jerk parameter and RMS errors by the proposed method compared with the augmented method.

8

x 104

Range (m)

Actual range Proposed method

6

CA method Augmented method

4 2 0

RMSE of Range (m)

0

150

10

20

30

40 50 time (sec)

60

70

80

90

30

40 50 time (sec)

60

70

80

90

Proposed method CA method

100

Augmented method

50 0 10

20

Fig. 5. The actual value and the estimation of target range and corresponding RMS errors estimation by the proposed method compared with the CA and the augmented methods.

ARTICLE IN PRESS A. Karsaz, H. Khaloozadeh / Signal Processing 89 (2009) 532–547

Azimuth (deg)

Fig. 4 displays the actual values and the estimation of jerk parameter which are evaluated by the proposed and the augmented methods. Note that the jerk value as the acceleration rate in the CA approach is equal to zero. The actual values and the estimations of the target range, azimuth and course are displayed in Figs. 5–7, respectively. Simulation results show a high performance of the proposed innovative algorithm especially in azimuth and

course tracking which is major deficiency of the most common algorithms for highly maneuvering targets. Due to the length constraints, a complete study of the angles tracking in all its aspects and target scenarios is beyond the scope of this paper. Since the performance of the proposed algorithm are presented as estimation of position, velocity, acceleration, jerk, range, azimuth, and course tracking, comparison based on the Monte Carlo simulation is also made to evaluate the performances

140 120 100 80 60 40 20

Actual azimuth Proposed method CA method Augmented method

10

RMSE of azimuth (deg)

543

20

30

40 50 time (sec)

60

70

80

90

1 Proposed method CA method Augmented method

0.8 0.6 0.4 0.2 0 0

10

20

30

40 50 time (sec)

60

70

80

90

Fig. 6. The actual value and the estimation of target azimuth and corresponding RMS errors estimation by the proposed, the CA and the augmented methods.

Course (deg)

165 Actual course Proposed method CA method Augmented method

160

155

150

RMSE of Course (deg)

0

10

20

30

40 50 time (sec)

60

70

80

90

3 Proposed method

2.5

CA method Augmented method

2 1.5 1 0.5 0 0

10

20

30

40 50 time (sec)

60

70

80

90

Fig. 7. The actual value and the estimation of target course and its RMS error estimation by the proposed, the CA and the augmented methods.

ARTICLE IN PRESS 544

A. Karsaz, H. Khaloozadeh / Signal Processing 89 (2009) 532–547

Table 2 The Monte Carlo simulation results for 300 runs on low, medium, and high jerks scenarios. Scenarios

Low jerk T ¼ 0.5 s jx(n) ¼ 0.007 m/s3 Jy(n) ¼ 0.004 m/s3

Medium jerk T ¼ 0.5 s jx(n) ¼ 0.07 m/s3 jy(n) ¼ 0.04 m/s3

High jerk T ¼ 0.01 s jx(n) ¼ 0.7 m/s3 jy(n) ¼ 0.4 m/s3

Performance (RMSE)

Proposed method

CA method

Augmented method

X-position (m) Y-position (m) X-velocity (m/s) Y-velocity (m/s) X-acceleration (m/s2) Y-acceleration (m/s2) Jerk (m/s3) Range (m) Azimuth (deg) Course (deg)

11.2021 7.5614 2.8382 1.8276 0.2550 0.1643 0.0086 13.6784 0.1915 3.4683

50.7160 26.3637 13.7643 11.8899 0.4420 0.4341 – 56.4084 1.6415 21.0823

11.2060 7.5558 2.8478 1.8206 0.2364 0.1601 0.0086 13.6631 0.2225 3.5666

X-position (m) Y-position (m) X-velocity (m/s) Y-velocity (m/s) X-acceleration (m/s2) Y-acceleration (m/s2) Jerk (m/s3) Range (m) Azimuth (deg) Course (deg)

9.2781 5.8344 5.5573 3.4733 0.1576 0.0950 0.0050 11.9841 0.1520 1.3886

120.1938 70.1533 54.2356 29.1381 1.8901 1.0750 – 138.3253 0.2924 1.4088

9.5022 5.9034 5.6573 3.4736 0.1770 0.1053 0.0050 12.0041 0.1691 1.4285

X-position (m) Y-position (m) X-velocity (m/s) Y-velocity (m/s) X-acceleration (m/s2) Y-acceleration (m/s2) Jerk (m/s3) Range (m) Azimuth (deg) Course (deg)

8.3623 4.2942 7.0096 8.0455 0.0622 0.0850 0.0003 10.2753 0.0747 8.3708

403.1173 270.1547 91.3607 101.0016 2.5068 3.7299 – 485.2701 1.1847 27.9692

8.3655 4.2949 7.0575 8.0485 0.0689 0.0816 0.0003 10.2754 0.0651 8.4094

pffiffiffiffiffiffiffi 0 sptp90 s, t ¼ 81 s, sx ¼ 10 10 m, sy ¼ 20 m.

Table 3 The Monte Carlo simulation results for 300 runs on medium jerk value with two sampling time T ¼ 2 s and T ¼ 1 s. Scenarios

Performance (RMSE)

Proposed method

CA method

Augmented method

Medium jerk value T ¼ 2s jx ðnÞ ¼ 0:07 m=s3 jy ðnÞ ¼ 0:04 m=s3

X-position (m) Y-position (m) X-velocity (m/s) Y-velocity (m/s) X-acceleration (m/s2)

12.0782 7.8339 1.9573 1.8723 0.1976

230.1278 150.1525 20.2302 18.1181 2.9871

12.0785 7.8544 1.9579 1.9015 0.1987

Y-acceleration (m/s2) Jerk (m/s3) Range (m) Azimuth (deg) Course (deg)

0.1410 0.0070 14.891 0.1967 3.3843

2.1730 – 274.5823 1.2014 4.4078

0.1347 0.0070 14.893 0.1977 3.4093

X-position (m) Y-position (m) X-velocity (m/s) Y-velocity (m/s) X-acceleration (m/s2) Y-acceleration (m/s2)

8.0782 5.833 1.6573 1.4723 0.1576 0.1310

180.1838 100.1345 20.2345 18.1671 3.0701 2.9826

8.0785 5.8344 1.6579 1.4715 0.1587 0.1357

Jerk (m/s3) Range (m) Azimuth (deg) Course (deg)

0.0050 9.8910 0.1567 2.3843

– 219.5073 1.2814 5.4074

0.0050 9.8920 0.1577 2.3993

Medium jerk value T ¼ 1s jx ðnÞ ¼ 0:07 m=s3 jy ðnÞ ¼ 0:04 m=s3

pffiffiffiffiffiffi 0 sptp180 s, t ¼ 171 s, sx ¼ 10 10 m, sy ¼ 20 m.

ARTICLE IN PRESS A. Karsaz, H. Khaloozadeh / Signal Processing 89 (2009) 532–547

Table 4 Computational cost of the proposed and the augmented methods for 1000 different scenarios. Algorithm

The proposed method

The augmented method

Averaged computation time (s)

0.0176

0.0254

between the proposed, the CA and the augmented methods. Performance evaluation of the methods is measured by the RMSE index over 300 runs simulations. This statistic is computed ignoring initial transitions. In Table 2, the results of the state estimations for the above mentioned methods are summarized. The RMSE is computed with 300 Monte Carlo runs on three different scenarios: low, medium, and high jerk values with different values of sampling time 0.01, 0.1, and 0.5 s. The Monte Carlo simulation results for 300 runs on medium jerk value with two sampling time T ¼ 2 s and T ¼ 1 s are also shown in Table 3. Tables 2 and 3 show that the proposed and the augmented methods can significantly improve the results of the target parameters estimations and they are superior to the CA method. Due to some random parameters in the algorithm, some little differences between the obtained results by the proposed and the augmented method would become. Computational loads of the proposed and the augmented methods for more than 1000 different scenarios are compared in Table 4. All results are obtained based on simulations using the tic/toc functions of Matlab software on a 2.4 GHz Intel P4 platform. It can be seen in this table that the proposed method requires smaller computational (about 30%) effort than the augmented method. In the highly maneuvering target tracking situations, the proposed and the augmented methods have a much better performance than the CA approach, because of the presence of jerk in maneuvering. This is true even when the jerk value is rather low. The performance of the proposed algorithm is as well as the results obtained by augmented method in the low, typical and high jerk value situations.

5. Conclusions For an arbitrary highly maneuvering target, the motion profile have higher order non-zero time-derivatives of acceleration which should be considered in target motion modeling for more accurate target tracking. Therefore, the performance in traditional IDE approaches which are based on the state only up to the second order of timederivatives is reduced in tracking the highly maneuvering targets. The proposed scheme in this paper intends to overcome the drawbacks of some similar works, such as computational expenses, based on the augmented methods using a new generalized dynamic modeling of the acceleration and the state vector dimensions reduction. The optimality of the proposed two-stage modified Kal-

545

man filter has been proved. Comparisons with two other estimators show that the proposed scheme has a computational advantage over the augmented algorithms and it has a significant improvement on the input estimation techniques. Appendix A In the following, we are going to derive a recursive P form of (47) to estimate UGam(n+1) such that U Gam ðn þ 1jn þ 1Þ is minimized. Based on (41) it is desired that obtaining the recursive algorithm for the acceleration state estimation in the form of the SKF as below: ^ þ 1jn þ 1Þ ¼ U ^ Gam ðn þ 1jnÞ Uðn ^˜ þ 1Þ ˜ þ 1Þ  Zðn þ K Gam ðn þ 1Þ½Zðn ^ Gam ðn þ 1jnÞ þ K Gam ðn þ 1Þ½Zðn ˜ þ 1Þ ¼U ^ Gam ðn þ 1jnÞ(A.1)  Hðn þ 1ÞM Gam ðn þ 1ÞU The residue of UGam(n+1) defined in (42), by using (A.1) and (41) gives, U˜ Gam ðn þ 1jn þ 1Þ ^ Gam ðn þ 1jn þ 1Þ ¼ U Gam ðn þ 1Þ  U ^ Gam ðn þ 1jnÞ  K Gam ðn þ 1Þ½Zðn ˜ þ 1Þ ¼ U Gam ðn þ 1Þ  U ^ Gam ðn þ 1jnÞ  Hðn þ 1ÞMGam ðn þ 1ÞU ¼ U˜ Gam ðn þ 1jnÞ  K Gam ðn þ 1Þ½Hðn þ 1ÞM Gam ðn þ 1ÞU Gam ðn þ 1Þ m ^ Gam ðn þ 1jnÞ þ Z˜ ðn þ 1jnÞ  Hðn þ 1ÞM Gam ðn þ 1ÞU ¼ U˜ Gam ðn þ 1jnÞ  K Gam ðn þ 1ÞHðn þ 1ÞMGam ðn þ 1ÞU˜ Gam ðn þ 1jnÞ m  K Gam ðn þ 1ÞZ˜ ðn þ 1jnÞ

(A.2)

Also the covariance matrix of U˜ Gam ðn þ 1jn þ 1Þ defined in (A.2), using (54) gives X ðn þ 1jn þ 1Þ U Gam

¼

X

ðn þ 1jnÞ

U Gam

þ K Gam ðn þ 1ÞHðn þ 1ÞMGam ðn þ 1Þ

X

ðn þ 1jnÞ

U Gam

 M Gam ðn þ 1ÞT Hðn þ 1ÞT K Gam ðn þ 1Þ X þ K Gam ðn þ 1Þ ðn þ 1jnÞK Gam ðn þ 1ÞT 

X

Zm

ðn þ 1jnÞM Gam ðn þ 1ÞT Hðn þ 1ÞT K Gam ðn þ 1ÞT

U Gam



X

ðn þ 1jnÞK Gam ðn þ 1ÞT

U Gam Z m

 K Gam ðn þ 1ÞHðn þ 1ÞMGam ðn þ 1Þ

X

ðn þ 1jnÞ

U Gam

þ K Gam ðn þ 1ÞHðn þ 1ÞMGam ðn þ 1Þ X  ðn þ 1jnÞK Gam ðn þ 1ÞT U Gam Z m

 K Gam ðn þ 1Þ

X

ðn þ 1jnÞ þ K Gam ðn þ 1Þ

m



X Z m U Gam

Z U Gam

ðn þ 1jnÞM Gam ðn þ 1ÞT Hðn þ 1ÞT K Gam ðn þ 1ÞT

ARTICLE IN PRESS 546

A. Karsaz, H. Khaloozadeh / Signal Processing 89 (2009) 532–547

¼

X

P The minimum UGam ðn þ 1jn þ 1Þ is obtained by setting KGam(n+1)W equal to D. Therefore,

ðn þ 1jnÞ

U Gam

þ K Gam ðn þ 1ÞHðn þ 1ÞM Gam ðn þ 1Þ

X

ðn þ 1jnÞ

U Gam T

T

K Gam ðn þ 1Þ ¼ 2

T

 M Gam ðn þ 1Þ Hðn þ 1Þ K Gam ðn þ 1Þ X þ K Gam ðn þ 1Þ ðn þ 1jnÞK Gam ðn þ 1ÞT 

X

Z

X

ðn þ 1jnÞM Gam ðn þ 1ÞT Hðn þ 1ÞT ðW T ÞW 1

U Gam

¼2

X

ðn þ 1jnÞM Gam ðn þ 1ÞT Hðn þ 1ÞT

U Gam

m

 ½3Hðn þ 1ÞM Gam ðn þ 1Þ X  ðn þ 1jnÞMGam ðn þ 1ÞT Hðn þ 1ÞT

ðn þ 1jnÞMGam ðn þ 1ÞT Hðn þ 1ÞT K Gam ðn þ 1ÞT

U Gam



X

ðn þ 1jnÞMGam ðn þ 1ÞT Hðn þ 1ÞT K Gam ðn þ 1ÞT

U Gam

 K Gam ðn þ 1ÞHðn þ 1ÞM Gam ðn þ 1Þ

X

U Gam

þ

X ðn þ 1jnÞ1

ðn þ 1jnÞ

U Gam

þ K Gam ðn þ 1ÞHðn þ 1ÞM Gam ðn þ 1Þ X  ðn þ 1jnÞMGam ðn þ 1ÞT Hðn þ 1ÞT K Gam ðn þ 1ÞT U Gam

 K Gam ðn þ 1ÞHðn þ 1ÞM Gam ðn þ 1Þ

X

The minimum error covariance obtained: X

ðn þ 1jnÞ

þ K Gam ðn þ 1ÞHðn þ 1ÞM Gam ðn þ 1Þ X  ðn þ 1jnÞMGam ðn þ 1ÞT Hðn þ 1ÞT K Gam ðn þ 1ÞT U Gam (A.3)

¼

X

¼

¼

X

ðn þ 1jnÞ

U Gam

ðn þ 1jnÞ þ 3K Gam ðn þ 1ÞHðn þ 1ÞM Gam ðn þ 1Þ

U Gam



T

þ K Gam ðn þ 1Þ 2

T

ðn þ 1jnÞMGam ðn þ 1Þ Hðn þ 1Þ K Gam ðn þ 1Þ

U Gam

X

X

ðn þ 1jnÞM Gam ðn þ 1ÞT Hðn þ 1ÞT ðWW T Þ1

U Gam

 ½2

X

þ 1jn þ 1Þ is

ðn þ 1jnÞ  DDT

 ½2 X

U Gam ðn

U Gam

ðn þ 1jn þ 1Þ

U Gam

P

ðn þ 1jn þ 1Þ

U Gam

U Gam

X

(A.8)

Zm

T

¼

X

ðn þ 1jnÞM Gam ðn þ 1ÞT Hðn þ 1ÞT

U Gam

U

ðn þ 1jnÞM Gam ðn þ 1ÞT Hðn þ 1ÞT K Gam ðn þ 1ÞT X

ðn þ 1jnÞ  4

 ½3Hðn þ 1ÞMGam ðn þ 1Þ X  ðn þ 1jnÞMGam ðn þ 1ÞT Hðn þ 1ÞT

Zm

U Gam

ðn þ 1jnÞM Gam ðn þ 1ÞT Hðn þ 1ÞT T

U Gam

U Gam

X ðn þ 1jnÞK Gam ðn þ 1ÞT

 2K Gam ðn þ 1ÞHðn þ 1ÞM Gam ðn þ 1Þ

X

X

þ

Gam X X ðn þ 1jnÞ1 Hðn þ 1ÞM Gam ðn þ 1Þ ðn þ 1jnÞ

Zm

ðn þ 1jnÞ

U Gam

(A.9)

U Gam

(A.4) The second and third terms of (A.4) can be written as T (n+1)WWTK Gam Gam(n+1) , where

K

WW T ¼ 3Hðn þ 1ÞM Gam ðn þ 1Þ X  ðn þ 1jnÞM Gam ðn þ 1ÞT Hðn þ 1ÞT X ðn þ 1jnÞ

X

ðn þ 1jnÞ ¼ AGam ðnÞ

U Gam

X

ðnjnÞATGam ðnÞ þ Q U Gam ðnÞ

(A.10)

U Gam

If the value of KGam(n+1) given by (A.8) substitute in ^ Gam ðn þ 1Þ leads to a minimum (A.1), the estimation of U error covariance matrix.

U Gam

þ

Based on (16), we have

(A.5)

Zm

References P

Since the covariance matrices U Gam ðn þ 1jnÞ and P Z m ðn þ 1jnÞ are symmetric, we can find a decomposition in the form of WWT, for appropriate matrix W, thus X

ðn þ 1jn þ 1Þ ¼

U Gam

X

ðn þ 1jnÞ

U Gam

þ ½K Gam ðn þ 1ÞW  D½K Gam ðn þ 1ÞW  DT  DDT

(A.6)

Comparing (A.4) and (A.6), the term D is obtained as: D¼2

X U Gam

ðn þ 1jnÞM Gam ðn þ 1ÞT Hðn þ 1ÞT ðW T Þ1

(A.7)

[1] S.S. Blackman, R.F. Popoli, Design and Analysis of Modern Tracking Systems, Artech House, Norwood, MA, 1999. [2] Y. Bar-Shalom, X.R. Li, Estimation and Tracking: Principles, Techniques, and Software, Artech House, Boston, MA, 1993. [3] Y. Bar-Shalom, X.R. Li, T. Kirubarajan, Estimation with Applications to Tracking and Navigation: Theory, Algorithms, and Software, Wiley, New York, 2001. [4] H. Khaloozadeh, A. Karsaz, A new state augmentation for maneuvering targets detection, in: Proceedings of the Third IEEE International Conference on Signal Processing & Communication (SPCOM), December 2005, pp. 65–69. [5] K. Mehrotra, P.R. Mahapatra, A jerk model for tracking highly maneuvering targets, IEEE Trans. Aerosp. Electron. Syst. 33 (4) (1997) 1094–1105. [6] P. Mookerjee, F. Reifler, Reduced state estimators for consistent tracking of maneuvering targets, IEEE Trans. Aerosp. Electron. Syst. 41 (2) (1999) 608–619.

ARTICLE IN PRESS A. Karsaz, H. Khaloozadeh / Signal Processing 89 (2009) 532–547

[7] X.R. Li, V.P. Jilkov, A survey of maneuvering target tracking–part IV: decision-based methods, in: Proceedings of the 2002 SPIE Conference on Signal and Data Processing of Small Targets, vol. 4728, Orlando, FL, April 2002, pp. 511–534. [8] X.R. Li, Y. Bar-Shalom, A recursive multiple model approach to noise identification, IEEE Trans. Aerosp. Electron. Syst. 30 (3) (1994) 671–684. [9] G.P. Cardillo, A.V. Mrstik, T. Plambeck, A track filter for reentry objects with uncertain drag, IEEE Trans. Aerosp. Electron. Syst. 35 (2) (1999) 395–409. [10] T.C. Wang, P.K. Varshney, A tracking algorithm for maneuvering targets, IEEE Trans. Aerosp. Electron. Syst. 29 (3) (1993) 910–924. [11] Y.T. Chan, A.G.C. Hu, J.B. Plant, A Kalman filter based tracking scheme with input estimation, IEEE Trans. Aerosp. Electron. Syst. (15) (1979) 237–244. [12] H. Lee, M.J. Tahk, Generalized input-estimation technique for tracking maneuvering targets, IEEE Trans. Aerosp. Electron. Syst. 35 (4) (1999) 1388–1402. [13] X.R. Li, Y.M. Zhang, X.R. Zhi, Multiple-model estimation with variable structure-part IV: design and evaluation of model-group switching algorithm, IEEE Trans. Aerosp. Electron. Syst. 35 (1) (1999) 242–254. [14] I.H. Whang, J. Lee, T. Sung, Modified input estimation technique using pseudo-residuals, IEEE Trans. Aerosp. Electron. Syst. 30 (1) (1994) 220–228. [15] P.L. Bogler, Tracking a maneuvering target using input estimation, IEEE Trans. Aerosp. Electron. Syst. 23 (3) (1987) 298–310. [16] Y.H. Park, J.H. Seo, J.G. Lee, Tracking using variable-dimension filter with input estimation, IEEE Trans. Aerosp. Electron. Syst. 31 (1) (1995) 399–408. [17] H.A.P. Blom, T. Bar-Shalom, The interacting multiple model algorithm for systems with Markovian switching coefficients, IEEE Trans. Autom. Control 33 (8) (1988) 780–783. [18] E. Mazor, A. Averbuch, Y. Bar-Shalom, J. Dayan, Interacting multiple model methods in target tracking: a survey, IEEE Trans. Aerosp. Electron. Syst. 34 (1) (1998) 103–123. [19] Y. Bar-Shalom, X. Li, Multitarget–Multisensor Tracking: Principles and Techniques, YBS Press, Storrs, CT, 1995. [20] R.R. Pitre, V.P. Jilkov, R.A. Li, Comparative study of multiple-model algorithms for maneuvering target tracking, in: Proceedings of the 2005 SPIE Conference on Signal Processing, Sensor Fusion, and Target Recognition vol. XIV, Orlando, FL, March–April 2005.

547

[21] F.B. Duh, C.T. Lin, Tracking a maneuvering target using neural fuzzy network, IEEE Trans. Syst. Man Cybern. Part B 34 (1) (2004) 16–33. [22] H. Bahari, A, Karsaz, H. Khaloozadeh, A new algorithm based on combined fuzzy logic and Kalman filter for target maneuver detection, in: International IEEE Conference ISSCAA, Harbin, China, 19-21 January 2006, pp. 520–524. [23] H. Bahari, A, Karsaz, H. Khaloozadeh, High maneuver target tracking based on combined Kalman filter and fuzzy logic, in: IEEE 2007 Conference Information, Decision and Control (IDC), Australia, 11–14 February 2007, pp. 1020–1027. [24] B. Friedland, Treatment of bias in recursive filtering, IEEE Trans. Autom. Control 14 (1969) 359–367. [25] W.D. Blair, Fixed-gain two-stage estimator for tracking maneuvering targets, IEEE Trans. Aerosp. Electron. Syst. 29 (3) (1993) 1004–1014. [26] C.S. Hsieh, F.C. Chen, Optimal solution of the two-stage Kalman estimator, IEEE Trans. Autom. Control 44 (1) (1999) 195–199. [27] C.S. Hsieh, F.C. Chen, General two-stage Kalman filters, IEEE Trans. Autom. Control 45 (4) (2000) 819–824. [28] T. Kawase, H. Tsurunosono, N. Ehara, I. Sasase, Two-stage Kalman estimator using an advanced circular prediction for tracking highly maneuvering targets, in: Proceedings of the ICASSP, 1998. [29] H.Z. Qiu, H.Y. Zhang, X.F. Sun, Solution of two-stage Kalman filter, IEE Proc. Control Theory Appl. 152 (2) (2005) 152–156. [30] F.C. Schweppe, Uncertain Dynamic Systems, Prentice-Hall, Englewood Cliffs, NJ, 1973. [31] F. Gustafsson, Adaptive Filtering and Change Detection, Wiley, New York, 2001. [32] F. Gustafsson, The marginalized likelihood test for detecting abrupt changes, IEEE Trans. Autom. Control 41 (1996) 66–78. [33] D.P. Malladi, J.L. Speyer, A generalized Shiryaev sequential probability ratio test for change detection and isolation, IEEE Trans. Autom. Control 44 (8) (1999) 1522–1534. [34] A. Karsaz, H. Khaloozadeh, N. Pariz, A new algorithm based on generalized target maneuver detection, in: Sixth IEEE International Conference on Control and Automation (ICCA), Guangzhou, China, 30 May 1 June. 2007, pp. 3034–3039. [35] R.A. Singer, Estimating optimal tracking filter performance for manned maneuvering targets, IEEE Trans. Aerosp. Electron. Syst. 6 (1970) 473–483.