An Improved Yaw Estimation Algorithm for Land Vehicles Using MARG Sensors

This paper presents a linear Kalman filter for yaw estimation of land vehicles using magnetic angular rate and gravity (MARG) sensors. A gyroscope measurement update depending on the vehicle status and constraining yaw estimation is introduced. To determine the vehicle status, the correlations between outputs from different sensors are analyzed based on the vehicle kinematic model and Coriolis theorem, and a vehicle status marker is constructed. In addition, a two-step measurement update method is designed. The method treats the magnetometer measurement update separately after the other updates and eliminates its impact on attitude estimation. The performances of the proposed algorithm are tested in experiments and the results show that: the introduced measurement update is an effective supplement to the magnetometer measurement update in magnetically disturbed environments; the two-step measurement update method makes attitude estimation immune to errors induced by magnetometer measurement update, and the proposed algorithm provides more reliable yaw estimation for land vehicles than the conventional algorithm.


Introduction
As a set of Euler angles, the yaw, pitch, and roll represent the orientation of a body frame with respect to a reference frame. The pitch and roll are also referred to as attitude. Yaw and attitude estimation are widely used in vehicular technologies including driver assistance [1][2][3], vehicle safety [4,5], etc. In recent years, magnetic angular rate and gravity (MARG) sensors [6] are widely used in orientation estimation. A MARG sensor consists of a triaxis magnetometer, a triaxis gyroscope, and a triaxis accelerometer. Reasonable installation and calibration make it acceptable to assume that the sensor frames are aligned with the body frame. Hence, a MARG sensor can measure the geomagnetic field, angular rate of the body frame, and the gravity resolved in the body frame in undisturbed environments.
In order to obtain an orientation estimation of the body frame, we can integrate the gyroscope output based on an initial value, but the result will drift away with time because of gyroscope measurement errors [7]. Alternatively, we can solve the Wahba problem [8] using magnetometer and accelerometer outputs, but the sensor outputs are apt to be interfered by motion accelerations and magnetic disturbances [9,10]. Therefore, to make the most of the information from MARG sensors and obtain robust orientation estimation, many fusion algorithms have been studied. These algorithms can be classified into two categories: one is based on complementary filters, which realize the fusion in frequency domain [11][12][13], and the other is based on Kalman filters, which employ a which can be used as supplement to the sensors information. We note that the yaw of a vehicle running along a straight road remains basically unchanged, and this character can be utilized to improve the yaw estimation as a supplement to the magnetometer measurement update. To achieve this, it should be known whether a vehicle is running straight. A straightforward way to determine the vehicle status is comparing the gyroscope output against some preset parameters, but the gyroscope output suffers from various measurement errors and vehicle bumps. Considering the turning motion cause changes of MARG outputs, we can determine the vehicle status by analyzing the correlation between outputs from different sensors. In addition, we can reform the measurement update process of a direct Kalman filter using the condition proved in [26] and thus make the attitude estimation immune to magnetic disturbances. In fact, more accurate attitude estimation is also helpful to improve the yaw estimation because they are coupled in the time propagation. Motivated by above discussion, we propose an improved yaw estimation algorithm, and the main contributions of this paper are as follows: (1) A new measurement update robust to magnetic disturbances is introduced, and its weight can be adjusted according to the vehicle status. The correlation coefficients between outputs from a MARG sensor are analyzed based on the vehicle kinematic model and Coriolis theorem, and a vehicle status marker is constructed. (2) A new two-step measurement update method is designed. The method implements measurement updates in two successive steps, and make a special processing of the magnetometer measurement update to eliminate its impact on attitude estimation.
In this paper, we construct a yaw estimation algorithm using a typical linear Kalman filter to implement basic quaternion estimation and applying the conventional strategy of reducing the measurement weight to deal with magnetic disturbances. This algorithm is referred to as conventional algorithm. Then, the conventional algorithm is improved with the new measurement update and the new two-step measurement update method. Finally, the performances of the improved algorithm are evaluated through comparing its results against that of the conventional algorithm in experiments. The rest of this paper is organized as follows: Section 2 describes a conventional yaw estimation algorithm. Section 3 details the improved algorithm. Section 4 provides the experiment results and discussion. Finally, the work is concluded in Section 5.

Preliminaries
Because of not suffering from the gimbal lock, low dimension, and offering a linear formulation of the orientation dynamics [9], quaternion is widely used for orientation representation. Any orientation of a body frame with respect to a reference frame can be represented by a unit quaternion q defined as: where q 0 is the scalar part; q 1 q 2 q 3 T is the vector part; α is the rotation angle; e x e y e z T is the unit vector that represents the rotation axis. In quaternion estimation algorithms, the gyroscope output is used to depict the quaternion dynamics; therefore, the gyroscope bias is an important factor that affects the estimation accuracy.
In this section, we describe a Kalman filter-based yaw estimation algorithm with the unit quaternion and gyroscope bias as states. This algorithm, referred to as conventional algorithm, provides the basement and benchmark for the improved algorithm presented in the next section. The conventional algorithm is based on a typical linear Kalman filter presented in [16]. In addition, adaptive measurement weights are applied to deal with measurement disturbances, and the yaw is computed with the quaternion estimation. The detailed derivations of the system model are presented in [16]; for conciseness, we only list the results and give necessary explanations in this section. Note that we use a different quaternion definition from that in [16] (The scalar part is the first component of a quaternion, whereas in [16], the scalar part is the last component); therefore, we rewrite the concerned equations accordingly, which will not affect the performances of the algorithm.
In this paper, the reference frame is the East, North, Up frame; the body frame is the Right, Forward, Above frame; the sensor frame is assumed aligned with the body frame; the ZXY sequence of Euler angles is chosen, and the yaw, pitch, and roll are respectively z-axis rotation angle, x-axis rotation angle, and y-axis rotation angle. For a clear writing, we define some notations that will be used throughout this paper as follows.
Vectors: x is the state vector; b is the gyroscope bias; m r and g r are respectively the geomagnetic filed and the gravity resolved in the reference frame; m, ω, and a are respectively the output of magnetometer, gyroscope, and accelerometer.
Matrices: C b r is the rotation matrix from the reference frame to the body frame; 0 and I are respectively null matrices and identity matrices, and their subscripts indicate dimensions, for example, the dimensions of 0 3×4 , 0 3 , and I 4 are respectively 3 × 4, 3 × 3, and 4 × 4; [v×] is the skew-symmetric matrix of vector v.
Subscripts: x, y, and z denote respectively the x-axis, y-axis, and z-axis component of a vector; k denotes the time step.

System Model
The state vector, composed of the unit quaternion and gyroscope bias, is defined as: Based on the well-known quaternion dynamics model and state augmentation, the process equation is written as: where: In (3), Ψ k is the transition matrix; Γ k is the process noise input matrix; ∆t is the sample interval; q v is the vector part of the quaternion state, and n k is the gyroscope noise with covariance matrix diag σ 2 ω , σ 2 ω , σ 2 ω . The linear pseudo-measurement equation of the accelerometer is: where: In (4), H a,k+1 is the measurement matrix, and − 1 2 Ξ k+1 δa k+1 is the measurement noise. The covariance matrix of δa k+1 is σ 2 a I 3 . Similarly, the linear pseudo-measurement equation of the magnetometer is: where: In (5), H m,k+1 is the measurement matrix, and − 1 2 Ξ k+1 δm k+1 is the measurement noise. The covariance matrix of δm k+1 is σ 2 m I 3 . From (4) and (5), the overall measurement equation can be written as: where: The noise covariance matrices of the process and measurement equations are written as: where ρ a and ρ m are adaptive weight coefficients, which can be adjusted in real time. The computation formulas for Q k , R a,k+1 , and R m,k+1 can be found in [16]. From (8) and (9), the noise covariance matrix of (6a) can be written as:

Adaptive Fusion Algorithm Based on Kalman Filter
Filter Initialization: Set initial values for state vector estimationx, i.e.,q andb, and error covariance matrix P.
Time Propagation: The process model, gyroscope output and posteriori estimation (or the initial estimation) at step k is used to compute the priori estimation at step k + 1 by: Measurement weight adjustment: To deal with motion accelerations and magnetic disturbances, the weight coefficients ρ a and ρ m are adjust in real time according to the disturbance intensity, which is equivalent to adjust the measurement weight. The adjusting expressions are: where the relative distances between the norms of measured vectors (a k+1 , m k+1 ) and reference vectors (g r , m r ) are used to represent the disturbance intensity, and the exponential functions are used to map the disturbance intensity to weight coefficients. The function of the parameters λ a and λ m is to adjust the mapping relations. If the weight coefficient should be more sensitive to the disturbance intensity, the corresponding parameter should be increased, otherwise, it should be reduced. The values used for λ a and λ m can be determined experimentally. The exponential function instead of a linear function is applied because the former can reduce the measurement weight more quickly when the relative distance increases [28]. Measurement Update: The measurement model, magnetometer and accelerometer outputs and the priori estimation at step k + 1 is used to compute the posteriori estimation at step k + 1 by: Unit Constraint: To preserve the unit-norm property of the quaternion estimation, the updated quaternion is normalized by:

Yaw Computation
The rotation matrix from the reference frame to the body frame can be represented as a function of either Euler angles or a unit quaternion, and the expressions are respectively: where ψ, φ, and γ are respectively the yaw, pitch, and roll; s and c denote sine and cosine function, respectively.
Using the first two elements of the second row of C(q), denoted as C 21 and C 22 respectively, and defining the range of the yaw as (−180 • , 180 • ], we obtain the following formulas for yaw computation:

Improved Algorithm
In this section, we keep the time propagation step of the conventional algorithm unchanged, and improve its measurement update step in two ways. Firstly, to improve the accuracy of the yaw estimation in the presence of magnetic disturbances, a measurement equation of the gyroscope bias is derived and its weight can be adjusted according to the vehicle status. Secondly, we design a two-step measurement update method to eliminate the impact of the magnetometer measurement update on the attitude estimation.

New Measurement Equation and Its Adaptive Weight
The gyroscope output can be written as: and its z-axis component is: where ω b is the angular rate resolved in the body frame, and n 1 is a white Gaussian noise vector. When vehicles are running straight, ω bz can be regarded as zero and (19b) can be rewritten as: where: Equation (20a) is a new measurement equation named as gyroscope measurement, and its characters are as follows. Firstly, it only updates the gyroscope bias estimation in measurement update, but it will affect the quaternion estimation in time propagation through (11); specifically speaking, it will constrain the yaw estimation from changing by correcting the z-axis angular rate to zero. Secondly, it is immune to magnetic disturbances because the latter cannot impact the gyroscope output. Thirdly, the gyroscope measurement only hold true in running straight status, in other words, its validity depends on the vehicle status.
Note that, the yaw is computed using quaternion estimation, and all the gyroscope bias components can impact the quaternion propagation. But, in the fusion algorithm, the quaternion estimation is not only based on propagation, but also on measurement update. The yaw derivative [12] can be written as: Obviously, the yaw derivative is related to pitch, roll, ω bx , and ω bz . In the fusion algorithm, the estimation of pitch and roll will be corrected by the accelerometer measurement (in the form of correcting quaternion), which is immune to magnetic disturbances. The two-step measurement update method which will be introduced in Section 3.2 makes pitch and roll estimation immune to the errors induced by magnetometer measurement update. In addition, generally speaking, the pitch and roll are small for land vehicles, hence the absolute value of the coefficient of ω bz is greater than that of ω bx . Therefore, the z-axis bias component is of more importance, especially in magnetically disturbed environments. Equation (21) is based on Euler angles, but it should be noted that both quaternion and Euler angles are representation of orientation, and orientation obeys a unique dynamic rule. In other words, the quaternion propagation models is different from that of Euler angles, but when they are used to represent a same orientation parameter, for example yaw, the parameter obeys a same dynamic rule no matter what representation is used.
The above characteristics suggest that the gyroscope measurement can be used to correct the yaw estimation in the presence of magnetic disturbances when vehicle is running straight. In practice, a vehicle cannot always run on a straight road, and even on a straight road the yaw may fluctuate. In order to employ this measurement correctly, we define the covariance of the measurement noise n 1z as ρ ω σ 2 ω , where ρ ω , similar to ρ a and ρ m , works as a adaptive weight coefficient that can be adjusted in real time. Obviously, the adjusting principle should be that the higher the extent of running straight is, the smaller the ρ ω is. Now, the problem is how to determine the extent of running straight. To address this problem, we analyze the correlations between different sensor outputs. The accelerometer output can be written as: where a m is the motion acceleration, and n 2 is a noise vector. According to the vehicle kinematic model [29], a m is given by: where v b is the velocity of the vehicle in the body frame.
The reference fame and body frame are shown in Figure 1. For a land vehicle, the running direction can be arbitrary in reference frame, but in the body frame, which is attached to the vehicle, its running direction is always forward, although the forward direction may change in reference frame. The "forward direction" is the y direction in the body frame defined in this paper. Dissanayake et al. [30] have pointed out that, when the vehicle does not jump off the ground and does not slide on the ground, velocity of the vehicle in the plane perpendicular to the forward direction is zero, hence v bx and v bz can be assumed as zero. Using this assumption and substituting (16) and (23) into (22), we simplify the x-axis component of (22) as: In practical situation, the assumption is somewhat violated due to the presence of side slip during cornering and vibrations caused by the engine and suspension system [30]. When the assumption is violated, the side slips and vibrations will cause corresponding accelerations, but these accelerations can be regarded as noise and a part of n 2x in (24). In practical situation, the assumption is somewhat violated due to the presence of side slip during cornering and vibrations caused by the engine and suspension system [30]. When the assumption is violated, the side slips and vibrations will cause corresponding accelerations, but these accelerations can be regarded as noise and a part of The magnetic field resolved in the body frame, denoted as b m , can be written as: where r d is the magnetic disturbances resolved in the reference frame. According to the Coriolis theorem, the change rate of b m in body frame can be written as: , the x-axis and y-axis components of b m  can be written as: The magnetic field resolved in the body frame, denoted as m b , can be written as: where d r is the magnetic disturbances resolved in the reference frame. According to the Coriolis theorem, the change rate of m b in body frame can be written as: Because . m r = 0, the x-axis and y-axis components of . m b can be written as: .
where C b r1· and C b r2· denote the first row and second row of C b r respectively. We rewrite (27) as: .
When the angle between the road plane and the horizontal plane remains constant and the vehicle do not vibrate, ω by and ω bx are zero. In practice, the angle between the two planes may change, however, in general, the rate of change will be small, and the duration of the vehicle attitude change is also very short. Hence, the corresponding ω by and ω bx will be small and close to zero for most of the time. In the case of vibration, which can be caused by bumps or ditches, ω by and ω bx may have large absolute value, but they will oscillate with high frequency. Hence, regarding the last item of the right side as noise and using the magnetometer output and difference to approximate the left side of (28), we obtain: where n 3 and n 4 are noises. Because m y,k and m x,k may be equal or close to zero and thus cause numerical instability, we define ∆m as: Equations (19b) and (24) show that the correlation between ω z and a x will be low if ω bz is stable; otherwise, the correlation will be high. Equations (19b), (29), and (30) show that ω z and ∆m have similar relationships; more importantly, d r will not affect the correlation between ω z and ∆m if it is constant and will not affect the correlation markedly except for an abrupt changing. Hence, the correlation between ω z and ∆m is robust to magnetic disturbances.
Considering ω bz is stable in running straight status and unstable in turning status, we construct a vehicle status marker c, whose expression is: where cor(·) estimates the correlation coefficient using the samples of the sensor outputs, and subscript N denotes the N points samples from k − N + 1 to k. In fact, in running straight status, MARG sensor outputs are stable except for uncorrelated noises; whereas in turning status, the outputs change, and the output from different sensor is correlative because their changes are caused by the same reason i.e., the vehicle is turning. The rationale of (31) is that use estimated correlation coefficients to distinguish status. In (31), the weighted sum combine the two correlation coefficient to make full use of the correlation between different sensors, and the second term of the sum is minus because ω z is negatively correlated to a x when ω bz is changing. Considering the sum is theoretically nonnegative, we set c to zero when correlation coefficient estimation errors result in a negative sum. Similar to ρ a and ρ m , the adjusting expression for ρ ω is: where λ ω is a coefficient that transforms the status marker to a proper range.

Two-Step Measurement Update Method
This part presents the two-step measurement update method, and the concerned parameters of the first and second update step are respectively denoted by subscript f and s in the following. First Step: The accelerometer and gyroscope measurement updates are implemented in this step. The measurement matrix and noise covariance matrix are respectively: The update expressions are: where: Second Step: In this step, the magnetometer measurement update is implemented and the ultimate quaternion estimation is computed. The measurement matrix and noise covariance matrix are respectively: The update expressions are: Extracting quaternion part ofx f ,k+1 andx s,k+1 , we obtainq f ,k+1 andq s,k+1 . Define q σ as: whereq * f ,k+1 is the conjugate quaternion ofq f ,k+1 and the symbol ⊗ is the quaternion multiplication operator. The conjugate quaternion of a unit quaternion represents the inverse rotation and a sequence of rotations can be represented by quaternion multiplication; therefore, q σ can be viewed as the correction quaternion induced by the magnetometer measurement update. In theory, q σ should only correct the yaw estimation; however, it also modify the attitude estimation in practice. The corresponding rotation matrix equation to Equation (39) can be written as: The condition for C(q s,k+1 ) and C(q f ,k+1 ) to share a same attitude [26] is: where q σv is the vector part of q σ and C(q f ,k+1 ) ·3 is the third column of C(q f ,k+1 ). The geometric meaning of (41) is that the rotation axis of q σ should be parallel to the z-axis of the reference frame. Based on this condition and (1), a new correction quaternionq σ is defined as: where: Note thatq σ is actually the projection of q σ on the direction depicted by C(q f ,k+1 ) ·3 , and this treatment makesq σ satisfy (41). Finally, the ultimate quaternion estimation is computed by:

Complete Improved Yaw Estimation Algorithm
Adding the adaptive gyroscope measurement update to the conventional algorithm and adopting the two-step measurement update method, we obtain the improved algorithm shown in Figure 2. Note that the improved algorithm preserves the linearity and all its measurement updates have an adaptive weight.

Complete Improved Yaw Estimation Algorithm
Adding the adaptive gyroscope measurement update to the conventional algorithm and adopting the two-step measurement update method, we obtain the improved algorithm shown in Figure 2. Note that the improved algorithm preserves the linearity and all its measurement updates have an adaptive weight.

Experiments
In this section, we evaluate the performances of the improved algorithm experimentally. The real data from a MARG sensor mounted on a test vehicle are processed by the conventional and improved algorithm, respectively. Then, we analyze the performances of the vehicle status marker, gyroscope measurement update, two-step measurement update method, and yaw estimation by comparing the results of the improved algorithm against that of the conventional algorithm and

Experiments
In this section, we evaluate the performances of the improved algorithm experimentally. The real data from a MARG sensor mounted on a test vehicle are processed by the conventional and improved algorithm, respectively. Then, we analyze the performances of the vehicle status marker, gyroscope measurement update, two-step measurement update method, and yaw estimation by comparing the results of the improved algorithm against that of the conventional algorithm and reference values.

Experimental Setup and Parameters
The test vehicle and experimental devices are shown in Figure 3. A Motion Tracker MTi-28A53G35 (Xsens, Enschede, The Netherlands) [31] was used as the MARG sensor. The Global Position System (GPS) unit (Unicoremm, Beijing, China) provided reference values of longitude, latitude, and forward velocity. 3-axis gyro module STIM210 (Sensonor AS, Horten, Norway) was used as the attitude and heading reference system (AHRS), which provided reference values of the yaw and attitude. The initial value of the AHRS was computed with the static accelerometer and magnetometer outputs [32]. The laptop (Lenovo, Beijing, China) supplied power for the MARG sensor, GPS unit, and the AHRS, and logged data from them at 100Hz. The GPS unit, MARG sensor, and AHRS were mounted on the test vehicle with the sensor and AHRS frames aligned with the vehicle body frame. In experiments, the vehicle was driven along test trajectories, which consisted of straight lines, corners with different angles and a circular line, involving a full range of yaw. Note that the magnetic sensor had been calibrated inline according to the sensor manual [20]; hence, we assumed that the impact of the hard and soft iron effects had been eliminated.
We analyzed the MARG sensor output in static condition, and set the noise parameters as follows: σ m = 0.0015, σ ω = 0.0056 rad/s, σ a = 0.008 m/s 2 . We found the proper values for λ m , λ ω , and λ a by trial and error, and set them to 20, 15, and 50 respectively. The estimation of q 0 was computed with the initial outputs of the magnetometer and accelerometer; the estimation of b 0 was set to the gyroscope output before the vehicle was started, and P 0 = 100I 7 . The initial values of Q, R a , and R m were computed with q 0 and P 0 [16]. inline according to the sensor manual [20]; hence, we assumed that the impact of the hard and soft iron effects had been eliminated. We analyzed the MARG sensor output in static condition, and set the noise parameters as follows:

Vehicle Status Marker
To analyze the experimental results more clearly, we intercepted a piece of the MARG sensor outputs corresponding to a segment of the test trajectories. Using these data, we computed the status marker c and the adaptive weight coefficient   , which are shown in Figure 4. In Figure 4a, the reference yaw indicates that the vehicle undergone a turning (about 90°) between 70 s and 80 s, and mainly run straight with small fluctuation of yaw in other times. Corresponding to the reference yaw, c increases markedly between 70 s and 80 s; and oscillates with small values in other times. One key to computing c effectively is to select a proper N, because too many sample points will cause severely lagged weight adjustment, whereas too few sample points will reduce the computation accuracy. In the experiments, we set

Vehicle Status Marker
To analyze the experimental results more clearly, we intercepted a piece of the MARG sensor outputs corresponding to a segment of the test trajectories. Using these data, we computed the status marker c and the adaptive weight coefficient ρ ω , which are shown in Figure 4. In Figure 4a, the reference yaw indicates that the vehicle undergone a turning (about 90 • ) between 70 s and 80 s, and mainly run straight with small fluctuation of yaw in other times. Corresponding to the reference yaw, c increases markedly between 70 s and 80 s; and oscillates with small values in other times.
One key to computing c effectively is to select a proper N, because too many sample points will cause severely lagged weight adjustment, whereas too few sample points will reduce the computation accuracy. In the experiments, we set N = 1/∆t.
To analyze the experimental results more clearly, we intercepted a piece of the MARG sensor outputs corresponding to a segment of the test trajectories. Using these data, we computed the status marker c and the adaptive weight coefficient   , which are shown in Figure 4. In Figure 4a, the reference yaw indicates that the vehicle undergone a turning (about 90°) between 70 s and 80 s, and mainly run straight with small fluctuation of yaw in other times. Corresponding to the reference yaw, c increases markedly between 70 s and 80 s; and oscillates with small values in other times. One key to computing c effectively is to select a proper N, because too many sample points will cause severely lagged weight adjustment, whereas too few sample points will reduce the computation accuracy. In the experiments, we set

Gyroscope Measurement Update
We processed the sensor outputs using the conventional and improved algorithm, respectively, and the yaw estimation results are shown in Figure 5. The results demonstrate that the improved algorithm has better estimation accuracy, and the reasons may be analyzed as follows. In the conventional algorithm, the magnetometer measurement update corrects the yaw estimation from drifting, but the measurement weight will be adjusted down in the presence of magnetic disturbances (Figure 4b shows obvious differences between the norm of magnetometer output and one from about 15 s to 50 s implying the presence of magnetic disturbances). As a result, the correction effect attenuates and thus the yaw estimation drifts. In contrast, the gyroscope measurement update provides another correction in the improved algorithm, and more importantly, the correction can hardly be attenuated by magnetic disturbances as analyzed in Section 3. The adaptive weight coefficients shown in Figure 4b verify the function of the gyroscope measurement update: ρ ω is much lower than ρ m in the presence of magnetic disturbances. In fact, ρ ω represent the validity of the gyroscope measurement, which is not affected by magnetic disturbances. Hence, the gyroscope measurement can be used in magnetically disturbed environments. The real-time adjusting ρ ω is critical for the gyroscope measurement to work effectively. Figure 5 also shows the yaw estimation result when ρ ω is set to constant 1. The result demonstrates that the constant ρ ω causes erroneous correction from the gyroscope measurement when the vehicle is turning, and hence results in significant estimation errors.

Gyroscope Measurement Update
We processed the sensor outputs using the conventional and improved algorithm, respectively, and the yaw estimation results are shown in Figure 5. The results demonstrate that the improved algorithm has better estimation accuracy, and the reasons may be analyzed as follows. In the conventional algorithm, the magnetometer measurement update corrects the yaw estimation from drifting, but the measurement weight will be adjusted down in the presence of magnetic disturbances (Figure 4b shows obvious differences between the norm of magnetometer output and one from about 15 s to 50 s implying the presence of magnetic disturbances). As a result, the correction effect attenuates and thus the yaw estimation drifts. In contrast, the gyroscope measurement update provides another correction in the improved algorithm, and more importantly, the correction can hardly be attenuated by magnetic disturbances as analyzed in Section 3. The adaptive weight coefficients shown in Figure 4b verify the function of the gyroscope measurement update:   is much lower than m  in the presence of magnetic disturbances. In fact,   represent the validity of the gyroscope measurement, which is not affected by magnetic disturbances. Hence, the gyroscope measurement can be used in magnetically disturbed environments. The real-time adjusting   is critical for the gyroscope measurement to work effectively. Figure 5 also shows the yaw estimation result when   is set to constant 1. The result demonstrates that the constant   causes erroneous correction from the gyroscope measurement when the vehicle is turning, and hence results in significant estimation errors. The gyroscope bias estimation results are shown in Figure 6. We estimated the reference value of z b by calculating the mean value of the gyroscope output in a not-moving interval, and the result was 0.00058 rad/s. Obviously, the z b estimation is updated more effectively in the improved algorithm, which is verified by its better yaw estimation accuracy. The gyroscope bias estimation results are shown in Figure 6. We estimated the reference value of b z by calculating the mean value of the gyroscope output in a not-moving interval, and the result was 0.00058 rad/s. Obviously, the b z estimation is updated more effectively in the improved algorithm, which is verified by its better yaw estimation accuracy. The gyroscope bias estimation results are shown in Figure 6. We estimated the reference value of z b by calculating the mean value of the gyroscope output in a not-moving interval, and the result was 0.00058 rad/s. Obviously, the z b estimation is updated more effectively in the improved algorithm, which is verified by its better yaw estimation accuracy.
(a) (b) In Section 3.1, we draw a conclusion based on Equation (21) that z b estimation is more important than x b and y b estimation in the fusion algorithm. To verify this conclusion based on experimental data, we constructed a conventional adaptive fusion algorithm using quaternion as In Section 3.1, we draw a conclusion based on Equation (21) that b z estimation is more important than b x and b y estimation in the fusion algorithm. To verify this conclusion based on experimental data, we constructed a conventional adaptive fusion algorithm using quaternion as state variables, and using Equation (6) as measurement equation. Note that, in this algorithm, the gyroscope bias is not estimated and compensated, and hence it will always impact the quaternion propagation. To examine the impact of the gyroscope bias, we added a constant bias on the reference angular rate provided by the AHRS and used the biased value to realize the quaternion propagation in the fusion algorithm.
We used this algorithm to estimate yaw for three times, and the constant bias were respectively set to [b, 0, 0], [0, b, 0], and [0, 0, b]. The results were denoted as bx, by, and bz respectively. Note that, in this process, all the parameters of the algorithm are same except for the constant bias. The root mean square (RMS) errors of the yaw estimation corresponding to different b are listed in Table 1. Obviously, b z has a more significant impact on yaw estimation. It should be noted that the gyroscope measurement updates only provide a "partial" correction, in other words, they cannot provide absolute yaw information but that the yaw is unchanged to some extent, and the "extent" is indicated by the adaptive weight. Therefore, it is reasonable to employ the gyroscope measurement updates in combination with the magnetometer measurement updates. In addition, the gyroscope measurement updates only work effectively in running straight status. In practice, turning status is inevitable, but the duration of turning is relatively short, and running along a straight road is a more usual status for most land vehicles.
In practice, a running vehicle cannot avoid bumps and ditches, which will cause oscillations of a x , ω x , and ω y . We assume these oscillations as part of the noise items in Equations (24) and (29) respectively, and construct the vehicle status marker c based on the correlations between different sensor outputs. To evaluate these noise assumptions, we analyzed a piece of the sensor outputs corresponding to a straight road with some bumps and ditches. The raw signals of a x , ω x , and ω y are shown in Figures 7 and 8, respectively. The term c computed by Equation (31) is shown in Figure 9.
In practice, a running vehicle cannot avoid bumps and ditches, which will cause oscillations of x a , x  , and y  . We assume these oscillations as part of the noise items in Equations (24) and (29) respectively, and construct the vehicle status marker c based on the correlations between different sensor outputs. To evaluate these noise assumptions, we analyzed a piece of the sensor outputs corresponding to a straight road with some bumps and ditches. The raw signals of x a , x  , and y  are shown in Figures 7 and 8, respectively. The term c computed by Equation (31) is shown in Figure 9.   The vehicle met a bump or ditch at about 7 s, 19 s, 27 s, 45 s, and 58 s. Obviously, the bumps and ditches cause oscillations of the sensor outputs. Figure 9 demonstrate that the bumps and ditches do not affect c significantly, which displays small values and is consistent with the running straight status.

Two-Step Measurement Update Method
In the experiments, the test road is basically level except for some speed breaks. To evaluate the performances of the two-step measurement update method, we computed attitudes using the quaternion [33] from the first step update, the conventional algorithm and the improved algorithm respectively, and the results are shown in Figure 10. Obviously, the improved algorithm and the first step update have the same attitude estimation, which verify that the second step update do not modify the attitude estimation. In addition, the attitude estimation of the improved algorithm is more consistent with the reference values than the conventional algorithm. The comparison demonstrates that the two-step measurement update method eliminates the attitude estimation errors induced by the magnetometer measurement update. Besides magnetic disturbances, the magnetometer measurement error can also cause attitude estimation errors. For example, the attitude estimation of the conventional algorithm change incorrectly during the turning of the vehicle (70-80 s), and similar phenomena can be found in other turning processes, which, we think, is due to the dynamic errors   The vehicle met a bump or ditch at about 7 s, 19 s, 27 s, 45 s, and 58 s. Obviously, the bumps and ditches cause oscillations of the sensor outputs. Figure 9 demonstrate that the bumps and ditches do not affect c significantly, which displays small values and is consistent with the running straight status.

Two-Step Measurement Update Method
In the experiments, the test road is basically level except for some speed breaks. To evaluate the performances of the two-step measurement update method, we computed attitudes using the quaternion [33] from the first step update, the conventional algorithm and the improved algorithm respectively, and the results are shown in Figure 10. Obviously, the improved algorithm and the first step update have the same attitude estimation, which verify that the second step update do not modify the attitude estimation. In addition, the attitude estimation of the improved algorithm is more consistent with the reference values than the conventional algorithm. The comparison demonstrates that the two-step measurement update method eliminates the attitude estimation errors induced by the magnetometer measurement update. Besides magnetic disturbances, the magnetometer measurement error can also cause attitude estimation errors. For example, the attitude estimation of the conventional algorithm change incorrectly during the turning of the vehicle (70-80 s), and similar phenomena can be found in other turning processes, which, we think, is due to the dynamic errors of the magnetometer. The vehicle met a bump or ditch at about 7 s, 19 s, 27 s, 45 s, and 58 s. Obviously, the bumps and ditches cause oscillations of the sensor outputs. Figure 9 demonstrate that the bumps and ditches do not affect c significantly, which displays small values and is consistent with the running straight status.

Two-Step Measurement Update Method
In the experiments, the test road is basically level except for some speed breaks. To evaluate the performances of the two-step measurement update method, we computed attitudes using the quaternion [33] from the first step update, the conventional algorithm and the improved algorithm respectively, and the results are shown in Figure 10. Obviously, the improved algorithm and the first step update have the same attitude estimation, which verify that the second step update do not modify the attitude estimation. In addition, the attitude estimation of the improved algorithm is more consistent with the reference values than the conventional algorithm. The comparison demonstrates that the two-step measurement update method eliminates the attitude estimation errors induced by the magnetometer measurement update. Besides magnetic disturbances, the magnetometer measurement error can also cause attitude estimation errors. For example, the attitude estimation of the conventional algorithm change incorrectly during the turning of the vehicle (70-80 s), and similar phenomena can be found in other turning processes, which, we think, is due to the dynamic errors of the magnetometer.
It should be noted that modifications of the attitude estimation made by magnetometer measurement update essentially arise from the disagreements between the accelerometer and magnetometer outputs, which can also be caused by non-gravitational acceleration and accelerometer measurement errors. The improved algorithm prevents the magnetometer output from modifying the attitude estimation because firstly the gravity resolved in the body frame can provide sufficient attitude information, and secondly the accelerometer output is assumed more reliable than the magnetometer output in magnetically disturbed environments as has been accepted and verified in many applications.

Yaw Estimation
We tested two trajectories, which are referred to as A and B. MARG sensor outputs corresponding to the trajectories were processed. To show the performances of the yaw estimation intuitively, we used dead reckoning [29] to reproduce the test trajectories based on the yaw estimations and the vehicle velocity data provided by the GPS unit. The results are shown in Figure 11, where the reference trajectories recorded by the GPS unit is also plotted. Clearly, the reproduced trajectories from the improved algorithm are closer to the reference than that from the conventional algorithm. The root mean square (RMS) errors of the improved algorithm in trajectory A and B are respectively 1.8° and 2.9°. The reason the improved algorithm performed better in A is, as pointed out earlier, the gyroscope measurement only work effectively in running straight status. A was approximately rectangular and has longer straight roads, whereas B consisted of more turning road with different angles and curvature. Obviously, A was more suitable for the gyroscope measurement to take effect. It should be noted that vehicle status cannot affect the work of the two-step measurement update method, which is also helpful for the yaw estimation.
The RMS errors of the yaw and attitude estimation are listed in Table 2. In addition, we used the

Yaw Estimation
We tested two trajectories, which are referred to as A and B. MARG sensor outputs corresponding to the trajectories were processed. To show the performances of the yaw estimation intuitively, we used dead reckoning [29] to reproduce the test trajectories based on the yaw estimations and the vehicle velocity data provided by the GPS unit. The results are shown in Figure 11, where the reference trajectories recorded by the GPS unit is also plotted. Clearly, the reproduced trajectories from the improved algorithm are closer to the reference than that from the conventional algorithm.

Yaw Estimation
We tested two trajectories, which are referred to as A and B. MARG sensor outputs corresponding to the trajectories were processed. To show the performances of the yaw estimation intuitively, we used dead reckoning [29] to reproduce the test trajectories based on the yaw estimations and the vehicle velocity data provided by the GPS unit. The results are shown in Figure 11, where the reference trajectories recorded by the GPS unit is also plotted. Clearly, the reproduced trajectories from the improved algorithm are closer to the reference than that from the conventional algorithm. The root mean square (RMS) errors of the improved algorithm in trajectory A and B are respectively 1.8° and 2.9°. The reason the improved algorithm performed better in A is, as pointed out earlier, the gyroscope measurement only work effectively in running straight status. A was approximately rectangular and has longer straight roads, whereas B consisted of more turning road with different angles and curvature. Obviously, A was more suitable for the gyroscope measurement to take effect. It should be noted that vehicle status cannot affect the work of the two-step measurement update method, which is also helpful for the yaw estimation.
The RMS errors of the yaw and attitude estimation are listed in Table 2. In addition, we used the Kalman filter with vector selection [7] and the complementary filter with varying gains [21], referred The root mean square (RMS) errors of the improved algorithm in trajectory A and B are respectively 1.8 • and 2.9 • . The reason the improved algorithm performed better in A is, as pointed out earlier, the gyroscope measurement only work effectively in running straight status. A was approximately rectangular and has longer straight roads, whereas B consisted of more turning road with different angles and curvature. Obviously, A was more suitable for the gyroscope measurement to take effect. It should be noted that vehicle status cannot affect the work of the two-step measurement update method, which is also helpful for the yaw estimation.
The RMS errors of the yaw and attitude estimation are listed in Table 2. In addition, we used the Kalman filter with vector selection [7] and the complementary filter with varying gains [21], referred to as VSKF and VGCF respectively, to process the experiment data, and the RMS errors are also listed in Table 2. Clearly, the improved algorithm performed best in the yaw and attitude estimation. In fact, as pointed out in the introduction, the VSKF and VGCF essentially use the same strategy as the conventional algorithm to handle magnetic disturbances: reducing the confidence in the magnetometer measurements, which will result in poor correction for the yaw estimation. Whereas the improved algorithm not only reduce the confidence in the disturbed measurements but also introduce new measurements based on vehicle status, which provide supplemental information for the yaw estimation. In addition, the two-step measurement update method eliminates the impacts of magnetometer errors and magnetic disturbances on attitude estimation, which is also helpful for the yaw estimation. Note that the better yaw estimation accuracy of the improved algorithm is due to its enhanced ability to deal with magnetic disturbances; and it will have the same level of yaw estimation accuracy as the conventional algorithm in magnetically homogeneous environments.

Conclusions
An improved yaw estimation algorithm for land vehicle using a MARG sensor was proposed in this paper. Under running straight assumption, we derived the gyroscope measurement equation, which update the gyroscope bias and thus constrain the yaw estimation from changing. The validity of the gyroscope measurement depends on the vehicle status; to determine the vehicle status, we analyzed the correlations between different sensors based on the vehicle kinematic model and Coriolis theorem, and constructed a vehicle status marker used to adjust the weight of the gyroscope measurement. In addition, we designed a two-step measurement update method, which implements the magnetometer measurement update separately and eliminates its impact on attitude estimation. Adopting the gyroscope measurement update and the two-step measurement update method, we improved the conventional yaw estimation algorithm. The improved algorithm was tested in experiments and compared against the conventional algorithm. Based on the experiment results, the performances and characters of the improved algorithm were discussed and the conclusion is as follows. The gyroscope measurement update is an effective supplement to the magnetometer measurement update in magnetically disturbed environments; the two-step measurement update methods make attitude estimation immune to the errors induced by magnetometer measurement update; and the improved algorithm provides more reliable yaw estimation for land vehicles than the conventional algorithm. Finally, it should be noted that the vehicle status marker is based on statistic characteristics between different sensors, which make it robust to disturbances, but on the other hand, insensitive to the varying status. The improvements of its real-time performance and ability to detect the turning status with small angular rates will be the topics of further work.