Application of Improved 5th-Cubature Kalman Filter in Initial Strapdown Inertial Navigation System Alignment for Large Misalignment Angles.

In view of the fact the accuracy of the third-degree Cubature Kalman Filter (CKF) used for initial alignment under large misalignment angle conditions is insufficient, an improved fifth-degree CKF algorithm is proposed in this paper. In order to make full use of the innovation on filtering, the innovation covariance matrix is calculated recursively by an innovative sequence with an exponent fading factor. Then a new adaptive error covariance matrix scaling algorithm is proposed. The Singular Value Decomposition (SVD) method is used for improving the numerical stability of the fifth-degree CKF in this paper. In order to avoid the overshoot caused by excessive scaling of error covariance matrix during the convergence stage, the scaling scheme is terminated when the gradient of azimuth reaches the maximum. The experimental results show that the improved algorithm has better alignment accuracy with large misalignment angles than the traditional algorithm.

1. Introduction

Initial alignment is the basis for obtaining accurate misalignment angles of Strapdown Inertial Navigation Systems (SINS), and its accuracy and efficiency will directly affect the SINS performance [1,2]. In an application environment with large misalignment angles and moving bases, the linear initial alignment algorithm is greatly restricted [3], while nonlinear algorithms have been increasingly used for initial alignment [4,5].

In order to apply the linear Kalman filter to a nonlinear system, the Extended Kalman Filter (EKF) [6] truncates the Taylor series expansion of the nonlinear function. The accuracy of EKF is low for strongly nonlinear systems, and there may be a problem in that the Jacobi matrix cannot be solved [7]. An Unscented Transformation (UT) is used in the Unscented Kalman Filter (UKF) [8] to replace the local linearization for obtaining the first two moments of the states. Although the precision of UKF is higher than EKF, there is a “dimension disaster” problem [9] in high dimensional applications of UKF, and parameter adjustment is needed when using UKF. Spherical Radial Cubature Rule (SRCR) is used in the third-degree CKF (hereinafter referred to as CKF3) to approximate the nonlinear integral and obtain the two moments of the states, which can achieve third-order accuracy [10]. CKF has a better numerical stability than UKF, and has been widely applied [11,12]. In order to further improve the accuracy of the CKF3, fifth-degree CKFs (hereinafter referred to as CKF5) [13,14] have been applied gradually [15]. CKF5 has problems of numerical stability and big calculation [16], and in this sense, CKF5 is more suitable for low dimensional systems. The adaptive Kalman filter is an effective method to restrain the divergence and improve the accuracy of filtering. The main idea of the adaptive Kalman filter is to scale the part of error covariance matrix or system noise matrix by calculating the ratio of matrix ranks of actual updating and filter updating [17,18]. The traditional method for estimating the updating covariance is the open window method, which uses the mean value of historical data and not takes into account the different weights of new and old updates. The accuracy and real-time performance of the traditional adaptive Kalman filter should be improved. Combining EKF and adaptive filtering [19], the update covariance matrix is directly introduced into the calculation of gain matrix, which improves the utilization rate of updates. However, EKF cannot effectively meet the accuracy requirement of the nonlinear initial alignment with large misalignment angles. Adaptive filtering is introduced into Unscented Particle Filter (UPF) [20], and the efficiency of convergence can be improved by partially scaling the error covariance matrix. Nevertheless, the precision of in-motion and nonlinear initial alignment with large misalignment angles should be further verified experimentally. The 5-th CKF is directly applied to the initial alignment [12], but experiments on initial alignment with large misalignment angles are lacking. The singular value decomposition can improve the robustness of the initial alignment algorithm [21], but there is no verification of actual data except for simulations.

In order to improve the efficiency and accuracy of the 5-th CKF in initial alignment situations, an entire scaling of the error covariance matrix method with the ratio of matrix ranks of actual updates and filtered updates is introduced in this paper. To solve the problem of oscillation or divergence of improved 5-th CKF caused by excessive scaling of the error covariance matrix, the scaling is terminated when the gradient of the azimuth error reaches the maximum. The improved algorithm needs to store the historical IMU data, GNSS data and the calculated azimuth misalignment angle data, and then recalculate it for alignment. To improve the stability of the 5-th CKF algorithm, the Singular Value Decomposition (SVD) algorithm is used instead of the Cholesky’s method in the matrix decomposition. The large misalignment angle model is simplified in this paper, which means that only the velocity errors and the misalignment angles are used as the states. The reduction of the states has the advantages of requiring less calculation and better calculation stability, which can reduce the chances of filter divergence.

In this paper, the nonlinear model of initial alignment with large misalignment angles is first analyzed and a simplified alignment model with 6-dimensional state variables is proposed. Secondly, the traditional CKF5 algorithm is briefly described in this paper, and then the improved algorithms based on fading factor for the updating covariance matrix, short-term scaling of the error covariance matrix and SVD decomposition are introduced. Finally, the effectiveness of the improved algorithms is verified by simulation experiments and the dynamic vehicle experiments.

2. Nonlinear Model of Initial Alignment with Large Misalignment Angle

The SINS updating algorithm includes attitude-updating, speed-updating and position-updating, and the corresponding differential equations are expressed with Equations (1) to (3) [22]: where is the direction cosine matrix of the attitude converting from the measurement reference frame (abbreviated as b frame) to the navigation reference frame (abbreviated as n frame), is the rotational angular rate of b frame relative to n frame, i is the inertial reference frame, ( is the antisymmetric matrix of , the vector expresses the eastern velocity, northern velocity and vertical velocity, respectively, is the output of the accelerometer, L is the latitude of the location of carrier, λ is the longitude of the location of carrier, h is the height of the location of carrier, R is the radius of meridian of the location of carrier and R is the radius of prime vertical circle of the location of carrier. Assuming n′ and n are the calculation reference frame and ideal navigation reference frame, respectively, the misalignment (expressed as Euler angle) between the two navigation reference frames can be expressed as: where is the relative rotational angular rate between the n reference frame and n′ reference frame, ε is the random constant drift of Fiber Optic Gyroscope (FOG) and is the random drift of FOG. and can be expressed as:

Under the reference frame, Equation (2) can be written as:

Subtracting Equation (2) from Equation (8) and neglecting the calculation error of gravitational acceleration and higher-order terms, we obtain: where , and we have:

The purpose of initial alignment is to quickly identify the misalignment angles. It is difficult to estimate the errors of inertial sensors accurately, and that is not the purpose of initial alignment, then the errors of ε and ∇ can be regarded as system noises. Therefore, the states used in initial alignment are misalignment angles φ and the velocity error δν, which can reduce the dimensions of the filter, and thus reduce the computation and improve the numerical stability of CKF5. In order to simplify the filtering error equation, a SINS nonlinear error model which can reflect the characteristics of SINS errors has important practical value in the initial alignment [20]. When the SINS and GNSS are integrated for initial alignment under moving base conditions, GNSS can provide relatively accurate velocity and location information, then Equations (4) and (9) can be simplified as: where is the random drift of accelerometer, is the angular velocity of the Earth determined by the information from GNSS. With the difference between the velocity of SINS measured by GNSS and the calculated velocity working as the measurement of z, the nonlinear filtering equations of initial alignment with large misalignment angles can be expressed as: where x = [ ()], , H = [03×3 3×3], is the error vector of measurement and the function of f can be obtained from Equation (12). The model of initial alignment (Equations (12) and (13)) takes into account of the moving base situation, so the estimation of the initial alignment angles will theoretically not be affected by the line movement and angle movement of SINS.

3. CKF5 and Improved CKF5

3.1. CKF5 Algorithm

For the integral of the global nonlinear function (f) and variable x = ry in the spherical-radial coordinate system, (f) can be written as: where δ = 1, r ∈ [0, ∞), U is n-dimensional unit hypersphere and σ(g) is the element of U. (f) can be divided into two parts, spherical integral S(r) and radial integral R:

The numerical integration is approximated with the summation of weighted sigma points in the algorithm of CKF, then Equation (16) can be expressed as:

Then:

The m-order approximation of S(r) can be written as: where d = [p1,p2,…,p], , , and Γ(g) is the gamma function. Note that when m = 5, the value of S is calculated with and . The points set on the U is defined as: where and are defined as: where is the ith column of identity matrix . Then the Equation (19) can be expressed as: where , and . For the symmetry, and are needed to calculated in the solution of 5-order integral of . Then

Let , then we have:

Then it follows from Equations (19) and (23) that the integrated result of Equation (18) can be obtained. The algorithm flow of CKF5 is shown as follows:

Prediction phase:

Let be the estimated state at time and decompose with Cholesky’s method:

Compute cubature points: where and .

Compute the predicted mean and predicted covariance : where , , and is system noise covariance matrix.

Update phase:

Decompose with Cholesky’s method:

Compute cubature points:

Compute the estimated measurement , the gain matrix and the covariance matrixes of , :

Update and : where is the measurement at time .

3.2. Improved CKF5 Algorithm

The filtering precision of CKF5 can reach the fourth order, and the accuracy of CKF3 can reach the second order. The amount of sigma points of CKF5 is , and that of CKF3 is , so the cost of computing of CKF3 is less than that of CKF5. In order to estimate the numerical stability of nonlinear filter, the stability factor can be defined as [23]: where is the weight of sigma point and is the amount of sigma points. If , the numerical integration has absolute stability theoretically, on the contrary, non-absolute stability for . When , the cumulative truncation error may leads to a non-positive error covariance matrix, even further to the error of calculation or divergence. For CKF3, and for CKF5, , which means that CKF3 has better numerical stability than CKF5. For CKF5, we have:

According to Equation (14), and , the robustness of CKF is poor in this situation. The non-positive covariance matrix caused by strong nonlinear or data errors may bring many problems in Kalman filters, such as filtering error or divergence, which also can lead to the error of Cholesky decomposition in CKF. In order to improve the numerical stability of CKF5, and considering that the singular value decomposition (SVD) has better numerical robustness than Cholesky [24], SVD is introduced into CKF5.

For s linear optimal Kalman filter, when the gain matrix is optimal, the innovation sequence should be orthogonal. This conclusion can be extended to nonlinear filtering. The adaptive Kalman (hereinafter referred to as AKF) scales the error covariance matrix or system noise matrix by calculating the ratio of matrix ranks of actual innovation and filtered innovation [17,18]. Essentially, the traditional AKF only partly scale the error covariance matrix in the next filtering period, and the utilization rate of innovation is low. For the innovation sequence, the contribution of the new innovation to the state estimation is large than the old one. Therefore, ICKF5 uses the exponent fading factor to weigh the innovation sequence, and then scale the entire error covariance matrix in the next filtering period. Assuming the current filter period is and the weight set of innovation sequence is , we have: where is fading factor and , empirically, . To illustrate the impact of on , was taken as 0.7, 0.8, and 0.9, respectively, for the convenience of observation, the innovations are weighted by the exponent fading factor from . The distribution of is shown in Figure 1.

Figure 1

The distribution of weights of innovations under different .

Figure 1 shows that if the is big, the weight of the current innovation is small, correspondingly, the weights of old innovations are big, which means that the old innovation has bigger influence than the new innovation on calculating the covariance matrix of innovation sequence. The opposite conclusion can be obtained if is small. A recursive formula for the covariance matrix of the innovation sequence is given below. When the forgetting factor is determined, then:

It follows from Equations (39) and (40) that:

The open window method is used in traditional AKF, the covariance matrix of innovation at time can be expressed as:

The weight of in Equation (42) can be replaced by the weights of in Equation (41), then a recursive formula between and can be expressed as:

In order to improve the feedback effect of innovation, we scale the error covariance matrix of the next filtering period times entirely. The increase of the error covariance matrix makes the filter gain increase and the error correction ability of the innovation in the next filtering period is strengthened. This method can reduce the dependence of the filter on the model, and then improve the filtering accuracy and reduce the risk of divergence. is defined as: where is the operator of obtaining the rank of a matrix, is the measurement error covariance matrix at time , is the measurement matrix at time . The error covariance matrix of next filtering period is scaled as:

The algorithm of ICKF5 is summarized as follows:

Prediction phase:

Input the initial values for filtering: .

Compute cubature points according to Equation (26).

Compute the predicted mean and predicted covariance according to Equations (27) and (28), respectively.

Scaling the entire : .

Update phase:

Decompose with the SVD method.

Compute cubature points according to Equation (30).

Compute the estimated measurement , the gain matrix and the covariance matrixes of , according to Equations (31)–(34), respectively.

Update , and according to Equations (35), (36) and (44), respectively.

The scaling of the error covariance matrix by ICKF5 can effectively improve the utilization of innovation and speed up the convergence of initial alignment. ICKF5 is a kind of suboptimal algorithm, if the scaling value of is too large, the system will be in a state of overshoot, which may lead to oscillation or divergence. For example, when getting the initial alignment in the convergence stage, the gain and the error covariance matrix should ideally remain unchanged or little changed. If the feedback for error covariance matrix is big at the stage of alignment convergence, it may cause the filtering convergence to be unstable or lead to the divergence. In order to reduce the possibility of divergence of initial alignment filtering and ensure the stationarity after convergence, we use a short-term feedback strategy that we stop the scaling of in ICKF5 before the filtering occurs oscillating or diverging, which can ensure the convergence of initial alignment. For convenience, we call the short-term feedback algorithm IICKF5. The first step of IICKF5 is that the ICKF5 is used to complete the initial alignment firstly, then the gradient sequence of the set of azimuth misalignment angle outputted by ICKF5 can be obtained, where is gradient operator. The second step, the time corresponding to the maximum value of the gradient is chosen as the time to terminate the feedback to the error covariance matrix . can be expressed as:

To be brief, ICKF5 is used in the period of and CKF5 is used in , where is the end time of initial alignment filtering. The flow chart of IICKF5 is shown in Figure 2.

Figure 2

The flow chart of IICKF5.

Compared to ICKF5, IICKF5 needs to recalculate the stored IMU data and the misalignment angle set {x} output by ICKF5. IICKF5 has a certain delay compared with ICKF5, and some storage space is needed.

3.3. Simulation and Its Analysis

In order to verify the effectiveness of ICKF5 and IICKF5, a simulation of initial alignment under static base conditions is performed. The simulation conditions are set as follows: the gyro constant drift is , random drift is ; the accelerometer zero bias is , random deviation is ; the local geographic latitude is , longitude is ; the misalignment angles are set as and the simulation period is . The general parameters of the initial alignment algorithms are set with the same values. The initial value of is set as . The matrix of is determined according to the random drift of IMU, and . The initial error covariance matrix is set as . Where is the operator of obtaining the diagonal matrix with diagonal elements. For UKF, , and . EKF, UKF [25], CKF3, CKF5, ICKF5 and IICKF5 algorithms are applied to the nonlinear initial alignment. In addition to IICKF5, the calculated misalignment angles are as shown in Figure 3.

Figure 3

The estimated misalignment angles (): (a) Pitching misalignment angle; (b) rolling misalignment angle; (c) azimuth misalignment angle.

Figure 3 shows that the pitching misalignment angles converge to 20° and the rolling misalignment angles converge to −40° with the exception of EKF, and only the azimuth misalignment angles of CKF5 and ICKF5 converge to 170°.

The gradient values of the azimuth angle errors of ICKF5 are shown in Figure 4.

Figure 4

The gradient values of azimuth angle error.

Figure 4 shows that there is an obvious mutation in the azimuth angle error at , which is the overshoot caused by ICKF5, so The estimated azimuth angles of CKF5, ICKF5 and IICKF5 are shown in Figure 5.

Figure 5

The estimated azimuth misalignment angles of CKF5, ICKF5 and IICKF5.

As seen in Figure 5, IICKF5 enters the convergence stage at about 400 s, while CKF5 and ICKF5 converge slower than IICKF5. Figure 3 also shows that the convergences of horizontal misalignment angles are faster than that of the azimuth misalignment angles, which is due to the fact that the observability of the horizontal misalignment angles is high under the velocity matching mode. There is a coupling relationship between the horizontal misalignment angles and the bias of accelerometer. When the variation of heading angle is small, the accuracies of the horizontal misalignment angles depend on the bias of the accelerometer. The simulation experiment in this paper is based on static base, so the variation of heading angle is very small, and the horizontal misalignment has a slow drift due to the bias of the accelerometer. From this point of view, the estimation period of the horizontal misalignment angles is not the longer the better. For the sake of accuracy, we select the estimations of horizontal misalignment angles at 50 s and select the estimation value of azimuth misalignment angle at the last moment (800 s) as the outputs of alignment, and the alignment errors are shown in Table 1.

Table 1

Comparison of alignment errors (°).

	EKF	UKF	CKF3	CKF5	ICKF5	IICKF5
error of pitch (°)	−7.0205	0.0899	0.1270	0.1590	0.0829	0.0291
error of roll (°)	7.8234	−0.1464	−0.1634	−0.2078	−0.2949	−0.2738
error of heading (°)	62.7783	157.7361	61.2902	13.6403	−1.2458	−0.0878

Table 1 shows that the pitch and the azimuth misalignment errors of IICKF5 are the smallest, and the azimuth misalignment error is close to the limitd alignment precision of 0.05°, which meets the initial alignment requirement.

4. Dynamic Vehicle Experiment

The experimental equipment is shown in Figure 6. The GNSS real-time differential system consisted of two model ProPak6 receivers (NovAtel, Calgary, AB, Canada). The attitude reference system is a NovAtel SPAN system (NovAtel, Calgary, AB, Canada), which can reach centimeter level precision. An LN-200 IMU named from the Novatel Company is used for the initial alignment experiment. The bias of the accelerometer is 0.3 mg and its scale factor is 300 ppm. The constant drift, the random drift and the scale factor of gyro are , and 100 ppm, respectively. The output rates of the LN-200, SPAN system and GNSS system are 200 Hz, 100 Hz and 5 Hz, respectively.

Figure 6

The field test setup.

The period of initial alignment is 700 s. In addition to IICKF5, the calculated azimuth angle errors are as shown in Figure 7.

Figure 7

The alignment errors of large azimuth misalignment angle of dynamic vehicle experiment: (a) Pitch errors; (b) roll errors; (c) heading errors.

For the azimuth angle errors, Figure 7 shows that EKF, UKF, CKF3 did not converge to the accurate azimuth angle, ICKF5 has faster convergence rate than CKF5, and ICKF5 has much more oscillations in the convergence process than CKF5 because of the overshoot. For the pitch angle and the roll angle errors, we can see that CKF5 and ICKF5 have greater oscillations than other algorithms, which is caused by the inaccuracy of model noises, such as the GPS’s measurement noise, IMU noise and so on. The nonlinear filter of high order precision amplifies the noises with high frequencies, which easily results in larger estimation errors if the estimated state values are small. When obtaining the gradient of the azimuth angles of ICKF5 in Figure 7c, the gradient values are as shown in Figure 8.

Figure 8

The gradient values of azimuth angle error (dynamic vehicle experiment).

Figure 8 shows that there is a mutation in the azimuth angle error at , so The estimated misalignment angles of all the algorithms are shown in Figure 9.

Figure 9

The alignment errors of different algorithms: (a) Pitch errors; (b) roll errors; (c) heading errors.

From Figure 9 we can find that IICKF5 has a faster convergence rate than CKF5 and ICKF5 except for the pitching misalignment angle, and the azimuth misalignment angle of IICKF5 has better convergence stability than ICKF5. The alignment errors at the end of alignment (at 700 s) are shown in Table 2.

Table 2

Comparison of alignment errors with in-motion base (°).

	EKF	UKF	CKF3	CKF5	ICKF5	IICKF5
error of pitch (°)	−3.2814	−1.2868	−2.6771	−3.7189	3.0270	5.9149
error of roll (°)	−1.2960	−0.2867	3.3751	5.4245	−2.9030	1.8357
error of heading (°)	−139.5622	−266.3212	−204.3457	7.9986	42.5004	−1.9540

Table 2 shows that IICKF5 has better accuravy than CKF5 in estimating the rolling misalignment angle and azimuth misalignment angle. When the alignment is finished, the alignment errors of ICKF5 and IICKF5 did not converge to small values (e.g., 1°), so more time is needed for fine alignment. In order to further compare the convergence abilities of different algorithms, the mean and standard deviation of azimuth misalignment errors are calculated in the range of [100 s, 700 s], which can be regarded as a relatively complete convergence period for ICKF5 and IICKF5, and the results are shown in Table 3.

Table 3

Mean and Standard Deviation (SD) of azimuth misalignment errors (°).

	EKF	UKF	CKF3	CKF5	ICKF5	IICKF5
Mean (°)	−155.2618	−264.1352	−201.1962	1.6920	−6.1625	−0.6852
SD (°)	38.8935	3.8062	11.5359	10.2307	41.2056	4.6229

From Table 3 we can see that IICKF5 has the smallest errors of mean and standard deviation. In this sense, the improved algorithm is suitable for quickly identifying a large azimuth misalignment angle in a complex alignment environment. ICKF5 is a kind of real-time algorithm that represents an improvement on CKF5, which does not need to process historical IMU data and no storage space for the IMU is needed. IICKF5 needs to obtain the gradient value of the historical azimuth angle, so the historical IMU data need to be stored for the second initial alignment filtering. For the 700 s period of in-motion initial alignment, the length of IMU data is 70,000 (in order to synchronize with the reference system), and the storage space for IMU data is 2869 KB. The data stored are double and floating types with 64 bits. The update rate of GNSS is 5 Hz, so the storage space it needs is small, and the size is 10 KB. The capacity of the azimuth misalignment angle data is 6 KB. The total storage space needed for IICKF5 is 2885 KB, and the additional burden on the navigation computer is small.

5. Conclusions

It is not only convenient to calculate the innovation variance matrix in a recursive and exponent fading factor way, but also to effectively use the current and old updates. The speed and precision of convergence can be effectively improved by scaling error covariance matrix with the ratio of matrix ranks of actual updates and filter updates. In order to solve the problem of oscillation or divergence caused by excessive scaling of the error covariance matrix, the scaling scheme is terminated when the azimuth gradient reaches the maximum. The improved CKF can improve the precision of initial alignment better than the traditional CKF. The length of filtering time is also dependent on the experience, and the robustness of the improved algorithm is also worthy of further study.