Quaternion Rotation: The Blocks and Plasticine Version

Posted on October 20, 2008 by Chris at 12:41 pm

As you know, I’m currently working on Mini-Game: Dwarf vs Zombie. As I’m working on the next section, I’ve come across some code that uses Quaternions to perform the rotation. In the past, I’ve kind of glossed over how they work, because they somehow seem beond my understanding, copied what I can from the available code and moved on. This time, however, I am determined to know exactly how they work, and why we use them as opposed to other rotation methods, such as rotation matrices.

I thought you said this was “Blocks and Plasticine” ?!!!

Before you go too far in, I have made a few assumptions about your maths skills. Firstly, I would expect that you know how to add, subtract and multiply. Secondly, you may want to know basic vector and matrix oprations like cross products and dot products.

If you can’t be bothered with any of that, I’ve included some code and samples at the end to speed you along.

What is a Quaternion?

A quaternion is a mathematical entity that can be used for a variety of operations. Most of their usefulness has been superseded by other areas of maths, but they still provide some significant performance boosts in the calculation of 3d rotations.

Quaternions are commonly expressed as a scalar (w) plus a 3 dimensional vector V, Sometimes, you will also see the notation [w,V].  They are also often just stored as a 4 dimensional vector <w,x,y,z>.

q = w + V
q = [w, V]
q = <w,x,y,z>

Supposing that these are the values, we could create a structure to hold these values like the one below. I’ve also included a few methods that we will cover later.

struct Quaternion {
public:
     float w,x,y,z;

public:
     // constructor
     Quaternion(float qw, float qx, float qy, float qz);

     // FromAxisAngle - creates a quaternion to do a rotation,
     //     from an <x,y,z> axis and an angle in degrees
     static Quaternion FromAxisAngle(float angle, float ax, float ay, float az); 

     // GetConjugate - returns the conjugate of the current quaternion
     Quaternion GetConjugate();

     // Multiply - Multiplies two quaternions togeather
     Quaternion Multiply(Quaternion q2);

     // RotatePoint - Rotates the given point by the current quaternion
     Quaternion RotatePoint(float x, float y, float z);

     // Slerp - Performs Spherical Linear IntERPolation
     Quaternion Slerp(Quaternion q2, float time);
};

A Quaternion for Rotation

The most common use for Quaternions that we find lately is to perform a short cut for rotation about an arbitraty axis.

To perform a rotation, we need a few initial values.

First, we need the <x,y,z> axis that we will be rotating around. This will be a standard vector, which we will call “Axis” in the calculations below.

Next, we will need the angle of rotation. For simplicity, I’ve called this “angle” in the calculations.

We also need to use the point that we will be rotating. This is another <x,y,z> vector, that is convinently named “Point”.

Now we have the basic values, we need to create the following quaternions.

q = [cos(angle/2), Axis * (sin(angle/2))]
q' = [cos(angle/2), -(Axis * (sin(angle/2)))]
P = [0, Point]

The (q) value defined above represents the rotation,

The (q’) value is known as the conjugate, which is simply calculated using one of the formulas below.

q' = [q.w, -q.V]
q' = [q.w, -q.x, -q.y, -q.z]

The (P) value here is simply the point that we want to rotate expressed as a quaternion. In this case, we use the Point vector without modification, and specify zero for the w value.

To make things a little easier, you could add a couple of functions to the Quaternion struct above.

Quaternion::Quaternion(float qw, float qx, float qy, float qz){
     w = qw; x = qx; y = qy; z = qz;
}

Quaternion Quaternion::FromAxisAngle(float angle, float ax, float ay, float az){
     float sinangle, cosangle;
     float length;

     // we need to normalize the axis values so that they are a unit vector
     length = (float) sqrt(ax*ax+ay*ay+az*az);
     ax /= length;
     ay /= length;
     az /= length;

     // precalculate the sin and cos of the angle so we dont need to keep using the sin function
     // NOTE: sin uses radians, and we want degrees as angles, so the formula is
     //     angle = angle * (PI / 180);
     //     sinangle = sin(angle/2);
     // which is the same as saying
     //     sinangle = sin(angle * PI/360);
     sinangle = (float) sin(angle * PI / 360.0);
     cosangle = (float) cos(angle * PI / 360.0);

     return Quaternion(cosangle, ax * sinangle, ay * sinangle, az * sinangle);
}

Quaternion Quaternion::GetConjugate(){
     return Quaternion(w,-x,-y,-z);
}

Now we have our preliminary quaternions defined, we can use the following formulas to calculate the new rotated point.

R = q * P * q'
Result = <R.x, R.y, R.z>

Multiplying the three quaternions q, P, and q’ will give us a third quaternion (R), which we can then extract the the Result vector from it’s x,y,z values. Multiplying two quaternions togeather requires a special formula, which we will cover next.

Quaternion Multiplication

As you can see in the example above, quaternions are multiplied togeather to give us the result that we need. Just as there is a special way of multiplying matrices togeather, there is also a way of multiplying quaternions with each other.

If you know matrix multiplication the formula below will give you the correct multiplication. If you let w1 = q1.w, and V1 = q1.V, and likewise for w2, and V2.

q1 * q2 = w1 * w2 - DotProduct(V1, V2) + w1 * V2 + w1 * V1 + CrossProduct(V1, V2)

If you don’t know matrix multiplication, the formula below is for you. I’ve written it as it would appear in code, so you should just be able to copy and paste :).

Quaternion Quaternion::Multiply(Quaternion q2){
     Quaternion Result;
     float w1, x1, y1, z1;
     float w2, x2, y2, z2;

     w1 = w; x1 = x; y1 = y; z1 = z;
     w2 = q2.w; x2 = q2.x; y2 = q2.y; z2 = q2.z;

     Result.w = w1 * w2 - x1 * x2 - y1 * y2 - z1 * z2;
     Result.x = w1 * x2 + x1 * w2 + y1 * z2 - z1 * y2;
     Result.y = w1 * y2 - x1 * z2 + y1 * w2 + z1 * x2;
     Result.z = w1 * z2 + x1 * y2 - y1 * x2 + z1 * w2;

     return Result;
}

So, now that we have multiplication defined, we are on easy street. We just plug in that nifty “q * P * q’” thingame to get our result… but what if we could be even lazier? How ’bout some code that makes you not need to remember that little rule?

Quaternion Quaternion::RotatePoint(float x, float y, float z){
     Quaternion p(0,x,y,z);
     Quaternion cq = GetConjugate();
     Quaternion result;

     // now for the magic (q * P * q')
     result = Multiply(p).Multiply(cq);
     return result;
}

It’s polite to Slerp your Interpolation

Interpolation is just a nifty word for “moving from one point to another at small intervals”. It is used to make the animation that you see between two keyframes as a smooth transition from one to the next.

To make an angle move smoothly from one quaternion to another, we use a method called Spherical Linear Interpolation (SLERP). The formula for Slerp is given by

q(time) = ((sin(angle * (1 - time))/(sin(angle))) * q1
          + ((sin(angle * time)) / sin(angle)) * q2

where
0.0 <= time <= 1.0  and
(q1 dot q2) >= 0
angle = cos-1(q1 dot q2)

NOTE q === -q in regards to angles, so we can simply change the sign of
     either q1 or q2 to make the rule fit.

Fortunately, I’ve come prepared with some code that you can use to do the task.

Quaternion Quaternion::Slerp(Quaternion q2, float time){
     float dotprod;
     float angle;
     double sinangle;
     float a, b;
     Quaternion result;
     float w1, x1, y1, z1;
     float w2, x2, y2, z2;

     w1 = w; x1 = x; y1 = y; z1 = z;
     w2 = q2.w; x2 = q2.x; y2 = q2.y; z2 = q2.z;

     // get the dot product, and ensure it's >= 0
     // to ensure the shortest path (< 180) is chosen
     dotprod = w1*w2 + x1*x2 + y1*y2 + z1*z2;
     if(dotprod < 0.0f){
          w1 = -w1; x1 = -x1; y1 = -y1; z1 = -z1;
          dotprod = w1*w2 + x1*x2 + y1*y2 + z1*z2;
     }

     // get the angle
     angle = (float) acos(dotprod);
     if(angle == 0){
          return Quaternion(w1,x1,y1,z1);
     }
     sinangle = sqrt(1 - (dotprod * dotprod));

     // calculate the multipliers
     a = (float)(sin(angle * (1.0f - time)) / sinangle);
     b = (float)(sin(angle * time) / sinangle);

     // create the result quaternion
     return Quaternion(a*w1+b*w2, a*x1+b*x2, a*y1+b*y2, a*z1+b*z2);
}

No Really, I Want my Blocks and Plasticine!!!!

Ok, so you’re not interested in all the theory, here’s the basics.

  • Quaternions are used to rotate a point around an axis
  • Download and use my source code [c++][java] (there’s your blocks.. of code)
  • If you are interested in how it all works, read the stuff above. It doesn’t get any easier, but it can be a whole lot harder.
  • Take a look at how to use the code below, and play with the samples (malleable, just like plasticine)

Example: Rotating a point <1,0,0> around the z-axis by 45º.

void main(){
     Quaternion angleQuaternion;
     Quaternion result;

     angleQuaternion = Quaternion::FromAxisAngle(45.0f, 0.0f, 0.0f, 1.0f);
     result = angleQuaternion.RotatePoint(1.0f,0.0f,0.0f);
     cout << "(" << result.x << ", " << result.y << ", " << result.z << ")" << endl;
}

Sample: Rotating point <1,0,0> around the z-axis

>

Example: Using SLERP to rotate a point from 45º to 135º around the z-axis.

void main(){
     Quaternion q1, q2, result;

     q1 = Quaternion::FromAxisAngle(45.0f, 0.0f, 0.0f, 1.0f);
     q2 = Quaternion::FromAxisAngle(135.0f, 0.0f, 0.0f, 1.0f);

     for(int i = 0; i < 10; i++){
          result = q1.Slerp(q2, 0.1f * i);
          cout << "(" << result.x << ", " << result.y << ", " << result.z << ")" << endl;
     }
}

Sample: SLERP around the z-axis


Why Use Quaternion Rotation?

The reasons for using quaternions in your rotations kind of play out like an add for a stocktake sale.

  • Less Storage Space - A quaternion holds 4 values, while a 3×3 matrix holds 9 values. “That’s a saving of more than 50%.. what a bargain!”
  • Less Operations - Each time you multiply two quaternions together, you need to perform 16 operations, while a 3×3 matrix will require 27, and a 4×4 will need 64… “that’s up to 75% off!!!”
  • Easier to Interpolate - Making a smooth transition between two quaternions is much easier than doing the same between two matrices. “open 24 hours.. hurry in while stocks last”
Filed under: Articles, Game Development Tags: , , ,

Ooh, a Wizard.. sounds like my kind of toon

Posted on October 16, 2008 by Chris at 6:31 pm

For some reason, I just cant keep myself away from playing the bam-bam caster types… in Diablo my first toon was a Mage, in Diablo 2 I had about 3 Sorceresses, in WoW my first toon was a Mage. It would seem that whatever dps caster Blizzard comes out with, I’m willing to lap it up.. and the new Wizard class that was recently announced for Diablo 3 will likley satisfy my zap-lust.

It does seem strange that Blizzard are deciding to choose Wizard for this class. It’s almost as if they are going to lengths to avoid using the same class name between different versions of Diablo.. with the exception of the Barbarian, which seems to have made his way from Diablo 2 to Diablo 3 without a name change. Me being stuck in my own ways, I will still call it a mage.. my Sourceresses are also called mages.. probably more to do with the fact that I like one syllable words more than anything else.

Perhaps I could look up a thesaurus and find a few more names that Blizzard can call a mage for other games.. like Diablo 4, 5, 6, etc… Conjurer, Thaumaturge, Warlock..

Chime in with suggestions if you can think of more synonims for mage :)

Filed under: General Tags:

Mini-Game: Dwarf vs Zombie (Part 2)

Posted on October 2, 2008 by Chris at 6:21 pm

Well, here I am again writing my next part of the Dwarf Vs Zombie game. Feels like months since I’ve actually been able to get to sit down and do some more work on this.

This one’s not going to be a big update, even though it took more work than Part 1, there isn’t really much to say.

Articles
   Part 1 - Stuff to download, getting into the development
   Part 2 - Setting the Scene, and Model Drawing

Downloads
   Source: Click Here

Setting the Scene

Before you look too far, you will notice that I wrote about a Scene Graph back in Building a Solar System in XNA (Part 2). It’s quite a simple system to implement, and is powerful enough to handle whatever I throw at it.

I’ve been at it again, this time writing one in C++ to handle the 3d environment we are working with. Here’s a little code that shows you how we can move recursively through the graph to draw the entire tree from bottom to top, without needing to have any additional code in the ISceneNode interface.

void SceneGraph::Render(){
   Render(m_BaseNode);
}

void SceneGraph::Render(ISceneNode *node){
  SceneNodeIterator it;
  node->Render();
  it = node->begin();
  while(it != node->end()){
    Render(*it);
    it++;
  }
}

You might have recognised the sneaky use of STL to implement the node list.. why re-invent the list, when you can simply inherit it :). Here’s the code for our SceneNode interface. Please DO remember the virtual destructor, unless you would prefer some wild and zany adventures with memory leaks later.

class ISceneNode : public std::list<ISceneNode *>{
public:
  virtual ~ISceneNode();
  virtual void Render() = 0;
  virtual void Update(float fStep) = 0;
};

typedef std::list<ISceneNode *>::iterator SceneNodeIterator;

Playing with Models (again)

So, what good is a Scene Graph without stuff to fill it?.. The next step is to get our models that we downloaded a couple of months ago showing up in the scene. To do this, I’ve written the ModelNode class that is responsible for the drawing of the models, and an MS3DLoader class that can read the data from a ‘.ms3d’ file.

I’ve basically inferred the MS3D file format from the MilkShape 3D Binary Model Viewer source code that you can download from the MilkShape site, and ignored the parts that I don’t intend to use. So, if you are after an MS3D loader, I’d suggest that you go there to find a full featured one.

Once the data is loaded, we can then use our ModelNode to render the models into our scene.

First, we will position and scale our models, so that they are about the right size and place. It is one of the problems with using art from multiple places, they are usually all different scales, and often face different directions when you get them.. this is just a simple hack to get them the same size and direction.

void ModelNode::Render(){
  int nTriangleIndex;
  int nVertIndex;
  int nMatIndex;

  glPushMatrix(); // push the matrix onto the stack
  glTranslatef(m_InitialPosition[0],m_InitialPosition[1],m_InitialPosition[2]);
  glRotatef(m_InitialRotation, 0.0f, 1.0, 0.0f);
  glScalef(m_Scale,m_Scale,m_Scale);

  // Do the rest of the rendering in here
  [...]

  glPopMatrix(); // pop the matrix we pushed earlier
}

Now we iterate through each of the model’s groups, and set the material information for each one. I’m almost certain that we use the glMaterial function for this.. but being that neither of the models do anything interesting with the materials, we are unlikely to find out if this is a mistake right now or not.

We also tell the system to use the appropriate texture here. If you look at the code, you will notice that I dont use the texture information from the ms3d file, but rather allocate it manually during the loading process.

  glColor3f( 1.0f, 1.0f, 1.0f);
  for(unsigned int i = 0; i < m_Data.Groups.size(); i++){
    // get the material information
    nMatIndex = m_Data.Groups[i].MaterialIndex;
    glMaterialfv(GL_FRONT_AND_BACK, GL_AMBIENT, m_Data.Materials[nMatIndex].Ambient);
    glMaterialfv(GL_FRONT_AND_BACK, GL_DIFFUSE, m_Data.Materials[nMatIndex].Diffuse);
    glMaterialfv(GL_FRONT_AND_BACK, GL_EMISSION, m_Data.Materials[nMatIndex].Emissive);
    glMaterialf(GL_FRONT_AND_BACK, GL_SHININESS, m_Data.Materials[nMatIndex].Shininess);
    glMaterialfv(GL_FRONT_AND_BACK, GL_SPECULAR, m_Data.Materials[nMatIndex].Specular);

    // set up the texture
    if(glIsTexture(m_Data.Materials[nMatIndex].TextureId)){
      glEnable(GL_TEXTURE_2D);
      glBindTexture(GL_TEXTURE_2D, m_Data.Materials[nMatIndex].TextureId);
    } else {
      glDisable(GL_TEXTURE_2D);
    }

    // we draw the triangles in here
    [...]
  }

Finally, within each group we draw the individual triangles, complete with normals, and texture co-ordinates from the loaded data. You can just ignore the funky business going on with the texture coordinates..

I thought I was being smart earlier when I wrote the loading process, but it loaded them in the order s1-s2-s3-t1-t2-t3 rather s1-t1-s2-t2-s3-t3, so I couldn’t use the glTexCoord3fv. Not that it really matters, they are all there.

    glBegin(GL_TRIANGLES);
    for(unsigned int j = 0; j < m_Data.Groups[i].Indices.size(); j++){
      nTriangleIndex = m_Data.Groups[i].Indices[j];

      // get the vertex information
      for(int k = 0; k < 3; k++){
        nVertIndex = m_Data.Triangles[nTriangleIndex].Indeces[k];

        glNormal3fv(&m_Data.Triangles[nTriangleIndex].Normals[k * 3]);
        glTexCoord2f(
          m_Data.Triangles[nTriangleIndex].TextureCoords[k],
          m_Data.Triangles[nTriangleIndex].TextureCoords[k + 3]);
        glVertex3fv(m_Data.Vertices[nVertIndex].verts);
      }
    }
    glEnd();