# Maneuvering Detection Using Multiple Parallel CUSUM Detector

Authors: Hongqiang Liu; Zhongliang Zhou; Chunguang Lu
Journal: Mathematical Problems in Engineering
 (2018)

---

## Abstract

The switching model tracking algorithm based on hard decisions is an important method to solve the maneuvering target tracking problem. The use of Doppler velocity not only helps shorten the delay time of maneuvering detection but also provides information about the target motion model. A novel target maneuvering detection method named Multiple Parallel Cumulative Sum (M-CUSUM) for target multiple maneuvering models is proposed in this paper based on Doppler velocity. The main scheme of the proposed approach consists of the following: firstly, the problem framework of multiple model maneuvering detection is put forward; secondly, the statistic of acceleration is obtained through modeling the mapping relationship between Doppler velocity and the normal/tangential acceleration according to the geometry and kinematics; thirdly, the joint empirical distribution of the normal/tangential acceleration is obtained by the statistical experiment method and then the approximate joint probability distribution function of the normal/tangential acceleration is acquired by use of Gaussian Mixture Model (GMM) with Expectation Maximization (EM) algorithm; fourthly, it is taken as the prior information of target maneuvering which is introduced to the likelihood ratio of prediction measurement residual by the marginalization method; finally, the standard Cumulative Sum (CUSUM) detector is extended as Multiple Parallel CUSUM detector. Simulation results show that M-CUSUM detector has a smaller maneuver onset detection delay time compared with similar detectors and has the ability of pattern recognition of target maneuvers.


---

## Body

## 1. Introduction

Tracking the manned maneuvering aircraft by use of radar is a maneuvering target tracking (MTT) problem [1–3]. It is hard to accurately estimate the motion state of maneuvering target by using only single motion model, because the current radar technology cannot obtain the observation about accelerations directly supporting the selection of the target motion model. The switching model tracking (SMT) algorithm [4] based on (the) hard decision and the interacting multiple model (IMM) tracking algorithm [5–7] based on the soft decision are two kinds of basic methods to solve MTT. The SMT algorithms use a maneuvering detector to judge the target maneuvering onset. When the target maneuvering is successfully detected, the motion model of tracking algorithm is immediately switched. The timely and correct maneuvering detection is the key to the SMT algorithm. The most maneuver detectors [8] only use the position measurement conditioned on the linear Gauss hypothesis. As the position measurement cannot directly reflect the maneuver state, however, radar can obtain not only the target position measurement but also Doppler velocity measurement directly affected by the maneuvering [9]. The maneuvering state of target can be obtained by use of Doppler velocity. Some scholars obtain the turn rate estimation from Doppler velocity [10–12]. As the acceleration is directly related to the target maneuvering type, other scholars obtain the statistics of the target acceleration directly through Doppler velocity. Bizup and Brown [13] propose a method of maneuvering target detection using Doppler velocity. By assuming the constant turn rate motion, a statistic cmin of acceleration is deduced from Doppler velocity. Ru et al. [14] combine the statistic cmin with the tangential acceleration to obtain the statistic of total acceleration, cmin2, and cmin3. The main difference between cmin2 and cmin3 is whether Doppler velocity measurement is applied in the filtering. It turns out that the statistic cmin2 or cmin3 using Cumulative Sum (CUSUM) detector has a shorter detection onset delay time than cmin [15]. Lu et al. [16] use Mahalanobis distance and Euclidean distance optimization method to give two statistics of the total acceleration Emin and Mmin, respectively, which overcome the problem that the acceleration has no solution when the noise is higher. Using the statistics in the above references to detect and track the target, there are the following problems.Assuming a Constant Turn Rate Motion.The statistics cmin, cmin2, and cmin3 are obtained by assuming a constant turn rate motion. However, in addition to the normal acceleration, the target maneuver usually has a tangential acceleration and the trajectory exhibits a curvilinear shape. There is an inherent system bias when the statistics cmin, cmin2, and cmin3 are used in the maneuver target detection.One-to-Many Mapping Relationship.There is a one-to-many mapping relationship between Doppler velocity and the normal acceleration. It is assumed that the target velocity direction change is less than 2π in the radar scanning period, so the one-to-many mapping is restricted to a pair of four mappings in the statistics cmin, cmin2, and cmin3. In practice, the tracking radar scanning period is less than 1 s and the maximum instantaneous turn rate of the manned maneuvering aircraft is not more than π/4 [17]. Therefore, if one-to-four mapping is applied to track such target, it will cause the dispersion of probability weight and affect the sensitivity of detector as some unlikely turns are considered.Unknown Motion Model.The detector of total acceleration based on Doppler velocity can detect target maneuver rapidly and distinguish the maneuvering and nonmaneuvering state. After the maneuvering is detected, the motion model with larger acceleration noise is usually chosen to filter. However, the motion model mismatch problem cannot be completely solved by increasing the acceleration noise. If a maneuver detector can not only detect whether the target is making a maneuvering or not but also recognize the target maneuvering model, we could select the matched motion model to improve the tracking accuracy with a lower acceleration noise.To solve the above problems, the mapping relationship between Doppler velocity and the normal/tangential acceleration is deduced without assuming a constant turn rate motion. Meanwhile, it is assumed that the target velocity direction change is less thanπ during the radar scanning period, and one-to-four mapping relationship is restricted to one-to-two which reduces the dispersion degree of probability weight. In this paper, two parameters, namely, normal/tangential acceleration, are used to describe the target maneuver, and the single total acceleration is not adopted. The joint empirical distribution of the normal/tangential acceleration is obtained by the statistical experiment method [15] and approximated by use of Gaussian Mixture Model (GMM) with Expectation Maximization (EM) algorithm [18]. Then the approximate joint probability distribution function (PDF) of the normal/tangential accelerations is obtained. The approximate joint PDF is taken as the a priori information of target maneuvering which is introduced to calculate the likelihood ratio of prediction measurement residual by the marginalization method. A Multiple Parallel Cumulative Sum (M-CUSUM) detector is proposed by extending the standard CUSUM detector. In addition, it has been proposed that the use of Doppler velocity measurement can improve tracking accuracy [19–21]. In this paper, Doppler velocity is used separately in the maneuvering detector and filtering. As the measurement equation with Doppler velocity is highly nonlinear, the traditional EKF algorithm [22] does not work well. We use the measurement conversion algorithm proposed by Jifeng and Huimin to construct the linear KF filter.The content of this paper is designed as follows. We describe the tracking and detection problem and propose a concept framework of multiple maneuver model detection and recognition in Section2. We present the method of using Doppler velocity to deduce the normal/tangential acceleration in Section 3. We present how to obtain the empirical distribution of the normal/tangential acceleration and the approximated GMM with EM algorithm in Section 4. We describe M-CUSUM detector for multiple maneuver model detection and recognition in Section 5. In Section 6, firstly the results of the empirical distribution of the normal/tangential acceleration and GMM fitting are shown in five scenes; secondly, M-CUSUM proposed in this paper is compared with other maneuver detection methods; finally, it is applied to the SMT algorithm and is compared with IMM algorithm and CSAF algorithm. We summarize the research content of this paper in Section 7. In Appendix A, the modified measurement conversion algorithm and filter algorithm are given. In Appendix B, The EM algorithm is given to learn the parameters in GMM. In Appendix C, the state space equation of CSAF algorithm is given in two-dimensional space.


## 2. Problem Descriptions

### 2.1. Target Maneuver Model and Measurement Model

An aircraft is moving in a two-dimensional plane whose acceleration can be decomposed into a normal acceleration perpendicular to the velocity direction and a tangential acceleration along the velocity direction. The normal/tangential accelerations are used as the input vector of the motion model and then the target maneuver model [18] is expressed as(1)Xk+1=FXk+ΓXkAk+Wk,where the state vector is Xk=[xk,vx,k,yk,vy,k]′ in Cartesian coordinate system; xk,yk are the position component; vx,k,vy,k are the velocity component; Ak=[al,k,an,k]′ is the acceleration vector; al,k is the tangential acceleration; an,k is the normal acceleration; Wk is a white Gaussian noise vector and its covariance matrix is Qk:(2)
Q
k
=
q
·
Q
′
O
O
Q
′
,
Q
′
=
1
3
T
3
1
2
T
2
1
2
T
2
T
,
O
=
0
0
0
0
.F and Γ(Xk) are as follows:(3)F=1T000100001T0001,ΓXk=T2ck2-T2sk2Tck-Tsk-T2sk2T2ck2Tsk-Tck,ck=cos⁡αk,sk=sin⁡αk,where T is the radar scanning period; αk is the target velocity direction.Assume that the radar is located at the origin of coordinateXOY. The measurement equation of radar which includes Doppler velocity is(4)Zk=HXk+Vk,where Zk=[rk,θk,vr,k]′ is the observation vector; rk is the range; θk is the azimuth; vr,k is Doppler velocity; Vk is a white Gaussian noise and its covariance matrix is Rk=diag⁡σr2,σθ2,σvr2; σr2,σθ2, and σvr2 are variances of the range, azimuth, and Doppler velocity, respectively; H(Xk) is the nonlinear measurement function vector:(5)HXk=xk2+yk2±arccos⁡xkxk2+yk2xkvx,k+ykvy,kxk2+yk2,where “±” is “+” if yk≥0 and “±” is “−” if yk<0.


### 2.2. The Problem Framework of Multiple Maneuvering Detection

We define the nonmaneuvering target model as a constant velocity motion model; that is,A∣A=0. Then the maneuvering target model is denoted as A∣A≠0. The traditional maneuvering detection problem is expressed as follows:(H0)
The target does not maneuver;Am∈A∣A=0 for m=1,…,k, where m is the radar sample time.(H1)
The target starts to maneuver at unknown timen; Am∈A∣A≠0 for m=n,…,k.For the above maneuvering detection problem, the traditional maneuver detector can accomplish the judgment between the maneuver and nonmaneuver, that is to say, to solve a simple binary hypothesis [23] test task. However, decision-makers often want to get the information of the target maneuver type. Then hypothesis H1 can be decomposed into multiple branches and the maneuver detection problem can be described as a new form:(H0)
The target does not maneuver;Am∈A∣A=0 for m=1,…,k.(Hi)
The target starts to maneuver at unknown timen; Am∈A∣A=Ai for m=n,…,k and i=1,…,M.A i represents the ith maneuvering model; M is the total number of maneuvering target models. In particular, the nonmaneuver A=0 is A=[0,0]′, which means that al=0 and an=0. The maneuvering A=Ai is Ai=[ali,ani]′, which can represent the different maneuvering model if [ali,ani] is the different nonzero vector.


## 2.1. Target Maneuver Model and Measurement Model

An aircraft is moving in a two-dimensional plane whose acceleration can be decomposed into a normal acceleration perpendicular to the velocity direction and a tangential acceleration along the velocity direction. The normal/tangential accelerations are used as the input vector of the motion model and then the target maneuver model [18] is expressed as(1)Xk+1=FXk+ΓXkAk+Wk,where the state vector is Xk=[xk,vx,k,yk,vy,k]′ in Cartesian coordinate system; xk,yk are the position component; vx,k,vy,k are the velocity component; Ak=[al,k,an,k]′ is the acceleration vector; al,k is the tangential acceleration; an,k is the normal acceleration; Wk is a white Gaussian noise vector and its covariance matrix is Qk:(2)
Q
k
=
q
·
Q
′
O
O
Q
′
,
Q
′
=
1
3
T
3
1
2
T
2
1
2
T
2
T
,
O
=
0
0
0
0
.F and Γ(Xk) are as follows:(3)F=1T000100001T0001,ΓXk=T2ck2-T2sk2Tck-Tsk-T2sk2T2ck2Tsk-Tck,ck=cos⁡αk,sk=sin⁡αk,where T is the radar scanning period; αk is the target velocity direction.Assume that the radar is located at the origin of coordinateXOY. The measurement equation of radar which includes Doppler velocity is(4)Zk=HXk+Vk,where Zk=[rk,θk,vr,k]′ is the observation vector; rk is the range; θk is the azimuth; vr,k is Doppler velocity; Vk is a white Gaussian noise and its covariance matrix is Rk=diag⁡σr2,σθ2,σvr2; σr2,σθ2, and σvr2 are variances of the range, azimuth, and Doppler velocity, respectively; H(Xk) is the nonlinear measurement function vector:(5)HXk=xk2+yk2±arccos⁡xkxk2+yk2xkvx,k+ykvy,kxk2+yk2,where “±” is “+” if yk≥0 and “±” is “−” if yk<0.


## 2.2. The Problem Framework of Multiple Maneuvering Detection

We define the nonmaneuvering target model as a constant velocity motion model; that is,A∣A=0. Then the maneuvering target model is denoted as A∣A≠0. The traditional maneuvering detection problem is expressed as follows:(H0)
The target does not maneuver;Am∈A∣A=0 for m=1,…,k, where m is the radar sample time.(H1)
The target starts to maneuver at unknown timen; Am∈A∣A≠0 for m=n,…,k.For the above maneuvering detection problem, the traditional maneuver detector can accomplish the judgment between the maneuver and nonmaneuver, that is to say, to solve a simple binary hypothesis [23] test task. However, decision-makers often want to get the information of the target maneuver type. Then hypothesis H1 can be decomposed into multiple branches and the maneuver detection problem can be described as a new form:(H0)
The target does not maneuver;Am∈A∣A=0 for m=1,…,k.(Hi)
The target starts to maneuver at unknown timen; Am∈A∣A=Ai for m=n,…,k and i=1,…,M.A i represents the ith maneuvering model; M is the total number of maneuvering target models. In particular, the nonmaneuver A=0 is A=[0,0]′, which means that al=0 and an=0. The maneuvering A=Ai is Ai=[ali,ani]′, which can represent the different maneuvering model if [ali,ani] is the different nonzero vector.


## 3. The Normal/Tangential Acceleration

The key of reasonable multiple maneuvering detector is to obtain the information of normal/tangential acceleration. In the aspect of radar measurements, Doppler velocity contains the target maneuvering information. And [13] establishes the mapping relationship between Doppler velocity and the normal acceleration on the assumption of a constant turn rate motion. This paper renews to build a mapping relationship between Doppler velocity and the normal/tangential acceleration without the assumption.We assume that the target moves in 2D plane of coordinateXOY, as shown in Figure 1.Figure 1
Geometric relationship among variables.In Figure1, β is the angle between the target velocity direction and range extending direction; v is the size of target velocity. The basic relationship can be obtained:(6)α=β+θ,vr=vcos⁡β.By formula (6), formula (7) can be obtained:(7)vr=vcos⁡α-θ.Asβ∈-π,π,(8)α=±arccos⁡vrv+θ.Assume that the target velocity direction denoted asαk-1 is known at time k-1. After the radar measurement (rk,θk,vr,k) coming at time k, αk may be as follows:(9)αk=αk1=arccos⁡vr,kvk+θkαk2=-arccos⁡vr,kvk+θk.The variable Δαk∈-π,π is denoted as the change of target velocity direction during T in the following formula:(10)Δαkαk,αk-1=αk-αk-1αk,αk-1∈0,π∪αk,αk-1∈-π,0αk-αk-1αk∈0,π,αk-1∈-π,0∪αk-αk-1≤παk-αk-1αk∈-π,0,αk-1∈0,π∪αk-αk-1≥-παk-αk-1-2παk∈0,π,αk-1∈-π,0∪αk-αk-1≥παk-αk-1+2παk∈-π,0,αk-1∈0,π∪αk-αk-1≤-π.If the direction fromαk-1 to αk is the clockwise direction, denote Δαk<0; otherwise, denote Δαk>0. In the radar scanning period T, there are an infinite number of ways from αk-1 to αk. It is assumed that the change of target velocity direction during T is less than 2π in [13]. Under this assumption, there are four turn ways from αk-1 to αk. The maximum instantaneous turn rate of the manned maneuvering aircraft is less than π/4, such as F-22. And the scanning period of the common Doppler radar is less than 1 s under the tracking pattern. Then, we assume that the change of target velocity direction during T is less than π. Under this assumption, there are two turn ways from αk-1 to αk, which are the turn ways from αk-1 to αk(1) with Δαk(αk(1),αk-1)≤π and from αk-1 to αk(2) with Δαk(αk(2),αk-1)≤π as shown in Figure 2. Thus, the possible turn ways may only be “1” and “2,” excluding “3” and “4” in Figure 2.Figure 2
The possible turning ways.The turn rate denoted asωk meets the left hand system:(11)ωk=ΔαkT.Then, ωk is obtained by formulas (9), (10), and (11):(12)ωk=Δαkαk1,αk-1T,Δαkαk2,αk-1T.The normal acceleration can be calculated byan,k=ωkvk. If the tangential acceleration al,k=0,vk=vk-1; if al,k≠0,  vk=vk-1+al,kT. This implies an assumption that an and al are constant in the period T. As T is very short, this assumption is satisfied in practice. Then, an,k is denoted as the following formula:(13)an,t=an,t1=Δαtarccos⁡vr,t/vt-1+al,tT+θt,αt-1Tvt-1+al,tT,an,t2=Δαt-arccos⁡vr,t/vt-1+al,tT+θt,αt-1Tvt-1+al,tT.The variableal,k is unknown at time k. Generally, the tangential acceleration of aircraft is provided by its engine thrust and the instantaneous change of thrust is slow. In a short period 2T, we can assume al,k=al,k-1. Then by formula (7), we can calculate(14)al,k=al,k-1=1Tvr,k-1cos⁡αk-1-θk-1-vr,k-2cos⁡αk-2-θk-2.If the exact values ofvr,k,vr,k-1,vr,k-2,θk,θk-1,θk-2,vk-1, and αk-1 are known, we could calculate an,k and al,k by formulas (13) and (14). Due to the random noise, it is assumed that the radar measurements vr,k,vr,k-1,vr,k-2,θk,θk-1,θk-2 are random variables with the Gauss distribution. In addition to the above radar measurements, other variables are unknown which are considered as random variables obeying some probability distributions. All variables in formulas (13) and (14) are treated as random variables. Then an,k and al,k are functions of random variables. The joint PDF f(an,k,al,t) describes the target maneuvering. However, there is a high degree of nonlinear relationship between an,k(al,k) and its arguments; thus PDFs of its arguments cannot directly pass to an,k(al,k). Here, we use the statistical simulation experiment to obtain the empirical distribution of f(an,t,al,t).


## 4. Empirical Distribution and Its Approximated PDF

The vector(an,al) has uncountable combinations that correspond to maneuvering models. Obviously, it is convenient for us to study maneuvering detection by determining some typical maneuvering models through some combinations. If we select M typical maneuvering models, we should build M joint PDFs denoted as fi(an,al)i=1M whose empirical distribution functions are obtained through statistical simulation experiments. In the experiment, the standard 2D curve motion model is adopted to obtain the true trajectory of the target as the following formula:(15)x˙=vcos⁡α+wxy˙=vsin⁡α+wyv˙=al+wvα˙=anv+wα,where wx,wy,wv, and wα are the independent Gauss white noise; σwx,σwy,σwv, and σwα are the corresponding standard deviation. In the experiment, multiple real trajectories given an and al could be obtained by formula (15) and by initializing different x,y,v,α. The measurement vector (r,θ,vr) from Doppler radar is obtained by the following formula:(16)rθvr=x2+y2±arccos⁡xx2+y2vcos⁡α∓arccos⁡xx2+y2+wrwθwvr,where “±” is “+” if y≥0 and “±” is “−” if y<0; “∓” is opposite to “±”; wr,wθ, and wvr are the independent Gauss white noise; σwr,σwθ, and σwvr are the corresponding standard deviation. The radar measurement sequence can be obtained through inputting a real trajectory into formula (16). In addition to radar measurement, the value of vk-1,αk-1 and αk-2 should be known to calculate an,k and al,k by formulas (13) and (14). We use a filter to estimate them. In order to solve the problem of the high nonlinearity of measurement equation (4), [14] proposes a measurement conversion algorithm that converts the measurement in Polar coordinate system to the pseudo-Gauss noise measurement in Cartesian coordinates system and then used KF for state estimation. For the better application in this paper, the filter has been modified as shown in Appendix A. X^k-1=[x^k-1,v^x,k-1,y^k-1,v^y,k-1]′ can be obtained by filtering and v^k-1,α^k-1, and α^k-2 can be calculated by the following formulas:(17)v^k-1=v^x,k-12+v^y,k-12,α^k-1=arccos⁡v^x,k-1v^x,k-12+v^y,k-12,v^y,k-1≥0,-arccos⁡v^x,k-1v^x,k-12+v^y,k-12,v^y,k-1<0.The variablesa^n,k and a^l,k can be calculated with those values by formulas (13) and (14). For the variable a^n,k, there exist two possible values, which are a^n,k(1) and a^n,k(2). In [13, 14], the statistic cmin=min⁡an,k(1),an,k(2) is employed. However, when |α-θ| is small, a^n(1) and a^n(2) are similar. They will be taken approximately with the same probability as an estimate of an. Then cmin will introduce a statistical error, since the smaller value between a^n(1) and a^n(2) is selected. A new statistic cR is proposed in this paper based on the statistic cmin; that is,(18)cR=min⁡a^n,k1,a^n,k2,arccos⁡vr,kv^k-1∣k-1≥η,rand⁡a^n,k1,a^n,k2,arccos⁡vr,kv^k-1∣k-1<η,where rand⁡a^n,k(1),a^n,k(2) is to select a value between a^n,t(1) and a^n,t(2) randomly; η≈π/18 is the experiential threshold value.We apply the Monte Carlo simulation method under every initial state given a maneuvering model(an,al). In each simulation, the target firstly moves at constant velocity in accordance with the initial state and then moves in accordance with the maneuvering model (an,al) after time tc. Record cR and al at time tc and the empirical distribution denoted as g(cR,al) could be acquired. We denote the approximate joint PDF f′(an,al) of f(an,al) which can be obtained by fitting g(cR,al). This is a fitting problem of multidimensional nonlinear PDF. GMM can be used to approximate an arbitrary joint PDF. Therefore, GMM is used to approximate g(cR,al) and EM algorithm is applied to learn the parameters of GMM as shown in Appendix B.


## 5. M-CUSUM Detector

The purpose of this paper is to design a multiple maneuvering model detector that can deal with the compound hypothesis test problem described in Section2.2. The maneuvering detector can be divided into two groups: one is a batch sliding window detector and the other is a sequential detector. Since the radar measurement is coming in order, the sequential detection is more suitable for the maneuvering detection. CUSUM as a kind of sequential detector can ensure the fastest detection at a certain decision error rate under the condition that the input vector A is known. For a simple binary hypothesis testing of the conventional maneuvering detection described in Section 2.2, the standard CUSUM detector is(19)Lk=max⁡Lk-1+log⁡fZ~k∣H1,Z1:k-1fZ~k∣H0,Z1:k-1,0,L0=0.Decision rules are as follows:(1)
AcceptH1, if Lk≥κ, and the maneuvering time m^=min⁡k:Lk≥λ.(2)
Continue to detect andk=k+1, if Lk<κ.Z 1 : k ≜ ( Z 1 , … , Z k ); κ is the detection threshold (decision error rate); Z~k is the measurement prediction residual at time k. The Multiple Parallel CUSUM detector, called M-CUSUM detector, is proposed in the paper for the detection and recognition of M maneuver patterns. The specific form is(20)Lki=max⁡Lk-1i+log⁡fZ~k∣Hi,Z1:k-1fZ~k∣H0,Z1:k-1,0,L0i=0i=1,2,…,M,and M detection thresholds are denoted as κ1,κ2,…,κM.Decision rules are follows:(1)
AcceptHi, if Lki≥κi and Lkj<κj(∀j≠i), and the maneuvering time m^=min⁡k:Lki≥κi.(2)
Continue to detect andk=k+1, if Lki<κi(∀i).By formula (20), it is the key of the detection and recognition of M maneuver patterns to calculate f(Z~k∣H0,Z1:k-1) and f(Z~k∣Hi,Z1:k-1). Assume that f(Z~k∣H0,Z1:k-1) is a zero mean Gaussian distribution in general:(21)fZ~k∣H0,Zk-1=NZ~k;0,Sk,where Sk is the measurement prediction residual error covariance at time k. f(Z~k∣Hi,Zk-1) is the measurement prediction residual likelihood function and is the Gaussian distribution with the mean H(Γ(X^k-1)Ai):(22)fZ~k∣Hi,Zk-1=NZ~k;HΓX^k-1Ai,Sk.The measurement conversion algorithm (see AppendixA) can approximately represent H(Γ(X^k-1)Ai) as a linear form HcΓ(X^k-1)Ai. Although the real value of Ai is unknown, the GMM approximation joint PDF f′(Ai) of Ai can be obtained by the statistical simulation experiment. Then f′(Ai) can be used as a priori information of target maneuvering and we adopt the marginal method to calculate f(Z~k∣Hi,Zk-1):(23)fZ~k∣Hi,Zk-1=ENZ~k;HcΓX^k-1Ai,Skf′Ai=∫NZ~k;HcΓX^k-1Ai,Skf′AidAi=∑j=kKλjiNZ~k;HcΓX^k-1μji,Sk+HcΓX^k-1ΣjiHcΓX^k-1′,f′Ai=∑j=1Kλji2π-1Σji-1/2exp⁡-12Ai-μjiΣji-1Ai-μji′,where f′(Ai) is the matrix form of GMM; μji is the mean vector; Σji is the covariance matrix; λji is the weight value.


## 6. Simulation Experiment

The main purpose of the simulation experiment is to verify M-CUSUM detector and compare it with other maneuvering detectors. The SMT algorithms using these detectors are also compared with state-of-the-art IMM algorithm and current statistical model adaptive filter (CSAF) algorithm [24, 25].


### 6.1. Simulation Settings

#### 6.1.1. Scenario Settings

The real trajectory and its radar measurement are generated by formulas (15) and (16), respectively. The standard deviations in formula (15) are set as σwx=2 m, σwy=2 m, σwv=0.1 m/s, and σwα=0.005 rad. The standard deviations in formula (16) are set as σwr=50 m, σwθ=0.01 rad, and σwvr=1 m/s. The initial state of target (x0,y0,v0,ϕ0) is set as (40 km, 60 km, 300 m/s, π/6). Doppler radar is located at the origin of coordinate XOY and its scanning period is set as T=0.5 s. Design five scenarios for performance evaluation. The target moves at a constant velocity (that is the nonmaneuver model A0=(0,0)) during 0 s~20 s in each scenario. The target starts to maneuver at tc=21 s and keeps 20 s. In five scenarios, the maneuver models of target are designed as A1=(0,10),A2=(0,-20),A3=(10,0),A4=(10,10), and A5=(5,20), as shown in Figure 3.Figure 3
The real trajectories of the five scenarios.


#### 6.1.2. Tracking Algorithm Settings

The filters are divided into two types according to the different state equations adopted. The first filter is one whose state equation [26] does not include the input vector Ut:(24)Xt+1=FXt+wt.Formula (24) is often taken as the state equation by the standard tracking algorithm where the common SMT algorithms could be compared. The state transform matrix F in formula (24) is the constant velocity model (CV) and we design two level noise covariance matrices Qlow and Qhigh, respectively, as the nonmaneuvering model and the maneuvering model. The SMT algorithm switches the motion model between Qlow and Qhigh based on the maneuvering detector. The model set of IMM algorithm consists of two CV motion models with Qlow and Qhigh, respectively. In the simulation, we especially set Qlow(q=10) and Qhigh(q=1000).The second filter is one whose state equation includes the input vectorUt, as described in the following formula:(25)Xk+1=FXk+Uk+Wk.For the STM algorithm with M-CUSUM detector and IMM algorithm, the input vectorUk is denoted as Γ(Xk)Ak. Then the state equation (25) is the same as the state equation (1). Given F and Qk, the different Ak represents the different maneuvering model. If Ak=A0, the state equation (25) describes the nonmaneuver model, that is, CV; if At∈{A1,A2,A3,A4,A5}, the state equation (25) describes the maneuver models corresponding to the five scenarios. In the simulation, the switching models of M-CUSUM detector are denoted as {Ai}i=15 and the model set of IMM algorithm is denoted as {Ai}i=05. The state transform matrix F in formula (25) also is CV and the noise covariance matrixes are set as QA0(q=10) and QAi(q=600),  i=1,2,3,4,5. Reference [25] proposes a direct method of estimating the tangential and normal accelerations of maneuvering targets in three-dimensional space. The method takes the modified Rayleigh distribution as the statistical model of the tangential and normal accelerations that directly are introduced into the state vector X and updates the state noise covariance matrix Qk online to adapt to the current maneuver of target. The method is called the current statistical model and adaptive Kalman filter (CSAF) algorithm in the paper. In order to compare with M-CUSUM, CSAF algorithm should be deduced in two-dimensional space and modified with the measurement from Doppler radar. The discrete Uk,Fk, and Qk in the state equation (25) can be deduced through (C.5) given in Appendix C for the current statistical model, and the adaptive Kalman filter algorithm can refer to [25] and Appendix A.Two types of filter are set up to provide a comparative platform for different tracking algorithms. The maneuver detector based on the binary hypothesis can only be used in the detection of nonmaneuver model and maneuver model but cannot give the specific information of maneuver model. Thus it can only be used in the first filter but not in the second filter. M-CUSUM is a multiple maneuver model detector, which can realize the hard decision among multiple maneuver models, and the IMM algorithm can realize the soft decision among multiple maneuver models. Both M-CUSUM and IMM can be used in the tracking problems with the first filter and second filter. CSAF algorithm can realize the adaptive modification ofUk and Qk which can be compared with M-CUSUM and IMM with the second filter only. For two types of filters, the measurement conversion algorithm in Appendix A is used. The difference is that the different state equation is used.The average onset detection delay time(n^-n) is used to measure the performance of the maneuver detector. The root mean square error (RMSE) of state estimation is used to measure the performance of the tracking algorithm:(26)RMSEx^k=1N∑i=1Nxk-x^k′xk-x^k,where N is Monte Carlo simulation time; xk is the real state vector; x^k is the state estimation vector.


### 6.2. Simulation Results

#### 6.2.1. Empirical Distribution and the Approximated GMM

The empirical distributions{gi(cR,al)}i=15 corresponding to the maneuver models {Ai}i=15 are acquired by the statistical experiment. We should assign many values to the initial state (x0,y0,v0,ϕ0) for experiments, because f(an,al) is affected by the velocity magnitude, velocity direction, and position. In each dimension, a number of typical levels are selected for uniform experiments, and the specific values of parameters are shown in Table 1.Table 1
Initial state parameter settings.
Parameters Maximum Minimum Grading x 0 80 km 40 km 20 km y 0 80 km 40 km 20 km v 0 400 m/s 200 m/s 100 m/s α 0 π - 5 π / 6 π / 6In Table1, there are 162 initial state values. We execute 30 Monte Carlo simulations for each maneuvering model in {Ai}i=15 and each initial state value in Table 1. Record the value of (cR,al) at time tc=21 s, and we can get 4860 groups of data for each maneuver model. Then the empirical distribution function gi(cR,al) corresponding to one maneuvering model Ai can be obtained. Two-dimensional GMM by EM algorithm in Appendix B is used to fit the empirical distribution, and then the approximate PDF fi′(an,al) can be acquired. Specifically, three two-dimensional Gauss distribution functions are used to fit the distributions. The scatter diagrams and the approximate GMM diagrams of each maneuver model in {Ai}i=15 are shown in Figure 4. And the parameters of GMM are recorded in Table 2. From Figure 4, different maneuvering models correspond to different distribution patterns, which proves that the probability distribution of variables in formulas (13) and (14) cannot be directly transmitted to an or al. It is also claimed that the normal acceleration is related to the tangential acceleration without the assumption of a constant turn rate motion. Thus, it will introduce an error if the joint PDF f(an,al) is decomposed into the product of marginal PDF f(an)f(al).Table 2
The parameters of GMM.
Maneuvering models μ Σ λ A 1 0.0850 7.1198 1.8693 −1.2086 0.1494 −1.2086 147.0940 0.2879 0.0207 75.3546 1.2346 0.3808 1.2346 1.2266 0.0057 8.8141 3.9931 0.1218 0.4699 0.1218 4.1001  A 2 0.2328 19.6408 91.6298 −2.8712 0.4534 −2.8712 1.3740 −0.1778 13.0432 1.6643 0.8671 0.1157 0.8671 154.5459 0.0761 17.5508 3.3611 0.0244 0.4310 0.0244 7.9094  A 3 8.4827 0.6988 5.1910 −0.0291 0.5962 −0.0291 16.2650 4.5230 0.6965 12.0374 0.0475 0.3000 0.0475 8.4115 9.8643 0.6820 1.9299 0.0232 0.1038 0.0232 83.8308  A 4 10.2197 4.2566 3.3473 0.6226 0.1979 0.6226 119.7231 11.5421 7.1672 4.3900 −7.5418 0.4087 −7.5418 25.3830 10.5928 9.4848 10.7989 −9.5082 0.3934 −9.5082 14.3750  A 5 3.9518 0.2270 59.0632 −1.2884 0.4022 −1.2884 2.8622 −19.1749 0.8473 36.5717 −13.2782 0.0195 −13.2782 75.4061 9.8141 18.4943 24.5112 0.9392 0.5783 0.9392 9.8745Figure 4
Scatter and GMM diagrams corresponding to five maneuver models.
(a)
Scatter diagram ofA1 (b)
f 1 ′ ( a n , a t ) (c)
Scatter diagram ofA2 (d)
f 2 ′ ( a n , a t ) (e)
Scatter diagram ofA3 (f)
f 3 ′ ( a n , a t ) (g)
Scatter diagram ofA4 (h)
f 4 ′ ( a n , a t ) (i)
Scatter diagram ofA5 (j)
f 5 ′ ( a n , a t )


#### 6.2.2. Detection Performance Comparison

In the simulation experiment, the M-CUSUM detector is compared with the measurement residual (MR), Input Estimate (IE), CUSUM, andcmin and cmin2 detectors. In the calculation of the average delay time (n^-n), the false alarm probability is set as pf=0.01 and the miss probability is set as pm=0.1. For MR and IE detectors, the window size is set as 5. For cmin detector, the threshold value of fading memory average (FMA) is set as 5000. For cmin2 and CUSUM detectors, the threshold value can be calculated by κ=log⁡((1-pf)/pm) [27], that is, κ=1.9542. For comparison, set the threshold values of M-CUSUM detector as κ1=κ2=⋯=κ5=κ.Table3 summarizes detection delay time of all above algorithms and the detectors have longer detection delay time for relatively small maneuver scenarios as A1 and A3. cmin, cmin2, and M-CUSUM detectors taking advantage of Doppler velocity measurement are more sensitive to the target maneuver compared with MR, IE, and CUSUM detectors. But when maneuver models (A4 and A5) have both the normal acceleration and the tangential acceleration, M-CUSUM detector has a shorter detection delay time. It is the main reason why M-CUSUM detector uses a more accurate prior knowledge of the normal/tangential acceleration.Table 3
The average delay time of maneuver detection (s).
n ^ - n MR IE CUSUM c m i n c m i n 2 M-CUSUM A 1 5.32 4.72 4.26 2.33 2.38 2.37 A 2 4.08 3.58 3.35 2.65 2.71 2.62 A 3 5.13 4.92 4.20 4.09 3.73 3.01 A 4 4.03 3.64 4.03 3.19 3.03 1.46 A 5 3.49 3.91 3.22 2.29 2.15 1.28


#### 6.2.3. Estimation Performance Comparison

In the first filter, the estimation performances of SMT algorithms with MR, IE, CUSUM,cmin, cmin2, and M-CUSUM detectors and IMM algorithm are compared. These tracking algorithms all use Doppler velocity measurement in the state estimation. The initial state vector X0 and its standard deviation are set as(27)X0=40km,245m/s,60km,150m/s,σx0=30m,σvx=5m/s,σy0=30m,σvy=5m/s.The initial estimation error covariance matrix denoted asP0, the initial model probability denoted as Pm0, and model transition probability of IMM algorithm denoted as Pm are set as(28)P0=1030000250000103000025,Pm0=0.90.1,Pm=0.90.10.10.9.We execute 100 Monte Carlo simulations for the five scenarios in the first filter. Table4 records the position and velocity of the average RMSE with different tracking algorithms for each scene. The velocity and position of RMSE of different tracking algorithms for each scene are shown at the sampling time in Figure 5, respectively. For the maneuvering target, the sooner it switches the motion model, the more conducive it is to the state estimation for the SMT algorithm. It is concluded that the smaller delay time leads to the lower RMSE, which means the better tracking accuracy according to Tables 3 and 4. Figure 5 shows that M-CUSUM has the smallest RMSE in target maneuver phase compared with MR, IE, CUSUM, cmin, cmin2, and CUSUM. The difference among the average position RMSE of each detector is small. However, the large difference is reflected in the average velocity RMSE in Table 4. When the target maneuvering is stronger, the average velocity RMSE is higher. This shows that the degree of mismatch between the real maneuvering model and the filter motion model is more directly reflected in the velocity than in the position, which is the reason why the efficiency of the detectors based on Doppler velocity measurement is higher than that based on the position measurement. The average position and velocity RMSE of IMM algorithm are higher than the SMT algorithms. The SMT algorithms using the matched target motion model Qlow have lower RMSE during the nonmaneuver phase 0~20 s as shown in Figure 5. As the model noise of IMM algorithm is larger than Qlow, it has a higher RMSE during the nonmaneuver phase 0~20 s. However, as IMM algorithm could adaptively adjust the weights of two models Qlow and Qhigh to match the target maneuver model, it has a lower RMSE during the maneuver phase 21 s~40 s. The difference in average RMSE between IMM algorithm and SMT algorithms mainly comes from the accumulation of nonmaneuver phase as shown in Table 4. The tracking algorithms based on M-CUSUM, cmin, and cmin2 detectors have the comparable estimation performance. Moreover, since M-CUSUM detector has the normal/tangential acceleration in the target maneuver, the corresponding algorithm obtains the better estimation accuracy. The CV motion model is more mismatched to the maneuver model with larger normal accelerations as A2 and A5 compared with other maneuver models as A1 and A3 as shown in Table 4 and Figure 5.Table 4
The average position RMSE (m) and average velocity RMSE (m/s).
Average RMSE MR IE CUSUM c m i n c m i n 2 M-CUSUM IMM A 1 Position 35.47 34.41 33.47 31.52 32.46 32.94 41.39 Velocity 184.01 180.97 170.37 139.58 138.84 129.35 190.37 A 2 Position 44.32 44.18 44.46 41.05 39.51 39.87 49.94 Velocity 236.47 241.97 238.69 193.63 180.59 182.88 207.73 A 3 Position 20.99 20.41 19.47 18.98 19.32 18.87 27.38 Velocity 116.54 117.70 109.31 110.07 113.13 115.36 166.03 A 4 Position 31.80 30.25 30.91 29.64 29.84 28.53 38.06 Velocity 134.03 134.72 126.22 129.89 131.33 126.38 174.31 A 5 Position 60.59 58.89 58.79 58.37 59.01 57.28 71.83 Velocity 270.96 246.01 248.49 220.72 228.27 213.27 272.90Figure 5
RMSE of the first filter.
(a)
Position RMSE ofA1 (b)
Velocity RMSE ofA1 (c)
Position RMSE ofA2 (d)
Velocity RMSE ofA2 (e)
Position RMSE ofA3 (f)
Velocity RMSE ofA3 (g)
Position RMSE ofA4 (h)
Velocity RMSE ofA4  (i)
Position RMSE ofA5 (j)
Velocity RMSE ofA5In the second filter, IMM algorithm adopts six models:{Ai}i=05. The initial state vector X0 and the initial estimation error covariance matrix P0 are the same as the first filter. The initial model probability denoted as Pm0 and model transition probability denoted as Pm of IMM algorithm are set as(29)Pm0=0.90.020.020.020.020.02,Pm=0.90.020.020.020.020.020.020.90.020.020.020.020.020.020.90.020.020.020.020.020.020.90.020.020.020.020.020.020.90.020.020.020.020.020.020.9.The parameters of CSAF algorithm are set as(30)almax=30m/s2,a-lmax=-30m/s2,anmax=70m/s2,a-nmax=-70m/s2,αl=0.1,αn=0.1,where almax and a-lmax are the positive and negative maximum of tangential acceleration; anmax and a-nmax are the positive and negative maximum of normal acceleration; αl and αn are the maneuvering frequency of the tangential and normal acceleration, respectively. The initial state vector X0 and its standard deviation of CSAF are set as(31)X0=40km,60km,245m/s,150m/s,2m/s2,2m/s2,σx0=30m,σy0=30m,σvx=5m/s,σvy=5m/s,σal=0.1m/s2,σan=0.2m/s2.The initial estimation error covariance matrix of CSAF denoted as P0 is set as(32)P0=103000000103000000250000002500000050000005.If the SMT algorithm with the M-CUSUM detector detects the target maneuver, it will adopt the maneuver modelAi,i=1,2,3,4,5, which firstly reaches the detection threshold to replace the term At=A0 in state equation (1). From the perspective of position and velocity, IMM algorithm and SMT algorithm in the second filter have higher estimation accuracy than in the first filter according to Figure 6. The more important reason is that the motion models {Ai}i=05 of state equation (1) are more suitable to the maneuver states of the target in the five scenes than the motion models Qlow and Qhigh. The SMT algorithm with M-CUSUM detector is significantly better than the IMM algorithm in two strong maneuver scenarios as A2 and A5 as shown in Figure 6.Figure 6
RMSE of the second filter.
(a)
Position RMSE ofA1 (b)
Velocity RMSE ofA1 (c)
Position RMSE ofA2 (d)
Velocity RMSE ofA2 (e)
Position RMSE ofA3 (f)
Velocity RMSE ofA3 (g)
Position RMSE ofA4 (h)
Velocity RMSE ofA4  (i)
Position RMSE ofA5 (j)
Velocity RMSE ofA5First of all, this is because the IMM algorithm has the competition problem among multiple motion models, which can cause the probability weights to be dispersed by the unmatched models, and the system error is introduced into the state fusion.Secondly, M-CUSUM detector can achieve higher recognition rates in strong maneuver scenes compared with other scenes, as shown in Table5. The recognition rate of A2 and A5 maneuver models based on M-CUSUM detector is 93% and 95%, respectively.Table 5
Recognition results of 100 Monte Carlo simulations based on M-CUSUM.
Recognition results A 1 A 2 A 3 A 4 A 5 A 1 89 3 3 5 0 A 2 2 93 1 1 3 A 3 4 3 84 7 2 A 4 3 1 3 91 2 A 5 2 2 0 1 95Then, the target state correction effect of the input vector is obvious in the strong maneuver scene, because the magnitude of input is higher than the noise level. However, as the magnitudes of input vector and noise are similar for weak maneuvering scenes, the target state correction of input vector will be interfered with the noise, which can lead to the performance degradation of the SMT algorithm with M-CUSUM detector.The variablesUk and Qk are modified online in CSAF algorithm as the current statistical model is used. Then the performance of CSAF is similar to IMM. CSAF algorithm has bigger position and velocity RMSE than SMT algorithm with M-CUSUM detector in the nonmaneuver stage, as shown in Figure 6. This is because the normal acceleration estimation obtained by the filter directly is more sensitive to the position measurement error than the target maneuvering in the nonmaneuvering stage. The tangential and normal accelerations in the iterative process are mutually independent with the modified Rayleigh distribution in CSAF algorithm. However, it has been proven in Section 6.2.1 that they are correlative. The error due to the independence assumption is constantly accumulated and then CSAF algorithm has significant tracking error around 30 s~40 s of the maneuvering scenarios A4 and A5 that have both the tangential and normal acceleration, as shown in Figures 6(g)–6(j). However, CSAF algorithm has the comparability with STM algorithm and IMM algorithm in the maneuvering scenarios A1, A2, and A3 that have only the normal acceleration or the tangential acceleration.


## 6.1. Simulation Settings

### 6.1.1. Scenario Settings

The real trajectory and its radar measurement are generated by formulas (15) and (16), respectively. The standard deviations in formula (15) are set as σwx=2 m, σwy=2 m, σwv=0.1 m/s, and σwα=0.005 rad. The standard deviations in formula (16) are set as σwr=50 m, σwθ=0.01 rad, and σwvr=1 m/s. The initial state of target (x0,y0,v0,ϕ0) is set as (40 km, 60 km, 300 m/s, π/6). Doppler radar is located at the origin of coordinate XOY and its scanning period is set as T=0.5 s. Design five scenarios for performance evaluation. The target moves at a constant velocity (that is the nonmaneuver model A0=(0,0)) during 0 s~20 s in each scenario. The target starts to maneuver at tc=21 s and keeps 20 s. In five scenarios, the maneuver models of target are designed as A1=(0,10),A2=(0,-20),A3=(10,0),A4=(10,10), and A5=(5,20), as shown in Figure 3.Figure 3
The real trajectories of the five scenarios.


### 6.1.2. Tracking Algorithm Settings

The filters are divided into two types according to the different state equations adopted. The first filter is one whose state equation [26] does not include the input vector Ut:(24)Xt+1=FXt+wt.Formula (24) is often taken as the state equation by the standard tracking algorithm where the common SMT algorithms could be compared. The state transform matrix F in formula (24) is the constant velocity model (CV) and we design two level noise covariance matrices Qlow and Qhigh, respectively, as the nonmaneuvering model and the maneuvering model. The SMT algorithm switches the motion model between Qlow and Qhigh based on the maneuvering detector. The model set of IMM algorithm consists of two CV motion models with Qlow and Qhigh, respectively. In the simulation, we especially set Qlow(q=10) and Qhigh(q=1000).The second filter is one whose state equation includes the input vectorUt, as described in the following formula:(25)Xk+1=FXk+Uk+Wk.For the STM algorithm with M-CUSUM detector and IMM algorithm, the input vectorUk is denoted as Γ(Xk)Ak. Then the state equation (25) is the same as the state equation (1). Given F and Qk, the different Ak represents the different maneuvering model. If Ak=A0, the state equation (25) describes the nonmaneuver model, that is, CV; if At∈{A1,A2,A3,A4,A5}, the state equation (25) describes the maneuver models corresponding to the five scenarios. In the simulation, the switching models of M-CUSUM detector are denoted as {Ai}i=15 and the model set of IMM algorithm is denoted as {Ai}i=05. The state transform matrix F in formula (25) also is CV and the noise covariance matrixes are set as QA0(q=10) and QAi(q=600),  i=1,2,3,4,5. Reference [25] proposes a direct method of estimating the tangential and normal accelerations of maneuvering targets in three-dimensional space. The method takes the modified Rayleigh distribution as the statistical model of the tangential and normal accelerations that directly are introduced into the state vector X and updates the state noise covariance matrix Qk online to adapt to the current maneuver of target. The method is called the current statistical model and adaptive Kalman filter (CSAF) algorithm in the paper. In order to compare with M-CUSUM, CSAF algorithm should be deduced in two-dimensional space and modified with the measurement from Doppler radar. The discrete Uk,Fk, and Qk in the state equation (25) can be deduced through (C.5) given in Appendix C for the current statistical model, and the adaptive Kalman filter algorithm can refer to [25] and Appendix A.Two types of filter are set up to provide a comparative platform for different tracking algorithms. The maneuver detector based on the binary hypothesis can only be used in the detection of nonmaneuver model and maneuver model but cannot give the specific information of maneuver model. Thus it can only be used in the first filter but not in the second filter. M-CUSUM is a multiple maneuver model detector, which can realize the hard decision among multiple maneuver models, and the IMM algorithm can realize the soft decision among multiple maneuver models. Both M-CUSUM and IMM can be used in the tracking problems with the first filter and second filter. CSAF algorithm can realize the adaptive modification ofUk and Qk which can be compared with M-CUSUM and IMM with the second filter only. For two types of filters, the measurement conversion algorithm in Appendix A is used. The difference is that the different state equation is used.The average onset detection delay time(n^-n) is used to measure the performance of the maneuver detector. The root mean square error (RMSE) of state estimation is used to measure the performance of the tracking algorithm:(26)RMSEx^k=1N∑i=1Nxk-x^k′xk-x^k,where N is Monte Carlo simulation time; xk is the real state vector; x^k is the state estimation vector.


## 6.1.1. Scenario Settings

The real trajectory and its radar measurement are generated by formulas (15) and (16), respectively. The standard deviations in formula (15) are set as σwx=2 m, σwy=2 m, σwv=0.1 m/s, and σwα=0.005 rad. The standard deviations in formula (16) are set as σwr=50 m, σwθ=0.01 rad, and σwvr=1 m/s. The initial state of target (x0,y0,v0,ϕ0) is set as (40 km, 60 km, 300 m/s, π/6). Doppler radar is located at the origin of coordinate XOY and its scanning period is set as T=0.5 s. Design five scenarios for performance evaluation. The target moves at a constant velocity (that is the nonmaneuver model A0=(0,0)) during 0 s~20 s in each scenario. The target starts to maneuver at tc=21 s and keeps 20 s. In five scenarios, the maneuver models of target are designed as A1=(0,10),A2=(0,-20),A3=(10,0),A4=(10,10), and A5=(5,20), as shown in Figure 3.Figure 3
The real trajectories of the five scenarios.


## 6.1.2. Tracking Algorithm Settings

The filters are divided into two types according to the different state equations adopted. The first filter is one whose state equation [26] does not include the input vector Ut:(24)Xt+1=FXt+wt.Formula (24) is often taken as the state equation by the standard tracking algorithm where the common SMT algorithms could be compared. The state transform matrix F in formula (24) is the constant velocity model (CV) and we design two level noise covariance matrices Qlow and Qhigh, respectively, as the nonmaneuvering model and the maneuvering model. The SMT algorithm switches the motion model between Qlow and Qhigh based on the maneuvering detector. The model set of IMM algorithm consists of two CV motion models with Qlow and Qhigh, respectively. In the simulation, we especially set Qlow(q=10) and Qhigh(q=1000).The second filter is one whose state equation includes the input vectorUt, as described in the following formula:(25)Xk+1=FXk+Uk+Wk.For the STM algorithm with M-CUSUM detector and IMM algorithm, the input vectorUk is denoted as Γ(Xk)Ak. Then the state equation (25) is the same as the state equation (1). Given F and Qk, the different Ak represents the different maneuvering model. If Ak=A0, the state equation (25) describes the nonmaneuver model, that is, CV; if At∈{A1,A2,A3,A4,A5}, the state equation (25) describes the maneuver models corresponding to the five scenarios. In the simulation, the switching models of M-CUSUM detector are denoted as {Ai}i=15 and the model set of IMM algorithm is denoted as {Ai}i=05. The state transform matrix F in formula (25) also is CV and the noise covariance matrixes are set as QA0(q=10) and QAi(q=600),  i=1,2,3,4,5. Reference [25] proposes a direct method of estimating the tangential and normal accelerations of maneuvering targets in three-dimensional space. The method takes the modified Rayleigh distribution as the statistical model of the tangential and normal accelerations that directly are introduced into the state vector X and updates the state noise covariance matrix Qk online to adapt to the current maneuver of target. The method is called the current statistical model and adaptive Kalman filter (CSAF) algorithm in the paper. In order to compare with M-CUSUM, CSAF algorithm should be deduced in two-dimensional space and modified with the measurement from Doppler radar. The discrete Uk,Fk, and Qk in the state equation (25) can be deduced through (C.5) given in Appendix C for the current statistical model, and the adaptive Kalman filter algorithm can refer to [25] and Appendix A.Two types of filter are set up to provide a comparative platform for different tracking algorithms. The maneuver detector based on the binary hypothesis can only be used in the detection of nonmaneuver model and maneuver model but cannot give the specific information of maneuver model. Thus it can only be used in the first filter but not in the second filter. M-CUSUM is a multiple maneuver model detector, which can realize the hard decision among multiple maneuver models, and the IMM algorithm can realize the soft decision among multiple maneuver models. Both M-CUSUM and IMM can be used in the tracking problems with the first filter and second filter. CSAF algorithm can realize the adaptive modification ofUk and Qk which can be compared with M-CUSUM and IMM with the second filter only. For two types of filters, the measurement conversion algorithm in Appendix A is used. The difference is that the different state equation is used.The average onset detection delay time(n^-n) is used to measure the performance of the maneuver detector. The root mean square error (RMSE) of state estimation is used to measure the performance of the tracking algorithm:(26)RMSEx^k=1N∑i=1Nxk-x^k′xk-x^k,where N is Monte Carlo simulation time; xk is the real state vector; x^k is the state estimation vector.


## 6.2. Simulation Results

### 6.2.1. Empirical Distribution and the Approximated GMM

The empirical distributions{gi(cR,al)}i=15 corresponding to the maneuver models {Ai}i=15 are acquired by the statistical experiment. We should assign many values to the initial state (x0,y0,v0,ϕ0) for experiments, because f(an,al) is affected by the velocity magnitude, velocity direction, and position. In each dimension, a number of typical levels are selected for uniform experiments, and the specific values of parameters are shown in Table 1.Table 1
Initial state parameter settings.
Parameters Maximum Minimum Grading x 0 80 km 40 km 20 km y 0 80 km 40 km 20 km v 0 400 m/s 200 m/s 100 m/s α 0 π - 5 π / 6 π / 6In Table1, there are 162 initial state values. We execute 30 Monte Carlo simulations for each maneuvering model in {Ai}i=15 and each initial state value in Table 1. Record the value of (cR,al) at time tc=21 s, and we can get 4860 groups of data for each maneuver model. Then the empirical distribution function gi(cR,al) corresponding to one maneuvering model Ai can be obtained. Two-dimensional GMM by EM algorithm in Appendix B is used to fit the empirical distribution, and then the approximate PDF fi′(an,al) can be acquired. Specifically, three two-dimensional Gauss distribution functions are used to fit the distributions. The scatter diagrams and the approximate GMM diagrams of each maneuver model in {Ai}i=15 are shown in Figure 4. And the parameters of GMM are recorded in Table 2. From Figure 4, different maneuvering models correspond to different distribution patterns, which proves that the probability distribution of variables in formulas (13) and (14) cannot be directly transmitted to an or al. It is also claimed that the normal acceleration is related to the tangential acceleration without the assumption of a constant turn rate motion. Thus, it will introduce an error if the joint PDF f(an,al) is decomposed into the product of marginal PDF f(an)f(al).Table 2
The parameters of GMM.
Maneuvering models μ Σ λ A 1 0.0850 7.1198 1.8693 −1.2086 0.1494 −1.2086 147.0940 0.2879 0.0207 75.3546 1.2346 0.3808 1.2346 1.2266 0.0057 8.8141 3.9931 0.1218 0.4699 0.1218 4.1001  A 2 0.2328 19.6408 91.6298 −2.8712 0.4534 −2.8712 1.3740 −0.1778 13.0432 1.6643 0.8671 0.1157 0.8671 154.5459 0.0761 17.5508 3.3611 0.0244 0.4310 0.0244 7.9094  A 3 8.4827 0.6988 5.1910 −0.0291 0.5962 −0.0291 16.2650 4.5230 0.6965 12.0374 0.0475 0.3000 0.0475 8.4115 9.8643 0.6820 1.9299 0.0232 0.1038 0.0232 83.8308  A 4 10.2197 4.2566 3.3473 0.6226 0.1979 0.6226 119.7231 11.5421 7.1672 4.3900 −7.5418 0.4087 −7.5418 25.3830 10.5928 9.4848 10.7989 −9.5082 0.3934 −9.5082 14.3750  A 5 3.9518 0.2270 59.0632 −1.2884 0.4022 −1.2884 2.8622 −19.1749 0.8473 36.5717 −13.2782 0.0195 −13.2782 75.4061 9.8141 18.4943 24.5112 0.9392 0.5783 0.9392 9.8745Figure 4
Scatter and GMM diagrams corresponding to five maneuver models.
(a)
Scatter diagram ofA1 (b)
f 1 ′ ( a n , a t ) (c)
Scatter diagram ofA2 (d)
f 2 ′ ( a n , a t ) (e)
Scatter diagram ofA3 (f)
f 3 ′ ( a n , a t ) (g)
Scatter diagram ofA4 (h)
f 4 ′ ( a n , a t ) (i)
Scatter diagram ofA5 (j)
f 5 ′ ( a n , a t )


### 6.2.2. Detection Performance Comparison

In the simulation experiment, the M-CUSUM detector is compared with the measurement residual (MR), Input Estimate (IE), CUSUM, andcmin and cmin2 detectors. In the calculation of the average delay time (n^-n), the false alarm probability is set as pf=0.01 and the miss probability is set as pm=0.1. For MR and IE detectors, the window size is set as 5. For cmin detector, the threshold value of fading memory average (FMA) is set as 5000. For cmin2 and CUSUM detectors, the threshold value can be calculated by κ=log⁡((1-pf)/pm) [27], that is, κ=1.9542. For comparison, set the threshold values of M-CUSUM detector as κ1=κ2=⋯=κ5=κ.Table3 summarizes detection delay time of all above algorithms and the detectors have longer detection delay time for relatively small maneuver scenarios as A1 and A3. cmin, cmin2, and M-CUSUM detectors taking advantage of Doppler velocity measurement are more sensitive to the target maneuver compared with MR, IE, and CUSUM detectors. But when maneuver models (A4 and A5) have both the normal acceleration and the tangential acceleration, M-CUSUM detector has a shorter detection delay time. It is the main reason why M-CUSUM detector uses a more accurate prior knowledge of the normal/tangential acceleration.Table 3
The average delay time of maneuver detection (s).
n ^ - n MR IE CUSUM c m i n c m i n 2 M-CUSUM A 1 5.32 4.72 4.26 2.33 2.38 2.37 A 2 4.08 3.58 3.35 2.65 2.71 2.62 A 3 5.13 4.92 4.20 4.09 3.73 3.01 A 4 4.03 3.64 4.03 3.19 3.03 1.46 A 5 3.49 3.91 3.22 2.29 2.15 1.28


### 6.2.3. Estimation Performance Comparison

In the first filter, the estimation performances of SMT algorithms with MR, IE, CUSUM,cmin, cmin2, and M-CUSUM detectors and IMM algorithm are compared. These tracking algorithms all use Doppler velocity measurement in the state estimation. The initial state vector X0 and its standard deviation are set as(27)X0=40km,245m/s,60km,150m/s,σx0=30m,σvx=5m/s,σy0=30m,σvy=5m/s.The initial estimation error covariance matrix denoted asP0, the initial model probability denoted as Pm0, and model transition probability of IMM algorithm denoted as Pm are set as(28)P0=1030000250000103000025,Pm0=0.90.1,Pm=0.90.10.10.9.We execute 100 Monte Carlo simulations for the five scenarios in the first filter. Table4 records the position and velocity of the average RMSE with different tracking algorithms for each scene. The velocity and position of RMSE of different tracking algorithms for each scene are shown at the sampling time in Figure 5, respectively. For the maneuvering target, the sooner it switches the motion model, the more conducive it is to the state estimation for the SMT algorithm. It is concluded that the smaller delay time leads to the lower RMSE, which means the better tracking accuracy according to Tables 3 and 4. Figure 5 shows that M-CUSUM has the smallest RMSE in target maneuver phase compared with MR, IE, CUSUM, cmin, cmin2, and CUSUM. The difference among the average position RMSE of each detector is small. However, the large difference is reflected in the average velocity RMSE in Table 4. When the target maneuvering is stronger, the average velocity RMSE is higher. This shows that the degree of mismatch between the real maneuvering model and the filter motion model is more directly reflected in the velocity than in the position, which is the reason why the efficiency of the detectors based on Doppler velocity measurement is higher than that based on the position measurement. The average position and velocity RMSE of IMM algorithm are higher than the SMT algorithms. The SMT algorithms using the matched target motion model Qlow have lower RMSE during the nonmaneuver phase 0~20 s as shown in Figure 5. As the model noise of IMM algorithm is larger than Qlow, it has a higher RMSE during the nonmaneuver phase 0~20 s. However, as IMM algorithm could adaptively adjust the weights of two models Qlow and Qhigh to match the target maneuver model, it has a lower RMSE during the maneuver phase 21 s~40 s. The difference in average RMSE between IMM algorithm and SMT algorithms mainly comes from the accumulation of nonmaneuver phase as shown in Table 4. The tracking algorithms based on M-CUSUM, cmin, and cmin2 detectors have the comparable estimation performance. Moreover, since M-CUSUM detector has the normal/tangential acceleration in the target maneuver, the corresponding algorithm obtains the better estimation accuracy. The CV motion model is more mismatched to the maneuver model with larger normal accelerations as A2 and A5 compared with other maneuver models as A1 and A3 as shown in Table 4 and Figure 5.Table 4
The average position RMSE (m) and average velocity RMSE (m/s).
Average RMSE MR IE CUSUM c m i n c m i n 2 M-CUSUM IMM A 1 Position 35.47 34.41 33.47 31.52 32.46 32.94 41.39 Velocity 184.01 180.97 170.37 139.58 138.84 129.35 190.37 A 2 Position 44.32 44.18 44.46 41.05 39.51 39.87 49.94 Velocity 236.47 241.97 238.69 193.63 180.59 182.88 207.73 A 3 Position 20.99 20.41 19.47 18.98 19.32 18.87 27.38 Velocity 116.54 117.70 109.31 110.07 113.13 115.36 166.03 A 4 Position 31.80 30.25 30.91 29.64 29.84 28.53 38.06 Velocity 134.03 134.72 126.22 129.89 131.33 126.38 174.31 A 5 Position 60.59 58.89 58.79 58.37 59.01 57.28 71.83 Velocity 270.96 246.01 248.49 220.72 228.27 213.27 272.90Figure 5
RMSE of the first filter.
(a)
Position RMSE ofA1 (b)
Velocity RMSE ofA1 (c)
Position RMSE ofA2 (d)
Velocity RMSE ofA2 (e)
Position RMSE ofA3 (f)
Velocity RMSE ofA3 (g)
Position RMSE ofA4 (h)
Velocity RMSE ofA4  (i)
Position RMSE ofA5 (j)
Velocity RMSE ofA5In the second filter, IMM algorithm adopts six models:{Ai}i=05. The initial state vector X0 and the initial estimation error covariance matrix P0 are the same as the first filter. The initial model probability denoted as Pm0 and model transition probability denoted as Pm of IMM algorithm are set as(29)Pm0=0.90.020.020.020.020.02,Pm=0.90.020.020.020.020.020.020.90.020.020.020.020.020.020.90.020.020.020.020.020.020.90.020.020.020.020.020.020.90.020.020.020.020.020.020.9.The parameters of CSAF algorithm are set as(30)almax=30m/s2,a-lmax=-30m/s2,anmax=70m/s2,a-nmax=-70m/s2,αl=0.1,αn=0.1,where almax and a-lmax are the positive and negative maximum of tangential acceleration; anmax and a-nmax are the positive and negative maximum of normal acceleration; αl and αn are the maneuvering frequency of the tangential and normal acceleration, respectively. The initial state vector X0 and its standard deviation of CSAF are set as(31)X0=40km,60km,245m/s,150m/s,2m/s2,2m/s2,σx0=30m,σy0=30m,σvx=5m/s,σvy=5m/s,σal=0.1m/s2,σan=0.2m/s2.The initial estimation error covariance matrix of CSAF denoted as P0 is set as(32)P0=103000000103000000250000002500000050000005.If the SMT algorithm with the M-CUSUM detector detects the target maneuver, it will adopt the maneuver modelAi,i=1,2,3,4,5, which firstly reaches the detection threshold to replace the term At=A0 in state equation (1). From the perspective of position and velocity, IMM algorithm and SMT algorithm in the second filter have higher estimation accuracy than in the first filter according to Figure 6. The more important reason is that the motion models {Ai}i=05 of state equation (1) are more suitable to the maneuver states of the target in the five scenes than the motion models Qlow and Qhigh. The SMT algorithm with M-CUSUM detector is significantly better than the IMM algorithm in two strong maneuver scenarios as A2 and A5 as shown in Figure 6.Figure 6
RMSE of the second filter.
(a)
Position RMSE ofA1 (b)
Velocity RMSE ofA1 (c)
Position RMSE ofA2 (d)
Velocity RMSE ofA2 (e)
Position RMSE ofA3 (f)
Velocity RMSE ofA3 (g)
Position RMSE ofA4 (h)
Velocity RMSE ofA4  (i)
Position RMSE ofA5 (j)
Velocity RMSE ofA5First of all, this is because the IMM algorithm has the competition problem among multiple motion models, which can cause the probability weights to be dispersed by the unmatched models, and the system error is introduced into the state fusion.Secondly, M-CUSUM detector can achieve higher recognition rates in strong maneuver scenes compared with other scenes, as shown in Table5. The recognition rate of A2 and A5 maneuver models based on M-CUSUM detector is 93% and 95%, respectively.Table 5
Recognition results of 100 Monte Carlo simulations based on M-CUSUM.
Recognition results A 1 A 2 A 3 A 4 A 5 A 1 89 3 3 5 0 A 2 2 93 1 1 3 A 3 4 3 84 7 2 A 4 3 1 3 91 2 A 5 2 2 0 1 95Then, the target state correction effect of the input vector is obvious in the strong maneuver scene, because the magnitude of input is higher than the noise level. However, as the magnitudes of input vector and noise are similar for weak maneuvering scenes, the target state correction of input vector will be interfered with the noise, which can lead to the performance degradation of the SMT algorithm with M-CUSUM detector.The variablesUk and Qk are modified online in CSAF algorithm as the current statistical model is used. Then the performance of CSAF is similar to IMM. CSAF algorithm has bigger position and velocity RMSE than SMT algorithm with M-CUSUM detector in the nonmaneuver stage, as shown in Figure 6. This is because the normal acceleration estimation obtained by the filter directly is more sensitive to the position measurement error than the target maneuvering in the nonmaneuvering stage. The tangential and normal accelerations in the iterative process are mutually independent with the modified Rayleigh distribution in CSAF algorithm. However, it has been proven in Section 6.2.1 that they are correlative. The error due to the independence assumption is constantly accumulated and then CSAF algorithm has significant tracking error around 30 s~40 s of the maneuvering scenarios A4 and A5 that have both the tangential and normal acceleration, as shown in Figures 6(g)–6(j). However, CSAF algorithm has the comparability with STM algorithm and IMM algorithm in the maneuvering scenarios A1, A2, and A3 that have only the normal acceleration or the tangential acceleration.


## 6.2.1. Empirical Distribution and the Approximated GMM

The empirical distributions{gi(cR,al)}i=15 corresponding to the maneuver models {Ai}i=15 are acquired by the statistical experiment. We should assign many values to the initial state (x0,y0,v0,ϕ0) for experiments, because f(an,al) is affected by the velocity magnitude, velocity direction, and position. In each dimension, a number of typical levels are selected for uniform experiments, and the specific values of parameters are shown in Table 1.Table 1
Initial state parameter settings.
Parameters Maximum Minimum Grading x 0 80 km 40 km 20 km y 0 80 km 40 km 20 km v 0 400 m/s 200 m/s 100 m/s α 0 π - 5 π / 6 π / 6In Table1, there are 162 initial state values. We execute 30 Monte Carlo simulations for each maneuvering model in {Ai}i=15 and each initial state value in Table 1. Record the value of (cR,al) at time tc=21 s, and we can get 4860 groups of data for each maneuver model. Then the empirical distribution function gi(cR,al) corresponding to one maneuvering model Ai can be obtained. Two-dimensional GMM by EM algorithm in Appendix B is used to fit the empirical distribution, and then the approximate PDF fi′(an,al) can be acquired. Specifically, three two-dimensional Gauss distribution functions are used to fit the distributions. The scatter diagrams and the approximate GMM diagrams of each maneuver model in {Ai}i=15 are shown in Figure 4. And the parameters of GMM are recorded in Table 2. From Figure 4, different maneuvering models correspond to different distribution patterns, which proves that the probability distribution of variables in formulas (13) and (14) cannot be directly transmitted to an or al. It is also claimed that the normal acceleration is related to the tangential acceleration without the assumption of a constant turn rate motion. Thus, it will introduce an error if the joint PDF f(an,al) is decomposed into the product of marginal PDF f(an)f(al).Table 2
The parameters of GMM.
Maneuvering models μ Σ λ A 1 0.0850 7.1198 1.8693 −1.2086 0.1494 −1.2086 147.0940 0.2879 0.0207 75.3546 1.2346 0.3808 1.2346 1.2266 0.0057 8.8141 3.9931 0.1218 0.4699 0.1218 4.1001  A 2 0.2328 19.6408 91.6298 −2.8712 0.4534 −2.8712 1.3740 −0.1778 13.0432 1.6643 0.8671 0.1157 0.8671 154.5459 0.0761 17.5508 3.3611 0.0244 0.4310 0.0244 7.9094  A 3 8.4827 0.6988 5.1910 −0.0291 0.5962 −0.0291 16.2650 4.5230 0.6965 12.0374 0.0475 0.3000 0.0475 8.4115 9.8643 0.6820 1.9299 0.0232 0.1038 0.0232 83.8308  A 4 10.2197 4.2566 3.3473 0.6226 0.1979 0.6226 119.7231 11.5421 7.1672 4.3900 −7.5418 0.4087 −7.5418 25.3830 10.5928 9.4848 10.7989 −9.5082 0.3934 −9.5082 14.3750  A 5 3.9518 0.2270 59.0632 −1.2884 0.4022 −1.2884 2.8622 −19.1749 0.8473 36.5717 −13.2782 0.0195 −13.2782 75.4061 9.8141 18.4943 24.5112 0.9392 0.5783 0.9392 9.8745Figure 4
Scatter and GMM diagrams corresponding to five maneuver models.
(a)
Scatter diagram ofA1 (b)
f 1 ′ ( a n , a t ) (c)
Scatter diagram ofA2 (d)
f 2 ′ ( a n , a t ) (e)
Scatter diagram ofA3 (f)
f 3 ′ ( a n , a t ) (g)
Scatter diagram ofA4 (h)
f 4 ′ ( a n , a t ) (i)
Scatter diagram ofA5 (j)
f 5 ′ ( a n , a t )


## 6.2.2. Detection Performance Comparison

In the simulation experiment, the M-CUSUM detector is compared with the measurement residual (MR), Input Estimate (IE), CUSUM, andcmin and cmin2 detectors. In the calculation of the average delay time (n^-n), the false alarm probability is set as pf=0.01 and the miss probability is set as pm=0.1. For MR and IE detectors, the window size is set as 5. For cmin detector, the threshold value of fading memory average (FMA) is set as 5000. For cmin2 and CUSUM detectors, the threshold value can be calculated by κ=log⁡((1-pf)/pm) [27], that is, κ=1.9542. For comparison, set the threshold values of M-CUSUM detector as κ1=κ2=⋯=κ5=κ.Table3 summarizes detection delay time of all above algorithms and the detectors have longer detection delay time for relatively small maneuver scenarios as A1 and A3. cmin, cmin2, and M-CUSUM detectors taking advantage of Doppler velocity measurement are more sensitive to the target maneuver compared with MR, IE, and CUSUM detectors. But when maneuver models (A4 and A5) have both the normal acceleration and the tangential acceleration, M-CUSUM detector has a shorter detection delay time. It is the main reason why M-CUSUM detector uses a more accurate prior knowledge of the normal/tangential acceleration.Table 3
The average delay time of maneuver detection (s).
n ^ - n MR IE CUSUM c m i n c m i n 2 M-CUSUM A 1 5.32 4.72 4.26 2.33 2.38 2.37 A 2 4.08 3.58 3.35 2.65 2.71 2.62 A 3 5.13 4.92 4.20 4.09 3.73 3.01 A 4 4.03 3.64 4.03 3.19 3.03 1.46 A 5 3.49 3.91 3.22 2.29 2.15 1.28


## 6.2.3. Estimation Performance Comparison

In the first filter, the estimation performances of SMT algorithms with MR, IE, CUSUM,cmin, cmin2, and M-CUSUM detectors and IMM algorithm are compared. These tracking algorithms all use Doppler velocity measurement in the state estimation. The initial state vector X0 and its standard deviation are set as(27)X0=40km,245m/s,60km,150m/s,σx0=30m,σvx=5m/s,σy0=30m,σvy=5m/s.The initial estimation error covariance matrix denoted asP0, the initial model probability denoted as Pm0, and model transition probability of IMM algorithm denoted as Pm are set as(28)P0=1030000250000103000025,Pm0=0.90.1,Pm=0.90.10.10.9.We execute 100 Monte Carlo simulations for the five scenarios in the first filter. Table4 records the position and velocity of the average RMSE with different tracking algorithms for each scene. The velocity and position of RMSE of different tracking algorithms for each scene are shown at the sampling time in Figure 5, respectively. For the maneuvering target, the sooner it switches the motion model, the more conducive it is to the state estimation for the SMT algorithm. It is concluded that the smaller delay time leads to the lower RMSE, which means the better tracking accuracy according to Tables 3 and 4. Figure 5 shows that M-CUSUM has the smallest RMSE in target maneuver phase compared with MR, IE, CUSUM, cmin, cmin2, and CUSUM. The difference among the average position RMSE of each detector is small. However, the large difference is reflected in the average velocity RMSE in Table 4. When the target maneuvering is stronger, the average velocity RMSE is higher. This shows that the degree of mismatch between the real maneuvering model and the filter motion model is more directly reflected in the velocity than in the position, which is the reason why the efficiency of the detectors based on Doppler velocity measurement is higher than that based on the position measurement. The average position and velocity RMSE of IMM algorithm are higher than the SMT algorithms. The SMT algorithms using the matched target motion model Qlow have lower RMSE during the nonmaneuver phase 0~20 s as shown in Figure 5. As the model noise of IMM algorithm is larger than Qlow, it has a higher RMSE during the nonmaneuver phase 0~20 s. However, as IMM algorithm could adaptively adjust the weights of two models Qlow and Qhigh to match the target maneuver model, it has a lower RMSE during the maneuver phase 21 s~40 s. The difference in average RMSE between IMM algorithm and SMT algorithms mainly comes from the accumulation of nonmaneuver phase as shown in Table 4. The tracking algorithms based on M-CUSUM, cmin, and cmin2 detectors have the comparable estimation performance. Moreover, since M-CUSUM detector has the normal/tangential acceleration in the target maneuver, the corresponding algorithm obtains the better estimation accuracy. The CV motion model is more mismatched to the maneuver model with larger normal accelerations as A2 and A5 compared with other maneuver models as A1 and A3 as shown in Table 4 and Figure 5.Table 4
The average position RMSE (m) and average velocity RMSE (m/s).
Average RMSE MR IE CUSUM c m i n c m i n 2 M-CUSUM IMM A 1 Position 35.47 34.41 33.47 31.52 32.46 32.94 41.39 Velocity 184.01 180.97 170.37 139.58 138.84 129.35 190.37 A 2 Position 44.32 44.18 44.46 41.05 39.51 39.87 49.94 Velocity 236.47 241.97 238.69 193.63 180.59 182.88 207.73 A 3 Position 20.99 20.41 19.47 18.98 19.32 18.87 27.38 Velocity 116.54 117.70 109.31 110.07 113.13 115.36 166.03 A 4 Position 31.80 30.25 30.91 29.64 29.84 28.53 38.06 Velocity 134.03 134.72 126.22 129.89 131.33 126.38 174.31 A 5 Position 60.59 58.89 58.79 58.37 59.01 57.28 71.83 Velocity 270.96 246.01 248.49 220.72 228.27 213.27 272.90Figure 5
RMSE of the first filter.
(a)
Position RMSE ofA1 (b)
Velocity RMSE ofA1 (c)
Position RMSE ofA2 (d)
Velocity RMSE ofA2 (e)
Position RMSE ofA3 (f)
Velocity RMSE ofA3 (g)
Position RMSE ofA4 (h)
Velocity RMSE ofA4  (i)
Position RMSE ofA5 (j)
Velocity RMSE ofA5In the second filter, IMM algorithm adopts six models:{Ai}i=05. The initial state vector X0 and the initial estimation error covariance matrix P0 are the same as the first filter. The initial model probability denoted as Pm0 and model transition probability denoted as Pm of IMM algorithm are set as(29)Pm0=0.90.020.020.020.020.02,Pm=0.90.020.020.020.020.020.020.90.020.020.020.020.020.020.90.020.020.020.020.020.020.90.020.020.020.020.020.020.90.020.020.020.020.020.020.9.The parameters of CSAF algorithm are set as(30)almax=30m/s2,a-lmax=-30m/s2,anmax=70m/s2,a-nmax=-70m/s2,αl=0.1,αn=0.1,where almax and a-lmax are the positive and negative maximum of tangential acceleration; anmax and a-nmax are the positive and negative maximum of normal acceleration; αl and αn are the maneuvering frequency of the tangential and normal acceleration, respectively. The initial state vector X0 and its standard deviation of CSAF are set as(31)X0=40km,60km,245m/s,150m/s,2m/s2,2m/s2,σx0=30m,σy0=30m,σvx=5m/s,σvy=5m/s,σal=0.1m/s2,σan=0.2m/s2.The initial estimation error covariance matrix of CSAF denoted as P0 is set as(32)P0=103000000103000000250000002500000050000005.If the SMT algorithm with the M-CUSUM detector detects the target maneuver, it will adopt the maneuver modelAi,i=1,2,3,4,5, which firstly reaches the detection threshold to replace the term At=A0 in state equation (1). From the perspective of position and velocity, IMM algorithm and SMT algorithm in the second filter have higher estimation accuracy than in the first filter according to Figure 6. The more important reason is that the motion models {Ai}i=05 of state equation (1) are more suitable to the maneuver states of the target in the five scenes than the motion models Qlow and Qhigh. The SMT algorithm with M-CUSUM detector is significantly better than the IMM algorithm in two strong maneuver scenarios as A2 and A5 as shown in Figure 6.Figure 6
RMSE of the second filter.
(a)
Position RMSE ofA1 (b)
Velocity RMSE ofA1 (c)
Position RMSE ofA2 (d)
Velocity RMSE ofA2 (e)
Position RMSE ofA3 (f)
Velocity RMSE ofA3 (g)
Position RMSE ofA4 (h)
Velocity RMSE ofA4  (i)
Position RMSE ofA5 (j)
Velocity RMSE ofA5First of all, this is because the IMM algorithm has the competition problem among multiple motion models, which can cause the probability weights to be dispersed by the unmatched models, and the system error is introduced into the state fusion.Secondly, M-CUSUM detector can achieve higher recognition rates in strong maneuver scenes compared with other scenes, as shown in Table5. The recognition rate of A2 and A5 maneuver models based on M-CUSUM detector is 93% and 95%, respectively.Table 5
Recognition results of 100 Monte Carlo simulations based on M-CUSUM.
Recognition results A 1 A 2 A 3 A 4 A 5 A 1 89 3 3 5 0 A 2 2 93 1 1 3 A 3 4 3 84 7 2 A 4 3 1 3 91 2 A 5 2 2 0 1 95Then, the target state correction effect of the input vector is obvious in the strong maneuver scene, because the magnitude of input is higher than the noise level. However, as the magnitudes of input vector and noise are similar for weak maneuvering scenes, the target state correction of input vector will be interfered with the noise, which can lead to the performance degradation of the SMT algorithm with M-CUSUM detector.The variablesUk and Qk are modified online in CSAF algorithm as the current statistical model is used. Then the performance of CSAF is similar to IMM. CSAF algorithm has bigger position and velocity RMSE than SMT algorithm with M-CUSUM detector in the nonmaneuver stage, as shown in Figure 6. This is because the normal acceleration estimation obtained by the filter directly is more sensitive to the position measurement error than the target maneuvering in the nonmaneuvering stage. The tangential and normal accelerations in the iterative process are mutually independent with the modified Rayleigh distribution in CSAF algorithm. However, it has been proven in Section 6.2.1 that they are correlative. The error due to the independence assumption is constantly accumulated and then CSAF algorithm has significant tracking error around 30 s~40 s of the maneuvering scenarios A4 and A5 that have both the tangential and normal acceleration, as shown in Figures 6(g)–6(j). However, CSAF algorithm has the comparability with STM algorithm and IMM algorithm in the maneuvering scenarios A1, A2, and A3 that have only the normal acceleration or the tangential acceleration.


## 7. Conclusions

A theoretical framework to support the detection and recognition of multiple maneuvering models is proposed in the paper. Multiple Parallel CUSUM detector named as M-CUSUM is proposed by extending the standard CUSUM detector. The normal and tangential accelerations based on Doppler velocity measurement are deduced without the assumption of a constant turn rate motion. A new statistic of normal acceleration is proposed to reduce the error based on the statisticcmin. The joint PDF of normal and tangential accelerations is used to describe the target maneuvering model. Its joint empirical PDF is obtained by the statistical experiment method and approximated by use of GMM with EM algorithm. The computer simulation results show that M-CUSUM detector not only has a good detection performance but also can recognize the target maneuver model to increase the matching degree of motion model.



---

*Source: 5062184-2018-03-12.xml*