Quaternions

Note that a mathematica ".nb" file and all necessary accompaniments are included in the file Quaternion.zip.

Extensions implemented by Peter Meulbroek.  All extensions are covered by copywrite, but are released under the GNU public license (*copyleft *).   For a copy of the license agreement, check out http://www.fsf.org/copyleft/gpl.html

We will extend the standard Mathematica package "Algebra`Quaternions`".  Note that this package is now on "MathSource", the Wolfram Research's repository of mathematica packages.  It can be found here.

Introduction

This package extends the Quaternion package of Jason Kastner and Matthew Markert to add the functionality to use quaternions as surrogates for rotation.  Quaternions can represent a rotation by defining the plane of the rotation and the angle of the rotation.  The plane is defined as normal to the quaternion vector, and the angle is defined as alpha, where the pure component of the quaternion is Cos[alpha/2].  The extensions add a rotation function, and the included documentation describes how to use quaternions to replace e.g., rotation matricies.

Several websites have been very helpful with understanding the graphical implications of quaternions . They are

Though the mathematics of the existing mathematica package is nice, it is not easy to use in the context of graphics.  I needed to play with quaternions as a replacement for rotation matrices (bad, bad , bad), and thus,  modified this package.  In this context, given a point {a,b,c} a quaternion {w,x,y,z} represents a rotation of α within the plane defined the point {a,b,c} normal to the vector {x,y,z}, where α is defined by w = [Graphics:Images/index_gr_2.gif].  Most people solving a rotation problem use principle rotation matrices (rotations around each major axes) or euler angles (rotations around z-, x-, z-axes).  The problem with these methods are twofold:

1:  the system is overdetermined.  For example, the rotation matrices have too many elements (some of which are identical).  This means computational and storage expense, as well as a more insidious problem: drift.  Due to roundoff errors in computation, the values that are supposed to be identical can drift apart.  When this occurs, the matrix changes from a simple rotation to a shearing rotation, resulting in distortion.

2:  Both systems suffer from "gimble lock", where a rotation maps one principle axis onto another, resulting in a loss of a degree of freedom.  

However, quaternions solve both these problems.  Any rotation of a point to another point occurs within a plane around a single axis, so no gimble lock.  All members of the quaternion are significant, so no drift.  Finally, as a bonus, it becomes easy to interpolate between quaternions, hence we can do smooth animations using them.

Setup

[Graphics:Images/index_gr_3.gif]

Note that we change the format of quaternions, so they print nicely.

[Graphics:Images/index_gr_4.gif]
<< 1 1 1 1 >>

Translations

To use quaternions in graphics, we need to be able to create them from a variety of sources.  

To/From Vectors

Test the new translation rules

[Graphics:Images/index_gr_5.gif]
<< 1 1 2 3 >>
[Graphics:Images/index_gr_6.gif]
[Graphics:Images/index_gr_7.gif]

First, from/to rotation matrices.  Note that these functions are inverse to each other modulo a constant

Matrix to Quaternions

  A rotation may be converted back to a quaternion through the use of  the following algorithm:

Produce a Quaternion for a known rotation.  A rotation around the z-axis has the following form:

[Graphics:Images/index_gr_8.gif]

We will test a rotation of π/2.  Note that this would tend to cause "gimble lock".

[Graphics:Images/index_gr_9.gif]
<< [Graphics:Images/index_gr_10.gif] 0 0 [Graphics:Images/index_gr_11.gif] >>

This is correct. [Graphics:Images/index_gr_12.gif] = [Graphics:Images/index_gr_13.gif] ⇒ α ⩵ [Graphics:Images/index_gr_14.gif];  Around the axis (0,0,z).  Note that the quaternion is normalized

Quaternions to Matrix

We now try to reverse this.  

[Graphics:Images/index_gr_15.gif]
[Graphics:Images/index_gr_16.gif]
[Graphics:Images/index_gr_17.gif]
<< [Graphics:Images/index_gr_18.gif] 0 0 [Graphics:Images/index_gr_19.gif] >>

Hm.  The problem is that the original Quaternion isn't normalized.  We get back the right FORM of the quaternion, but no longer symmetric.

Normalising Quaternions

Test the normalization function.  We can also assure ourselves that the To and From functions work

[Graphics:Images/index_gr_20.gif]
<< [Graphics:Images/index_gr_21.gif] 0 0 [Graphics:Images/index_gr_22.gif] >>
[Graphics:Images/index_gr_23.gif]
[Graphics:Images/index_gr_24.gif]
[Graphics:Images/index_gr_25.gif]
<< [Graphics:Images/index_gr_26.gif] 0 0 [Graphics:Images/index_gr_27.gif] >>

Bingo!

Quaternion Inverse

[Graphics:Images/index_gr_28.gif]
<< 1 [Graphics:Images/index_gr_29.gif] 0 [Graphics:Images/index_gr_30.gif] >>
[Graphics:Images/index_gr_31.gif]
<< 3 0 0 0 >>

The conjugate is actually the inverse, modulo a constant.  we define a  new implementation of Inverse that takes this into account

[Graphics:Images/index_gr_32.gif]
<< [Graphics:Images/index_gr_33.gif] [Graphics:Images/index_gr_34.gif] 0 [Graphics:Images/index_gr_35.gif] >>
[Graphics:Images/index_gr_36.gif]
<< 1 0 0 0 >>

We get back a unit Quaternion.

[Graphics:Images/index_gr_37.gif]
[Graphics:Images/index_gr_38.gif]

Rotations

The point of the extensions are to use them to predict and implement rotations.  Lets experiment on a geometric figure.

[Graphics:Images/index_gr_39.gif]

we define a set of axes to make the rotations clearer.  The axes are: red=x, green=y, blue=z.

[Graphics:Images/index_gr_40.gif]

[Graphics:Images/index_gr_41.gif]

[Graphics:Images/index_gr_42.gif]

Define a function to plot any points generated by the functions that follow.  To help visualize the points, we'll connect them with polygons.  To make clear what's going on, we print the list of points that are graphed.

[Graphics:Images/index_gr_43.gif]

Rotations

Helper functions

Quaternions are defined (for some, odd reason) using 1/2 angles.  Thus, to easily implement them, and to avoid inverse trig functions, we'll use the 1/2 angle formulas In the functions below, cα  is the cosine of the angle alpha.  the functions return cos[α/2] and sin[α/2], resp.

[Graphics:Images/index_gr_44.gif]
[Graphics:Images/index_gr_45.gif]
[Graphics:Images/index_gr_46.gif]

CreateQuaternion

We are given two points, and want to create the quaternion that defines the least rotation of one point to get to the other. To create the quaternion, we must define the rotation plane (using a normal vector), and the angle of rotation.  If we assume that the points are defined by unit vectors, then the axes of rotation is the cross product, and the angle of rotation is the dot product (otherwise, it gets tricky).  Note that we define rotations around the origin.  

A special note:  in the case the defining points form a line, the cross product fails, since an infinite number of planes are defined by a line running through a point.  In that case, we arbitrarily choose to rotate through an easy plane to visualize, the yz plane, defined by the x-axis

[Graphics:Images/index_gr_47.gif]
[Graphics:Images/index_gr_48.gif]
[Graphics:Images/index_gr_49.gif]

Note that the default tetrahedron is vertical.  This leads to all sorts of issues, since spherical coordinates have a pole (no pun intended) at z = 0.

[Graphics:Images/index_gr_50.gif]
[Graphics:Images/index_gr_51.gif]

[Graphics:Images/index_gr_52.gif]

[Graphics:Images/index_gr_53.gif]

Example 1:  rotation around the z axis

Create a quaternion that rotates {1,0,0} to the vector defined by {1,1,0}.  Note that we implicitly normalize the points, so that the rotation is a simple one.

[Graphics:Images/index_gr_54.gif]
<< [Graphics:Images/index_gr_55.gif] 0 0 [Graphics:Images/index_gr_56.gif] >>

Rotate the tetrahedra using this  quaternion

[Graphics:Images/index_gr_57.gif]
[Graphics:Images/index_gr_58.gif]

Get back a set of vectors to plot

[Graphics:Images/index_gr_59.gif]
[Graphics:Images/index_gr_60.gif]

Plot both the original and the rotated version

[Graphics:Images/index_gr_61.gif]
[Graphics:Images/index_gr_62.gif]
[Graphics:Images/index_gr_63.gif]

[Graphics:Images/index_gr_64.gif]

[Graphics:Images/index_gr_65.gif]

Example 2:  rotation around an arbitrary axis

Create a quaternion that rotates {1,0,0} to the vector defined by {1,1,1}.  Note that we implicitly normalize the points, so that the rotation is a simple one.

[Graphics:Images/index_gr_66.gif]
<< [Graphics:Images/index_gr_67.gif] 0 [Graphics:Images/index_gr_68.gif] [Graphics:Images/index_gr_69.gif] >>

Rotate the tetrahedron using this  quaternion

[Graphics:Images/index_gr_70.gif]
[Graphics:Images/index_gr_71.gif]

Get back a set of vectors to plot

[Graphics:Images/index_gr_72.gif]
[Graphics:Images/index_gr_73.gif]

Plot both the original and the rotated version

[Graphics:Images/index_gr_74.gif]
[Graphics:Images/index_gr_75.gif]
[Graphics:Images/index_gr_76.gif]

[Graphics:Images/index_gr_77.gif]

[Graphics:Images/index_gr_78.gif]