How to integrate 3D rotation with quaternions

How to integrate 3D rotation with quaternions

Post by John Hencke » Tue, 20 Apr 1999 04:00:00



I am writing a 3D rigid body simulation. I use a quaternion for the
orientation, but I use a "spin vector" for the angular velocity. The
direction of the spin vector (omega) is the spin axis and the magnitude is
the rate (radians per second). I think that this is a common practice.

When I step forward in time, for example using Euler integration, how do I
add the spin vector to the quaternion?

Given:
    dt = duration of timestep
    q = orientation of the rigid body (a quaternion)
    omega = angular velocity of the body (a 3d vector)

How do I compute Euler integration:

    q' = q + dt*omega

I have come up with four alternatives,

(1)
   a. Convert q to a spin vector, s
   b. compute r = s + dt*omega using vector arithmetic
   c. q' = convert r back to a quaternion

(2)
   a. convert dt*omega to a quaternion, q2.
   b. use quaternion composition q' = q * q2.

   (You can see definition of quaternion conversion and composition
    at http://www.cs.berkeley.edu/~laura/cs184/quat/quaternion.html )

(3)
   a. convert omega to a quaternion, q2.
   b. q' = SLERP(dt from q to q2)

(4)
   a. compute dq using the matrix...

                |-x -y -z |
     dq = 1/2 * | w -z  y | * omega
                | z  w -x |
                |-y  x  w |

   b. compute q' = q + dt*dq using element-wise arithmetic,
      i.e. ( q.w + dt*dq.w, q.x + dt*dq.x, ....)
   c. normalize q' to the unit hypersphere

My opinion is that (1) is dead wrong, (2) and (3) might be equivalent and
are the most accurate, (4) is ok for very small dt*omega, but is bad as q
+ dt*dq leaves the surface of the hypersphere.

Any help or pointers are GREATLY APPRECIATED!!!

If you post a reply, please cc: my email address.


 
 
 

How to integrate 3D rotation with quaternions

Post by John Blackbur » Tue, 20 Apr 1999 04:00:00




> I am writing a 3D rigid body simulation. I use a quaternion for the
> orientation, but I use a "spin vector" for the angular velocity. The
> direction of the spin vector (omega) is the spin axis and the magnitude is
> the rate (radians per second). I think that this is a common practice.

> When I step forward in time, for example using Euler integration, how do I
> add the spin vector to the quaternion?

> Given:
>     dt = duration of timestep
>     q = orientation of the rigid body (a quaternion)
>     omega = angular velocity of the body (a 3d vector)

> How do I compute Euler integration:

>     q' = q + dt*omega

The correct formula (for the relationship between q, the rotation
quaternion, and w, the angular velcoity vector) is:

      dq/dt =1/2 w * q

where 'w * q' is shorthand for [0,w] * q, i.e. the quaternion product of a
quaternion derived from w and q. From this you get:

      dq = dt * 1/2 w * q

giving:

      q' = q + dt * 1/2 w * q

for the calculation step.

Note that the sign of this/order of the product depends on your
conventions for the way you measure angular velocity (e.g. which direction
is +ve) and the way you use quaternions to describe rotations.

After doing this you may want to re-normalise q, as what you're in effect
doing is using a straight tangent to approximate a curve, and this will
steadliy diverge from true: although a non-unit quaternion can be used for
rotations it's often easisest to stick with unit quaternions at all times,
as then many calcuations (in paticular inverses) are much simpler. You can
also improve accuracy by decreasing the step size, i.e. dt, and doing
additional calculations per frame.

For small enough dt and over short periods this produces pretty good
results: often friction is more significant than the inaccuracies of this
approach, or something happens (such as a collision) before you notice
problems.

Also remember that for general rigid body motion, even with no external
forces/torques, w is non-constant and needs to be re-calculated each
iteration.

John