1. INTRODUCTION
The working principle of transfer alignment is to utilise relatively high-accuracy reference information to calibrate the Slave Inertial Navigation System (SINS) (Groves, Reference Groves2008; Titterton and Weston, Reference Titterton and Weston2004; Chattaraj et al., Reference Chattaraj, Mukherjee and Chaudhuri2013). The reference information is provided by the Master Inertial Navigation System (MINS). This alignment technique can fulfil three goals as follows: in-flight alignment of SINS; online calibration of main error parameters of SINS and related aiding parameters estimation, such as lever arm and time delay. Once the aiding parameters are compensated, the other two goals will be presented together, namely that, in order to improve the alignment accuracy, the main error parameters of inertial sensors should be calibrated; if the main error parameters of inertial sensors are calibrated, the alignment accuracy is certain to be improved.
The traditional transfer alignment methods are generally Kalman filtering techniques. For Kalman filtering-based methods, different measurement matching patterns can be selected, which process the difference between MINS and SINS. The inertial measurement matching pattern converges faster and has less computation, compared with velocity and attitude matching patterns. However, simultaneously more prior information is needed for initialisation which is hard to obtain; more theoretical derivations and practical engineering tests are required to guarantee the correctness of the Kalman filtering-based models. If the prior information is not accurate enough, the filtering will converge to a false “true value” and may even diverge, caused by the incorrectness of the measurement model and the inconformity between sensor error models and practical tests.
The implementation of linear Kalman filtering techniques should satisfy certain conditions, such as small attitude error. Other than filtering techniques, Wu et al. (Reference Wu, Wu, Hu and Hu2011) proposed the Optimisation-Based Alignment (OBA) method, which transforms the alignment problem into an attitude optimisation problem by transforming the velocity differential equation, and later the OBA method is applied to initial alignment and parameter identification (Wu et al., Reference Wu, Wu, Hu and Hu2011; Silson, Reference Silson2011; Kang et al., Reference Kang, Fang and Wang2012; Wu and Pan, Reference Wu and Pan2013; Wu et al., Reference Wu, Wang and Hu2014). Inspired by the OBA method and aiming to restrain the influences of both the random noise and bias, fixed integral interval sliding methods for optimisation-based alignment were proposed (Li et al., Reference Li, Xu, Chang and Zha2014; Chang et al., Reference Chang, Li and Chen2015a; Reference Chang, Hu and Li2015b). The length of the integrated data is supposed to be a key factor due to the ambivalent effect of integrating the sensor noise.
Inspired by the idea of the aforementioned optimisation-based methods, a new optimisation-based transfer alignment and calibration method is proposed for the inertial measurement vector integration matching in this paper, which has not been previously researched. In the proposed method, an inertial measurement vector integration is used instead of inertial measurement vector since using the inertial measurement vector directly may cause divergence of estimations due to noise coverage. In addition, proper selection of the data integration segment is critical for guaranteeing the error observability. If the integration time interval is so long that the dynamic trajectory characteristics will be lost, the estimations will also diverge due to the decrease of degree of observability. For example, during level flight, including a short time of pitching manoeuvres, if long-time integration matching is selected, the characteristics of pitching will be lost in the level flight data, resulting in the decrease of observability degree.
Compared with filtering techniques, no prior information of initial attitude and sensor errors is needed. Not only the initial attitude but also main sensor error parameters are estimated simultaneously by solving a nonlinear optimisation problem. The validity of the proposed method has been evaluated through simulation.
2. PRINCIPLE OF THE PROPOSED METHOD
We denote the MINS body frame with m, the SINS body frame with s, and the inertial frame with i. The inertial measurements of MINS and SINS are obtained in their own body frames, namely m and s frames. The difference of the corresponding inertial measurements results from the misalignment between m and s frames. Simultaneously, SINS has less accuracy compared with MINS, thus the sensor errors are also one main factor for the difference of the corresponding inertial measurements. In view of the difference between MINS and SINS, the misalignment and sensor error parameters can be fused into the inertial measurement vector integration matching equation and the corresponding algorithm is proposed as follows.
2.1. Relationship of Inertial Measurement Vectors between MINS and SINS
The relationship between accelerometer and gyroscope outputs of MINS and SINS is expressed as
where
$C_s^m $
is the transformation matrix from s frame to m frame, f
m
and
$\omega _{im}^m $
the accelerometer and gyroscope outputs of MINS in m frame,
$\tilde f^s $
and
$\tilde \omega _{is}^s $
the corresponding ideal outputs of SINS in the s frame. Compared with SINS, MINS has higher measurement accuracy by at least one order of magnitude and inertial sensor errors of MINS are ignored here. For the pre-calibrated Inertial Measurement Unit (IMU), the lever arm and installation errors change slowly compared with other error terms. They are considered as constant and are compensated before the alignment. Therefore the error models of SINS can be expressed as
where f
s
and
$\omega _{is}^s $
are the error-contaminated outputs of SINS, aSF and gSF are the scale factor of accelerometer and gyroscope of SINS, ba and bg are the corresponding constant bias. aSF and gSF are matrices as follows
$$aSF = \left[ {\matrix{ {aSF_x} & 0 & 0 \cr 0 & {aSF_y} & 0 \cr 0 & 0 & {aSF_z} \cr}} \right] \quad gSF = \left[ {\matrix{ {gSF_x} & 0 & 0 \cr 0 & {gSF_y} & 0 \cr 0 & 0 & {gSF_z} \cr}} \right]$$
In practical situations, the ideal outputs
$\tilde f^s $
and
$\tilde \omega _{is}^s $
are not available. However, their coefficients aSF and gSF are small quantities. They can be substituted by f
s
and
$\omega _{is}^s $
$$\eqalign{\tilde f^s \times \tilde \omega _{is}^s & = (\,f^s - aSF \times \tilde f^s - ba) \times (\omega _{is}^s - gSF \times \tilde \omega _{is}^s - bg) \cr & \approx (\,f^s - aSF \times f^s - ba) \times (\omega _{is}^s - gSF \times \omega _{is}^s - bg)} $$
Expanding Equation (5) further
$$\eqalign{& (\,f^s - aSF \times f^s - ba) \times (\omega _{is}^s - gSF \times \omega _{is}^s - bg) \cr &\qquad = f^s \times \omega _{is}^s - f^s \times gSF \times \omega _{is}^s - f^s \times bg - aSF \times f^s \times \omega _{is}^s + ba \times bg \cr & \quad\qquad{\rm +} aSF \times f^s \times gSF \times \omega _{is}^s - aSF \times f^s \times bg + \omega _{is}^s \times ba + ba \times gSF \times \omega _{is}^s} $$
Turning to the terms on the right side of Equation (6), ba × bg, aSF × f s × gSF × ω s is ,aSF × f s × bg and ba × gSF × ω s is are second order small quantities and can be omitted. Equation (6) is rewritten as
Then Equation (1) is rewritten
The misalignment and sensor error parameters of SINS are fused into the inertial measurement vector matching Equation (8). Considering the effect of noise, integration is done for the matching equation Equation (8).
2.2. Integration and Problem Transformation
The current time t N = NT (N = 0, 1, 2…), where T is the sampling time interval. Integrating Equation (8) in time interval [0, t N ], we can obtain the equalities
$$\eqalign{ \rho (t_N ) = C_s^m (\mu (t_N ) + \gamma (t_N )bg + \chi (t_N )ba + S_1 (t_N )\mathop {aSF}\limits^ \to + S_2 (t_N )\mathop {gSF}\limits^ \to ) \cr & N = 0,1,2,3...} $$
where
$$S_1 (t_N ) = \;\int_0^{t_N} {(\omega _{is}^s (\tau ) \times )\left[ {\matrix{ {\,f_x^s (\tau )} & 0 & 0 \cr 0 & {\,f_y^s (\tau )} & 0 \cr 0 & 0 & {\,f_z^s (\tau )} \cr}} \right]} d\tau $$
$$S_2 (t_N ) = - \!\int_0^{t_N} {(\,f^s (\tau ) \times )\left[ {\matrix{ {\omega _{isx}^s (\tau )} & 0 & 0 \cr 0 & {\omega _{isy}^s (\tau )} & 0 \cr 0 & 0 & {\omega _{isz}^s (\tau )} \cr}} \right]} d\tau $$
$$\mathop {aSF}\limits^ \to = \left[ {\matrix{ {aSF_x} \cr {aSF_y} \cr {aSF_z} \cr}} \right] \quad \mathop {gSF}\limits^ \to = \left[ {\matrix{ {gSF_x} \cr {gSF_y} \cr {gSF_z} \cr}} \right] $$
where
$\left( {a \times} \right)$
is a skew symmetric matrix formed from the elements of vector a.
In order to eliminate the different error characteristics caused by the different integration time interval, we construct new equalities by even-interval time difference
$$\eqalign{ & \rho (t_N ) - \rho (t_{N - \nabla} ) = C_s^m [(\mu (t_N ) - \mu (t_{N - \nabla} )) + (\gamma (t_N ) - \gamma (t_{N - \nabla} ))bg + (\chi (t_N ) - \chi (t_{N - \nabla} ))ba \cr & \quad\qquad\qquad\qquad\qquad\qquad\qquad+ (S_1 (t_N ) - S_1 (t_{N - \nabla} ))\mathop {aSF}\limits^ \to + (S_2 (t_N ) - S_2 (t_{N - \nabla} ))\mathop {gSF}\limits^ \to ) \cr & \qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad N = \nabla, \nabla + 1,\nabla + 2,\nabla + 3...} $$
Define
$$\eqalign{& \rho \;\triangleq \;\rho (t_N ) - \rho (t_{N - \nabla} ) \quad \mu \;\triangleq \;\mu (t_N ) - \mu (t_{N - \nabla} ) \quad \gamma \;\triangleq \;\gamma (t_N ) - \gamma (t_{N - \nabla} ) \cr & \chi \;\triangleq \;\chi (t_N ) - \chi (t_{N - \nabla} ) \quad S_1 \;\triangleq \;S_1 (t_N ) - S_1 (t_{N - \nabla} ) \quad S_2 \;\triangleq \;S_2 (t_N ) - S_2 (t_{N - \nabla} )} $$
Then Equation (17) can be expressed as
Equation (19) is used to form the objective function of the optimisation problem. The time difference interval
$(\nabla + 1)T$
should be adjusted according to the practical conditions. Proper selection of the time difference interval is critical for guaranteeing the observability degree of the system and inertial sensor error parameters. If the integration time interval is so long that the dynamic trajectory characteristics will be lost, the estimations will also diverge due to the decrease of degree of observability.
We aim to estimate the transformation matrix
$C_s^m $
, the accelerometer scale factor aSF and bias ba, the gyroscope scale factor gSF and bias bg, which are treated as constant during the alignment process. They are connected together by Equation (19). We adopt the Newton-Lagrange method, in Wu et al. (Reference Wu, Wang and Hu2014) to obtain the optimal estimations. The solution procedure is same as that in Wu et al. (Reference Wu, Wang and Hu2014), and is listed below for convenience.
Let
$q = [\matrix{ s & \eta \cr} ]^T $
be the unit quaternion corresponding to the transformation matrix
$C_s^m $
, where s is the scale part and η is the vector part. Replace
$C_s^m $
with q, then Equation (19) can be expressed as
where the superscript * denotes the quaternion conjugate operator, and ° denotes the quaternion multiplication operator. Define operators as follows
Therefore, we can obtain from Equation (19)
The problem of alignment and calibration, namely estimating the constant parameters, is transformed to a typical nonlinear optimisation, subject to an equality constraint.
The Newton-Lagrange method is applied to solve the problem. The Lagrange function for Equation (22) is
where
$X = [q^T \quad{\rm } bg^T \quad{\rm } ba^T {\rm } \quad\mathop {aSF}\limits^ \to {\rm } \quad\mathop {gSF}\limits^ \to ]^T $
, λ is the Lagrange multiplier. The iterative step starting from
$(X_k, \lambda _k )$
is given by
and
$(\Delta X,\Delta \lambda )$
is calculated from the following linear equations
Where
Where J is the Jacobian matrix of K T K and H is the Hessian matrix of K T K. When the starting point X 0 is close enough to the true value, the iteration is quadratically convergent. In order to get a relatively good initial value q 0 for X 0, namely attitude, we ignored the sensor errors in Equation (22). Equation (23) is reduced to the pure attitude alignment (Wu et al., Reference Wu, Wu, Hu and Hu2011).
$$\mathop {\min} \limits_q \sum\limits_N {q^T \left( {\left[ \matrix{ - \hfill \cr \mu \hfill} \right] - \left[ \matrix{ + \hfill \cr \rho \hfill} \right]} \right)^T \left( {\left[ \matrix{ - \hfill \cr \mu \hfill} \right] - \left[ \matrix{ + \hfill \cr \rho \hfill} \right]} \right)q} \quad {\rm s}{\rm. t}. \quad q^T q = 1$$
It is a typical Wahba problem which can be solved by q-method and the estimation result is good enough to serve as the initial attitude.
The corresponding optimisation procedure is presented as follows:
For the current time t N ,
-
1) compute ρ, μ, γ, χ, S 1, S 2 according to Equations (10)–(15);
-
2) set k = 0 and start with a good initial point
$\left( {X_k, \lambda _k} \right)$
using Equation (31); -
3) compute
$\nabla _X L(X_k, \lambda _k )$
,
$\nabla h(X_k )$
and
$\nabla _{XX}^2 L(X_k, \lambda _k )$
according to Equations (27)–(30); -
4) solve
$\left( {\Delta X,\Delta \lambda} \right)$
according to Equation (26); -
5) update
$\left( {X_k, \lambda _k} \right)$
with Equation (25); -
6) set k = k + 1 and go to step 3) until convergence;
-
7) obtain the current attitude.
2.3. The Compared Kalman Filter
In order to validate the performance of the proposed method, the Kalman method is selected for comparison. For the compared Kalman method, choose velocity error δv x , δv y , δv z , attitude error φ x , φ y , φ z , misalignment ξ x , ξ y , ξ z , gyroscope bias bg x , bg y , bg z , gyroscope scale factor δKg x , δKg y , δKg z , accelerometer bias ba x , ba y , ba z , accelerometer scale factor δKa x , δKa y , δKa z as the state vector of the system:
$$\eqalign{& X = [\delta v_x \quad \delta v_y \quad \delta v_z \quad \phi _x \quad \phi _y \quad \phi _z \quad bg_x \quad bg_y \quad bg_z \quad \delta Kg_x \quad \delta Kg_y \quad \delta Kg_z \cr & \qquad\quad\qquad ba_x \quad ba_y \quad ba_z \quad \delta Ka_x \quad \delta Ka_y \quad \delta Ka_z \quad \xi _x \quad \xi _y \quad \xi _z ]^T} $$
Then the system state equation is constructed as
Where W is the system noise, and F is the system transfer matrix.
Take the differences of velocity and attitude between MINS and SINS as the measurement vector. The measurement matrix can be expressed as
where C is the direct cosine matrix (DCM) of SINS.
3. SIMULATION
In order to validate the performance of the proposed method, we simulate the trajectory of a missile vertically launched from the ground lasting for 350 s. 0–60 s is the overload period and 60–350 s is the free flight period. We choose the overload period as the alignment process and compare the proposed method with a standard Kalman filter technique with velocity and attitude matching pattern, the filter time lasting for 350 s. For the body frame, x directs the right of the missile body, y the front, z the upward and vertical. The accelerometer and gyroscope outputs of the MINS are shown in Figures 1 and 2 and the sampling rate is 50 Hz. The time difference interval is set as 1 s according to the dynamic trajectory characteristics. The misalignment between the MINS and SINS is set as (2°, 3°, 4°) and initial position (110·32° 47·80° 953·73 m). The inertial sensor error parameters of SINS are: 1) gyroscope bias 2°/h, scale factor 1000 ppm, noise 0·2°/h/
$\sqrt Hz $
; 2) accelerometer bias 1000 ug, scale factor 1000 ppm, noise 100 ug/
$\sqrt Hz $
. Both the methods are evaluated by 50 Monte Carlo runs. The results are shown in Figures 3–7 and Tables 1 and 2.
-
1) For the misalignment, both the methods converge before 60 s and the proposed method converges faster. The estimations of the proposed method and Kalman filtering technique stabilise at (2·0002° ± 0·0004°, 3·0001° ± 0·0001°, 4·0001° ± 0·0032°) and (2·0143° ± 0·0026°, 3·0376° ± 0·0098°, 3·9743° ± 0·0024°) at 60 s. The proposed method is more accurate.
-
2) For the accelerometer bias, both the methods are effective and the estimations of the proposed method and the Kalman filter technique are (979·734 ± 46·421, 972·958 ± 47·621, 990·414 ± 29·744) ug and (899·952 ± 6·039, 907·917 ± 7·450, 1075·860 ± 4·125) ug respectively. However, the proposed method converges at 45 s and the Kalman filter technique at 280 s. Obviously, the proposed method has higher accuracy and faster convergence for the misalignment and accelerometer bias.
-
3) For the accelerometer scale factors, the proposed method converges at 45 s and the estimation is (1005·997 ± 96·424, 1043·771 ± 77·556, 913·253 ± 176·391) ppm. The Kalman filter technique dose not converge until 280 s and is only effective on the y-direction scale factor error, and the estimation is (−5·7696 ± 3·909, 1005·2471 ± 5·7982, −40·9304 ± 26·6239) ppm.
-
4) For the gyroscope scale factors, the estimations of the proposed method and Kalman filtering technique are (965·892 ± 61·332, 1047·036 ± 69·557, 950·047 ± 80·881) ppm and (294·3267 ± 260·5873, 275·5503 ± 187·2138, 305·8104 ± 242·1307) ppm respectively. Similar to the accelerometer scale factors, the proposed method also converges at 45 s; although the Kalman filter technique converges, not all the estimations are meaningful.
-
5) For the gyroscope bias, it is not effectively estimated by either of the two methods. Only two direction biases converge for the proposed method and the accuracy is not good enough (2·2623° ± 1·2293°, −0·3487° ± 13·1622°, 2·1416° ± 0·6595°); although all directions converge for the Kalman filter technique, the accuracy is not satisfactory (1·2488° ± 0·3276°, 1·4892° ± 0·2636°, 0·9718 ± 0·1649°).
Accelerometer outputs of the MINS.
Gyroscope outputs of the MINS.
Misalignment estimate across 50 Monte Carlo runs. a) optimisation-based method b) Kalman filter method.
Accelerometer scale factor estimate across 50 Monte Carlo runs. a) optimisation-based method b) Kalman filter method.
Gyroscope scale factor estimate across 50 Monte Carlo runs. a) optimisation-based method b) Kalman filter method.
Accelerometer bias estimate across 50 Monte Carlo runs. a) optimisation-based method b) Kalman filter method.
Gyroscope bias estimate across 50 Monte Carlo runs. a) optimisation-based method b) Kalman filter method.
Estimation results.
Estimation results.
In conclusion, compared with other error terms, the gyroscope bias is not effectively estimated by either of the two methods. This phenomenon is predictable because during the alignment process the observability of the gyroscope bias is the weakest compared with other parameters. Taking the overload period as an example, for the Kalman filter technique, the accelerometer scale factors along the y-direction and bias mainly affect the velocity errors; attitude errors are mainly affected by the gyroscope scale factor along the y-direction and accelerometer errors. Therefore the observability of the gyroscope bias is weak and prone to be disturbed when the velocity or attitude is the observation for Kalman filter techniques. For the proposed optimisation-based method, all the inertial sensor errors are observed as “equivalent bias”. Compared with gyroscope scale factors, the gyroscope bias contributes less to the “equivalent bias” and is more likely to be disturbed.
The results of the pure attitude alignment in Equation (31) and the proposed method are compared in Table 3. It can be seen that the resulting accuracy of the proposed method is higher than that of the pure attitude alignment. When the main error terms of inertial sensors are calibrated, the alignment accuracy is improved.
Results at 60 s.
4. CONCLUSION
The Kalman filtering-based transfer alignment method is limited by its dependence on prior information. In this paper, a new optimisation-based method is proposed with the capability of self-initialisation and no prior information is needed. The misalignment and inertial sensor error parameters are estimated simultaneously. For the proposed method, accelerometer and gyroscope outputs of the MINS and SINS are used directly and navigation rate equations are not needed to derive the corresponding error models. In the proposed method, inertial measurement vector integration is used instead of the inertial measurement vector to eliminate the effect of the noise. In addition, proper selection of the time difference interval is made which is critical for guaranteeing the observability degree of the system and inertial sensor errors. Simulation results demonstrate that the proposed method has faster convergence and higher accuracy than Kalman filter techniques. It provides a new way to solve the alignment and calibration problem.