Angular rates from quaternions

Hello,

I have a matrix containing quaternions in each row. And, each row represents a different time instant.

I would like to find the angular velocities along x, y and z ( ) using the quaternions.

I noticed that there is a function called angvel in matlab that should do it. However, I do not have this function in my matlab license.

Is there another way to find the angular velocities?

Thank you

[James Tursa]

You can calculate these angular velocities directly from the quaternions, but you need to know two things:

1) The times of the quaternions

2) The convention of the quaternions

I am guessing you already have (1). For (2), you need to know whether these quaternions are left chain vs right chain in order to back out the angular rate vector. (Basically, you start with the qdot derivative expression and solve for w). Do you know this? Where are these quaternions coming from? What are the coordinate systems involved?

E.g., depending on the quaternion convention and coordinate systems involved, there are four different qdot equations:

qdot = (1/2) q * w

qdot = -(1/2) q * w

qdot = (1/2) w * q

qdot = -(1/2) w * q

Which one of these needs to be solved for w depends on the quaternion convention and how you want w coordinatized.

[Elie Hatem]

Thank you for response, I managed to find the rates, however, I think the way I found qdot could not be ideal.

I have a roll angle that is changing with time, while the pitch and yaw remain at 0 degrees.

So, I found the quaternions by stacking all the angles at each time instant in their respective vectors and using the eul2quat function provided by matlab:

euler = [ phi' , theta , psi ];

quat = eul2quat(euler,'XYZ');

Then, I normalized the quaternion in each row of the obtained matrix:

%% extracting number of rows and columns from phi
[r, c] = size(phi');

%% normalizing quat
for i = 1:r
    qi = quat(i);
    qi_norm = sqrt(sum(qi)^2);
    quat(i) = 1 / qi_norm * qi;
end

However, I am not sure about qdot. Since I do not actually have it's values, I calculated them from quat as follows:

% computing q_dot
q_dot = diff(quat) / t_step;

With t_step being the time step between each angle change.

Then, I calculated the rates as follows:

%% angular velocities
w = zeros(r-1,4);

for i=1:r-1
    qi = quat(i,:);
    q_dot_i = q_dot(i,:);
    temp = 2 * quatmultiply(q_dot_i, quatinv(qi));
    w(i,:) = temp;
end

% removing the unwanted column
w = w(:,2:end);

So, I think I got the rates, but I'm afraid at some point in the computation, they reach values that are around 40rad/s, and in my case this is huge. Do you think it could be caused by the way I computed qdot? If so, do you think there's a better way to find it?

Thank you