## Abstract

This paper focuses on a heterogeneous redundant inertial navigation system (INS) that consists of one three-axis rotation-modulation fiber optic gyroscope (FOG) INS and one strapdown ring laser gyroscope (RLG) INS. A novel co-calibration method is proposed to estimate the FOG scale factor instability and RLG bias instability error in real time, in which the geometric constraints of the FOG-INS rotating axis are deployed to establish the observation equation without external reference information. In addition, the proposed method has a global navigation capability, which is achieved by using the normal vector of the earth ellipsoid and the state transformation extended Kalman filter to predict positioning errors. An output position error compensation method is used without disturbing the original INSs. Simulation and physical experiment results show that the positioning error of the heterogeneous redundant marine INS is reduced by more than 30% over 300 h of navigation.

## 1 INTRODUCTION

Nowadays, system-level redundant inertial navigation systems (INSs) are adopted in modern strategic navigation applications to improve reliability. System-level redundancy in nautical navigation is normally achieved by applying more than one set of rotation-modulation INSs (RINSs) to improve navigation accuracy (Levinson & Majure, 1987). Although single-axis and two- or three-axis RINSs are widely used, these systems are hindered by the fact that the INS reliability can be undermined by mechanical gimbals. Moreover, the RINS navigation information contains a saw-tooth pattern of the velocity and attitude error, which is not expected in high-accuracy transfer alignment applications or other high-accuracy attitude measurement applications (Levinson & Majure, 1987). Although an RINS can provide smoother output information during a halt (Levinson & Majure, 1987), interrupting the rotation-modulation process inevitably reduces the navigation accuracy, while increasing the relative complexity of control. For this reason, a heterogeneous redundant configuration consisting of one RINS with a strapdown INS has been developed to reach a trade-off in some marine applications (Liu et al., 2014, 2016).

Rotation-modulation technology was initially applied and then developed for a ring laser gyroscope (RLG)-INS. For example, a 16-position rotation-modulation method for a double-axis INS was established for MARLIN (marine ring laser inertial navigator) by Sperry Marine (Levinson & Majure, 1987). All constant error sources are successfully averaged out, including the gyroscope and accelerometer bias, misalignment angles, and scale factors over a modulation cycle. Subsequently, a single-axis four-position rotation-modulation method was developed and used in the AN/WSN-7B RLG-INS produced by Litton Marine Systems, which averages out the horizontal gyroscope and accelerometer constant bias in one rotation period of 20 min, achieving 1 n mile within 24 h (Tucker & Levinson, 2000). The success of rotation-modulation technology in RLG-INSs has motivated interest in fiber optic gyroscope (FOG)-INSs, the performance of which has been greatly improved in recent years. Today, strategic-grade FOGs surpass RLGs in terms of a lower angular random walk noise and better cost-effectiveness (Lefèvre, 2013).

However, an engineering challenge still exists for the application of rotation-modulation technology in FOG-INSs. Although FOG-INSs utilize the same theoretical principle as RLGs, the FOG features greater sensitivity to the thermal environment, which results in a scale factor instability of 10 ppm or even larger for a typical high-grade FOG over a full temperature range (Korkishko et al., 2017; Lefevre, 2012), whereas the RLG instability after compensation can be maintained at less than 1 ppm. When the FOG-INS is in a start-up process or the environmental temperature varies dramatically, the FOG scale factor instability can be even greater (Shen & Chen, 2012; Zhang, 2018). Consequently, although most of the constant errors of a rotation-modulation FOG-INS are averaged out in one multi-position rotation-modulation cycle, the saw-tooth pattern of the velocity and attitude error caused by the FOG scale factor residual error is inevitable. Therefore, continuous-rotation modulation methods have been proposed for double-axis FOG-INSs to mitigate this drawback (Du et al., 2017; Liu et al., 2017), and self-calibration technologies have been advanced to calibrate the turn-on error during or after the docking alignment stage (Gao et al., 2016). However, the residual gyroscope scale factor error coupled with the earth’s rotation and vehicle angular motion can generate an equivalent gyroscope drift in the northward and downward directions, leading to a linear build-up of longitude positioning error (Levinson & Majure, 1987) in long-term navigation.

Effectual approaches to realizing the full potential of heterogeneous RINSs have been explored for decades. The information fusion methods can be roughly divided into two categories according to the fusion level. One category includes methods based on loosely coupled fusion in which the fault detection methodology is adopted to directly weight or select navigation results for each independent INS (He et al., 2014; Li et al., 2011). The other category comprises methods based on tightly coupled fusion in which a parameter estimation methodology, such as the least squares or Kalman filter method, is utilized to estimate and compensate for the error source of each INS (L. Wang et al., 2018; Wu & Li, 2019).

Polar region navigation ability plays an essential part in global navigation. As the traditional method of representing latitude and longitude position is not suitable for polar region navigation because of the meridian convergence problem, the transverse latitude and longitude position representation method has been adapted (Li et al., 2016; Qin et al., 2018; Yao et al., 2016). Accordingly, an algorithm with the capability of switching between polar and non-polar regions is required for global navigation. A new position representation method based on the normal vector of the surface of the earth ellipsoid model has been proposed for global navigation, which is more concise and does not require such switches (Gade, 2010; Liu et al., 2020). In addition, with this method, it is simple to convert the position representation from the normal vector to latitude and longitude as well as transverse latitude and longitude.

This paper proposes a novel co-calibration method for a specific type of dual heterogeneous redundant marine INS. This dual heterogeneous redundant marine INS consists of a three-axis rotation-modulation FOG-INS and a strapdown RLG-INS. The gimbal rotation axes of the FOG-INS can be predetermined as known vectors in the inertial measurement unit (IMU) body frame, and the direction cosine matrix (DCM) between the IMU body frame of the FOG-INS and of the strapdown RLG-INS can also be predetermined when all three gimbal rotation axes of the FOG-INS are at their initial position. The geometrical constraints offer relative attitude observations to estimate the predominant error sources of each INS (scale factor error for the FOG and bias for the RLG). High-accuracy FOGs potentially have a smaller angular random walk error than RLGs (Mead & Mosor, 2020) whereas RLGs have a smaller scale factor error than FOGs (Lèfevre et al., 2020); thus, there are complementary advantages for the dual heterogeneous redundant marine INS consisting of both a three-axis rotation-modulation FOG-INS and a strapdown RLG-INS.

The remainder of this paper is organized as follows: Section 2 describes the proposed three-axis continuous-rotation modulation method for the FOG-INS. Section 3 derives detailed Kalman filter formulations of a co-calibration method for the dual heterogeneous redundant marine INS consisting of one three-axis rotation-modulation FOG-INS and one strapdown RLG-INS. Section 4 presents an output compensation architecture to correct the position errors of each INS by using the estimated error source without disturbing the original inertial navigation process. Polar region navigation applicability is achieved by using the normal vector model. Simulations of the heterogeneous RINS are presented in Section 5, followed by physical experiments in Section 6. Finally, conclusions are presented in Section 7.

## 2 THREE-AXIS ROTATION MODULATION

### 2.1 Three-Axis Rotation Modulation Method

The IMU frame adopts right-hand coordination, and the x, y, and z axes are defined as the right, forward, and upward directions. To focus on the topic of this paper, the nonorthogonal angles between three gimbals are preliminarily considered to be well calibrated. That is to say, the inner gimbal axis is identical to the IMU y-axis, and the middle gimbal axis is identical to the x-axis when the inner gimbal is at the initial position. The outer gimbal axis is identical to the z-axis when the other gimbals are at their initial positions.

The original purpose of the rotation-modulation technique is to restrain the systematic errors in a marine INS. In this paper, a novel rotation-modulation technique is proposed that can also provide certain observation information to estimate the scale factor error for the FOG and the bias for the RLG in a dual heterogeneous redundant marine INS. Continuous rotation modulation is adopted in the FOG-INS. The three-axis rotation-modulation method proposed in this paper includes two processes, namely, a rotation process and a stopping process. During the rotation process, at any given moment, two axes are continuously rotating: the outer gimbal is rotating relatively rapidly by 360° clockwise or counterclockwise, and the inner or middle gimbal is rotating slowly by 360° clockwise or counterclockwise. The period of the inner or middle gimbal is an integral multiple of the period of the outer gimbal. During the stopping process, when the outer gimbal returns to the initial position after completing a 360° continuous-rotation cycle, it stops briefly while the inner or middle gimbal is still rotating. At this time, the vector of the rotation axis in the body frame can naturally provide an observation for INS error estimation. Figure 1 visualizes the rotation sequence of the physical turntable.

A detailed mathematical description of the rotation modulation is given below. The three-axis rotation-modulation scheme is designed as follows. The continuous-rotation modulation period of the IMU should be chosen to avoid the Schuler period, which is 84.4 min. Let Ω be a constant angular rate and *T*_{0} - 2*π* / Ω. The IMU rotates 360° about the inner gimbal axis or the middle gimbal axis at a speed of ±Ω in a period of *T*_{0}. The outer gimbal adopts a shorter rotation-modulation cycle. Let us suppose that *T _{Out}* =

*T*

_{0}/ 2

*m*, where ,

*m*≫ 1. Let

*T*be a short period; then, the corresponding rotation rate of the outer gimbal is as follows:

_{s}1

The IMU rotates 360° about the outer gimbal axis at a speed of ±Ω_{Out} in a period of *T _{Out}* −

*T*).

_{t}For *T _{s}* = 1

*s*, with the initial angular position of the three gimbals at the zero point, the angles of the three axes and the measurement moment are illustrated in Figure 2.

According to Figures 1 and 2, the rotation rates of the three gimbals can be written as piecewise functions. The rotation rate of the inner gimbal is as follows:

2

The rotation rate of the middle gimbal is given by the following:

3

The rotation rate of the outer gimbal is as follows:

4

where *k* = 0,1,2,…, *p* = 0,1,2…, *q* = 0,1,2… (*m* − 1).

Compared with the rotation speed of the outer gimbal, the rotation speed of the inner and middle gimbals is much lower. Thus, the outer gimbal completes several consecutive reciprocating rotation cycles before the inner (or middle) gimbal completes one rotation cycle. When the outer gimbal is rotated to return to its initial position, the inner (or middle) gimbal axis that is still rotating can be predetermined as a fixed vector in the coordinate frame of the RLG-INS. This vector provides a relative attitude measurement, which enables one to estimate the scale factor error of the FOG and the bias of the RLG. The stopping rotation period *T _{s}* of the outer axis is very short; thus, its negative effect on rotation modulation can be ignored.

### 2.2 Error Source in Long-Term Navigation

The proposed method can be decomposed into two two-axis rotation-modulation processes in which the inner and middle gimbal axes rotate in turn respectively with the outer axis. Owing to the contributions of Levinson (Levinson & Majure, 1987), the two-axis rotation-modulation error source has been investigated. The constant bias, scale factor error, and misalignment can be averaged out over a complete modulation period (also called an indexing cycle). However, the influence of the earth’s rotation accumulates during the long-term inertial navigation process. The earth’s rotation coupled with the scale factor error generates an equivalent gyroscope bias that becomes the main error source of this three-axis rotation FOG-INS, with the equivalent gyroscope bias given below (a detailed derivation can be found in Appendix A):

5

where , , and represent the equivalent eastern, northern, and upward gyroscope bias in the local navigation frame, *δk _{gx}*,

*δk*and

_{gy}*δk*denote the gyroscope scale factor errors of the

_{gz}*x, y,*and

*z*axes,

*ω*is the earth’s rotation rate, which is 7.292115 × 10

_{ie}^{−5}

*rad*/

*s*, and

*L*is the local latitude. The * symbol indicates the periodic trigonometric function elements that are averaged out in one rotation-modulation period.

## 3 CO-CALIBRATION METHOD FOR A DUAL HETEROGENEOUS INS

The main error of the rotation-modulated FOG-INS generated during long-term navigation is caused by the residual gyroscope scale factor error coupled with the earth’s rotation and the vehicle’s angular motion, whereas the residual gyroscope bias error is the main error source for the fixed strapdown RLG-INS.

### 3.1 System Equation Modeling of the Dual INS

In the earth-centered earth-fixed frame (ECEF, *e* frame), the attitude error is defined by , where the tilde represents the DCM containing the attitude error. The differential equations of the attitude error of the two INSs are directly given by the following:

6

where represents the earth’s rotation rate vector with respect to the inertial frame (*i* frame) represented in the *e* frame, *b*_{1} and *b*_{2} are the IMU body frame (*b* frame) of the rotation-modulated FOG-INS and fixed RLG-INS, respectively, , are the DCMs between the *b* frame and *e* frame, and , are the gyroscope errors.

The joint attitude error and its differential equation are defined as follows:

7

The gyroscope errors of the two INSs are modeled as the residual of the scale factor error and the bias error (considering that the gyroscope bias of the FOG-INS is averaged out by rotation modulation and the gyroscope scale factor of the RLG-INS is well calibrated), with the corresponding angular random walk terms and modeled as zero-mean white noise:

8

Here, *diag*(·) is a diagonal matrix formed from a certain vector, and is the gyroscope output of the FOG-INS. The approach for is given in Section 4.3.

and are modeled as first-order Markov processes, with the following differential equations:

9

where *τ*_{1}, *τ*_{2} are correlation times of the first-order Markov process, which are typically determined by a stationary test (El-Diasty & Pagiatakis, n.d., 2009), and *w*_{δk}, ** w_{δb}** are white noise.

We select the joint attitude error and the error sources of each INS as states of the Kalman filter:

10

and establish the system equation as follows:

11

12

13

14

15

16

Here, × denotes the transformation of a vector to a skew-symmetric matrix as well as the cross product.

### 3.2 Measurement Equation of the Dual INS

The inner and middle gimbal axes of the rotation-modulation FOG-INS are represented as vectors and , respectively, in the IMU body frame of the FOG-INS, which can be acquired by calibration and simple computation by rotating the corresponding axis. Therefore, a measurement of the rotating vector can be established as follows:

17

where is the DCM between the FOG-INS and RLG-INS and the initial (when all three gimbal rotation axes of the FOG-INS are at the initial position) can be obtained beforehand by precise mechanism installation and calibration. The tilde represents the DCM containing the attitude error. Note that Equation (17) relies on the hypothesis of a perfect three-axis turntable.

When the outer gimbal of the rotation-modulation FOG-INS is at the initial position and *ω _{Out}*(

*t*) = 0,

*ω*(

_{Mid}*t*) = 0, then holds. When the outer gimbal of the rotation-modulation FOG-INS is at the initial position and

*ω*(

_{Out}*t*) = 0,

*(*

**ω**_{In}*t*) = 0, then holds.

The DCM containing the attitude error is decomposed as follows:

18

where . The measurement vector can be derived as follows:

19

Thus, the measurement equation can be given as the following:

20

where:

21

*ν _{inner}* and

*ν*are measurement noise vectors modeled as Gaussian white noise vectors. The system equation and measurement equation are given in Equations (10), (11), and (20), respectively. As filter states, the gyroscope scale factor error of the FOG and the bias error of the RLG can be estimated by the Kalman filter. A brief Kalman filter introduction is given in Appendix C.

_{middle}## 4 OUTPUT COMPENSATION METHOD

### 4.1 Output Compensation of the Dual INS

Directly compensating the scale factor error for the FOG and the bias for the RLG by feeding back the estimation into the navigation calculation process would result in changes in the dual INS internal procedure, whereas incorrect estimations of the FOG scale factor error and RLG bias would cause irreparable INS error. It is recommended that the independence of the original dual INSs be maintained for system reliability. Therefore, an output compensation method is proposed by which the original navigation results of the INSs are preserved along with the compensated navigation results. The real-time output switches between these two sets of navigation results, as schematically illustrated in Figure 3. In addition, the normal vector of the surface of the earth ellipsoid model is used for global positioning representation, which is described in detail in Appendix B.

and are the original normal vector position outputs of each INS represented in the *e* frame, and and are the normal vector position compensation results of each INS calculated by the error prediction model that takes and as the input. In addition, and are the normal vector positions after the correction of each INS. The process for transforming the normal vector position representation to local latitude and longitude or transversal latitude and longitude is given in Section 4.2.

Switch 1 is turned to the right before the Kalman filter stabilizes; hence, the positioning result is the original output of each INS. Then, switch 1 turns left to output the compensated result of each INS. Generally, the switching time can be configured at the end of a complete modulation period of the FOG-INS. The error prediction model adopts a linear system model based on the state transformation extended Kalman filter (ST-EKF) (M. Wang et al., 2018, 2019, 2021), which reduces the complexity of the output compensation algorithm compared with the traditional EKF in that the ST-EKF uses gravity rather than the dynamic specific force of each INS in the EKF. Note that the prediction only contains the propagation step of the standard Kalman filter, and detailed equations are given in Section 4.2.

### 4.2 Global Capability of the Compensation Method

To enhance the capacity of the proposed method for global navigation, the normal vector of the surface of the earth ellipsoid model is adopted to represent the vehicle position on the earth.

The prediction states of the FOG-INS and RLG-INS are as follows:

22

where , , and represent the attitude error, horizontal state-transformed velocity error, and normal vector position error represented in the *e* frame, respectively. The state-transformed velocity error and horizontal state-transformed velocity error are defined as follows:

23

where , are the velocity output and velocity error of the INS represented in the *e* frame, respectively.

The differential equation of the linear prediction model can be represented as the following:

24

25

where *g* is the norm of the local gravity, , *j* = 1,2, and * K_{R}* is a scale matrix:

26

with *R _{E}*,

*R*,

_{N}*h*as the transverse radius of curvature, meridian radius of curvature, and geodetic altitude, respectively. A detailed derivation can be found in Appendix B.3.

The inputs of each linear prediction model are as follows:

27

with the following driving matrix:

28

We assign •_{j}(*t _{k}*) = (•)

_{j,k}. for brevity, and the discrete prediction models are as follows:

29

30

where and represent the predicted error and estimated gyroscope error source of each INS at time *t _{k}* with initial values of and . Δ

*t*=

*t*–

_{k}*t*

_{k–1}is the discrete time.

The predicted normal vector error can be easily converted into the local longitude error *δλ* and latitude error *δL* at low and middle latitudes or the transversal longitude error *δλ _{t}* and transversal latitude error

*δL*at high latitudes. Therefore, the output is designed to adapt the global navigation in this paper, as illustrated in Figure 4.

_{t}The transformation from the normal vector position representation to local longitude and latitude in non-polar regions as well as transversal longitude and latitude in polar regions is as follows:

31

32

The above result was acquired from Appendix B. Note that the vertical velocity error and the height error are neglected, as they can be calibrated by other sensors, such as a depthometer.

### 4.3 Recovering Gyroscope Measurements from the Attitude Matrix Output

Note that Equation (13) contains the angular velocity of the INS. Yet, some INSs may not provide original gyroscope data. To solve this problem, the angular velocity of the INS is recovered through the relationship between the equivalent rotation vector and the DCM. Let us suppose that the output frequency of the INS DCM is 100 Hz, then the angular velocity can be calculated as follows:

33

where is the equivalent rotation vector of the INS rotation, is its norm, and *T* is the period of computation. can be acquired by the chain law:

34

where:

35

## 5 SIMULATIONS

A semi-physical simulation was performed to verify the effectiveness of the proposed method. The simulation process is illustrated in Figure 5. The theoretical gyroscope and accelerometer sampling data of the INS were calculated based on the IMU attitude, angular velocity, local gravity, and earth’s rotation. The actual gyroscope and accelerometer static noise and slow drift error were obtained based on 15 days (360 h) of 200-Hz sampling data from an actual FOG-INS and an actual RLG-INS. Allan variance analysis was applied to determine the parameters of the noise model and the Kalman filter. Constant bias and scale factor errors of the gyroscopes were also added to the simulated sampling data based on the IMU error model. The error parameters are shown in Table 1.

During the initial alignment stage of 6 h, the equivalent north and upward gyroscope biases of the RLG-INS are estimated, whereas the gyroscope bias estimation of the FOG-INS is unnecessary because it has been averaged out. After the initial alignment has been completed, the estimated gyroscope bias of the fixed RLG-INS is −0.0033, −0.00003, and 0.0006°/h. Most of the north and upward gyroscope bias has been estimated whereas the eastern bias has not yet been estimated because of its unobservability. Fortunately, the eastern gyroscope bias has no contribution to the build-up of longitudinal positioning error. The residual gyroscope scale factor error of the FOG-INS can exceed 10 ppm, considering the repeatability; thus, the simulated gyroscope scale factor errors of the three FOGs are set as 15, −25, and −35 ppm, respectively, taking the instability into account. Note that the constant bias errors and the scale factor errors of the accelerometers are ignored in the simulation process because their contribution to the positioning error is of less significance than those of the gyroscopes. The simulation results are displayed in Figures 6 and 7.

In the simulation, the measurement noise of the positioning mechanism of the rotated FOG-INS is 2 arcseconds, that is, the positioning repeatability of the outer gimbal is 2 arcseconds. After a complete 5-h rotation-modulation cycle, the estimated scale factor errors of the three FOGs are 14.82, −25.18, and −34.22 ppm; these values are very close to the simulated value, with an estimation accuracy better than 1 ppm. The residual gyroscope bias of the RLG-INS with slow drift is also estimated with an accuracy better than 0.0001°/h. The results of the positioning error of eachINS are normalized by dividing by the respective maximum of the total position error and are compared in Figure 7. Based on the proposed co-calibration method, the total positioning errors are reduced by 57% for the FOG-INS and 32% for the RLG-INS. Note that the error source is selected as a constant, but the time-varying error can also be estimated by using proper parameters of correlation time and noise strength in Equation (9).

To further explore the possibility of estimating the time-varying errors, for example, to estimate the time-varying errors with a 24-h period, a simulation was performed in which sine waves with a 24-h period and an amplitude of 10 ppm are added to the FOG measurements while sine waves with a 24-h period and an amplitude of 0.001°/h are added to the RLG measurements.

The results are displayed in Figures 8 and 9. It can be seen that the sine-wave time-varying error can be estimated in real time and that the position errors of the FOG-INS and RLG-INS are reduced by 66% and 27%, respectively, after compensation.

## 6 PHYSICAL EXPERIMENT

A high-precision FOG-INS, a high-precision RLG-INS, and a three-axis flight simulation turntable were utilized to conduct a dual-INS co-calibration experiment. The main specifications of the FOG-INS and RLG-INS are listed in Table 2. The physical experiment results are presented in Figures 10 and 11.

It can be seen from Figures 10 and 11 that the proposed co-calibration method can compensate for most of the longitudinal position error, which linearly increases with time. During the 300-h inertial navigation process, the normalized maximum position error of the FOG-INS decreases from 1 to 0.63, corresponding to a reduction of 37%, while the normalized maximum position error of the RLG-INS declines from 1 to 0.33, corresponding to a decrease of 67%.

## 7 CONCLUSIONS

For long-term navigation of a three-axis rotation-modulation FOG-INS, the linear accumulated positioning error is primarily attributed to the instability of the gyroscope scale factor whereas linear accumulated positioning error of strapdown RLG-INS mainly results from gyroscope bias. The co-calibration method proposed in this paper for a dual INS can estimate and compensate the main error sources of each INS. The accuracy and effectiveness of this method were verified by simulations and physical experiments. During a 300-h physical experiment, the position errors were reduced by 37% and 67% for the FOG-INS and RLG-INS, respectively. In addition, this method adopts an output compensation approach to maintain the independence of each INS without comprehensive reformation of the original dual INSs, which facilitates upgrades of the original system.

## HOW TO CITE THIS ARTICLE

Cui, J., Wu, W., Ma, T., Wang, M., & Ji, C. (2023). A novel co-calibration method for a dual heterogeneous redundant marine INS. *NAVIGATION, 70*(4). https://doi.org/10.33012/navi.620

## APPENDIX A EQUIVALENT GYROSCOPE BIAS GENERATED BY THE SCALE FACTOR ERROR

A derivation of the equivalent gyroscope bias generated by the scale factor error coupled with earth rotation is given here. The earth rotation vector is expressed in the local *n* frame (in eastern, northern, upward sequence):

1

where *ω _{ie}* is the earth’s rotation rate, which is 7.292115 × 10

^{−5}

*rad*/

*s,*and

*L*is the local latitude.

Without a loss of generality, the three-axis rotation rates are assigned as *ω _{x}*,

*ω*, and

_{y}*ω*. The DCM between the

_{z}*n*frame and

*b*frame is given for two segments according to Figure 2:

2

where corresponds to the period during which the middle gimbal is at the initial position while the inner and outer gimbal axes rotate. corresponds to the period during which the inner gimbal is at the initial position while the middle and outer gimbal axes rotate.

The gyroscope scale factor error and misalignment error are presented in a matrix:

3

where *δk _{gj}* (

*j*=

*x, y, z*) is the gyroscope scale factor error and

*δk*(

_{jm}*j*=

*x, y m*=

*y, z*

*j*≠

*m*) denotes the misalignment error.

The equivalent gyroscope bias can be derived as follows:

4

which can be calculated via symbolic computation software such as Mathematica or MATLAB.

The result is expressed as follows:

5

where , , and represent the equivalent eastern, northern, and upward gyroscope bias in the *n* frame. The * symbol denotes the periodic trigonometric function elements that can be averaged out over a considerable long-term rotation-modulation process. It can be determined from Equation (5) that the three gyroscope scale factor errors generate equivalent gyroscope biases in the northern and upward directions; hence, these terms must be estimated and compensated for during long-term navigation.

## APPENDIX B INTRODUCTION OF THE NORMAL VECTOR POSITION REPRESENTATION

The normal vector (n-vector) is a non-singular two-dimensional representation of horizontal position on the surface of the earth ellipsoid model with three elements (Gade, 2010). Figure 12 illustrates the geographical definition of the normal vector, earth frame, and transverse earth frame. The n-vector can represent all positions on earth uniquely, avoiding singularities or discontinuities near the polar area.

### B.1 From n-Vector Representation to the Navigation Frame and Transverse Navigation Frame

According to the geometry of the system, the following relation is established:

6

where *L*, *λ* are the local latitude and longitude and *L _{t}*,

*λ*are the transverse latitude and longitude. The DCMs between the ECEF frame and the local navigation frame as well as the transverse navigation frame are as follows:

_{t}7

By taking the differential of both sides of Equation (6), we can express the relationship between the n-vector error and navigation frame and the transversal navigation frame as follows:

8

9

10

Hence, Equations (31) and (32) in the main text can be proved based on the above mathematical derivations.

### B.2 From n-Vector Representation to the ECEF Frame

The n-vector can represent horizontal positioning on the reference ellipsoid. The ECEF frame position can be represented by the n-vector and the geodetic height, while the positioning error can be represented by the n-vector error and the geodetic height error as follows:

11

where *R _{E}* is the transverse radius of curvature,

*h*is the geodetic height,

*e*is the eccentricity of the earth ellipse model, and and

*δh*are the horizontal and vertical position errors, respectively.

### B.3 Derivation of Navigation Equations via n-Vector Positioning Representation

We can rewrite the differential equations of the attitude error and state transform velocity error (M. Wang et al., 2018):

12

The horizontal and vertical elements in can be isolated by the n-vector as follows:

13

in which is the horizontal element of and is the vertical element.

The differential equation of the n-vector is as follows:

14

where *R _{E}*,

*R*,

_{N}*h*are the transverse radius of curvature, meridian radius of curvature, and altitude, respectively, and .

*is a scale matrix:*

**K**_{R}15

The differential equation of geodetic altitude is used:

16

Note that the vertical velocity error and height error are neglected because they can be calibrated by other sensors, such as a depthometer. We have the following:

17

where the approximation of Equation (18) is used and the vertical error is regarded as zero:

18

is the norm of the local gravity, which can be calculated as follows (Groves, 2012):

19

with:

20

Here, *R*_{0} = 6378137*m*, , *R _{P}* = 6356752.31425

*m*, and

*μ*= 3.986004418×10

^{14}

*m*

^{3}·

*s*

^{−2}.

The n-vector differential equation represented by the state transform velocity error is as follows:

21

Thus, the derivation of the error prediction model is accomplished.

## APPENDIX C INTRODUCTION OF THE KALMAN FILTER

For a given set of state variables * x_{k}* =

*x*_{k/k}:=

*(*

**x***t*), system equation

_{k}*(*

**F***t*), and measurement

_{k}*:=*

**z**_{k}**(**

*z**t*), we initialize the initial state and initial covariance:

_{k}22

The hat denotes an estimated variable. To apply the Kalman filter in discrete form, five essential equations are required, corresponding to two steps, namely, the propagation step and the update step. For the propagation step, we apply the following:

23

We apply an update when a measurement is received:

24

where *Φ*_{k/k−1} and *Γ*_{k−1} are a discretized system matrix and a noise driving matrix with discrete time Δ*T*, as given in Equation (25).

25

** Q** and

**are the spectrum density of the device random noise and measurement noise:**

*R*26

This is an open access article under the terms of the Creative Commons Attribution License, which permits use, distribution and reproduction in any medium, provided the original work is properly cited.