A cross product calibration method for micro-electro mechanical system gyroscope in unmanned aerial vehicle attitude determination system

Abstract

This paper presents a novel calibration method for micro-electro mechanical system gyroscope in attitude measurement system of small rotor unmanned aerial vehicles. This method is based on an observation vector and its cross product, which is especially valuable for the in-field calibration without the aid of external equipment. By analysing the error model of the tri-axial gyroscope, the principle of calibration is proposed. Compared with other algorithms, numerical simulations are performed to evaluate the effectiveness of integral form of the cross product calibration method. Experiment on the hex-rotor unmanned aerial vehicle platform shows that the proposed method has great advantages in low-cost integrated navigation system.

Keywords

Micro-electro mechanical system gyroscope in-field calibration cross product calibration observation vector integrated navigation system

Introduction

Attitude determination of unmanned aerial vehicles (UAVs) is essentially an integrated navigation problem.¹ The theory and technique of the calibration for conventional inertial devices (especially in integrated navigation systems (INS)) are well established.^2–5 These calibration methods all rely on high-precision equipment, which can provide accurate attitude reference and angular velocity.^6–10 For an INS equipped with low-end micro-electro mechanical system (MEMS) sensors, using specific and precision equipment to calibrate tri-axial gyroscopes will raise the cost and is usually impractical. Therefore, the above-mentioned calibration methods are probably unsuitable, and using attitude information for calibration is the best strategy.

Liu and Fang proposed a six-position rotation calibration method for field calibration of a fibre optic gyroscope (FOG) inertial measurement unit (IMU).¹¹ The reported experimental results show that the calibration accuracy is at the same level as the conventional method using a turntable. However, it is noteworthy that the three axes of the IMU are required to point to the sky and the earth in turn and rotate 360 degrees in both directions. If there are errors in the actual direction or rotation angle, the calibration results will be directly affected. Solutions to this problem are not described in detail in their research reports.

Wan and Huang proposed a gyroscope field calibration algorithm based on heading angle, pitch angle and roll angle information. The pitch angle and the roll angle were measured by an accelerometer, and the heading angle was obtained by a magnetic compass or the Global Positioning System.¹²

W. T. Fong et al. proposed a method for IMU calibration without external equipment. This method directly refers to gravitational acceleration and realizes the calibration of tri-axial gyroscope.¹³ Cheuk et al.¹⁴ used the W. T. Fong’s method for IMU self-calibration. E. Dorveaux et al. proposed another self-calibration method of a tri-axial gyroscope with reference to gravitational acceleration.^15–20

According to the research of X. Li and Z. Li,²¹ a novel calibration method that employs the invariance of the dot product of two constant vectors is presented for tri-axial field sensors in strap-down navigation systems, and this method has been used to calibrate three-axis magnetometer. But we cannot find an auxiliary vector that has constant dot product with the angular velocity vector, so this method is not suitable for tri-axial gyroscope.

Based on the above analysis, a vector cross product calibration method for a tri-axial gyroscope is proposed, which is expressed as the cross product of the reference vector and the angular velocity. This method is derived from the derivative of the reference vector with respect to time, which can be used in INS.

Methodology

Attitude algorithms for INS rely on three 3D vectors, namely the geomagnetic vector, gravity vector and angular velocity. In ideal cases, the sensor output should be proportional to the measurand, that is, $v = k u$ , in which k is the sensitivity (scaling factor), while $u = (\begin{matrix} u_{1} & u_{2} & u_{3} \end{matrix})^{T}$ and $v = (\begin{matrix} v_{1} & v_{2} & v_{3} \end{matrix})^{T}$ indicate the vector being measured and the sensor output, respectively. However, the sensor usually has different error sources that come from manufacture and fabrication, such as bias, scale factor, misalignment and noise.^21,22 A generalized error model of the tri-axial vector sensor can be described by equation (1)

\begin{matrix} (\begin{matrix} v_{1} \\ v_{2} \\ v_{3} \end{matrix}) = (\begin{matrix} k_{11} & k_{12} & k_{13} \\ k_{21} & k_{22} & k_{23} \\ k_{31} & k_{32} & k_{33} \end{matrix}) (\begin{matrix} u_{1} \\ u_{2} \\ u_{3} \end{matrix}) + (\begin{matrix} b_{1} \\ b_{2} \\ b_{3} \end{matrix}) + (\begin{matrix} ε_{1} \\ ε_{2} \\ ε_{3} \end{matrix}) \\ = Ku + b + ε_{0} + H . O . T . \end{matrix}

(1)

In equation (1), the vector $ε = (\begin{matrix} ε_{1} & ε_{2} & ε_{3} \end{matrix})^{T}$ represents sensor noise and higher order error terms, and it can be further divided into the sensor noise $ε_{0}$ and the higher order terms (i.e. the non-linear part) H.O.T. The calibration of low-end three-axis gyroscope is usually aimed at constant error and linear error, that is, K and b in equation (1), while the terms of ε ₀ and higher order error are neglected. Thus, equation (1) will become a linear error model as $v = Ku + b$ , in which all the constant errors can be included in $b = (\begin{matrix} b_{1} & b_{2} & b_{3} \end{matrix})^{T}$ , and all the errors that are proportional to vector u can be described by the matrix $K = (k_{ij})_{3 \times 3}$ .

For a low-end tri-axial gyroscope, if the angular velocity in the aircraft coordinate frame is ω , while the output of gyroscope is $\tilde{ω}$ , the error model can be simplified as equation (2)

\begin{matrix} (\begin{matrix} {\tilde{ω}}_{1} \\ {\tilde{ω}}_{2} \\ {\tilde{ω}}_{3} \end{matrix}) = (\begin{matrix} k_{11} & k_{12} & k_{13} \\ k_{21} & k_{22} & k_{23} \\ k_{31} & k_{32} & k_{33} \end{matrix}) (\begin{matrix} ω_{1} \\ ω_{2} \\ ω_{3} \end{matrix}) + (\begin{matrix} b_{1} \\ b_{2} \\ b_{3} \end{matrix}) \\ \tilde{ω} = K ω + b \end{matrix}

(2)

For a time-invariant vector u in the geographic coordinate frame, its time derivative can be calculated by its cross product with ω , as shown in equation (3)

\overset{\cdot}{u} = u \times ω

(3)

Equation (3) can also be written in the form of matrix multiplication, as shown in equation (4)

(\begin{matrix} {\overset{\cdot}{u}}_{1} \\ {\overset{\cdot}{u}}_{2} \\ {\overset{\cdot}{u}}_{3} \end{matrix}) = (\begin{matrix} 0 & - ω_{3} & ω_{2} \\ ω_{3} & 0 & - ω_{1} \\ - ω_{2} & ω_{1} & 0 \end{matrix}) (\begin{matrix} u_{1} \\ u_{2} \\ u_{3} \end{matrix}) \equiv u [ω \times]

(4)

In equation (4), the matrix $[ω \times]$ on the right is called ‘cross product matrix’, which corresponds to $ω$ . For any non-zero vector u , $u [ω \times] = u \times ω$ if u is written as a row vector, and $[ω \times] u = ω \times u$ if u is written as a column vector. If K and b are determined, we denote $L = K^{- 1}$ and $d = K^{- 1} b$ , then substitute the error model equation (2) into equation (4), and it leads to equation (5)

\overset{\cdot}{u} = u [(L \tilde{ω} - d) \times]

(5)

According to equation (5), if u and $\tilde{ω}$ are both available, L and d can be determined, and the calibration of gyroscope can be implemented. Equation (5) is the basis of the vector-aided calibration method.

In MEMS-based INS, the vector u can be provided by calibrated tri-axial magnetometer or accelerometer, that is, the geomagnetic or gravity vector can serve as u . In the process of data acquisition, the INS needs to be rotated to provide an excitation signal, that is, a non-zero angular velocity.

Equation (5) is in a differential form. Naturally we can transform equation (5) into an integral form. Assuming that the starting and ending time of the rotation process are $t_{0}$ and T, respectively, the change of vector u (denoted as $Δ u$ ) can be calculated by equation (6)

Δ u = \int_{t_{0}}^{T} u [(L \tilde{ω} - d) \times] dt

(6)

To get a closed-form solution of the elements in L and d , the symbol $I (u_{i} {\tilde{ω}}_{j})$ is introduced to represent the integral $I (u_{i} {\tilde{ω}}_{j}) \equiv \int_{t_{0}}^{T} u_{i} {\tilde{ω}}_{j} dt$ , and equation (6) can also be rewritten in a component-wise form, which is shown in equation (7)

\begin{matrix} (\begin{matrix} Δ u_{1} \\ Δ u_{2} \\ Δ u_{3} \end{matrix}) = \\ [\begin{matrix} 0 & 0 & 0 & I (- u_{3} {\tilde{ω}}_{1}) & I (- u_{3} {\tilde{ω}}_{2}) & I (- u_{3} {\tilde{ω}}_{3}) & I (u_{2} {\tilde{ω}}_{1}) & I (u_{2} {\tilde{ω}}_{2}) & I (u_{2} {\tilde{ω}}_{3}) & 0 & I (u_{3}) & I (- u_{2}) \\ I (u_{3} {\tilde{ω}}_{1}) & I (u_{3} {\tilde{ω}}_{2}) & I (u_{3} {\tilde{ω}}_{3}) & 0 & 0 & 0 & I (- u_{1} {\tilde{ω}}_{1}) & I (- u_{1} {\tilde{ω}}_{2}) & I (- u_{1} {\tilde{ω}}_{3}) & I (- u_{3}) & 0 & I (u_{1}) \\ I (- u_{2} {\tilde{ω}}_{1}) & I (- u_{2} {\tilde{ω}}_{2}) & I (- u_{2} {\tilde{ω}}_{3}) & I (u_{1} {\tilde{ω}}_{1}) & I (u_{1} {\tilde{ω}}_{2}) & I (u_{1} {\tilde{ω}}_{3}) & 0 & 0 & 0 & I (u_{2}) & I (- u_{1}) & 0 \end{matrix}] \\ {[\begin{matrix} l_{11} & l_{12} & l_{13} & l_{21} & l_{22} & l_{23} & l_{31} & l_{32} & l_{33} & d_{1} & d_{2} & d_{3} \end{matrix}]}^{T} \end{matrix}

(7)

Equation (7) is a linear expression of all the elements in L and d , and thus it is a linear regression problem.

Equations (5) and (6) are both derived from the cross product in equation (3), and the only difference is that equation (5) is in differential form, while equation (6) is in integral form. We will focus on the integral form in the following discussion.

It can be seen that the above-mentioned cross product calibration method shares the same principle with the calibration method described in previous studies,^13–21 but the proposed method greatly simplifies the calculation in the calibration process. Meanwhile, the proposed calibration method only needs a time-invariant vector u (gravity vector g or geomagnetic vector h ) without using turntable or other precision equipment. Moreover, when using the magnetometer or accelerometer to provide the reference vector u in calibration, the misalignment between gyroscope and magnetometer or accelerometer can be eliminated simultaneously.

Numerical simulation and experiment

In the following simulations, we use the North-East-Down (NED) coordinate frame and use the gravity vector g ₀ = (0 0 9.8)^T m s⁻² as the reference vector u .

We use the trapezoidal method (described in Särkkä et al.²³) to calculate the integration in equation (7) numerically. The sampling rate is 100 Hz, and the corresponding time interval is 0.01 s.

Rotation arrangement in calibration

To implement the gyro calibration, we first collect raw data during several rotations of the UAV. First, we point the x, y and z axis of the IMU vertically upwards and downwards, respectively, and thus there are six different cases. Then, we rotate the IMU both clockwise and counter-clockwise around x, y and z axis, respectively, so that a total of 36 different rotations can be obtained. Furthermore, in each rotation, the angular velocity increases linearly from 0 to 90°/s within 1 s and then decreases linearly to zero again within 1 s. So the rotation angle is exactly 90.

According to the above arrangement, there will be 36 rotations corresponding to different combinations of the angular velocity ω and the reference vector u . However, among these 36 rotations, there are 12 rotations in which the angular velocity is parallel or anti-parallel to the reference vector, and thus the reference vector u is unchanged during the rotation. In such cases, both sides of equations (5) and (6) become zero. In other words, these 12 rotations do not provide effective data for gyroscope calibration. In fact, this kind of rotation is called ‘useless rotation’ in Li et al.²⁴ In contrast, the reference vector and the angular velocity are orthogonal to each other in the other 24 rotations, and such rotations can be referred to as ‘effective rotations’, as shown in Figure 1.

Figure 1.

Useless rotation versus effective rotation for gyroscope calibration.

As stated in section ‘Introduction,’ the principle of the proposed cross product calibration method is quite similar to that of the Fong’s method.¹³ In fact, if the reference vector u on the right side of equation (7) is calculated by differential equation method, equation (7) can be regarded as a first-order approximation of Fong’s method. However, in Fong’s method, the reference vector u is calculated by INS updating algorithm (such as the fourth-order Runge-Kutta algorithm in the following simulation), rather than provided by accelerometer. As a result, the elements of L and d will be accumulated in the integrals at the right side of equation (7) and that will inevitably lead to deviation of the calibration results. What is more, it can be expected that such deviation will increase with larger gyroscope error.

To compare Fong’s calibration method with the proposed algorithm, we carry out the following numerical simulations. For Fong’s method, we use the fourth-order Runge-Kutta (R-K) algorithm to calculate the reference vector (i.e. the gravity vector), as shown in equation (8)

{\begin{matrix} g_{k + 1} = g_{k} + \frac{h}{6} (k_{1} + 2 k_{2} + 2 k_{3} + k_{4}) \\ k_{1} = g_{k} \times {\tilde{ω}}_{k} \\ k_{2} = (g_{k} + \frac{h}{2} k_{1}) \times {\tilde{ω}}_{k + 0.5} \\ k_{3} = (g_{k} + \frac{h}{2} k_{2}) \times {\tilde{ω}}_{k + 0.5} \\ k_{4} = (g_{k} + h k_{3}) \times {\tilde{ω}}_{k + 1} \end{matrix}

(8)

Equation (8) describes the renewal process of gravity vector from the moment k to k + 1. In equation (8), h is the time interval of the numerical algorithm. It should be noted that ${\tilde{ω}}_{k + 0.5}$ in equation (8) represents the measured angular velocity at h/2 after the moment k. In other words, the sampling rate of the gyroscope must be twice of that of the reference vector. Therefore, for Fong’s method, the sampling rate of gyroscope is set to 200 Hz.

We assume that the gyro has the following errors

\begin{matrix} L = (\begin{matrix} 1.070 & 0.005 & - 0.020 \\ - 0.015 & 0.990 & 0.030 \\ 0.025 & - 0.010 & 0.960 \end{matrix}), \\ b = (\begin{matrix} 3.000 \\ - 1.000 \\ - 2.000 \end{matrix}) (° / s) \end{matrix}

We first compare the proposed cross product calibration method with Fong’s method with different rotation arrangements, that is, the above-mentioned 36 and 24 rotations. The simulation results are presented in Table 1.

Table 1.

Calibration results under different rotation arrangements.

Rotation arrangement		36 rotations, range 90°	24 rotations (without the useless rotations), range 90°
Cross product calibration (integral form)	L	$(\begin{matrix} 1.0700 & 0.0050 & - 0.0200 \\ - 0.0150 & 0.9900 & 0.0300 \\ 0.0250 & - 0.0100 & 0.9600 \end{matrix})$	$(\begin{matrix} 1.0700 & 0.0050 & - 0.0200 \\ - 0.0150 & 0.9900 & 0.0300 \\ 0.0250 & - 0.0100 & 0.9600 \end{matrix})$
Cross product calibration (integral form)	b (°/s)	${(\begin{matrix} 3.0000 & - 1.0000 & - 2.0000 \end{matrix})}^{T}$	${(\begin{matrix} 3.0000 & - 1.0000 & - 2.0000 \end{matrix})}^{T}$
Imitating Fong’s method	L	$(\begin{matrix} 1.0526 & 0.0190 & - 0.0236 \\ - 0.0261 & 0.9918 & 0.0023 \\ 0.0274 & - 0.0136 & 0.9674 \end{matrix})$	$(\begin{matrix} 1.0535 & 0.0047 & - 0.0174 \\ - 0.0141 & 0.9933 & 0.0297 \\ 0.0271 & - 0.0097 & 0.9686 \end{matrix})$
Imitating Fong’s method	b (°/s)	${(\begin{matrix} 2.5986 & - 0.5179 & - 1.4620 \end{matrix})}^{T}$	${(\begin{matrix} 2.7606 & - 0.9423 & - 1.8104 \end{matrix})}^{T}$

According to Table 1, the performance of the cross product calibration method is good with either 36 or 24 rotations, but Fong’s method shows unsatisfactory results, especially with 36 rotations. That means Fong’s method is more significantly affected by the useless rotations.

Next, we evaluate the impacts of the rotation angle, which is also an important issue. It can be expected that the cumulative error and its influence on Fong’s calibration method will increase with larger rotation angle. We set the rotation angle to 180°, 270° and 360°, respectively, and the simulation results are presented in Table 2.

Table 2.

Simulation results for different rotation ranges.

Rotation arrangement		24 rotations, range 180°	24 rotations, range 270°	24 rotations, range 360°
Cross product calibration (integral form)	L	$(\begin{matrix} 1.0700 & 0.0050 & - 0.0200 \\ - 0.0150 & 0.9900 & 0.0300 \\ 0.0250 & - 0.0100 & 0.9600 \end{matrix})$	$(\begin{matrix} 1.0699 & 0.0050 & - 0.0200 \\ - 0.0150 & 0.9899 & 0.0300 \\ 0.0250 & - 0.0100 & 0.9599 \end{matrix})$	$(\begin{matrix} - 0.1909 & - 0.0060 & 0.0107 \\ 0.0038 & - 0.1675 & - 0.0030 \\ - 0.0049 & 0.0023 & - 0.1643 \end{matrix})$
Cross product calibration (integral form)	b (°/s)	${(\begin{matrix} 3.0000 & - 1.0000 & - 2.0000 \end{matrix})}^{T}$	${(\begin{matrix} 3.0000 & - 1.0000 & - 2.0000 \end{matrix})}^{T}$	${(\begin{matrix} 3.0000 & - 1.0000 & - 2.0000 \end{matrix})}^{T}$
Imitating Fong’s method	L	$(\begin{matrix} 0.9988 & - 0.0011 & - 0.0149 \\ - 0.0137 & 1.0029 & 0.0306 \\ 0.0286 & - 0.0045 & 1.0004 \end{matrix})$	$(\begin{matrix} 0.8795 & - 0.0244 & - 0.0213 \\ - 0.0285 & 1.0261 & 0.0380 \\ 0.0036 & 0.0303 & 1.0948 \end{matrix})$	$(\begin{matrix} 0.0000 & 0.0000 & 0.0000 \\ 0.0000 & 0.0001 & 0.0000 \\ 0.0000 & 0.0000 & 0.0001 \end{matrix})$
Imitating Fong’s method	b (°/s)	${(\begin{matrix} 1.4202 & - 0.6092 & - 0.8712 \end{matrix})}^{T}$	${(\begin{matrix} - 10.9139 & 2.3808 & 1.9687 \end{matrix})}^{T}$	${(\begin{matrix} 8.6989 & - 2.6782 & - 6.0390 \end{matrix})}^{T}$

According to Table 2, the result of Fong’s method is still poorer than that of the proposed method when the rotation angle is 180° or 270°. However, it is particularly noteworthy that both methods cannot calculate the matrix L correctly when the rotation angle is 360°. This is mainly because the reference vector has no change after a 360° rotation and thus it can be viewed as another kind of ‘useless rotation’.

According to the above simulation results, the proposed cross product calibration method can give accurate results in most cases, but the rotation range should be less than 360°. Therefore, when using the proposed cross product calibration method, a sufficient and feasible rotation arrangement is as follows:

Point the three axes of the IMU vertically upwards and downwards, respectively (six different cases in total).

In each of the above cases, rotate the IMU around the two horizontal axes both clockwise and counter-clockwise (see Figure 1, 24 rotations in total).

The angle of each rotation is 90°. We will use this arrangement in the following numerical simulations.

Impacts of sensor noise and error magnitude

To make the test more comprehensive, we randomly select 1000 error samples for the numerical simulation. The distortion matrix is generated as $L = I_{3 \times 3} + D$ , where $I_{3 \times 3}$ is an identity matrix, and each element of D is zero-mean and normally distributed with the standard deviation $σ_{D} = 0.01$ . Besides, the biases are also assumed to be normally distributed with the standard deviation $σ_{b} = 3^{°} / s$ .

To evaluate the calibration results of the matrix L , we use an evaluation factor based on the Frobenius matrix norm. If the presumed distortion matrix is L _c, and the simulation result is L _p, the evaluation factor will be defined as

Ev (L) = \frac{{| I - L_{c} L_{p}^{- 1} |}_{F}}{{| I - L_{p}^{- 1} |}_{F}}

(9)

In equation (9), the symbol|·|_F indicates the Frobenius norm of a matrix. For any 3×3 matrix A , | I − A |_F reflects the difference between A and Identity matrix I . Obviously, if L _p is close to L _c, Ev( L ) will be close to zero. On the contrary, if Ev( L ) is close to (or even larger than) 1, it means that the calibration method cannot eliminate the error of gyroscope.

On the contrary, we evaluate the calibration result of the bias by the relative error Ev( b ), as shown in equation (10)

Ev (b) = \frac{| b_{c} - b_{p} |}{| b_{p} |} \times 100 %

(10)

In equation (10), b _p and b _c denote the presumed bias and the simulation result, respectively.

First, we add a Gaussian white noise item to the measurement of gyroscope, with the standard deviation $σ_{ω}$ varying from 0.2 to 2.0°/s. The simulation results are shown in Figure 2. As can be seen from Figure 2, when the noise of gyroscope enlarges, the evaluation factors Ev( L ) and Ev( b ) of the two calibration methods both increase, but the results of the proposed method are much better than those of Fong’s method.

Figure 2.

Simulation results at different gyroscope noise levels (PM = proposed method, FM = Fong’s method).

Second, we evaluate the calibration results when the reference vector has noise. We add a Gaussian white noise item to the gravity vector, and the standard deviation varies from 1 to 10 mg. The calibration effect is shown in Figure 3. Once again, the proposed method shows better performance than that of Fong’s method. Nevertheless, we can also see that the noise of reference vector has greater impacts to the proposed method, since this noise item will affect both sides of equation (7).

Figure 3.

Simulation results at different reference vector noise levels (PM = proposed method, FM = Fong’s method).

Moreover, we test the effectiveness of both calibration methods when dealing with different magnitudes of sensor errors. In this test, σ_D varies from 0.02 to 0.1, and σ_b changes from 2.5 to 10°/s, while the noise level $σ_{ω}$ is set to 0.1°/s. The simulation results are presented in Figures 4 and 5. It is clear that the proposed method always shows better performance compared to Fong’s method, even with large σ_D and σ_b.

Figure 4.

Simulation results of Ev( L ) versus σ_D and σ_b: (a) the proposed method and (b) Fong’s method.

Figure 5.

Simulation results of Ev( b ) versus σ_D and σ_b: (a) the proposed method and (b) Fong’s method.

Experiment on UAV

We apply all the above-mentioned methods to a hex-rotor UAV, which consists of a tri-axial magnetometer (HMC5883), a tri-axial accelerometer (ADXL345) and a MEMS gyroscope (L3G4200D), as shown in Figure 6.

Figure 6.

Flight controller with MEMS gyro.

In the experiment, the UAV is rotated in hand, and different calibration methods of raw data are collected. Data are collected in 24 rotations following the recommended procedure, that is, approximately 90° rotations around the horizontal axes. The accelerometer ADXL345 is calibrated before the experiment, in order to serve as the reference for the calibration of the gyroscope, and then L and b calculated by different calibration methods are compensated.

Finally, we installed the UAV on the test platform for static live test and observe hovering test, as shown in Figure 7.

Figure 7.

Hex-rotor UAV.

Figure 8 shows the comparison of gyroscope yaw angular deviations before and after calibration in the static live test. In Figure 8, the blue line represents the yaw angle calculated according to uncalibrated gyroscope data, which shows the drift up to 25.6° in 1 min. The red line represents the yaw angle according to gyroscope data calibrated by Fong’s method, and the drift error can keep within 8.0°. The green line represents the error of yaw angle after calibration with cross product algorithm, and the effect is better than that of Fong’s method. It can be seen from Figure 8 that the attitude angle would accumulate drift after the solution, because there still have been some errors after the calibration of the UAV gyroscope. The corrected sampling value can be combined with the data fusion algorithm,^25,26 for example, Kalman filter and complementary filter, to suppress the accumulated errors and ensure the attitude stability for long-time flight.

Figure 8.

Angular deviation for static live test.

Next, the UAV remote control throttle remains 50% for hovering test. Figure 9 shows the comparison of yaw angular deviations after calibration with two methods. As we can see in Figure 9, the red line represents the yaw angle obtained by Fong’s method, and the drift error of yaw angle can remain within 1.5°. The green line represents the error of yaw angle after calibration with the cross product algorithm, and the effect is obviously better than that of Fong’s method.

Figure 9.

Angular deviation for hovering test.

Conclusion

We introduce a gravity-aided calibration method for low-end MEMS gyros in this paper. The proposed method essentially unifies the calibration method based on gravity vector, which makes its expression more concise and greatly simplifies the calculation in the calibration process. This method requires no external equipment and is particularly suitable for the in-field calibration of a MEMS gyroscope in a UAV. Although its precision may be limited in practice, this method is efficient to provide a good initial value for the online parameter estimation.

The proposed method is based on the time derivative of the gravity vector, and the integral form cross product is used. Numerical simulations are carried out to test the performance of the cross product calibration method and to imitate Fong’s method. In the experiment on a UAV platform, the results show that the proposed method can offer higher accuracy in most cases.

The basis of the proposed method, that is, the time derivative of the gravity vector, can be replaced by the time derivative of any constant vector in the aircraft coordinate frame. Therefore, such a vector can be utilized instead of the gravity, and the proposed method will then be evolved to a vector-aided scheme.

Footnotes

Declaration of conflicting interests

The author(s) declared no potential conflicts of interest with respect to the research, authorship and/or publication of this article.

Funding

The author(s) disclosed receipt of the following financial support for the research, authorship and/or publication of this article: This work was supported by National Natural Science Foundation of China Grant #61361006 to Z.L., Guangxi Natural Science Foundation Grant #2015GXNSFBA139251 to Y.W. and Guangxi Key Laboratory of Automatic Detecting Technology and Instruments Grant #YQ15107 to X.L.

ORCID iDs

Yongjun Wang

Xiang Li

Data availability

The datasets used or analysed during the current study are available from the corresponding author on reasonable request. The data used to support the findings of this study are included within the article.

References

Dorveaux

Boudot

Hillion

, et al. Combining inertial measurements and distributed magnetometry for motion estimation. In: American control conference, San Francisco, CA, 29 June–1 July 2011, pp. 4249–4256. New York: IEEE.

ANSI/IEEE STD 293-1969. IEEE test procedure for single-degree-of-freedom spring-restrained rate gyros.

IEEE STD 529-1980. Supplement for strapdown applications to IEEE standard specification format guide and test procedure for single-degree-of-freedom rate-integrating gyros.

IEEE STD 952-1997. IEEE standard specification format guide and test procedure for single-axis interferometric fiber optic gyros.

IEEE STD 647-2006. IEEE standard specification format guide and test procedure for single-axis laser gyros.

Chatfield

. Notation, coordinate systems, and units. In: Chatfield

(ed.) Fundamentals of high accuracy inertial navigation [Progress in Astronautics & Aeronautics]. Reston, VA: AIAA, 1997, pp. 15–32.

Titterton

Weston

. Strapdown inertial navigation technology. 2nd ed. Reston, VA: AIAA; Stevenage: IET, 2004, pp. 148–173.

Syed

Aggarwal

Goodall

, et al. A new multi-position calibration method for MEMS inertial navigation systems. Meas Sci Technol 2007; 18(7): 1897–1907.

Frosio

Pedersini

Borghese

. Autocalibration of MEMS accelerometers. IEEE T Instrum Meas 2009; 58(6): 2034–2041.

10.

Won

Golnaraghi

. A triaxial accelerometer calibration method using a mathematical model. IEEE T Instrum Meas 2010; 59(8): 2144–2153.

11.

Liu

Fang

. Novel field calibration through rotation in six-position for FOG-IMU. Opto-Electron Eng 2008; 35(1): 60–65.

12.

Liang

Huang

. A field calibration approach for gyro errors. Comput Simul 2008; 25(6): 50–53.

13.

Fong

Ong

Nee

AYC

. Methods for in-field user calibration of an inertial measurement unit without external equipment. Meas Sci Technol 2008; 19(8): 817–822.

14.

Cheuk

Lau

Lin

, et al. Automatic calibration for inertial measurement unit. In: International conference on control automation robotics & vision, Guangzhou, China, 5–7 December 2012, pp. 1341–1346. New York: IEEE.

15.

Dorveaux

. Magneto-inertial navigation: principles and application to an indoor pedometer. Automatic. PhD Thesis, École Nationale Supérieure des Mines de Paris, Paris, 2012.

16.

Chow

Hol

Luinge

. Tightly-coupled joint user self-calibration of accelerometers, gyroscopes, and magneto-meters. Drones 2018; 6(2): 1–20.

17.

. Vector-aided in-field calibration method for low-end MEMS gyros in attitude and heading reference systems. IEEE T Instrum Meas 2014; 63(11): 2675–2681.

18.

Rosario

Lovell

Redmond

. Quaternion-based complementary filter for attitude determination of a smartphone. IEEE Sens J 2016; 16(15): 6008–6017.

19.

Ghasemi-Moghadam

Homaeinezhad

. Attitude determination by combining arrays of MEMS accelerometers, gyros, and magnetometers via quaternion-based complementary filter. Int J Numer Model 2017; 31(3): 2282.

20.

Qureshi

Golnaraghi

. An algorithm for the in-field calibration of a MEMS IMU. IEEE Sens J 2017; 17: 7479–7486.

21.

. A new calibration method for tri-axial field sensors in strap-down navigation systems. Meas Sci Technol 2012; 23(10): 2852–2855.

22.

Keighobadi

Vosoughi

Faraji

. Design and implementation of a model predictive observer for AHRS. GPS Solut 2018; 22(1): 29.

23.

Särkkä

Nieminen

Suuriniemi

, et al. A multi-position calibration method for consumer-grade accelerometers, gyroscopes, and magneto-meters to field conditions. IEEE Sens J 2017; 17(11): 3470–3481.

24.

Niu

Zhang

, et al. An in situ hand calibration method using a pseudo-observation scheme for low-end inertial measurement units. Meas Sci Technol 2012; 23(10): 105104.

25.

Tian

Chao

Flanagan

, et al. Design and evaluation of UAV flow angle estimation filters. IEEE T Aero Elec Sys 2019; 55(1): 371–383.

26.

Zhang

, et al. Attitude anti-interference federal filtering algorithm for MEMS-SINS/GPS/magnetometer/SV integrated navigation system. Meas Control 2020; 53(1): 46–60.