A Three-Stage Accelerometer Self-Calibration Technique for Space-Stable Inertial Navigation Systems

As a specific force sensor, the tri-axis accelerometer is one of the core instruments in an inertial navigation system (INS). During navigation, its measurement error directly induces constant or alternating navigation errors of the same order of magnitude. Moreover, it also affects the estimation accuracy of gyro drift coefficients during the initial alignment and calibration, which will indirectly result in navigation errors accumulating over time. Calibration can effectively improve measurement accuracy of the accelerometer. Device-level calibration can identify all of the parameters in the error model, and the system-level calibration can accurately estimate part of these parameters. Combining the advantages of both the methods and making full use of the precise angulation of the space-stabilized platform, this paper proposes a three-stage accelerometer self-calibration technique that can be implemented directly in the space-stable INS. The device-level calibration is divided into two steps considering the large amount of parameters. The first step is coarse calibration, which identifies parameters except for the nonlinear terms, and the second step is fine calibration, which not only identifies the nonlinear parameters, but also improves the accuracy of the parameters identified in the first step. The follow-on system-level calibration is carried out on part of the parameters using specific force error and attitude error to further improve the calibration accuracy. Simulation result shows that by using the proposed three-stage calibration technique in the space-stable INS, the estimation accuracy of accelerometer error can reach 1×10−6 g order of magnitude. Experiment results show that after the three-stage calibration, the accuracy of latitude, longitude, and attitude angles has increased by over 45% and the accuracy of velocity has increased by over 22% during navigation.


Introduction
Inertial navigation systems (INSs) tend to be the first choice for marine navigation because they are able to output motion parameters completely, autonomously, and continuously [1][2][3]. In order to fulfill the requirement for long-term and high-precision marine navigation, it is essential for INSs to employ high-precision inertial instruments. The gyroscope should be able to work stably for a long time, of which the drift error needs to be small and the constant drift can be precisely compensated. And as for the accelerometer, high stability of the bias and the scale factor is required [4].
All of the inertial instruments must be precisely tested and calibrated before they can work properly to provide accurate measurement in the INS [5,6]. According to the observed quantities, accelerometer calibration can be divided into two types, device-level calibration and system-level calibration. In device-level calibration, the accelerometer is set to measure the gravity of Earth at where g is the scalar gravity. Transforming it to the p-frame yields: where the direction cosine matrix from the n-frame to the b-frame is: cos ψ cos θ sin ψ cos θ − sin θ − sin ψ cos φ + cos ψ sin θ sin φ cos ψ cos φ + sin ψ sin θ sin φ cos θ sin φ sin ψ sin φ + cos ψ sin θ cos φ − cos ψ sin φ + sin ψ sin θ cos φ cos θ cos φ and the direction cosine matrix from the b-frame to the p-frame is: − cos S t sin h cos q − sin S t sin q cos S t sin h sin q − sin S t cos q − cos S t cos h − sin S t sin h cos q + cos S t sin q sin S t sin h sin q + cos S t cos q − sin S t cos h cos h cos q − cos h sin q − sin h   . (4) In Equation (3), ψ, θ, and φ are yaw, pitch, and roll, respectively. In Equation (4), S t is the rotation angle of the platform's revolution axis relative to the revolution axis of the inner gimbal, h is the rotation angle of the revolution axis of the middle gimbal relative to that of the outer gimbal, and q is the rotation angle of the revolution axis of the outer gimbal relative to the base. The rotation angle of the revolution axis of the inner gimbal relative to that of the middle gimbal is zero, for the outer gimbal is tracking the rotation angle of the inner gimbal, and for the revolution axes of the platform, the inner gimbal and the middle gimbal are orthogonal.
The geometric relationship between the orthogonal p-frame and the non-orthogonal a-frame is shown in Figure 1. It is assumed that installation error angles are all small quantities. β x , β y and β z denote the three installation error angles in the x-direction, y-direction, and z-direction, respectively. Each of these three small angles can be decomposed into two rotation angles around the axes of the p-frame, namely β xy and β xz , β yz and β yx , and β zx and β zy . According to the geometric relationship shown in Figure 1, the direction cosine matrix from the p-frame to the a-frame can be calculated as: Thus the specific force measured by the tri-axis accelerometer satisfies an equation, as: where SF is the matrix of scale factors, δb is the bias vector, and K nl is the matrix of nonlinear coefficients of the scale factors. In this paper, a three-stage calibration technique is proposed, to estimate and compensate all the parameters shown in Equation (6). The first two stages will lay a foundation of full-parameter and full-value estimations of model coefficients. The third stage will then correct part of the coefficients and further improve the calibration accuracy. Although the final calibration accuracy will mainly depend on the novel system-level calibration method, it must work on the basis of the device-level calibration. Therefore, the first two stages of the device-level calibration designed for the space-stable INS will also be introduced to form the integrated threestage calibration technique.

Device-Level Coarse Calibration
For the tri-axis accelerometer, there are totally 18 error coefficients that need to be calibrated: three scale factors, SF x , SF y and SF z ; three biases, δb x , δb y and δb z ; three nonlinear coefficients of the scale factors, K nlx , K nly and K nlz ; six installation error angles β xy , β xz , β yz , β yx , β zx and β zy ; and three zero-point offsets of the gimbals' angular sensors, θ 10 of the platform axis, θ 20 of the inner gimbal axis, and θ 30 of the middle gimbal axis. The platform needs to be rotated to 18 positions. Table 1 lists the order of rotation, and the corresponding rotation angles of the gimbals and components of the pframe gravity vector [3], where the symbol ↓ denotes that the angle remains unchanged. Figure 2 shows the directions of the accelerometer's sensitive axes in detail, where N, W, and U denote the local northern, western and up directions, respectively, and x, y, and z denote the orientations of the tri-axis accelerometer. According to the geometric relationship shown in Figure 1, the direction cosine matrix from the p-frame to the a-frame can be calculated as: Thus the specific force measured by the tri-axis accelerometer satisfies an equation, as: where SF is the matrix of scale factors, δb is the bias vector, and K nl is the matrix of nonlinear coefficients of the scale factors. In this paper, a three-stage calibration technique is proposed, to estimate and compensate all the parameters shown in Equation (6). The first two stages will lay a foundation of full-parameter and full-value estimations of model coefficients. The third stage will then correct part of the coefficients and further improve the calibration accuracy. Although the final calibration accuracy will mainly depend on the novel system-level calibration method, it must work on the basis of the device-level calibration. Therefore, the first two stages of the device-level calibration designed for the space-stable INS will also be introduced to form the integrated three-stage calibration technique.

Device-Level Coarse Calibration
For the tri-axis accelerometer, there are totally 18 error coefficients that need to be calibrated: three scale factors, SF x , SF y and SF z ; three biases, δb x , δb y and δb z ; three nonlinear coefficients of the scale factors, K nlx , K nly and K nlz ; six installation error angles β xy , β xz , β yz , β yx , β zx and β zy ; and three zero-point offsets of the gimbals' angular sensors, θ 10 of the platform axis, θ 20 of the inner gimbal axis, and θ 30 of the middle gimbal axis. The platform needs to be rotated to 18 positions. Table 1 lists the order of rotation, and the corresponding rotation angles of the gimbals and components of the p-frame gravity vector [3], where the symbol ↓ denotes that the angle remains unchanged. Figure 2 shows the directions of the accelerometer's sensitive axes in detail, where N, W, and U denote the local northern, western and up directions, respectively, and x, y, and z denote the orientations of the tri-axis accelerometer. Table 1. The order of rotation for the 18-position calibration, and the corresponding rotation angles of the gimbals and components of the p-frame gravity vector.

Rotation Angles of Gimbals
Components of the p-Frame Gravity Vector Table 1. The order of rotation for the 18-position calibration, and the corresponding rotation angles of the gimbals and components of the p-frame gravity vector. 30 10 30 11 During coarse calibration, the nonlinearity of scale factors is ignored. Thus the accelerometer measurement equation shown in Equation (6) can be simplified as: During coarse calibration, the nonlinearity of scale factors is ignored. Thus the accelerometer measurement equation shown in Equation (6) can be simplified as:

The Order of Rotation Rotation Angles of Gimbals Components of the p-Frame Gravity Vector
where C a p g p is the gravity vector sensed by the tri-axis accelerometer. Each position takes 20 s for the rapid measurement. Based on the outputs of the tri-axis accelerometer at different positions and Table 1, 18 equations can be obtained according to Equation (7) Sensors 2018, 18, 2888 6 of 16 to preliminarily estimate 15 parameters except the nonlinear coefficients. Scale factors and biases can be coarsely estimated by solving the equation set corresponding to the six rotation positions (No. 1, 3, 5, 7, 9, and 11) where each sensitive axis of the tri-axis accelerometer points upward or downward along the z n axis in turn. The other equations corresponding to the rest of the rotation positions (No. 2, 4, 6, 8, 10, and 12~18), where one of the sensitive axes is parallel to the x n axis or the y n axis, and other two sensitive axes have included angles of 45 • with the vertical direction, are used to estimate zero-point offsets of the angular sensors and installation error angles.
Estimated values of zero-point offsets,θ 20 andθ 30 , are used to coarsely correct the orientation of the p-frame, and to make the directions of the middle gimbal axis and the inner gimbal axis approach the directions of axes of the p-frame. Concretely, the inner gimbal axis and the middle gimbal axis are inversely rotatedθ 20 andθ 30 , respectively, and then coarse calibration is re-executed. This procedure will repeat until the re-estimatedθ 20 andθ 30 decrease to negligibly small quantities.

Device-Level Fine Calibration
Fine calibration will be executed on the basis of coarse calibration. The orientation of the platform has been corrected and new data will be measured in more accurate rotation positions. In order to improve the data collection accuracy, the measurement time of each position is extended to 10 min. Moreover, to improve the data processing accuracy, the estimation error of coarse calibration is estimated instead of the full values of error parameters.
According to Equation (6), the ideal specific force in the p-frame satisfies: If the result of coarse calibration and outputs of the tri-axis accelerometer during fine calibration ( f a ) have been obtained, the coarse estimation of f p can be calculated as: where "hatted" symbols denote the estimated values.
In order to reduce the requirement for position control accuracy of the platform during fine calibration, the observed variable is not the direct difference between the coarsely estimated value and the ideal value of the specific force at one certain position, but the difference between the squared amplitudes of them. Subtracting the squared amplitude of Equation (9) from that of Equation (8), and then taking half of the result yields the observed variable, as: where: According to Equations (10)-(15), the observation equation at position j can be written as: where the unknown error vector is: Based on the data collected at the 18 positions, an 18 × 12-dimension observation equation can be formed based on Equation (16), and the unknown variables shown in Equation (17) can be solved by the least squares method. Finally, using scale factors to convert the solved dx and adding the result to the coarse estimation will yield a full-value fine estimation of 12 error coefficients. It should be noted that β xz and β yz , β zy and β xy , and β yx and β zx are linearly correlated in pairs and cannot be decomposed from each other. Therefore the fine calibration can only improve the accuracy of f p , but it cannot calibrate the rotation characteristics of the a-frame.

System-Level Calibration Based on Specific Force Error and Attitude Error
Errors of inertial instruments will induce the navigational errors of INSs, including the positional error, the velocity error and the attitude error. Inversely, the independently observed navigation errors such as the specific force error, the attitude error and the position error can be used to calibrate the main error parameters in the system. As the position error is affected by both the gyroscope error and accelerometer error, it is generally difficult to directly extract the accelerometer-induced term from the position error. On the other hand, the specific force error and the horizontal attitude error are mainly caused by the accelerometer error. Therefore this paper employs the specific force and horizontal attitude angles to correct part of the accelerometer's error coefficients. The prerequisite for accurate observation in such a system-level calibration is that the system should be installed on a static and horizontal base.

Calibration Using Specific Force Error
Under the ideal condition, that there is no installation error between the base and the stabilized platform, and that the measurement of gimbals' rotation angle is accurate, there is an approximate relationship: where L denotes the latitude. On the static and horizontal base, the local scalar gravity g can be accurately calculated. Then, according to Equations (1) and (18), the ideal specific force vector in the p-frame can be written as: Similar to the device-level calibration method presented in the previous section, the difference of the squared amplitudes of the measured and the ideal specific force vectors is chosen as the observed variable, as: A 5-dimensional state vector is constructed as where: Thus the observation equation is: Considering Equation (19), the observation matrix is: The five components of x are actually the mean value, the sine and the cosine fundamental waves, and the sine and the cosine second harmonics of the Earth's rotation rate in the specific force error. Assuming that the observation noise is a zero-mean Gaussian white noise, x can be estimated by a Kalman filter.

Calibration Using Attitude Error
On the static and horizontal base, the attitude error of the INS with both vertical damping and horizontal velocity damping mainly consists of a constant value and a 24 h long period component, which is induced by the gyroscope error, the accelerometer error, and the rotation angle errors of platform's gimbals. Using the space-stable mechanization, the azimuth error will accumulate over time due to the gyroscope error, while the horizontal attitude errors are mainly caused by the accelerometer error and have no such accumulation. Therefore, the horizontal attitude errors are chosen as the observed quantities for accelerometer calibration.
Neglecting the short-period misalignment angle of gyro case rotation, the attitude error on the static and horizontal base has an expression as: where δκ b x , δκ b y , and δκ b z are the transverse tilt error, the longitudinal tilt error, and the azimuth error of the base, respectively, and the superscript e denotes the earth-fixed frame (e-frame).
Substituting the e-frame specific force error into Equation (28) yields the expressions of pitch error and roll error. Both of them consist of the mean value, the fundamental wave, and the second harmonic of the Earth's rotation rate. The five Fourier coefficients of the pitch error are: The other five Fourier coefficients of the roll error are: The subscripts indicate which component each expression is (1 for the fundamental waves, 2 for the second harmonics, s for the sine functions, and c for the cosine functions). After the error curves of pitch and roll are obtained, these coefficients can be extracted by Fourier expansion.

Parameter Separation Algorithm for Combined Error Coefficients of Accerelometer
In Sections 5.1 and 5.2, 15 combined coefficients (Equations (21)-(25), (29)- (33), and (34)-(38)) are estimated by static Kalman filter and Fourier expansion. Twelve of them are selected and divided into three groups to extract the residual error parameters of the tri-axis accelerometer. Since only the fundamental waves and second harmonics of pitch and roll errors are concerned, the requirement for the horizontality of the base need not be very strict.
The first group consists of three equations, which are the coefficients of the sine and cosine second Similarly, the second group consists of Equations (25), (33), and (38), and ∆β yz −∆β xz 2 can be solved as:∆ β yz −∆β xz 2 = 1 3 Another six parameters ( x g , y g , ∆β yx , ∆β zx , ∆β xy , ∆β zy ) can be solved using the least square method based on the remaining six combined coefficients, according to: The six parameters in Equation (41)  and ∆SF x 2SF x − ∆SF y 2SF y will be indirectly compensated by deducting the error terms about them from the n-frame specific force expression.

Simulation Result
As the final accuracy of the three-stage calibration mainly depends on the performance of the last stage, the simulation of the system-level calibration is conducted first. In the simulation, the location of the static base is set at 40 • N. The non-orthogonal installation error of the tri-axis accelerometer is quite a stable mechanical error, and after the 18-position two-stage device-level calibration, the residuals are generally less than 10 . Thus, the residuals of accelerometer's error parameters are assigned as Tables 2 and 3. In addition, a zero-mean random noise with an amplitude of 5 × 10 −3 m/s 2 is added to each sensitive axis of the tri-axis accelerometer.
Truth values of the system-level combined coefficients are calculated according to the error parameters shown in Tables 2 and 3, and are listed in Table 4 as reference values. The combined coefficients are estimated by system-level calibration based on specific force error and attitude error, and are also listed in Table 4. Comparison shows that the estimation accuracy of each combined coefficient is better than 10 −6 rad, which is equivalent to an accelerometer calibration accuracy that is better than 10 −6 g. It indicates that the method can be used for the accelerometer calibration of high-precision INSs.

Experiment Results
In order to validate the proposed three-stage accelerometer calibration method, a static-base navigation test has been conducted using a high-precision space-stable INS. The INS is in the navigation grade class and has been introduced in reference [36]. Two dual-axis free-rotor gyros and a tri-axis quartz accelerometer are mounted on a gimbaled platform. According to the outputs of the free-rotor gyros, the platform is servo-controlled to be stabilized in the inertial space. The gyros have an ultra-low drift of less than 0.001 • /h. The measurement error of the accelerometer is of a 1 × 10 −5 g order of magnitude, and its stability is better than 1 × 10 −5 g per year. The static measurement accuracy of gimbals' rotation angles can reach the level of 1 . The position error of this navigation system can reach the order of magnitude of 1 n mile per day. The data used in the experiment are collected when the system is installed on a static and horizontal base.
Navigation errors before calibration, after the original two-stage 18-position coarse and fine self-calibration, and after the entire three-stage calibration are compared in Figures 3-5 and Table 5.
Since the system employs vertical damping, height error and vertical velocity error can be neglected and are not presented. The output position of the system is compared with the reference position of the base to obtain the position error. The output velocity and horizontal attitude angles are actually the velocity error and horizontal attitude error, for the base is static and horizontal. All experiment results here have been normalized.
The system has employed gyro case rotation techniques to decrease the gyro drift. Although its influence on the platform's motion has been deducted during navigation calculating, there are still 8-minute residuals in the velocity and attitude error curves. As the cycle of these residuals are much shorter than the test time, they are severely squeezed in the error curves shown in Figures 4 and 5.
Accelerometer parameters will change after the high-precision INS works for a period of time, which will consequently lead to the increase of navigation errors. Therefore, on-site and online calibration is necessary to improve the navigation precision. Latitude error, longitude error, and attitude errors are reduced by over 30% after the two-stage calibration, and by over 45% after the three-stage calibration. Since the velocity error mainly consists of Schuler oscillation, and the error term is related to the gyro case rotation residual, the improvement of velocity accuracy is not conspicuous. The velocity error has decreased by over 10% after the two-stage calibration, and by over 22% after the three-stage calibration.
here have been normalized.
The system has employed gyro case rotation techniques to decrease the gyro drift. Although its influence on the platform's motion has been deducted during navigation calculating, there are still 8minute residuals in the velocity and attitude error curves. As the cycle of these residuals are much shorter than the test time, they are severely squeezed in the error curves shown in Figures 4 and 5.   here have been normalized. The system has employed gyro case rotation techniques to decrease the gyro drift. Although its influence on the platform's motion has been deducted during navigation calculating, there are still 8minute residuals in the velocity and attitude error curves. As the cycle of these residuals are much shorter than the test time, they are severely squeezed in the error curves shown in Figures 4 and 5.

Conclusions
A three-stage tri-axis accelerometer self-calibration technique is proposed with a thorough theoretical and mathematical analysis. Compared with existing calibration methods, the proposed three-stage calibration method has the following features and advantages: (1) Taking full advantage of high-precision angulation of the gimbaled platform in the space-stable INS, the two-stage device-level calibration is implemented directly on the platform of the space-stable INS. Thus the entire three-stage self-calibration can be implemented in the INS after the tri-axis accelerometer has been mounted on the platform. Besides, during the device-level calibration, zero-point offsets of the gimbals' angular sensors are simultaneously estimated and used to correct the platform's orientation. This will help improve the system's navigation performance. (2) Different from the existing calibration methods which usually estimate error parameters of both accelerometer and gyroscope simultaneously, the proposed system-level calibration is designed based on error propagation to choose observed quantities only resulting from accelerometer measurement error. This can get rid of the influence of gyro drifts.
(3) Combining device-level calibration and system-level calibration, the proposed method can achieve both full-parameter estimation of the accelerometer measurement model and high-precision correction for part of the parameters, while existing calibration methods normally concentrate on one of the two aspects. Both simulation and experiment results validate the proposed method. The simulation result shows that using the proposed three-stage calibration technique, the estimation accuracy of accelerometer error is better than 1 × 10 −6 g. Experiment results show that after implementing the three-stage calibration, during navigation, the accuracy of latitude, longitude, and attitude angles has increased by over 45%, and the accuracy of velocity has increased by over 22%.