Movable Surface Rotation Angle Measurement System Using IMU

Abstract

In this paper, we describe a rotation angle measurement system (RAMS) based on an inertial measurement unit (IMU) developed to measure the rotation angle of a movable surface. The existing IMU-based attitude (tilt) sensor can only accurately measure the rotation angle when the rotation axis of the movable surface is perfectly aligned with the X axis or Y axis of the sensor, which is always not possible in real-world engineering. To overcome the difficulty of sensor installation and ensure measurement accuracy, first, we build a model to describe the relationship between the rotation axis and the IMU. Then, based on the built model, we propose a simple online method to estimate the direction of the rotation axis without using a complicated apparatus and a method to estimate the rotation angle using the known rotation axis based on the extended Kalman filter (EKF). Using the estimated rotation axis direction, we can effectively eliminate the influence of the mounting position on the measurement results. In addition, the zero-velocity detection (ZVD) technique is used to ensure the reliability of the rotation axis direction estimation and is used in combination with the EKF as the switching signal to adaptively adjust the noise covariance matrix. Finally, the experimental results show that the developed RAMS has a static measurement error of less than 0.05° and a dynamic measurement error of less than 1° in the range of ±180°.

1. Introduction

The movable surfaces of an aircraft (such as elevators, flaps, and slats) are aerodynamic devices allowing a pilot to adjust and control the attitude of an aircraft. Movable surfaces are attached to the airframe on hinges or tracks (see Figure 1) and must be rotated to the angle assigned by the pilot to ensure the stability of the aircraft. Therefore, it is essential to accurately measure the rotation angles of the movable surfaces during an aircraft’s ground test.

The principle of measuring the rotation angle of a movable surface is equivalent to measuring the rotation angle of a rigid body rotating around a fixed axis, which has a specific direction. Optical motion capture systems are suitable for this task, however, these are expensive and require high computational costs [1]. Inertial measurement units (IMUs) containing tri-axis accelerometers and tri-axis gyroscopes are also suitable for this task. In recent years, microelectromechanical system (MEMS) IMUs have become widely available in many areas, such as state estimation and angle measurement, due to their small size and low cost [2,3].

Lian Hu et al. [4] proposed a Kalman filter-based algorithm that uses an IMU to estimate the roll angle of an agricultural machine in paddy fields. Milad Ghanbari et al. [5] proposed a tilt measurement algorithm consisting of a modified Kalman filter (KF) as a postfilter and a complementary filter as a prefilter to enhance the accelerometer bandwidth and eliminate the gyroscope drift. In ref. [6], the attitude (pitch and roll) was obtained by fusing the acceleration and angular velocity using a nonlinear complementary filter. In refs. [7,8,9], a KF-based algorithm using an IMU was developed to estimate the attitude during dynamic conditions. Although these methods can solve the attitude estimation problem, the pitch or roll is not equal to the rotation angle around the fixed rotation axis when neither the X nor Y axis of the IMU are aligned with the rotation axis of the movable surface.

Combining an IMU and a tri-axis magnetometer can form an attitude and heading reference system (AHRS) that can provide a complete orientation measurement relative to the direction of gravity and the Earth’s magnetic field [10,11,12]. We can parameterize the orientation measured by AHRS as two quantities: a unit vector indicating the direction of a rotation axis and an angle $\theta$ describing the magnitude of the rotation about the axis, where $\theta$ is the rotation angle of the movable surface. However, the performance of the AHRS degrades substantially due to the presence of significant magnetic disturbances in the aircraft ground test environment.

Aligning the measurement axis of the IMU with the rotation axis is hard. Estimating the mounting error (the misalignment of the IMU relative to the rotation axis) and correcting the measured pitch or roll is an effective method to measure the rotation angle. Although mounting error estimation methods can be found in numerous studies [13,14,15], these methods require an apparatus such as turntables and have to be performed in the laboratory. In ref. [16], a novel method to estimate mounting error was proposed, which requires only a tri-axis accelerometer. However, the rotation angle is obtained by solving a model equation, so the sensor noises directly affect the measurement results. Furthermore, it is impossible to measure the rotation angle while the movable surface is rotating without a gyroscope.

The primary objective of this research was to develop an IMU-based rotation angle measurement system (RAMS) that can accurately measure the rotation angle without using high-cost instruments and the RAMS can be randomly mounted on a movable surface. The measurement algorithm used in RAMS consists of two parts: rotation axis direction estimation and rotation angle estimation. The first part only requires the movable surface to be stationary for a few seconds in three different rotational positions after the RAMS is mounted on the surface while recording the accelerometer data. We can calculate the rotation axis direction in the RAMS (IMU) frame of reference from the recorded data. In the second part, an adaptive extended Kalman filter (EKF) is used to solve the rotation angle estimation problem. We use acceleration data to correct the rotation angle obtained from a gyroscope integration, where the measurement model in the EKF is built using the estimated rotation axis direction. Furthermore, the zero-velocity detection (ZVD) technique is used in both the rotation axis direction estimation and rotation angle estimation, which aims to detect whether the RAMS is stationary. The ZVD determines the moment of recording acceleration data during the rotation axis direction estimation and the moment of switching the EKF noise covariance matrix during the rotation angle estimation. The proposed algorithm is compared with the method in ref. [16], and the measurement accuracy of the RAMS is compared with a commercial high-precision AHRS (Ellipse2-N from SBG). The results show that the RAMS is able to measure the rotation angle within ±180°, and the maximum measurement error is less than 0.05°, even with it randomly mounted on a movable surface.

The rest of this article is organized as follows. Section 2 develops a model to describe the relationship between the IMU output and the rotation axis of the movable surface. Section 3 gives details of the proposed algorithm used in RAMS. Section 4 introduces the design for the RAMS. Section 5 performs a simulation to evaluate the performance of the proposed algorithm and the effect of major error sources. Section 6 presents the experimental results. Section 7 concludes this article with a summary of the designed measurement system.

2. Modeling

2.1. Problem Definition

In this section, we used IMU instead of RAMS since the IMU is the core of measuring the rotation angle of a movable surface. It is required to mount the IMU on the surface to measure the rotation angle of the aircraft’s movable surface. The following notation for frames of reference is used (see Figure 2):

• A denotes the frame fixed on the Earth.
• B denotes the IMU frame of reference.
• C denotes the movable surface frame of reference.

We assume that the X–Y plane of frame A is the horizontal plane, the X–Y plane of frame B is the measurement plane, and the Z axis of frame C is the rotation axis of the movable surface.

The coordinate transformation of the vector $\mathbf{e}\in {\mathbb{R}}^3$ between two different frames (e.g., A and B) is

$${}^B\mathbf{e}={}_B^A{\mathbf{R}}^T{}^A\mathbf{e}$$

(1)

where the left superscript of $\mathbf{e}$ implies that the corresponding vectors are expressed in different frames, and ${}_B^A\mathbf{R}\in {\mathbb{R}}^{3\times 3}$ is the rotation matrix of the frame B with respect to frame A.

The vector ${}^B\mathbf{v}={\left[\begin{array}{ccc}{v}_x& v_y& v_z\end{array}\right]}^T$ in frame B is used to represent the direction of the rotation axis with respect to the IMU, which satisfies the following conditions:

• $\left|\right|{}^B\mathbf{v}\left|\right|=\sqrt{v_x^2+v_y^2+v_z^2}=1$.
• $v_x\ge 0$, and if $v_x=0$, $v_y>0$.

where $\left|\right|·\left|\right|$ denotes the norm of the vector.

The rotation matrix $\mathbf{R}$ is expressed using a unit vector $\top ={\left[\begin{array}{ccc}x& y& z\end{array}\right]}^T$ and a rotation angle $\theta$ as [17]

$$\begin{array}{cc}\mathbf{R}\hfill & =R(\theta ,\top )\hfill \\ \hfill & =\left(\cos \theta \right)\mathbf{I}+\left(\sin \theta \right){\left[\top \right]}_{\times }+\left(1-\cos \theta \right)\left(\top {\top }^T\right)\hfill \\ \hfill & \begin{array}{c}\hfill =\left[\begin{array}{cc}\cos \theta +\left(1-\cos \theta \right)x^2& \left(1-\cos \theta \right)xy-\left(\sin \theta \right)z\\ \left(1-\cos \theta \right)xy+\left(\sin \theta \right)z& \cos \theta +\left(1-\cos \theta \right)y^2\\ \left(1-\cos \theta \right)zx-\left(\sin \theta \right)y& \left(1-\cos \theta \right)zy+\left(\sin \theta \right)x\end{array}\right.\\ \hfill \left.\begin{array}{c}\left(1-\cos \theta \right)zx+\left(\sin \theta \right)y\\ \left(1-\cos \theta \right)zy-\left(\sin \theta \right)x\\ \cos \theta +\left(1-\cos \theta \right)z^2\end{array}\right]\end{array}\hfill \end{array}$$

(2)

where $R(·)$ denotes the conversion of an axis-angle to a rotation matrix, and ${[·]}_{\times }$ denotes conversion of a vector to a skew-symmetric matrix.

The purpose of this paper was to estimate the rotation axis ${}^B\mathbf{v}$ and rotation angle $\theta$ of a movable surface.

2.2. Relationship between Accelerometer, Rotation Axis, and Rotation Angle

The tri-axis accelerometer can measure the external specific force acting on the sensor. When the accelerometer is stationary, the measured acceleration is the gravitational acceleration ${}^B\mathbf{g}={\left[\begin{array}{ccc}{a}_x& a_y& a_z\end{array}\right]}^T$. The distribution of gravitational acceleration ${}^A\mathbf{g}={\left[\begin{array}{ccc}0& 0& g\end{array}\right]}^T$ on each axis as depicted in Figure 3, and the following formulas are valid [18]

$$a_x=g\cos \alpha$$

(3)

$$a_y=g\cos \beta$$

(4)

$$a_z=g\cos \phi$$

(5)

where $\alpha$, $\beta$, and $\phi$ denote the angles included between gravitational acceleration and X, Y, and Z axes, respectively.

If both the measurement plane and the rotation axis are horizontal, the tilt angle can be measured in the stationary state by Equations (3)–(5). In most practical cases, however, the sensor cannot be mounted in an ideal position. The measured gravitational acceleration ${}^B\mathbf{g}$ can be expressed in terms of ${}_B^A\mathbf{R}$ as

$${}^B\mathbf{g}={}_B^A{\mathbf{R}}^T{}^A\mathbf{g}$$

(6)

Then, only consider the stationary state of the accelerometer, the rotation axis and the measured acceleration with respect to the frame C can be expressed as

$${}^C\mathbf{v}={}_C^B{\mathbf{R}}^T{}^B\mathbf{v}={\left[\begin{array}{ccc}0& 0& 1\end{array}\right]}^T$$

(7)

$${}^C\mathbf{g}={}_C^B{\mathbf{R}}^T{}^B\mathbf{g}$$

(8)

where ${}^C\mathbf{v}$ and ${}^C\mathbf{g}$ denote the rotation axis direction and the measured gravitational acceleration in frame C, respectively.

We assume that ${}^C{\mathbf{g}}_0={\left[\begin{array}{ccc}{}^C{a}_{x_0}& {}^C{a}_{y_0}& {}^C{a}_{z_0}\end{array}\right]}^T$ denotes the measured acceleration at the initial moment when the movable surface has not yet rotated, and ${\theta }_t$ denotes the rotation angle of the movable surface at time t. The rotation matrix with respect to the initial moment is

$$R\left({\theta }_t,{}^C\mathbf{v}\right)=\left[\begin{array}{ccc}\cos {\theta }_t& -\sin {\theta }_t& 0\\ \sin {\theta }_t& \cos {\theta }_t& 0\\ 0& 0& 1\end{array}\right]$$

(9)

At this moment, the measured acceleration in the frame C can be expressed as:

$$\begin{array}{cc}{}^C{\mathbf{g}}_t\hfill & =R{\left({\theta }_t,{}^C\mathbf{v}\right)}^T{}^C{\mathbf{g}}_0\hfill \\ \hfill & =\left[\begin{array}{c}{}^C{a}_{x_0}\cos {\theta }_t+{}^C{a}_{y_0}\sin {\theta }_t\\ -{}^C{a}_{x_0}\sin {\theta }_t+{}^C{a}_{y_0}\cos {\theta }_t\\ {}^C{a}_{z_0}\end{array}\right]\hfill \end{array}$$

(10)

We can obtain that the trajectory of ${}^C{\mathbf{g}}_t$ is a circle perpendicular to ${}^C\mathbf{v}$ as the rotation angle ${\theta }_t$ changes (see Figure 4).

For three different times $t_1$, $t_2$, and $t_3$, we can obtain

$${}^C{\mathbf{g}}_{t_1}=\left[\begin{array}{c}{}^C{a}_{x_0}\cos {\theta }_{t_1}+{}^C{a}_{y_0}\sin {\theta }_{t_1}\\ -{}^C{a}_{x_0}\sin {\theta }_{t_1}+{}^C{a}_{y_0}\cos {\theta }_{t_1}\\ {}^C{a}_{z_0}\end{array}\right]$$

(11)

$${}^C{\mathbf{g}}_{t_2}=\left[\begin{array}{c}{}^C{a}_{x_0}\cos {\theta }_{t_2}+{}^C{a}_{y_0}\sin {\theta }_{t_2}\\ -{}^C{a}_{x_0}\sin {\theta }_{t_2}+{}^C{a}_{y_0}\cos {\theta }_{t_2}\\ {}^C{a}_{z_0}\end{array}\right]$$

(12)

$${}^C{\mathbf{g}}_{t_3}=\left[\begin{array}{c}{}^C{a}_{x_0}\cos {\theta }_{t_3}+{}^C{a}_{y_0}\sin {\theta }_{t_3}\\ -{}^C{a}_{x_0}\sin {\theta }_{t_3}+{}^C{a}_{y_0}\cos {\theta }_{t_3}\\ {}^C{a}_{z_0}\end{array}\right]$$

(13)

Considering Equations (11)–(13) together, we can yield the following relationship:

$${}^C\mathbf{v}=\pm \frac{{\left[{}^C{\mathbf{g}}_{t_2}-{}^C{\mathbf{g}}_{t_1}\right]}_{\times }\left[{}^C{\mathbf{g}}_{t_3}-{}^C{\mathbf{g}}_{t_2}\right]}{\left|\right|{\left[{}^C{\mathbf{g}}_{t_2}-{}^C{\mathbf{g}}_{t_1}\right]}_{\times }\left[{}^C{\mathbf{g}}_{t_3}-{}^C{\mathbf{g}}_{t_2}\right]\left|\right|}$$

(14)

Multiplying both sides of Equation (14) by ${}_C^B\mathbf{R}$ yields

$${}^B\mathbf{v}=\pm \frac{{\left[{}^B{\mathbf{g}}_{t_2}-{}^B{\mathbf{g}}_{t_1}\right]}_{\times }\left[{}^B{\mathbf{g}}_{t_3}-{}^B{\mathbf{g}}_{t_2}\right]}{\left|\right|{\left[{}^B{\mathbf{g}}_{t_2}-{}^B{\mathbf{g}}_{t_1}\right]}_{\times }\left[{}^B{\mathbf{g}}_{t_3}-{}^B{\mathbf{g}}_{t_2}\right]\left|\right|}$$

(15)

Therefore, we can conclude that it is only needed to measure the gravitational acceleration at three different rotation positions to calculate the rotation axis direction ${}^B\mathbf{v}$.

When the external acceleration is not considered (the accelerometer is stationary), the relationship between the measured acceleration ${}^B{\mathbf{a}}_m$ and the rotation angle $\theta$ can be expressed as

$${}^B{\mathbf{a}}_m=R{\left(\theta ,{}^B\mathbf{v}\right)}^T{}^B{\overline{\mathbf{a}}}_0+{\mathbf{n}}_a$$

(16)

where ${}^B{\overline{\mathbf{a}}}_0={\left[\begin{array}{ccc}{\overline{a}}_{x_0}& {\overline{a}}_{y_0}& {\overline{a}}_{z_0}\end{array}\right]}^T$ is the mean value of the acceleration measurements over an initialization period with no motion, and ${\mathbf{n}}_a$ denotes the measurement noise, assumed to be Gaussian white noise, ${\mathbf{n}}_a\sim \mathcal{N}\left(0,{\sigma }_a^2{\mathbf{I}}_{3\times 3}\right)$.

2.3. Relationship between Gyroscope, Rotation Axis, and Rotation Angle

As part of the IMU, the tri-axis gyroscope provides measurements of the angular velocity about the X, Y, and Z axes. Gyroscopes are known to be affected by different error terms, such as a measurement noise error and a bias [19]. The relationship between the measured angular velocity ${}^B{\mathbf{\omega }}_m$ and the real angular velocity ${}^B{\widehat{\mathbf{\omega }}}_m$ can be simply modeled as [8,9]

$${}^B{\mathbf{\omega }}_m={}^B\widehat{\mathbf{\omega }}+\mathbf{b}+{\mathbf{n}}_{\omega }$$

(17)

$$\dot{\mathbf{b}}=\mathbf{w}$$

(18)

where $\mathbf{b}$ denotes the gyroscope bias, ${\mathbf{n}}_{\omega }$ denotes the measurement noise, and $\mathbf{w}$ denotes a random walk process. Both ${\mathbf{n}}_{\omega }$ and $\mathbf{w}$ are assumed to be Gaussian white noise, ${\mathbf{n}}_{\omega }\sim \mathcal{N}\left(0,{\sigma }_{\theta }^2{\mathbf{I}}_{3\times 3}\right)$, $\mathbf{w}\sim \mathcal{N}\left(0,{\sigma }_b^2{\mathbf{I}}_{3\times 3}\right)$.

With a fixed rotation axis, the relationship between rotation velocity, rotation axis, and rotation angle can be expressed as

$${}^B\widehat{\mathbf{\omega }}=\pm \left|\right|{}^B\widehat{\mathbf{\omega }}\left|\right|{}^B\mathbf{v}$$

(19)

$$\dot{\theta }={}^B{\mathbf{v}}^T{}^B\widehat{\mathbf{\omega }}={}^B{\mathbf{v}}^T({}^B{\mathbf{\omega }}_m-\mathbf{b}+{\mathbf{n}}_{\omega })$$

(20)

3. Algorithm Description

The whole angle measurement process is performed in two steps. The first step is to estimate the rotation axis direction ${}^B\mathbf{v}$. The second step is using ${}^B\mathbf{v}$ to construct the prediction model and the measurement model in the EKF, which is used to estimate the rotation angle. In addition, the ZVD is required for both rotation axis direction estimation and rotation angle estimation. Figure 5 shows the block diagram of the angle measurement algorithm. More details about the algorithm are introduced as follows.

3.1. Zero-Velocity Detection

The objective of ZVD is to decide whether the IMU is stationary or moving during a time epoch consisting of N measurements between time instants n and $n+N-1$ [20]. The ZVD is widely used in inertial navigation systems, providing the required information to reset the velocity error, and preventing the velocity error linearly increasing with time [21,22]. When the IMU is stationary, the following two conditions should be satisfied [23]:

• Acceleration condition: $\left|\right|{}^B{\mathbf{a}}_{\mathbf{m}}\left|\right|=g$.
• Angular velocity condition: $\left|\right|{}^B{\mathbf{\omega }}_{\mathbf{m}}\left|\right|=0$.

Isaac Skog et al. [20] confirmed that a simple threshold on the angular velocity magnitude works well. Antonio R. Jiménez et al. [24] applied thresholds to the acceleration magnitude, angular velocity magnitude, and local acceleration variance. In ref. [25], a time duration threshold was added. In this paper, the ZVD conditions are:

• The acceleration magnitude needs to satisfy a threshold:

  $$t{h}_{a_{\min }}<\left|\right|{}^B{\mathbf{a}}_m\left|\right|<t{h}_{a_{\max }}.$$

  (21)
• The local acceleration variance is less than a threshold:

  $${\zeta }_a^2={\mathrm{var}}_{t_w}\left(\left|\right|{}^B{\mathbf{a}}_{m_t}\left|\right|\right)<t{h}_{{\zeta }_a}$$

  (22)

  where ${\mathrm{var}}_{t_w}\left(\right||{}^B{\mathbf{a}}_{m_t}\left|\right|)$ is an operator that computes the variance of acceleration magnitude $\left|\right|{}^B{\mathbf{a}}_{m_t}\left|\right|$ measured in a time interval of length $t_w$ seconds.
• Both of the above conditions need to be satisfied, and the duration is longer than $t{h}_t$ seconds.

Figure 6 describes an example of how the ZVD works on actual acceleration data, and

Table 1. Parameter values of zero-velocity detection.

Parameter	Value
$t{h}_{a_{\max }}$	9.84 $\mathrm{m}/{\mathrm{s}}^2$
$t{h}_{a_{\min }}$	9.76 $\mathrm{m}/{\mathrm{s}}^2$
$t{h}_{{\zeta }_a}$	0.001 ${\mathrm{m}}^2/{\mathrm{s}}^4$
$t_w$	0.01 s
$t{h}_t$	0.01 s

shows the parameter settings at this point.

3.2. Estimation of Rotation Axis Direction

After mounting the RAMS on a movable surface, the first step is to estimate the rotation axis direction in the frame B. As introduced in Section 2, the key of the rotation axis direction estimation is to measure the gravitational acceleration at three different rotational positions, and then the rotation axis direction ${}^B\mathbf{v}$ can be calculated from Equation (15).

We average the acceleration data when the movable surface is stationary after being rotated to mitigate the measurement noise effect on the result of Equation (15), and rewrite Equation (15) as

$${}^B\mathbf{v}=\pm \frac{{\left[{}^B{\overline{\mathbf{a}}}_1-{}^B{\overline{\mathbf{a}}}_0\right]}_{\times }\left[{}^B{\overline{\mathbf{a}}}_2-{}^B{\overline{\mathbf{a}}}_1\right]}{\left|\right|{\left[{}^B{\overline{\mathbf{a}}}_1-{}^B{\overline{\mathbf{a}}}_0\right]}_{\times }\left[{}^B{\overline{\mathbf{a}}}_2-{}^B{\overline{\mathbf{a}}}_1\right]\left|\right|}$$

(23)

where ${}^B{\overline{\mathbf{a}}}_0$, ${}^B{\overline{\mathbf{a}}}_1$, and ${}^B{\overline{\mathbf{a}}}_2$ denote the average acceleration.

The ZVD ensures that stationary and motion intervals are correctly classified. Once the gravitational acceleration measurements at three different positions are obtained, the rotation axis direction ${}^B\mathbf{v}$ can be estimated according to Equation (23).

The operation procedure rotates the movable surface and then keeps it stationary for at least 1 s, while recording the acceleration data and inputting them into Algorithm A1 (see Appendix A). This is repeated operation until Algorithm A1 successfully obtains the rotation axis direction ${}^B\mathbf{v}$.

3.3. Estimation of Rotation Angle

An EKF, which consists of two stages (prediction and correction) [7,8,9,10,12], is used to estimate the rotation angle after the rotation axis direction is estimated. We define the rotation angle ${\theta }_k$ and the gyroscope bias ${\mathbf{b}}_k$ as state vectors

$${\mathbf{x}}_k=\left[\begin{array}{c}{\theta }_k\\ {\mathbf{b}}_k\end{array}\right].$$

(24)

3.3.1. Prediction

Assume that at time step k, we have the angular velocity ${}^B{\mathbf{\omega }}_m$ measured by gyroscope and the rotation axis direction ${}^B\mathbf{v}$. According to Equations (17)–(20), we can obtain the prediction equation of the EKF as

$$\begin{array}{cc}\hfill {\widehat{\mathbf{x}}}_k^{-}& =f({\widehat{\mathbf{x}}}_{k-1},{}^B{\mathbf{\omega }}_m)\hfill \\ \hfill & =\left[\begin{array}{c}{\widehat{\theta }}_{k-1}+{}^B{\mathbf{v}}^T\widehat{\mathbf{\omega }}\Delta t\\ {\widehat{\mathbf{b}}}_{k-1}\end{array}\right]\hfill \end{array}$$

(25)

$${\mathbf{P}}_k^{-}={\mathbf{F}}_{\mathbf{x}}{\mathbf{P}}_{k-1}{\mathbf{F}}_{\mathbf{x}}^T+\mathbf{Q}$$

(26)

where $\widehat{\mathbf{\omega }}={}^B{\mathbf{\omega }}_m-{\widehat{\mathbf{b}}}_{k-1}$, the minus superscript denotes the a priori (or predicted) estimate, the hat represents that the real system state is estimated by the EKF, ${\mathbf{P}}_{k-1}$ is the covariance matrix, ${\mathbf{F}}_{\mathbf{x}}$ is the Jacobi matrix of $f({\widehat{\mathbf{x}}}_{k|k},{\mathbf{\omega }}_m)$, and $\mathbf{Q}$ is the process noise covariance matrix.

The process noise covariance matrix $\mathbf{Q}$ is

$$\mathbf{Q}=\mathbf{B}\left[\begin{array}{cc}{\sigma }_{\theta }^2{\mathbf{I}}_{3\times 3}& \mathbf{0}\\ \mathbf{0}& {\sigma }_b^2{\mathbf{I}}_{3\times 3}\end{array}\right]{\mathbf{B}}^T$$

(27)

$$\mathbf{B}=\left[\begin{array}{cc}{}^B{\mathbf{v}}^T\Delta t& \mathbf{0}\\ \mathbf{0}& \Delta t{\mathbf{I}}_{3\times 3}\end{array}\right]$$

(28)

where ${\sigma }_{\theta }^2{\mathbf{I}}_{3\times 3}$ and ${\sigma }_b^2{\mathbf{I}}_{3\times 3}$ are the covariance matrix of ${\mathbf{n}}_{\omega }$ and $\mathbf{w}$.

3.3.2. Correction

We use the measured acceleration ${}^B{\mathbf{a}}_m$ to correct the prior estimate. According to Equation (16), the correction model is

$${}^B{\mathbf{a}}_m=h\left({\widehat{\mathbf{x}}}_k^{-}\right)+{\mathbf{n}}_a=R{\left(\theta ,{}^B\mathbf{v}\right)}^T{}^B{\overline{\mathbf{a}}}_0+{\mathbf{n}}_a$$

(29)

where ${}^B{\overline{\mathbf{a}}}_0$ and ${\mathbf{n}}_a$ are defined in the same way as in Equation (16).

The prediction equation can be written as

$${\widehat{\mathbf{x}}}_k={\widehat{\mathbf{x}}}_k^{-}+\mathbf{K}\left({}^B{\mathbf{a}}_m-h\left({\widehat{\mathbf{x}}}_k^{-}\right)\right)$$

(30)

$$\mathbf{K}={\mathbf{P}}_k^{-}{{\mathbf{H}}_{\mathbf{x}}}^T{\left({\mathbf{H}}_{\mathbf{x}}{\mathbf{P}}_k^{-}{{\mathbf{H}}_{\mathbf{x}}}^T+\mathbf{V}\right)}^{-1}$$

(31)

$${\mathbf{P}}_k=(\mathbf{I}-{\mathbf{KH}}_{\mathbf{x}}){\mathbf{P}}_k^{-}$$

(32)

where ${\mathbf{H}}_{\mathbf{x}}$ is the Jacobi matrix of $h\left({\widehat{\mathbf{x}}}_k^{-}\right)$, $\mathbf{K}$ is Kalman gain, and $\mathbf{V}$ is the measurement noise covariance matrix.

The measurement noise covariance matrix $\mathbf{V}$ is

$$\mathbf{V}={\sigma }_a^2{\mathbf{I}}_{3\times 3}$$

(33)

where ${\sigma }_a^2{\mathbf{I}}_{3\times 3}$ is the covariance matrix of ${\mathbf{n}}_a$.

The procedure of the EKF to estimate the rotation angle is shown in Figure 7.

3.3.3. Dynamic Adjustment of Noise Covariance

The above measurement model is only valid when the RAMS is stationary due to the external acceleration. Methods to resolve the covariance uncertainty due to external acceleration by adaptively adjusting the measurement noise covariance matrix can be found in many studies [7,8,9,26]. We chose a switching EKF structure to eliminate the effect of external acceleration.

Based on the standard EKF, we define two state covariance matrices $\mathbf{P}$ and ${\mathbf{P}}_z$, and two measurement noise covariance matrices $\mathbf{V}$ and ${\mathbf{V}}_z$. ${\mathbf{P}}_z$ and ${\mathbf{V}}_z$ are used to construct the EKF for stationary intervals. $\mathbf{P}$ and $\mathbf{V}$ are used to construct the EKF for motion intervals.

The ZVD is still used to distinguish whether the sensor is stationary or in motion. Figure 8 shows the complete rotation angle estimation procedure.

Figure 9 shows an example of how the EFK with ZVD and the standard EKF work on data recorded by the RAMS, respectively. The reference angle in Figure 9a was obtained by a precision turntable. The error in Figure 9b is defined as the difference between the angle ${\theta }_k$ estimated using different algorithms and the reference angle ${\widehat{\theta }}_k$:

$$\mathrm{error}={\theta }_k-{\widehat{\theta }}_k.$$

(34)

The complete motion process in Figure 9 is divided into three parts: stationary period 1 (0–1 s), moving period (1–2 s), and stationary period 2 (2–3 s). In stationary period 1, the two algorithms produced almost the same results (see Figure 9b). In the moving period, $\left|\right|{}^B{\mathbf{a}}_m\left|\right|\ne g$ and $\left|\right|{}^B{\mathit{\omega }}_m\left|\right|\ne 0$, the angle estimation error of the standard EKF was higher than that of the EKF with ZVD. In stationary period 2, the RAMS was stationary again, and the estimation error of the EKF with ZVD converged rapidly. Therefore, we can conclude that the improved EKF algorithm in this paper has significant advantages in dynamic processes.

4. System Design

The RAMS that was designed in this paper consists of a rotation angle measurement module and a 5G mobile communication module. The block diagram for the system is shown in Figure 10. The RAMS transfers the measurement results to the cloud server via 5G network. The data visualization application obtains the measurement data from the cloud server and visualizes them. Furthermore, the data visualization application can access RAMS measurements directly via USB. More details of the rotation angle measurement module are as follows.

The rotation angle measurement module (see Figure 10) is the core of the measurement in RAMS, which includes a tri-axis accelerometer (ADXL355), an IMU (BMI088) containing a tri-axis gyroscope, and a microcontroller (STM32F411). The ADXL355 [27] is a low noise density, low offset drift, low power, selectable measurement range tri-axis accelerometer with industry-leading noise, minimal offset drift over temperature, and long-term stability for precision applications with minimal calibration. The BMI088 [28] is a high-performance, low-cost inertial sensor consisting of a 16-bit digital tri-axis accelerometer and a 16-bit digital tri-axis gyroscope with high-vibration robustness and excellent temperature stability. The STM32F411 [29] is a microcontroller based on the ARM Cortex-M4 32-bit core operating frequency of up to 100 MHz. It features a single-precision floating-point unit, a full set of DSP instructions, up to five SPI interfaces (up to 50 Mbit/s), and a full-speed USB 2.0 controller.

The reason to choose two different sensors instead of having one combined one is that we need a low-noise accelerometer to achieve the measurement accuracy while saving costs.

Table 2. Accelerometer noise density prices for different sensor ICs.

IC	Noise Density $(\mathsf{\mu }\mathbf{g}/\sqrt{\mathbf{Hz}})$	Unit Price (USD)
ADXL355	22.5	52.14
ADIS16465	23	798.14
BMI088	175	28.75

 gives the prices and accelerometer noise densities for ADXL355, ADIS16465 (a high-precision MEMS IMU IC), and BMI088, where the prices are from Digi-Key. As can be seen from the information in Table 2, the noise density of ADXL355 is much less than that of the accelerometer in BMI088, however, the price is much lower than ADIS16465. This solution can achieve the desired measurement accuracy while keeping the cost of sensors under USD 100.

The STM32F411 uses SPI to read acceleration data from ADXL355 and angular velocity data from BMI088 and uses USB to send the measurement results to the communication module. The rotation angle measurement algorithm described in Section 3 is executed on the STM32F411.

Two sensors (ADXL355 and BMI088) were separately calibrated using the method in ref. [30] to reduce the sensor measurement errors.

5. Simulation and Discussion

In this section, we use simulations to verify the proposed algorithm and compare it with the method in ref. [16]. The properties of the IMU used in the simulations are presented in

Table 3. Properties of the IMU used in the simulations.

	Accelerometer	Gyroscope
Bias	50 mg	0.7°$/\mathrm{s}$
Noise Density	50 $\mathsf{\mu }\mathrm{g}/\sqrt{\mathrm{Hz}}$	0.02°$/\mathrm{s}/\sqrt{\mathrm{Hz}}$
Random Walk	10 $\mathsf{\mu }\mathrm{g}/\sqrt{\mathrm{Hz}}$	0.014°$\sqrt{\mathrm{Hz}/\mathrm{s}}$

.

Simulation

In the simulation, we first rotated the IMU around a rotation axis and estimated the misalignment between the rotation axis and the IMU using the method in ref. [16] and the method proposed in this paper, respectively. Then, the IMU was rotated around the rotation axis from −180° to 180° and 1000 measurements calculated by the two methods were recorded at 5° intervals. Figure 11 shows the maximum measurement errors of the two methods for different rotation angles with the rotation axis ${}^B\mathbf{v}={\left[\begin{array}{ccc}0.1734& 0.9835& -0.0523\end{array}\right]}^T$. Clearly, the measurement accuracy of the proposed rotation angle estimation method is superior against the method in ref. [16].

6. Discussion

Since the rotation angle measurement model in this paper relies on the rotation axis direction, the accuracy of the rotation axis direction estimation directly affects the accuracy of the measurement results.

For accelerometers, the axis misalignment, scaling factor, and fixed bias can be obtained and compensated by the method in ref. [30], while the measurement noise ${\mathbf{n}}_a$ directly affects the accuracy of the rotational axis estimation. To reduce the effect of ${\mathbf{n}}_a$, we should extend the stationary time as long as possible and average the acceleration data recorded during the stationary period.

In addition, the angle between the rotation axis and the horizontal plane also has an impact on the accuracy of the rotation axis direction estimation. As the rotation axis changes from horizontal to vertical, the noise weight of ${}^B{\overline{\mathbf{a}}}_1-{}^B{\overline{\mathbf{a}}}_0$ in Equation (23) becomes larger. Assume that $\gamma$ is the angle between the rotational axis ${}^B\mathbf{v}$ and the horizontal plane, ${\mathbf{g}}_0$ and ${\mathbf{g}}_1$ are the ideal values of the measured acceleration of gravity for two different rotational positions (${\mathbf{g}}_1=R{(30°,{}^B\mathbf{v})}^T{\mathbf{g}}_0$). The relationship between the variation of $\gamma$ and $\left|\right|{\mathbf{g}}_1-{\mathbf{g}}_0\left|\right|$ is shown in Figure 12. When the rotation axis is vertical, rotation about this axis does not cause a change in the measured acceleration (regardless of noise). When $\gamma =60°$, the magnitude of $\left|\right|{\mathbf{g}}_1-{\mathbf{g}}_0\left|\right|$ is half of that at $\gamma =0°$, so in the actual measurement process, the angle between the rotation axis and the horizontal plane should preferably be less than 60°.

7. Experiments Result

We used a precision three-axis turntable (SGT320E from China Aviation Industry Corporation) as an experimental platform to provide the actual reference data of rotation angles. The angular position accuracy of the SGT320E is $\pm 5^{″}$. The RAMS was compared with the method in ref. [16] and a high-precision commercial AHRS (Ellipse2-N from SBG). The estimates of the rotation angles obtained by the method in ref. [16] were calculated using raw data from the RAMS. The Ellipse2-N is a small, high-precision AHRS that contains a tri-axis gyroscope, accelerometer, and magnetometer, with a measurement accuracy of 0.1° for roll and pitch, and 0.8° for yaw.

To study the RAMS’s performance in different mounting positions, we adjusted the RAMS to different positions for several tests, taking the rotation axis direction close to the RAMS’s X axis direction (Test 1), Y axis direction (Test 2), and a tilt direction (Test 3) as examples.

All tests were divided into two steps. The first step controlled the specified rotation axis of the SGT320E to rotate three times to estimate the rotation axis direction. The second step was to control the rotation axis from −180° to 180° with a 1 s pause at 10° intervals.

The rotation axis direction estimated by the RAMS is shown in

Table 4. Rotation axis directions estimated by RAMS.

	Rotation Axis Direction
Test 1	${}^B{\mathbf{v}}_1={\left[\begin{array}{ccc}0.9932& -0.1044& 0.0523\end{array}\right]}^T$
Test 2	${}^B{\mathbf{v}}_2={\left[\begin{array}{ccc}0.1043& 0.9921& -0.0698\end{array}\right]}^T$
Test 3	${}^B{\mathbf{v}}_3={\left[\begin{array}{ccc}0.8070& -0.5864& 0.0698\end{array}\right]}^T$

. The results were retained to four decimal places. A total of four sets of rotation angle data were obtained, corresponding to the reference rotation angles obtained by the SGT320E, the RAMS measured angles, the method in ref. [16] measured angles, and the Ellipse2-N measured angles. The Ellipse2-N measurement results were obtained by converting the measured orientations into axis-angle form. The measurement errors were computed as the difference between measured values and the reference values.

It is common to quantify sensor performance as the static and dynamic root-mean-square error (RMSEs). The ZVD determined whether the RAMS was stationary or moving. The RMSEs of the angular measurements at different positions are summarized in

Table 5. Static and dynamic RMSE of Ellipse2-N and RAMS.

		RMSE Static (°)	RMSE Dynamic (°)
${}^B{\mathbf{v}}_1$	RAMS	0.0145	0.4138
	Ellipse2-N	0.1091	0.3545
	Method in ref. [16]	0.0672	3.6057
${}^B{\mathbf{v}}_2$	RAMS	0.0182	0.4133
	Ellipse2-N	0.1734	0.2624
	Method in ref. [16]	0.0804	2.7650
${}^B{\mathbf{v}}_3$	RAMS	0.0118	0.4611
	Ellipse2-N	0.0827	0.1299
	Method in ref. [16]	0.0709	3.6634

, where each value is retained to four decimal places. In all three sets of tests, the static RMSEs of the RAMS are less than 0.02°, and the moving RMSEs are less than 0.5°. The static RMSEs of the method in ref. [16] are less than 0.1°, however, this method cannot obtain accurate measurements under dynamic conditions. Both the static and moving RMSEs of Ellipse2-N exceed 0.1° because the magnetometer could not effectively suppress the cumulative error of the gyroscope at the heading angle due to the frequent changes in the magnetic field in the experimental environment. After running Ellipse2-N and the RAMS for 10 min, they were left to stand still for 1 min. During this minute, the RAMS’s measurement results were held constant while the Ellipse2-N’s measurement results shifted (see Figure 13), which shows that the designed RAMS is not affected by the gyroscope integration error.

The maximum measurement error for different angles is also an important indicator of system performance. As shown in Figure 14, the maximum measurement error of RAMS for various angles did not exceed 0.05°.

8. Conclusions

It is an important part of aircraft ground testing to accurately measure the rotation angle of movable surfaces. For existing IMU-based attitude sensors, the only way to obtain the desired measurement accuracy is to align the sensitive axis of the sensor with the movable surface rotation axis or calibrate the mounting error using precision equipment. Considering the fact that mounting errors are inevitable in the actual measurement process, this paper focuses on proposing a new measurement algorithm for rotational angle measurement without any additional calibration equipment.

The proposed algorithm uses a unit vector in the IMU frame of reference to represent the direction of the rotation axis, which can be estimated using only the stationary acceleration measurements for three different rotational positions. The process of the rotation axis direction estimation can be performed before or during the angle measurement (estimating the rotation axis direction during the measurement process requires saving the raw data from the RAMS and outputting the measurement results in real-time only after the rotation axis is estimated). The rotation angle estimation part of the proposed algorithm uses an adaptive EKF, whose correction model is built relying on the estimated rotation axis direction. The results of simulations and real experiments show that the static measurement accuracy of RAMS in the range of −180°–180° is better than that of a commercial AHRS and another measurement method in ref. [16], and the maximum static measurement error is less than 0.05°.

In addition, the RAMS designed in this paper is not only applicable to aircraft movable surface measurement but can also be applied to any rotation angle measurement project where the rotation axis is fixed, such as the measurement of a robot arm joint rotation angle and the calibration of a low-cost rotary encoder.