Quaternions angle error - correct?

[bkenwright]

Hi All,
I think this is right, but not sure after doing some of the maths.

If I have two rotated objects...lets say two sticks...and each has a rotation, in quaternions q0 and q1.

Now the difference, can be calculated as

qdiff = q0 * Conjugate( q1 )

Okay?

Of course both my object quaternions are normalized rotation quaternions.

But if I calculated the axis angle for my qdiff, I would have guess it to give the shortest angle error to rotate one to the other?...and the axis would be that by which its rotated to get it?

But I found after some debugging it doesn't seem right...if I work out the rotation angle on the two sticks, the rotation axis is different, and if I do the dot product on the two sticks to get the minimum angle to rotate one to the other, its different than what my qdiff axis angle value gives.

Is there something I' missing?

Many Thanx,

Ben.

 

[D H]

When computing q_diff = q₀ q₁^*, the eigen axis of q_diff doesn't have much meaning. It is the eigen angle that is meaningful.

 

[bkenwright]

Hi D H,
thanks for your reply.

What would you say is the best method to get the minimum axis angle between two objects then, if you have two quaternions representing there rotations, and you want to know the axis and angle to rotate on to match the other?

Thanx,

Ben.

 

[D H]

One way to get what you want is to use the quaternion "difference" (obviously not really a difference, but it is called that).

I'm assuming here that you have placed your quaternion into canonical form (unit quaternion with a non-negative real part). You can use either the scalar real part of the quaternion or the vector magnitude of the vector imaginary part of the quaternion to give you the desired angle. So, which should you use? The answer is the smaller of the two. The inverse sine is more accurate for angles between 0 and pi/4 while inverse cosine is better for angles between pi/4 and pi/2.

Note that either acos(real part) and asin(||vector part||) will give an answer between 0 and pi/2. If the angle is between pi/2 and pi you will need to subtract the computed angle from pi to yield the correct angle.

 

[bkenwright]

Hi D H,
thanks for your help. I've been fiddling around with what you said, but I think I must be missing something. I did a code snippet so you could see...in a few cases the returned angle difference doesn't seem to be correct.

And example case:

// Some normalized rotation quaternions from bodies
Quaternion qrot0 = new Quaternion( new Vector3(-0.15f, 0.68f, 0.15f), 0.68f );
Quaternion qrot1 = new Quaternion( new Vector3(-0.19f, -0.67f, 0.67f), 0.19f );

// Work out the difference
Quaternion qdiff = qrot0 * Quaternion.Conjugate(qrot1);


// Get Axis angle and see what we get (angle is in the w component)
// just acos, then x2
float angle;
Vector3 axis;
ToAxisAngle(ref qdiff, out axis, out angle);

// We get, angle = 3.53900266 ? should be 2.4ish?

// Logic checks to see if our angles are right!
Matrix mat0 = Matrix.CreateFromQuaternion( qrot0 );
Matrix mat1 = Matrix.CreateFromQuaternion( qrot1 );

Vector3 checkDir0 = new Vector3(mat0.M31, mat0.M32, mat0.M33);
Vector3 checkDir1 = new Vector3(mat1.M31, mat1.M32, mat1.M33);

float checkAngle = acos( Vector3.Dot( checkDir0, checkDir1 ) );
// testAngle => 2.471772 /// expected value

 

[D H]

Hi D H,
thanks for your help. I've been fiddling around with what you said, but I think I must be missing something. I did a code snippet so you could see...in a few cases the returned angle difference doesn't seem to be correct.

And example case:

// Some normalized rotation quaternions from bodies
Quaternion qrot0 = new Quaternion( new Vector3(-0.15f, 0.68f, 0.15f), 0.68f );
Quaternion qrot1 = new Quaternion( new Vector3(-0.19f, -0.67f, 0.67f), 0.19f );

Those are *not* unit quaternions, not even in single precision land.
Try normalizing your quaternions.

// Work out the difference
Quaternion qdiff = qrot0 * Quaternion.Conjugate(qrot1);// Get Axis angle and see what we get (angle is in the w component)
// just acos, then x2
float angle;
Vector3 axis;
ToAxisAngle(ref qdiff, out axis, out angle);

// We get, angle = 3.53900266 ? should be 2.4ish?

Should be 2.7316764521066457 ish. I would expect some error from using single precision arithmetic. However, you have done something here that has made you lose all precision.

Try normalizing your quaternions.

// Logic checks to see if our angles are right!
Matrix mat0 = Matrix.CreateFromQuaternion( qrot0 );
Matrix mat1 = Matrix.CreateFromQuaternion( qrot1 );

Vector3 checkDir0 = new Vector3(mat0.M31, mat0.M32, mat0.M33);
Vector3 checkDir1 = new Vector3(mat1.M31, mat1.M32, mat1.M33);

float checkAngle = acos( Vector3.Dot( checkDir0, checkDir1 ) );
// testAngle => 2.471772 /// expected value

That is just the angle between the z axes. It is not the eigen rotation angle.

Think of it this way. Consider the matrices A and B:

$$A = \bmatrix<br /> \phantom{-}1 & \phantom{-}0 & \phantom{-}0 \\<br /> \phantom{-}0 & \phantom{-}1 & \phantom{-}0 \\<br /> \phantom{-}0 & \phantom{-}0 & \phantom{-}1<br /> \endbmatrix$$

$$B = \bmatrix<br /> -1 & \phantom{-}0 & \phantom{-}0 \\<br /> \phantom{-}0 & -1 & \phantom{-}0 \\<br /> \phantom{-}0 & \phantom{-}0 & \phantom{-}1<br /> \endbmatrix$$

By your logic, the rotation between these matrices is zero. That obviously is not the case.

 

[bkenwright]

I'd cut the decimal places a few when I dumped the values out, the full values are:
(which are normalized)
Quaternion qrot0 = new Quaternion( new Vector3(-0.1564176f, 0.6895894f, 0.1564176f), 0.6895894f );
Quaternion qrot1 = new Quaternion( new Vector3(-0.1965395f, -0.6792439f, 0.6792439f), 0.1965395f );

All the quaternion rotations are normalized and I check for that in the code.


I draw the two bodies, and extract the extruded body length which I used to check for the angle (as can be seen in the screenshot).

https://www.physicsforums.com/attachment.php?attachmentid=31983&stc=1&d=1297024503

Not sure why it works for smaller angles, less than 90. I'll have a poke around and see why.

Thanx,

Ben

 

[D H]

It works for all angles.

In your example, the quaternion product q₀q₁^* is -0.2035261 - 0.4695329 i + 0.5287143 j - 0.6771832 k. Note that the magnitude of the vector part is 0.9790695. The angle can be computed using 2*acos(0.2035261) or by 2*asin(0.9790695), both of which yield 2.731676.

The matrices corresponding to q0 and q1 are

$$T_0 = \bmatrix<br /> -1.457167719820518\text{e-16} & -0.42070530006186835 & 0.90719736028047004 \\<br /> 0 & 0.90719736028046993 & 0.42070530006186835 \\<br /> -0.99999999999999989 & 0 & -1.457167719820518\text{e-16}<br /> \endbmatrix$$

$$T_1 = \bmatrix<br /> -0.8511340206185567 & 0 & -0.52494845360824727 \\<br /> 0.52494845360824727 & -2.2204460492503131\text{e-16} & -0.85113402061855647 \\<br /> 0 & -0.99999999999999978 & -2.2204460492503131\text{e-16}<br /> \endbmatrix$$

(Note: I use the left quaternion convention, v -> qvq*. If you use the right quaternion convention, v-> q*vq the matrices will be the transpose of the above.)

The product T₀T₁^T is

$$T_0T_1^T = \bmatrix<br /> -0.47623185139671659 & -0.77214653675005762 &0.42070530006186801 \\<br /> -0.22084859669227144 & -0.35807659353719445 & -0.90719736028046982 \\<br /> 0.8511340206185567 & -0.52494845360824705 & 3.2355623065705562\text{e-32}<br /> \endbmatrix$$

There are a number of ways to compute the single axis rotation from a transformation matrix. One simple way is to use the trace of the matrix: tr(T) = 1 + 2*cos(theta). The trace of the above matrix is -0.8343084449339111. Thus cos(theta) = -0.91715422246695555, or theta = 2.731676, which is the same as the value calculated using the quaternions.

(Note: The trace is not the most accurate technique in this case.)