Rotate accelerometer data with quaternion

Question

I have an inertial measurement unit consisting of an accelerometer and a gyroscope. My goal is to obtain the absolute acceleration of the x-axis inside the earth frame. Since my device is tilted I cannot use the raw data. Refering to the image, I would like to get x‘ and not x.

Fortunately I can retrieve the rotation in form of a quaternion for every measurement. Unfortunately I don’t know much about quaternions and many papers seem very advanced to me. So how do I need to proceed to rotate a vector of my acceleration like A = (ax, ay, az) to the global frame? Or is it a more intuitive way to convert the quaternion and work with euler angles?

Answer 1

As @neak mentioned, you can use AHRS to do quaternion calculations in python. I assume that you have arrays of acceleration data (acc_data) and gyroscope data (gyr_data) of shape (n, 3).

1. Calculate the quaternions using an Extended Kalman Filter:

   from ahrs.filters.ekf import EKF
   quaternion_arr = EKF(acc=acc_data, gyr=gyr_data).Q
   

   Since you have stated that you already have quaternions, simply replace the quaternion_arr with your data. The quaternions should be in Hamilton representation (w, x, y, z).

2. Rotate the acceleration data based on the quaternion array:

   from ahrs.common.quaternion import Quaternion
   rotated_acc = [Quaternion(q_np).rotate(acc_data[q_idx]) for q_idx, q_np in enumerate(quaternion_arr)]
   

Answer 2

You can check out AHRS toolbox

import numpy as np
from ahrs.filters import EKF
from ahrs.common.orientation import acc2q
ekf = EKF()
num_samples = 1000              # Assuming sensors have 1000 samples each
Q = np.zeros((num_samples, 4))  # Allocate array for quaternions
Q[0] = acc2q(acc_data[0])       # First sample of tri-axial accelerometer
for t in range(1, num_samples):
     Q[t] = ekf.update(Q[t-1], gyr_data[t], acc_data[t])

Once you got the quaternion, you can use the formula mentioned in the previous answer.

Answer 3

The answer depends on what quaternion convention is used by your device and what quaternion arithmetic routines you have available to you for manipulation. Start with these knowns:

q = your quaternion
A = your acceleration vector

Then your desired result will be one of these two calculations depending on the convention:

Arotated = q * A * inv(q)

or

Arotated = inv(q) * A * q

Where * is understood to be a quaternion multiply, inv(q) is the inverse of q, and A is understood to be used as a quaternion with 0 scalar part and A vector part. Then you simply pick off the vector part of Arotated for your answer.

Which brings up the big questions: Do you know the quaternion convention of your device? Is it active or passive? Is the scalar part first or last? Do you have any documentation that describes this? Do you have routines available for calculating quaternion multiply and inverses or will you have to write these from scratch?

-- EDIT --

If you don't know the quaternion convention that your device uses, then you will have to try and discover it through trial and error. That is, you will have to make assumptions and try the calculations to see which results make sense. E.g., assume the scalar part is the first element and the vector part is the last three elements. Then run the calculations above and see if one of them works. If not, then try vector part first and scalar part last and try the above calculations again.

Another trial and error technique is to rotate your device about a known axis and then see how the quaternion elements change in response to this rotation. Doing this about three different axes can pin down which elements are scalar vs vector in the quaternion.

For the quaternion arithmetic code, you can just Google "quaternion multiplication" to get the algorithms. E.g., Wikipedia has the basic formulas listed. Quaternion multiplication is listed here under "Hamilton product", and quaternion inverse for unit quaternions is simply the "conjugate" which negates the vector part. These would be fairly trivial for you to code up yourself without resorting to a library.

Answer 4

Think of a quaternion as a representation of axis-angle rotation.

If your quaternion rotates x' to x, then the inverse of this quaternion rotates x to x'.

So if you already have this quaternion, then invert it and apply it to the rotated axis x, so you could get the un-rotated version of x, which is x'.

=== EDIT ===

A quaternion has 4 values

quat=quaternion(x, y, z, w);

3 of these values are imaginary part, associated with the rotation axis.

and the remaining value is the real part of the quaternion.

Let me assume that the imaginary 3 values are x, y, z, (I'm just assuming because it depends on the implementation of the quaternion in your code)

then

x = sin(theta/2.0) * rx;
y = sin(theta/2.0) * ry;
z = sin(theta/2.0) * rz;

where theta is the rotation angle the quaternion represents, and r=(rx, ry, rz) is the rotation axis vector.

and

w = cos(theta/2.0);

so, to invert an axis-angle rotation, what you have to do is rotate the same amount of angle but in the opposite direction.

In this case, you just invert the sign of x, y, z and you have the inverse of the quaternion.

As for the multiplication between a quaternion and a vector, you treat the vector as a quaternion (with its real part coordinate being zero) so basically you're just multiplying two quaternions. (using Hamilton product).