Application of Interactive Multiple Model Adaptive Five-Degree Cubature Kalman Algorithm Based on Fuzzy Logic in Target Tracking

: Aiming at the shortcomings of low precision, hysteresis, and poor robustness of the general interactive multimodel algorithm in the “snake-like” maneuver tracking of anti-ship missiles, an interactive multimodel adaptive ﬁve-degree cubature Kalman algorithm based on fuzzy logic (FLIMM5ACKF) is proposed. The algorithm mainly includes adaptive ﬁve-degree cubature Kalman algorithm (A5CKF) and fuzzy logic algorithm (FL). A5CKF uses the Sage–Husa noise estimation principle to propose a state error covariance adaptive ﬁve-degree cubature Kalman algorithm to improve the performance of state estimation. Then, the fuzzy logic algorithm (FL) is added to the model probability update module to control the model probability update module. Finally, by setting the same tracking model simulation analysis, the algorithm has better convergence speed, tracking e ﬀ ect and robustness than the interactive multimodel cubature Kalman algorithm (IMMCKF), the interactive multimodel ﬁve-degree cubature Kalman algorithm (IMM5CKF) and the interactive multimodel adaptive ﬁve-degree cubature Kalman (IMMA5CKF).


Introduction
With the development of modern military, anti-ship missiles have various forms of maneuvering, and the most typical one is the "snake-like" maneuver [1], because the maneuvering method of anti-ship missiles is relatively variable, a model can no longer meet the needs of anti-ship missile tracking.Therefore, the multiple model algorithm is proposed, for example, the multiple model (MM) algorithm proposed by Magill [2], the GPB algorithm proposed by Ackerson and Fu [3], and the interactive multimodel (IMM) algorithm proposed by Blom [4,5], among these algorithms, the most widely used is IMM algorithm [6][7][8], which takes into account the characteristics of the model and considers that only one model matches the motion state at a certain time, it reduces the error of a single model, improves the effect of target tracking.However, IMM algorithm also has some shortcomings, in order to solve the problems that the traditional IMM algorithm has, such as low filtering precision, slow convergence speed, and linear filtering in target tracking, this article will improve the IMM algorithm.
Nowadays, the research of nonlinear filtering algorithms has attracted more and more researchers' attention.In the aspect of nonlinear filtering algorithm, the most widely used is the extended Kalman algorithm (EKF) [9,10], its core idea is to use linearization to approximate nonlinearization according to the first-order Taylor series expansion, however, the Jacobian matrix of nonlinear functions is difficult to find in many practical problems.In order to better improve the estimation performance, Julier et al. proposed the Unscented Kalman Filter (UKF) [11,12], which has better performance than EKF.In order to seek better filtering effects, the cubature Kalman algorithm (CKF) based on the third-order spherical-radial cubature rule was proposed by Haykin et al. [13,14]; it is different from UKF, CKF has a strict mathematical formula to prove, it has been shown that CKF has better performance than UKF when the dimension of the system is greater than three [15].In order to pursue better performance of CKF, 5CKF was proposed in the literature [16][17][18], and the simulation shows that 5CKF performance is better than CKF.For pursuing a better estimation performance of 5CKF filtering algorithm, combined with the improvement of Sage-Husa noise estimation in the literature [19], and the combination with the corresponding filtering algorithm in the literature [20][21][22][23], the adaptive five-degree cubature Kalman algorithm based on error covariance estimation is proposed, and it is used as a filtering algorithm under the IMM framework.
In the IMM algorithm, the transition probability between models is performed by the Markov transition matrix.When the target is tracked at a certain time, it works by a specific model, and its probability is approximately 1, other approximations are 0.However, this method requires a certain amount of time, and there is hysteresis, resulting in a decrease in tracking effect.For solving the problem of slow convergence caused by hysteresis, a fuzzy logic (FL) is proposed [24], and combine the FL algorithm with the IMM algorithm, and the corresponding fuzzy rules are formulated, when the model is transformed, the FL algorithm is used to judge whether the probability of the model is 1 or 0 to accelerate the update of the probability, this will speed up the convergence rate, the related references also proves the effectiveness of the combined algorithm [20,25].
Although the IMMCKF algorithm and the IMM5CKF algorithm have already achieved good results [14,18], they still cannot effectively solve the problem of low filtering precision and slow convergence in the tracking process, therefore, an interactive multimodel adaptive five-degree cubature Kalman algorithm based on fuzzy logic is proposed in this paper, it uses the maximum likelihood function obtained by the improved A5CKF algorithm in the parallel filtering, updates the probability through the FL algorithm, and finally obtains the result through the output data fusion.Finally, by setting the same simulation model analysis, compared with IMMCKF [14], IMMA5CKF, and IMM5CKF [18], FLIMMA5CKF has better tracking effect and robustness, and the hysteresis is also improved.The rest of the sections are arranged as follows.Section 2 introduces fuzzy logic algorithms, Section 3 introduces adaptive five-degree cubature Kalman algorithm, Section 4 introduces FLIMMA5CKF algorithm, Section 5 performs missile dynamics modeling, Section 6 provides simulation experiments with the analysis of the results, and the conclusions are provided in Section 7.

Fuzzy Logic Algorithm
The fuzzy logic algorithm is used to update the model probability [24], so that the model probability is quickly converted to accelerate the response speed of the filtering system.Figure 1 shows the structure of the fuzzy logic algorithm in the model probability update module.[11,12], which has better performance than EKF.In order to seek better filtering effects, the cubature Kalman algorithm (CKF) based on the third-order spherical-radial cubature rule was proposed by Haykin et al. [13,14]; it is different from UKF, CKF has a strict mathematical formula to prove, it has been shown that CKF has better performance than UKF when the dimension of the system is greater than three [15].In order to pursue better performance of CKF, 5CKF was proposed in the literature [16][17][18], and the simulation shows that 5CKF performance is better than CKF.For pursuing a better estimation performance of 5CKF filtering algorithm, combined with the improvement of Sage-Husa noise estimation in the literature [19], and the combination with the corresponding filtering algorithm in the literature [20][21][22][23], the adaptive five-degree cubature Kalman algorithm based on error covariance estimation is proposed, and it is used as a filtering algorithm under the IMM framework.
In the IMM algorithm, the transition probability between models is performed by the Markov transition matrix.When the target is tracked at a certain time, it works by a specific model, and its probability is approximately 1, other approximations are 0.However, this method requires a certain amount of time, and there is hysteresis, resulting in a decrease in tracking effect.For solving the problem of slow convergence caused by hysteresis, a fuzzy logic (FL) is proposed [24], and combine the FL algorithm with the IMM algorithm, and the corresponding fuzzy rules are formulated, when the model is transformed, the FL algorithm is used to judge whether the probability of the model is 1 or 0 to accelerate the update of the probability, this will speed up the convergence rate, the related references also proves the effectiveness of the combined algorithm [20,25].
Although the IMMCKF algorithm and the IMM5CKF algorithm have already achieved good results [14,18], they still cannot effectively solve the problem of low filtering precision and slow convergence in the tracking process, therefore, an interactive multimodel adaptive five-degree cubature Kalman algorithm based on fuzzy logic is proposed in this paper, it uses the maximum likelihood function obtained by the improved A5CKF algorithm in the parallel filtering, updates the probability through the FL algorithm, and finally obtains the result through the output data fusion.Finally, by setting the same simulation model analysis, compared with IMMCKF [14], IMMA5CKF, and IMM5CKF [18], FLIMMA5CKF has better tracking effect and robustness, and the hysteresis is also improved.The rest of the sections are arranged as follows.Section 2 introduces fuzzy logic algorithms, Section 3 introduces adaptive five-degree cubature Kalman algorithm, Section 4 introduces FLIMMA5CKF algorithm, Section 5 performs missile dynamics modeling, Section 6 provides simulation experiments with the analysis of the results, and the conclusions are provided in Section 7.

Fuzzy Logic Algorithm
The fuzzy logic algorithm is used to update the model probability [24], so that the model probability is quickly converted to accelerate the response speed of the filtering system.Figure 1 shows the structure of the fuzzy logic algorithm in the model probability update module.In this paper, we use two kinds of maneuver models of constant velocity and "snake-like" to simulate, so the input of the filter of the model probability update module is Λ 1 and Λ 2 .Firstly, the model probabilities ̂1 and ̂2, corresponding to the two models, are calculated according to the IMM algorithm to calculate the model probability.In this paper, we use two kinds of maneuver models of constant velocity and "snake-like" to simulate, so the input of the filter of the model probability update module is Λ 1 and Λ 2 .Firstly, the model probabilities μ1 and μ2 , corresponding to the two models, are calculated according to the IMM algorithm to calculate the model probability.
Symmetry 2019, 11, 767 3 of 16 It can be derived from Formulas (1)-(3): In the formula where p ji is the transition probability of the model j to i, and µ i is the probability of the first model, µ j (l) is the probability of IMM update, and c is the normalization constant.
Analyze model 1 (uniform linear motion), assuming that the input variables obtained by the model probability update module are I 1 and I 2 , and the resulting output variable is µ, let (1) The domain of input variables and output variables: (2) Fuzzy set of input variables and output variables: (3) Determining the membership function: Refer to Figures 2-4 to formulate fuzzy rules, I 1 represents the probability of model 1 at (l − 1) moment, I 2 represents the difference between the updated probability and the pre-update probability after updating at the moment of l. µ is the probability of time l, then, analyze these three variables: if I 2 is a negative value and is denoted as N, then the probability µ at k moment must be smaller than the model probability I 1 at (l − 1) moment; if I 2 is equal to 0, it is recorded as Z, then the probability µ at k moment should be equal to the model probability I 1 at (l − 1) moment; if I 2 is a positive value and is denoted as P, then the probability µ at l moment must be greater than the model probability I 1 at (l − 1) moment.Their fuzzy logic rules are shown in Table 1:          When the output of the fuzzy logic system is defuzzified, the median method is used to defuzzify and obtain the actual probability of the model.

Adaptive Five-Degree Cubature Kalman
In order to pursue the better performance of the 5CKF filtering algorithm [16,17], the error covariance is improved based on the principle of Sage-Husa noise estimator, and a five-degree cubature Kalman algorithm based on error covariance adaptive is proposed.
Assume random variable x, its mean is x = E[x], P = cov[x] is its covariance matrix, the A5CKF algorithm consists of two parts: time update and measurement update.

Time Update
(1) Performing Cholesky decomposition on the covariance matrix P l−1|l−1 : (2) Calculate the cubature point: ) where the dimension of the system is n; the weight of the cubature point is ω 0 , ω 1 , and ω 2 ; and the unit vector is e i .
(3) Propagating the cubature point: Substituting the cubature point into the equation of state to obtain the cubature point after propagation.
Symmetry 2019, 11, 767 (4) Prediction of state values: Calculate the predicted value by combining the transmitted sample points with the corresponding weighted values. xl|l l|l−1 (30) (5) Prediction of the state error covariance matrix: l|l−1 (34)

Measurement Update
(1) Performing Cholesky decomposition on the prediction covariance matrix P l|l−1 : (2) Calculate the cubature point: (3) Propagating the cubature point: Substituting the cubature point into the observation equation to obtain the cubature point after propagation.
(4) Predicted observation: Calculating the observed values by weighted summation of the propagated cubature points.
(5) Prediction of the error covariance matrix: Calculating the predicted error covariance matrix by weighted summation.
(6) Calculation of the cross-covariance matrix: Calculating the cross-covariance matrix by weighted summation.

Algorithm Overall Framework
In order to obtain good filtering accuracy and response speed, A5CKF and FL algorithms are added in IMM, and an interactive multimodel adaptive five-degree cubature Kalman algorithm based on FL algorithm is proposed, which has the advantages of three algorithms.In the FLIMMA5CKF algorithm, the A5CKF filter is used for parallel filtering, the obtained estimation values are weighted and fused, and finally the state estimation is output.The execution steps of the FLIMMA5CKF algorithm and the IMM algorithm are basically the same: it is divided into four steps: input interaction, parallel filtering, update probability, and output data fusion.

Input Interaction
Estimating the current time measurement value   according to the system state estimation and covariance estimation at the previous moment, and then performing a new initialization calculation on the model, wherein the new initial value is obtained by the Markov matrix between different models.Set the state optimal estimate and the estimated covariance matrix at model time ( − 1):

Input Interaction
Estimating the current time measurement value Z l according to the system state estimation and covariance estimation at the previous moment, and then performing a new initialization calculation on the model, wherein the new initial value is obtained by the Markov matrix between different models.Set the state optimal estimate and the estimated covariance matrix at model time (l − 1): Symmetry 2019, 11, 767 The following is the corresponding meaning: at time (l − 1), for model j, its conditional probability is µ j|i l−1|l−1 , its probability is µ j l−1 , its mathematical expectation is X0i l−1|l−1 , and its error covariance is p 0i l−1|l−1 .

Parallel Filtering
The choice of filtering algorithm can also directly determine the tracking effect.Therefore, the parallel filtering module selects the A5CKF algorithm proposed in Section 3 of this paper, which has better estimation effect.Then, based on the previously obtained a priori information combined with the execution steps of A5CKF, finally, we can get updated state variables and recursive estimation P l|l statistical characteristics.

Update Probability
Model updating is also the key factor to determine the effectiveness of the algorithm.It uses maximum likelihood function combined with FL algorithm in this paper, and uses the principle of fuzzy logic to control the updating of model probability.
The likelihood function: where the filter residual is represented by ε i l = z − z i l|l−1 , and the covariance matrix of the residual is S i,l

Output Data Fusion
According to the results of each model and their weights, the final output results are obtained.

Analysis of The Motion Characteristics of Missiles
The missile's motion can be seen as a six-degrees-of-freedom rigid body motion.The general situation is to describe the missile's motion by solving the missile's six-degree-of-freedom full-quantity equation [20].Since the computational complexity of the high-order solution is large, and in most cases, the aerodynamic data of the missile cannot be obtained.It is necessary to adopt a simplified dynamic model.Therefore, according to the law of missile movement, the lift, thrust and side force of the missile are planned.The position and attitude of the missile are obtained from the equation of motion of the missile's centroid, and then the missile's track is generated.The centroid motion equation of the missile is as follows: where v is the missile velocity; θ is the ballistic inclination of the missile; ψ v is the ballistic angle of the missile, x, y, and z are components in the inertial coordinate system; g x2 , g y2 , and g z2 are the components of gravity acceleration in the ballistic coordinate system; and a x2 , a y2 , and a z2 are the components of the ballistic overload in the inertial coordinate system.
Firstly, calculating the actual movement position of the missile using Runge-Kutta algorithm, and then the coordinate transformation is performed to convert the speed coordinate system to the projectile coordinate system for simulation.

Snake-Like Maneuvering
The "snake-like" is a maneuvering model: it is a special turning motion, when moving half a circle, its angular velocity direction will be opposite, and its radius is also changing over time [25].Suppose the target moves in a sinusoidal motion in the horizontal plane.The target makes a uniform linear motion in the X axis, the speed is V, the speed of the Y axis is 0, and the sine motion is performed in the Z-axis.
The discretization form of the discretization available state equation is: Deriving the discretized process noise input matrix G and the state transition matrix F: Symmetry 2019, 11, 767 11 of 16

Results and Discussion
Setting the same simulation conditions, the missile performs uniform linear motion in 0-15 s and snake motion in 16-60 s.Comparative analysis was performed using FLIMMA5CKF with IMMCKF, IMM5CKF, and IMMA5CKF algorithms.
The nonlinear observation equation of the target is described as follows: where Z(t) is the observation vector, V(t) is the measurement noise, and h[ * ] is the observation function.
The value is: The performance of the four algorithms can be measured by the displacement, velocity, acceleration Error and RMSE at time l: where m is the number of Monte Carlo simulations, n is the total number of sampling points, X m,k is the true value, and Xm,l is the estimated value.Assume that the target's maneuvering frequency is ω = 0.2π, the filtering period is 0.1 s, the Runge-Kutta method step size is 0.01 s, the initial probability of both models is 0.5, the initial state of the target is X(0) = 90, 000 2000 50 −1450 0 0 0 0 −20 T , the initial state error covariance matrix is P(0) = diag 100 100 100 25 25 25 1 1 1 , measuring noise standard deviation is 1 mrad, the standard noise standard deviation is 10, the Monte Carlo simulations are 50 times, and the transfer matrix between models is: 0.999 0.001 0.001 0.999 (80) Figures 6-8 are comparisons of error estimates for tracking targets using the improved FLIMMA5CKF algorithm, IMMA5CKF algorithm, IMM5CKF algorithm [18], and IMMCKF algorithm [14], respectively.
Figures 6-8 are comparisons of error estimates for tracking targets using the improved FLIMMA5CKF algorithm, IMMA5CKF algorithm, IMM5CKF algorithm [18], and IMMCKF algorithm [14], respectively.From Figures 6-8, combined with the data in Table 2, it can be seen that the estimation error of the IMMA5CKF algorithm on the X-axis and the Y-axis is significantly smaller than the IMM5CKF algorithm of Zhu [18], and the IMMCKF algorithm of Li [14], which proves the high-precision characteristics of the A5CKF algorithm proposed in this paper, it not only adapts the error covariance, but also adapts the noise.Comparing the FLIMMA5CKF algorithm with IMMA5CKF, we can find that the estimation error of the FLIMMA5CKF algorithm is smaller than that of the IMMA5CKF algorithm.It proves that the FL algorithm proposed in this paper not only improves the accuracy of the algorithm, but also accelerates the convergence speed.In the Z-axis direction, the FLIMMA5CKF algorithm is slightly worse than the other three algorithms before 35 s, but after 35 s, the error reduction is obvious and more stable, in contrast, the other three algorithms have varying degrees of mutation or jump, in Figure 6, the late estimation error is approximately ± 1 m, in Figure 7, the late speed estimation error is ~0.4 m/s, in Figure 8, and the late acceleration estimation error is ~0.2 m/s .From Figures 6-8, combined with the data in Table 2, it can be seen that the estimation error of the IMMA5CKF algorithm on the X-axis and the Y-axis is significantly smaller than the IMM5CKF algorithm of Zhu [18], and the IMMCKF algorithm of Li [14], which proves the high-precision characteristics of the A5CKF algorithm proposed in this paper, it not only adapts the error covariance, but also adapts the noise.Comparing the FLIMMA5CKF algorithm with IMMA5CKF, we can find that the estimation error of the FLIMMA5CKF algorithm is smaller than that of the IMMA5CKF algorithm.It proves that the FL algorithm proposed in this paper not only improves the accuracy of the algorithm, but also accelerates the convergence speed.In the Z-axis direction, the FLIMMA5CKF algorithm is slightly worse than the other three algorithms before 35 s, but after 35 s, the error reduction is obvious and more stable, in contrast, the other three algorithms have varying degrees of mutation or jump, in Figure 6, the late estimation error is approximately ±1 m, in Figure 7, the late speed estimation error is ~0.4 m/s, in Figure 8, and the late acceleration estimation error is ~0.2 m/s 2 .
Figures 9-11 shows the RMSE of the four algorithms.As can be seen from Figures 9-11, when the missile moves linearly at a constant speed, there is no significant difference in the tracking effect of several algorithms, but after 16 s, the missile begins to accelerate in a snake-like manner, as can be seen from Figures 9-11.Compared with the IMMCKF algorithm [14], the IMM5CKF algorithm [18], and the IMMA5CKF algorithm, the FLIMMA5CKF algorithm has a great improvement in tracking performance and stability.It can be seen from the Z-axis of Figures 9-11, FLIMMA5CKF has a better convergence speed than the other three algorithms, no peaks appear after 40 s, and the robustness is better.Figure 12 shows the model probability update graphs of the four algorithms, FLIMMA5CKF algorithm can better control model probability update due to the addition of fuzzy logic algorithm.It can be seen from the data in Table 3 that the performance of the FLIMMA5CKF algorithm in this paper is significantly higher than that of the IMMCKF algorithm [14], IMM5CKF algorithm [18], and IMMA5CKF algorithm in the X-axis and Y-axis directions.In the Z-axis direction, the algorithm of this paper is better than the other three, but the advantages of the algorithm are not particularly large.The shortcoming is that the FLIMMA5CKF algorithm in this paper is more complex than the other three algorithms, and the simulation time of the algorithm is relatively long.It can be seen from the data in Table 3 that the performance of the FLIMMA5CKF algorithm in this paper is significantly higher than that of the IMMCKF algorithm [14], IMM5CKF algorithm [18], and IMMA5CKF algorithm in the X-axis and Y-axis directions.In the Z-axis direction, the algorithm of this paper is better than the other three, but the advantages of the algorithm are not particularly large.The shortcoming is that the FLIMMA5CKF algorithm in this paper is more complex than the other three algorithms, and the simulation time of the algorithm is relatively long.

Figure 1 .
Figure 1.Fuzzy probabilistic model probability update module structure.

Figure 1 .
Figure 1.Fuzzy probabilistic model probability update module structure.

Table 2 .
X, Y, and Z direction maximum error.

Table 2 .
X, Y, and Z direction maximum error.

Table 3 .
Average RMSE comparison of four algorithms.

Table 3 .
Average RMSE comparison of four algorithms.