2. School of Automation and Electrical Engineering, University of Science and Technology of Beijing, Beijing 100083, China;
3. Social Robotics Laboratory, Interactive Digital Media Institute, and the Department of Electrical and Computer Engineering, National University of Singapore, Singapore 117576, Singapore
The requirement for accuracy is high in this area of autonomous navigation, especially with rapid developments in information technology. However, obtaining the motion state of a ballistic missile with high precision is very difficult when using only a single navigation technology. Strapdown inertial navigation system (SINS) has attracted more and more attention from transportation department and the space shuttle, in view of a lot of significant advantages of strong autonomy, high accuracy, and low computational cost. However, the errors of SINS increase with time, hence for long durations, it is difficult to use SINS independently [1]. Celestial navigation system (CNS) is also an autonomous navigation system where errors do not accumulate over time. Furthermore, CNS can provide the advantages of high location precision and good concealment. While it cannot independently accomplish positioning through navigation with the constraint of external conditions. Thus, there is a complementary relationship between SINS and CNS in terms of error characteristics, and we may increase the performance of navigation through their combination [2]. No matter domestic or abroad, integrated navigation has been widely used with a series of advantages and obtained important achievements. For instance, utilizing each subsystem's information can increase precision of the system positioning and raise the system reliability [3]. Accordingly, SINS/CNS integrated navigation can be applied to the ballistic missile movement.
For this navigation system, filtering algorithm is a significant factor affecting navigation performance, which is frequently used to estimate state and filter measurement. In the information fusion, Kalman filter (KF) is a basic algorithm aiming at acquiring the attitude, velocity, and position information of the target. Improved filtering methods, such as EKF and UKF, have been developed based on the standard Kalman filter [4].
It is difficult to gain the optimal navigation estimation solely relying on KF due to the nonlinearity of the missile's model, and both the EKF and UKF are suitable for nonlinear system [5]. EKF works by utilizing a firstorder Taylor series expansion to approximate linear system, then processing equations according to linear equations. However, linearization may cause the filter's divergence and poor performance when the system has strong nonlinearity. Julier and Uhlman provide the theory of UKF, which effectively solves the problem related to EKF based on the assumption that the state variables' distribution approximates a Gaussian distribution [6], [7]. Essentially, UKF adopts an unscented transformation (UT) strategy and provides a higher order moment of function without the need to linearize the nonlinear equations or compute the Jacobian matrix, as required in EKF. Therefore, UKF has no higherorder terms truncation error and can achieve higher precision compared with EKF for nonlinear system [8]. Under the hypothesis that the process and measurement noise of EKF and UKF are a Gaussian distribution with zeromean, we end up with the problem of filter tuning. This problem is timeconsuming and taxing. Fortunately, there are many available tools to aid the filter design. A popular approach to resolve this problem is to adopt an adaptive filter.
Mehra offered four distinct categories such as Bayes estimate, covariance matching, correlation, and maximum likelihood (ML) for the adaptive filter. The procedures used to implement the adaptive Kalman filter are classified into three main methods: innovationbased adaptive estimator (IAE), multiple model adaptive estimate (MMAE), and adaptive fading Kalman filter (AFKF) in [9]. MMAE is an important aspect in multimodel adaptive control, which is appropriate for systems with uncertain parameters or variable structure [10], [11]. The central idea of multiple model (MM) was put forward by Magill in 1956. This algorithm has become a mainstream approach to study estimation problems in hybrid systems. Meanwhile, the interactive multiple model (IMM) algorithm and static multiple model algorithm are most widely used in what. In MMAE method, the state estimation of the system is obtained by the weighted sum of each subfilter's estimation with the likelihood function of unknown elements. Thus, this method can adapt to the systems with uncertain or unknown parameters [12].
The contributions of this paper are to propose the MMAEUKF method and apply this proposed filtering approach to a SINS/CNS combined navigation system. At the same time, we carry out a simulation with motion of ballistic missile and make comparisons among EKF, MMAEEKF, UKF, and MMAEUKF algorithms. Simulation results are provided to indicate that MMAEUKF can obtain the highest navigation precision.
The organization of the remainder of this paper is shown below. Multiple model adaptive estimation method is presented in Section Ⅱ. Section Ⅲ describes the MMAEUKF algorithm in detail. In Section Ⅳ, the ballistic missile motion and integrated navigation is specified. Then, Section Ⅴ gives the simulation conditions and results. In the end, conclusions are drawn in Section Ⅵ.
Ⅱ. MULTIPLE MODEL ADAPTIVE ESTIMATIONThis section primarily researches the multiple model adaptive estimate algorithm, which is composed of a set of parallel filters and a hypothesis testing. In this filter bank, each filter has a specific model of the system which is expressed by independent vector parameters. Then, each filter forms estimated value of the current system state according to the input vector and selfmodel independently. Making use of the current state estimate to calculate predicted measurement, the residual error is obtained through subtraction of predicted measurement value from true value [13]. Then, based on the true vector parameters and measurement, hypothesis testing applies the residual error to compute the conditional probability of each model, and the obtained probability is applied to measure the estimation accuracy of each filter state. In the end, we obtain the optimal estimate of the system state with the state estimation of each model multiplied by the corresponding weights [14]. The MMAE architecture is demonstrated in Fig. 1.
Download:


Fig. 1 The MMAE architecture 
Multiple model adaptive estimation is a recursive estimator, which
utilizes a group of filters with unknown elements. The
elements are the parameters of process noise covariance matrix, indicated by the vector
$ \begin{align} p\left({\pmb p}^{(l)}\widetilde{{\pmb Y}}_{k}\right)= \dfrac{p\left(\widetilde{{\pmb Y}}_{k}{\pmb p}^{(l)}\right) p\left({\pmb p}^{(l)}\right)}{\sum\limits_{j=1}^{M}p\left(\widetilde{{\pmb Y}}_{k}{\pmb p}^{(j)}\right) p\left({\pmb p}^{(j)}\right)} \end{align} $  (1) 
where
And the posteriori probability density function calculated by
$ \begin{align}\label{eq.2} p\left({\pmb p}^{(l)}\widetilde{{\pmb Y}}_{k}\right) &= \frac{p\left(\widetilde{{\pmb y}}_{k}, {\pmb p}^{(l)}\widetilde{{\pmb Y}}_{k1}\right)} {p\left(\widetilde{{\pmb y}}_{k}\widetilde{{\pmb Y}}_{k1}\right)} \nonumber \\[2mm] &= \frac{p(\widetilde{{\pmb Y}}_{k} \widehat{{\pmb x}}^{(l)}_{k})p\left({\pmb p}^{(l)}\widetilde{{\pmb Y}}_{k1}\right)} {\sum\limits_{j=1}^{M}\left[p\left(\widetilde{{\pmb Y}}_{k}\widehat{{\pmb x}}^{(j)}_{k}\right) p\left({\pmb p}^{(j)}\widetilde{{\pmb Y}}_{k1}\right)\right]} \end{align} $  (2) 
where
$ \begin{align}\label{eq.3} \varpi^{(l)}_{k} \,&=\varpi^{(l)}_{k1}p\left(\widetilde{{\pmb y}}_{k} \widehat{{\pmb x}}^{(l)}_{k}\right) \nonumber\\[2mm] \varpi^{(l)}_{k}&\leftarrow \frac{\varpi^{(l)}_{k}}{\sum\limits_{j=1}^{M} \varpi^{(j)}_{k}} \end{align} $  (3) 
where
The conditional mean estimation is the sum of parallel filters estimations by weighting as below:
$ \begin{align} \widehat{{\pmb x}}^{}_{k} = \sum\limits_{l=1}^{M}\varpi^{(l)}_{k}\widehat{ {\pmb x}}^{(l)}_{k}. \end{align} $  (4) 
Meanwhile, the state estimations' error covariance is calculated through
$ \begin{align} {pmb P}^{}_{k} = \sum\limits_{l=1}^{M}\varpi^{(l)}_{k} \left(\widehat{{\pmb x}}^{(l)}_{k}  \widehat{{\pmb x}}^{}_{k}\right) \left(\widehat{{\pmb x}}^{(l)}_{k}  \widehat{{\pmb x}}^{}_{k}\right)^{T}. \end{align} $  (5) 
The likelihood function for the residual of measurement is constructed in this part. Primarily, the residual error is given by
$ \begin{align} \epsilon_{i} \equiv \left[ \begin{array}{c} {\pmb e}_{k}\\ {\pmb e}_{k1}\\ \vdots\\ {\pmb e}_{ki} \end{array}\right] \end{align} $  (6) 
where
$ \begin{align} L_{i} = \frac{1}{\left[\det(2\pi {{\pmb C}}_{i})\right]^{\frac{1}{2}}} \exp \left(\frac{1}{2}\epsilon_{i} ^{T}{{\pmb C}}_{i}^{1} \epsilon_{i}\right) \end{align} $  (7) 
where
When
$ \begin{align} L_{0} = &\ \frac{1}{\left\{\det\left[2\pi({{\pmb H}}_{k} {{\pmb P}}_{k}^{}{{\pmb H}}_{k}^{T} + {{\pmb R}}_{k})\right]\right\}^{\frac{1}{2}}} \nonumber\\[2mm] &\, \times \exp \left[\frac{1}{2}{\pmb e}_{k}^{T}(H_{k}{{\pmb P}}_{k}^{}H_{k}^{T} + {{\pmb R}}_{k})^{1}{\pmb e}_{k}\right]. \end{align} $  (8) 
This part presents the adaptive laws. The current measurements are just applied to update law given by (3) in standard MMAE method. Thus, the update law on account of the likelihood function is obtained through
$ \begin{align} & \varpi^{(l)}_{k} = \varpi^{(l)}_{k1}L^{(l)}_{0} \nonumber\\[2mm] & \varpi^{(l)}_{k} \leftarrow \frac{\varpi^{(l)}_{k}} {\sum\limits_{j=1}^{M}\varpi^{(j)}_{k}} \end{align} $  (9) 
since
$ \begin{align} L_{0}^{(l)} = &\ \frac{1}{\left\{\det\left[2\pi({{\pmb H}}_{k} {{\pmb P}}_{k}^{(l)}{{\pmb H}}_{k}^{T} + {{\pmb R}}_{k})\right]\right\}^{\frac{1}{2}}} \nonumber\\[2mm] &\, \times \exp \left[\frac{1}{2}{\pmb e}_{k}^{(l)T}({{\pmb H}}_{k} {{\pmb P}}_{k}^{(l)}{{\pmb H}}_{k}^{T} + {{\pmb R}}_{k})^{1}{\pmb e}_{k}^{(l)}\right]. \end{align} $  (10) 
In what follows, we consider the form of the discrete nonlinear system model
$ \begin{align}\label{eq.11} \begin{cases} {{\pmb x}}_{k} = f({{\pmb x}}_{k1}, {{\pmb u}}_{k1}, {{\pmb w}}_{k1})\\[1mm] {{\pmb z}}_{k} = h({{\pmb x}}_{k}, {{\pmb v}}_k) \end{cases} \end{align} $  (11) 
where
In order to resolve this problem of nonlinearity and nonGaussian, the UKF algorithm based on unscented transformation (UT) is provided
to improve the EKF algorithm. UT adopts sigma symmetric sampling
method, which needs
$ \begin{align} &\chi_{0, k} = \hat{{{\pmb x}}}_k \end{align} $  (12) 
$ \chi_{i, k} = \hat{{{\pmb x}}}_k +\sqrt{(\lambda+n)}(\sqrt{{{\pmb P}}(kk)})_i, \\ \qquad\qquad\qquad\qquad\qquad\quad i=1, 2, \ldots, n $  (13) 
$ \chi_{i, k} =\hat{{{\pmb x}}}_k \sqrt{(\lambda+n)}(\sqrt{({{\pmb P}}(kk)})_i, \\ \qquad\qquad\qquad\qquad\qquad\quad i=n+1, 2, \ldots, 2n $  (14) 
where
First, the basic UKF algorithm is detailed as follows:
1) Initialize:
$ \begin{align*} \hat{{{\pmb x}}}_{0}=E[{{\pmb x}}_{0}], \quad {{\pmb P}}_{0}=E[({{\pmb x}}_{0}\hat{{{\pmb x}}}_{0})({{\pmb x}}_{0}\hat{{{\pmb x}}}_{0})^{T}]. \end{align*} $ 
2) Calculate sigma points:
$ \begin{align} \chi_{k1} = [\chi_{0, k1} ~~\chi_{1, k1}~~ \ldots ~~\chi_{2n, k1}]. \end{align} $  (15) 
3) Time update:
a) Sigma points spread by the system equation are
$ \chi_{kk1} = \ f(\chi_{k1}, k1) $  (16) 
$ \hat{{{\pmb x}}}_{kk1} = \ \sum\limits^{2n}_{i=0} w_{i}^{(m)}\chi_{i, kk1} $  (17) 
$ {{\pmb P}}_{kk1} = \ \sum\limits^{2n}_{i=0}w_{i}^{(c)}[\chi_{i, kk1}\hat{{{\pmb x}}}_{kk1}]\nonumber\\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\ \left[\chi_{i, kk1}\hat{{{\pmb x}}}_{kk1}\right]^{T}+{{\pmb Q}}_{k}. $  (18) 
b) Sigma points propagated by the measurement equation are
$ {{\pmb z}}_{kk1} = h(\chi_{kk1}, k) $  (19) 
$ \hat{{{\pmb z}}}_{kk1} = \sum\limits^{2n}_{i=0} w_{i}^{(m)}{{\pmb z}}_{i, kk1}. $  (20) 
4) Measurement update:
a) Estimate covariance of the prediction measure
$ {{\pmb P}}_{\hat{z}_{k}\hat{z}_{k}} = \ \sum\limits^{2n}_{i=0} w_{i}^{(c)}\left[{{\pmb z}}_{i, kk1}\hat{{{\pmb z}}}_{kk1}\right]\\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\; \, \times \left[{{\pmb z}}_{i, kk1}\hat{{{\pmb z}}}_{kk1}\right]^{T}+{{\pmb R}}_{k} $  (21) 
$ {{\pmb P}}_{\hat{x}_{k}\hat{z}_{k}} = \ \sum\limits^{2n}_{i=0} w_{i}^{(c)}\left[\chi_{i, kk1}\hat{{{\pmb x}}}_{kk1}\right]\\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\; \, \times\left[{{\pmb z}}_{i, kk1}\hat{{{\pmb z}}}_{kk1}\right]^{T} $  (22) 
where
b) Compute filter gain and covariance correction
$ {{\pmb K}}_{k} = {{\pmb P}}_{\hat{x}_{k}\hat{z}_{k}}{{\pmb P}}_{\hat{z}_{k}\hat{z}_{k}}^{1} $  (23) 
$ \hat{{{\pmb x}}}_{k} = \hat{{{\pmb x}}}_{kk1} + {{\pmb K}}_{k}(\textbf{z}_{k}\hat{\textbf{z}}_{kk1}) $  (24) 
$ {{\pmb P}}_{k} = {{\pmb P}}_{kk1}  {{\pmb K}}_{k}{{\pmb P}}_{\hat{z}_{k}\hat{z}_{k}} {{\pmb K}}_{k}^{T}. $  (25) 
Then, the MMAEUKF algorithm is summarized as follows:
1) Initialization:
$ \begin{align*} \hat{{{\pmb x}}}_{0}=E[{{\pmb x}}_{0}], ~~~{{\pmb P}}_{0}=E\left[({{\pmb x}}_{0}\hat{{{\pmb x}}}_{0})({{\pmb x}}_{0}\hat{{{\pmb x}}}_{0})^{T}\right]. \end{align*} $ 
2) Estimate with UKF:
a) Calculation of sigma points;
b) Time update;
c) Measurement update.
3) Calculate likelihood equation using multiple model adaptive estimation
$ \begin{align} L_{0}^{(l)} = &\ \frac{1}{\left\{\det\left[2\pi({{\pmb H}}_{k} \hat{{{\pmb P}}}_{k}^{(l)}{{\pmb H}}_{k}^{T} + {{\pmb R}}_{k})\right]\right\}^{\frac{1}{2}}} \nonumber\\[2mm] &\, \times \exp \left[\frac{1}{2}{\pmb e}_{k}^{(l)T}({{\pmb H}}_{k} \hat{{{\pmb P}}}_{k}^{(l)}{{\pmb H}}_{k}^{T} + {{\pmb R}}_{k})^{1}{\pmb e}_{k}^{(l)}\right]. \end{align} $  (26) 
Evaluate weights and normalize the importance weight to be a normalization constant
$ \varpi^{(l)}_{k} = \, \varpi^{(l)}_{k1}L^{(l)}_{0} $  (27) 
$ \varpi^{(l)}_{k} \leftarrow \frac{\varpi^{(l)}_{k}} {\sum_{j=1}^{M}\varpi^{(j)}_{k}} $  (28) 
for
$ \hat{{\pmb x}}_{k} = \sum\limits_{l=1}^{M}\varpi^{(l)}_{k}\hat{{{\pmb x}}}^{(l)}_{k} $  (29) 
4) Output
$ {\pmb p}_{k} = \sum\limits_{l=1}^{M}\varpi^{(l)}_{k}(\hat{{{\pmb x}}}^{(l)}_{k}  \hat{{\pmb x}}_{k})(\hat{{{\pmb x}}}^{(l)}_{k}  \hat{{\pmb x}}_{k})^{T}. $  (30) 
The movement of the missile can be divided into three phases, including the vertical rising part of the boost stage, the turning part of the boost stage and the final flameout stage [17]. Under the assumption that ballistic missile is a particle and the earth is a fixed sphere, then the nonlinear dynamics model of the missile is shown as
$ \begin{align}\begin{cases} \dot{v} = \frac{PX_{d}}{m}g\sin\theta \\[1mm] \dot{x} = v\cos\theta \\ \dot{y} = v\sin\theta \end{cases} \end{align} $  (31) 
where
$ \begin{align}\label{eq31} \theta = \begin{cases} \theta_{0},&t \leq t_{1} \\[1mm] \frac{(\theta_{0}\theta_{k})(t_{2}t)^{2}}{(t_{2}t_{1})^{2}+ \theta_{k}},&t_{1} < t \leq t_{2} \\[1mm] \theta_{k},&t_{2} < t \leq t_{k}. \end{cases} \end{align} $  (32) 
where
To avert the coordinate transform as a result of the earth's
rotation, the gravity coordinate system at the launch site is turned
into inertial coordinate system. Then, the position of the ballistic
missile is obtained by making use of the velocity's firstorder
integration under the inertial coordinate system at the launch site, and the derivative of velocity is taken to obtain the absolute
acceleration of ballistic missile. The accelerometer's error is
modeled as the firstorder Markov process [18]. Because of the
relative equation (33) between absolute acceleration and apparent
acceleration
$ \begin{align}\label{eq.32} \dot{v}=\alpha + g. \end{align} $  (33) 
The accelerometer output is:
$ \begin{align} \alpha i_{ib}= {{\pmb B}}^{b}_{ib}(\alpha_{ib})={{\pmb B}}^{b}_{ib}(\dot{v}g) \end{align} $  (34) 
where
After the firstorder calculation, the gyroscope gains the attitude
angular rate denoted by
$ \begin{align} wi_{ib}=\left[ \begin{array}{ccc} \sin\psi&0&1 \\ \sin\gamma \cos\psi&\cos\gamma&0 \\ \cos\gamma \sin\psi &\sin\gamma&0 \end{array} \right] \left[ \begin{array}{c} \dot{\varphi} \\ \dot{\psi} \\ \dot{\gamma} \end{array} \right] \end{align} $  (35) 
where
Based on this assumption that the system has additive noise, the mathematical model can be indicated by
$ {{\pmb x}}_{k}= f({{\pmb x}}_{k1})+{{\pmb w}}_{k} $  (36) 
$ {{\pmb y}}_{k}= h({{\pmb x}}_{k})+{{\pmb v}}_{k} $  (37) 
High accuracy and high availability of navigation information are necessary for many nonlinear systems. The integrated navigation system can make the most of each subsystem's information, and further enhances the precision of navigation information. Thus, we combine SINS with CNS to build integrated navigation system that utilizes SINS as the main body. Next, we correct the attitude angle, velocity, and position information in SINS by the information obtained from CNS [19]. Employing improved filter technology, the errors of accelerometer, gyro or other navigation parameters are obtained in the system of combined navigation.
The SINS/CNS integrated navigation system uses the scheme of strapdown that is connecting the carrier with star sensor, which has the advantages of flexible application, low cost, and a simple structure, but the disadvantage of requiring a larger field angle due to the star sensors only working in a dynamic environment [20]. Fig. 2 shows the work process of integrated navigation system. In view of specific force and angular speed that is obtained from missile trajectory generator, the SINS can get the inertial instrument output. Then, the SINS attitude information turns up. Meanwhile, this CNS subsystem obtains the attitude of vehicle through the star sensors' measured information and directly outputs the angular error information. After filtering the angular error, we gain the SINS's attitudeangle error estimate and the final output of this integrated navigation system.
Download:


Fig. 2 The flow diagram of SINS/CNS integrated navigation 
Integrated navigation systems are more accurate than simple system
of inertial navigation because the fixed star's orientation in
the inertial space is essentially constant. Though the precession
of the earth's polar axis, parallax, and the star sensors' aberration
cause small changes in the direction of the stars, they result in
an attitude error of less than 1
In this integrated navigation system, the CNS usually applies the
star sensors to seek starry sky, capturing the stars, recognize star
charts, and then drawing the fixed stars' imaging location. After image
recognition, we eventually calculate the attitude matrix of star
sensor to gain the attitude information [21]. The star sensors
output the attitude angle's observed value under the inertial
coordinate system. In addition,
$ \begin{align} {{\pmb M}}=\left[ \begin{array}{ccc} 0&\cos\varphi &\cos\psi \cos\varphi \\[1mm] 0&\sin\varphi&\cos\psi \cos\varphi \\[1mm] 1&0&\sin\psi \end{array} \right]. \end{align} $  (38) 
And then, the attitude angle is indicated by
$ \begin{align} (\psi, \varphi, \gamma)=F(x, y, \alpha, \beta, f) \end{align} $  (39) 
where (
In this part, we carry out simulations for SINS/CNS integrated navigation system of the missile to verify the filtering performance of this improved algorithm. The performance of four filtering methods (EKF, MMAEEKF, UKF, and MMAEUKF) are also analyzed and compared.
A. Simulation ConditionsIn general, we assume a set of timevarying parameters of process
noise covariance, but in practice, we do not know these
values of the system in advance. Therefore, due to the
unpredictability of parameters, we need to adopt a method to
generate an optimal random distribution. Here, we select two methods
and contrast the simulation result of these two methods that are
producing random numbers. The first method is linear congruence
method to generate uniform distribution; the second is a Hammersley
[22] quasirandom number generation method. Fig. 3 shows a contrast
between the first method and Hammersley quasirandom method, for
randomly producing 250 elements between
Download:


Fig. 3 Comparison of uniform distribution and Hammersley sequence 
In this simulation, the filtering methods include the EKF, MMAEEKF, UKF, and MMAEUKF algorithms. Table Ⅰ shows the initial simulation parameters.
In the system of integrated navigation, the true values of the
parameters
$ \begin{align} {{\pmb Q}}_{k}=\left[ \begin{array}{cc} wgt^{2}\left[ \begin{array}{ccc} 1&0&0 \\ 0&1&0 \\ 0&0&1 \end{array} \right] &O_{3\times3} \\ O_{3\times3}&wga^{2}\left[ \begin{array}{ccc} 1&0&0 \\ 0&1&0 \\ 0&0&1 \end{array} \right] \\ \end{array} \right]. \end{align} $  (40) 
To highlight the advantage of filtering algorithm in this actual timevarying system, the Hammersley point set is used to distribute the parameters for each sub filter when setting the filter's parameters. In MMAE, the filters' number is 10.
For the parameter setting of system, we assume the set of parameters
for the process noise covariance matrix changes over time. That is to
say, the total flight time of this ballistic missile is 300 seconds;
during the first 100 seconds, the parameters of process noise
covariance
This part provides simulation result with filtering performance of navigation system by the simulations, and includes some comparisons of the four filtering algorithms. From Figs. 46, we can compare and analyze the error of yaw angle, pitch angle, and roll angle of the different filtering methods. According to the simulation results, the MMAE method indicates better performance in aspect of accuracy than the traditional EKF and UKF algorithms in information fusion. Furthermore, the MMAE method with UKF has obviously higher accuracy than the one with EKF algorithm in the state estimation process, and the MMAEUKF algorithm gains the minimal attitude angle error.
Download:


Fig. 4 The error of yaw angle 
Download:


Fig. 5 The error of pitch angle 
Download:


Fig. 6 The error of roll angle 
To demonstrate the advantages and disadvantages of these four algorithm in navigation system more clearly, we assess the attitude error navigation performance of these filtering methods by the root mean square error (RMSE). The RMSE is defined by
$ \begin{align} RMSE =\sqrt{\frac{1}{n}\sum\limits_{k=1}^{n}(x_{k}\hat{x}_{k})^{2}} \end{align} $  (41) 
where
Table Ⅱ shows the simulation result with root mean square error of three estimated attitude angle for different algorithms. Thus, we can see that the optimal estimate of the following four algorithms is to employ this improved MMAE method with unscented Kalman filter.
This paper has derived a new filtering approach based on multiple model adaptive estimate method, which is suited for nonlinear systemBallistic Missile Motion. In this method, we design UKF as sub filter in MMAE rather than EKF to obtain the estimate of system state. To indicate the navigation performance of the improved filtering algorithm more clearly, we carry out the simulations, and the behaviors of EKF, MMAEEKF, UKF, and MMAEUKF have also been compared and analyzed in the SINS/CNS integrated navigation system. Based on the simulation results, compared with conventional EKF and UKF methods, MMAEEKF and MMAEUKF algorithms obtain better navigation performance. MMAEUKF method can gain the highest navigation accuracy at the cost of the highest computational cost, and the improved methods can meet the demand of the ballistic missile's movement. With a view that integrated navigation system has some uncertain factors, filtering methods with adaptive parameters are to be studied and some real experiments are to be done in the future.
[1]  H. Y. He, J. N. Xu, F. J. Qin, and F. Li, "Platform bench strapdown algorithminertial navigation system, " in Proc. Chinese Automation Congr., Wuhan, China, 2015, pp. 15501554. http://ieeexplore.ieee.org/document/7382747/ 
[2]  W. R. Wu, D. Y. Wang, and X. L. Ning, Deep Space Probes Autonomous Navigation Principles and Techniques. Beijing, China: China Astronautic Publishing House, 2011. 
[3]  S. Y. Cho, "IMfilter for INS/GPSintegrated navigation system containing lowcost gyros, " IEEE Trans. Aerosp. Electron. Syst., vol. 50, no. 4, pp. 26192629, Oct. 2014. http://dx.doi.org/10.1109/TAES.2014.130128 
[4]  Y. G. Zhang, Y. L. Huang, Z. M. Wu, and N. Li, "A high order unscented Kalman filtering method, " Acta Autom. Sinica, vol. 40, no. 5, pp. 838848, May 2014. http://en.cnki.com.cn/Article_en/CJFDTotalMOTO201405006.htm 
[5]  Q. Z. Yan and M. X. Sun, "Suboptimal learning control for nonlinear systems with both parametric and nonparametric uncertainties, " Acta Autom. Sinica, vol. 41, no. 9, pp. 16591668, Sep. 2015. Suboptimal learning control for nonlinear systems with both parametric and nonparametric uncertainties 
[6]  H. M. T. Menegaz, J. Y. Ishihara, G. A. Borges, and A. N. Vargas, "A systematization of the unscented Kalman filter theory, " IEEE Trans. Autom. Control, vol. 60, no. 10, pp. 25832598, Oct. 2015. http://www.mendeley.com/research/systematizationunscentedkalmanfiltertheory 
[7]  F. Deng, J. Chen, and C. Chen, "Adaptive unscented Kalman filter for parameter and state estimation of nonlinear highspeed objects, " J. Syst. Eng. Electron., vol. 24, no. 4, pp. 655665, Aug. 2013. http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=xtgcydzjse201304013 
[8]  Z. R. Xing and Y. Q. Xia, "Comparison of centralised scaled unscented Kalman filter and extended Kalman filter for multisensor data fusion architectures, " IET Signal Process., vol. 10, no. 4, pp. 359365, Jun. 2016. http://www.growkudos.com/publications/10.1049%252Fietspr.2015.0205 
[9]  R. K. Mehra, "Approaches to adaptive filtering, " IEEE Trans. Autom. Control, vol. 17, no. 5, pp. 693698, Oct. 1972. http://www.mendeley.com/catalog/approachesadaptivefiltering/ 
[10]  W. He, S. Zhang, and S. S. Ge, "Adaptive boundary control of a nonlinear flexible string system, " IEEE Trans. Control Syst. Technol., vol. 22, no. 3, pp. 10881093, May 2014. http://www.scititles.com/articles/995e125ca3b94c1db5211bba9b8afd1a 
[11]  M. Huang, X. Wang, and Z. L. Wang, "Multiple model adaptive control for a class of nonlinear multivariable systems with zeroorder proximity boundedness, " Acta Autom. Sinica, vol. 40, no. 9, pp. 20572065, Sep. 2014. http://en.cnki.com.cn/Article_en/CJFDTotalMOTO201409028.htm 
[12]  J. C. Fang and X. L. Ning, Autonomous Celestial Navigation Method for Deep Space Detector. Shanxi, China: Northwestern Polytechnical University Press, 2010. 
[13]  W. L. Li and Y. M. Jia, "An information theoretic approach to interacting multiple model estimation, " IEEE Trans. Aerosp. Electron. Syst. , vol. 51, no. 3, pp. 18111825, Jul. 2015. http://dx.doi.org/10.1109/TAES.2015.140542 
[14]  H. Peng, F. F. Zhao, S. F. Fan, Z. L. Tang, and W. He, "Multiple model adaptive estimation for the celestial navigation system, " in Proc. 34th Chinese Control Conf. , Hangzhou, China, 2015, pp. 53035308. http://www.mendeley.com/research/multiplemodeladaptiveestimationcelestialnavigationsystem/ 
[15]  F. F. Zhao, G. Q. Zhao, S. F. Fan, Z. L. Tang, and W. He, "Multiple model adaptive estimation algorithm for SINS/CNS integrated navigation system, " in Proc. 34th Chinese Control Conf. , Hangzhou, China, 2015, pp. 52865291. http://www.en.cnki.com.cn/Article_en/CJFDTOTALSKTC201404005.htm 
[16]  V. Hassani, A. M. Pascoal, A. P. Aguiar, and M. Athans, "Multiple model adaptive estimation for open loop unstable plants, " in Proc. 2013 European Control Conf. , Zurich, Switzerland, pp. 16211626. http://www.sintef.no/home/Publications/Publication/?pubid=CRIStin+1057958 
[17]  J. Kim, S. S. Vaddi, P. K. Menon, and E. J. Ohlmeyer, "Comparison between nonlinear filtering techniques for spiraling ballistic missile state estimation, " IEEE Trans. Aerosp. Electron. Syst. , vol. 48, no. 1, pp. 313328, Jan. 2012. Comparison between nonlinear filtering techniques for spiraling ballistic missile state estimation 
[18]  D. Zhang, K. Fu, S. S. Ge, Z. L. Tang, and W. He, "Analysis of filtering methods for the SINS/CNS integrated navigation system of missile motion, " in Proc. 11th World Congr. Intelligent Control and Automation, Shenyang, China, 2014, pp. 38543859. http://dx.doi.org/10.1109/WCICA.2014.7053360 
[19]  Q. T. Gu and J. Fang, "Global optimality for generalized federated filter, " Acta Autom. Sinica, vol. 35, no. 10, pp. 13101316, Oct. 2009. http://dx.doi.org/10.1016/S18741029(08)601123 
[20]  W. Quan, J. C. Fang, F. Xu, and W. Sheng, "Hybrid simulation system study of SINS/CNS integrated navigation, " IEEE Aerosp. Electron. Syst. Mag., vol. 23, no. 2, pp. 1724, Feb. 2008. http://www.mendeley.com/research/hybridsimulationsystemstudysinscnsintegratednavigation/ 
[21]  K. Badshah, Y. Y. Qin, and J. L. Zhang, "SINS/CNS integration algorithm and simulations for extended time flights using linearized Kalman filtering, " in Proc. 2015 IEEE Int. Conf. Communication Software and Networks, Chengdu, China, pp. 3337. http://dx.doi.org/10.1109/ICCSN.2015.7296122 
[22]  J. M. Hammersley, "Monte Carlo methods for solving multivariable problems, " Ann. N. Y. Acad. Sci. , vol. 86, no. 3, pp. 844874, May 1960. http://dx.doi.org/10.1111/j.17496632.1960.tb42846.x 