3-axes quaternion rotation in OpenGL

Question

I'm trying to create an OpenGL program where the model of a bird is supposed to follow a defined path along the surface of a sphere described by Seiffert's spherical spiral. However, I've been stuck on getting the rotations right for quite some time now.

As a first step, I make the bird just follow a circular path in the x-z-plane:

// 1. Circle in x-z plane
float phi =  TWO_PI * t; // t = [0..1]

float x = boundingSphereRadius * cos(phi);
float y = 0.0f;
float z = boundingSphereRadius * sin(phi);

float rotationAngle = glm::orientedAngle(glm::vec3(0.0f, 0.0f, 1.0f),    
                                         glm::normalize(glm::vec3(x, 0, z)),
                                         glm::vec3(0.0f, 1.0f, 0.0f)) - HALF_PI;
glm::fquat rotation = glm::angleAxis(rotationAngle, glm::vec3(0.0f, 1.0f, 0.0f));

The fixed -HALF_PI is necessary so that the bird is correctly aligned. This works perfectly fine, and similarly I could achieve a circular rotation in the x-y- and the y-z-plane.

The problem occurs when I try to accumulate all the different rotations. The path I'm trying to follow looks like this:

As a requirement, the bird's belly is supposed to always face the surface of the sphere and the bird should fly in the forward direction.

My current approach looks like this, which consists of just combining the three orientation quaternions:

glm::fquat rotationX  = glm::angleAxis(glm::orientedAngle(glm::normalize(glm::vec3(0.0f, 0.0f, 1.0f)), glm::normalize(glm::vec3(x, 0, z)), glm::vec3(0.0f, 1.0f, 0.0f)) - HALF_PI, glm::vec3(0.0f, 1.0f, 0.0f));
glm::fquat rotationY1 = glm::angleAxis(-HALF_PI, glm::vec3(0.0f, 1.0f, 0.0f));
glm::fquat rotationY2 = glm::angleAxis(glm::orientedAngle(glm::vec3(0.0f, 1.0f, 0.0f), glm::normalize(glm::vec3(x, y, 0)), glm::vec3(0.0f, 0.0f, 1.0f)), glm::vec3(0.0f, 0.0f, 1.0f));
glm::fquat rotationY  = rotationY2 * rotationY1;
glm::fquat rotationZ  = glm::angleAxis(glm::orientedAngle(glm::vec3(0.0f, 0.0f, 1.0f), glm::normalize(glm::vec3(0, y, z)), glm::vec3(1.0f, 0.0f, 0.0f)) + HALF_PI, glm::vec3(1.0f, 0.0f, 0.0f));
glm::fquat rotation   = rotationZ * rotationY * rotationX;

However, the orientation changes are completely wrong and at some angles there are jumps occurring.

EDIT:

I'm trying different circles on the sphere now where more than one rotation is necessary. For beta = gamma = 0.0f and alpha = HALF_PI the circle is again in the x-z-plane and the value of rotationAngleXZ is changing, while rotationAngleXY is either -HALF_PI of HALF_PI and rotationAngleYZ is either 0.0f or PI. I guess I'm running into a Gimbal Lock here and I've read a multitude of articles about it, however I'm still not sure how I can prevent it in this case.

// 10. `Arbitrary` circles on sphere surface
// http://math.stackexchange.com/questions/643130/circle-on-sphere
//
// Parameters:
//      alpha = 0...HALF_PI - For alpha = 0, the circle is just a point - For alpha = HALF_PI, the circle is a Great Circle
//      (beta, gamma) = center of circle in spherical coordinates
float phi =  TWO_PI * t;

float x = boundingSphereRadius * ( (sin(alpha) * cos(beta) * cos(gamma)) * cos(phi) + (sin(alpha) * sin(gamma)) * sin(phi) - (cos(alpha) * sin(beta) * cos(gamma)));
float y = boundingSphereRadius * ( (sin(alpha) * sin(beta)) * cos(phi) + cos(alpha) * cos(beta));
float z = boundingSphereRadius * (-(sin(alpha) * cos(beta) * sin(gamma)) * cos(phi) + (sin(alpha) * cos(gamma)) * sin(phi) + (cos(alpha) * sin(beta) * sin(gamma)));

float rotationAngleXZ = glm::orientedAngle(glm::normalize(glm::vec3(0.0f, 0.0f, 1.0f)), glm::normalize(glm::vec3(x, 0, z)), glm::vec3(0.0f, 1.0f, 0.0f));
std::cout << "Rotation Angle XZ = " << rotationAngleXZ << std::endl;
glm::fquat rotationXZ = glm::angleAxis(rotationAngleXZ - HALF_PI, glm::vec3(0.0f, 1.0f, 0.0f));

float rotationAngleXY = glm::orientedAngle(glm::vec3(0.0f, 1.0f, 0.0f), glm::normalize(glm::vec3(x, y, 0)), glm::vec3(0.0f, 0.0f, 1.0f));
std::cout << "Rotation Angle XY = " << rotationAngleXY << std::endl;
glm::fquat rotationXY_Y = glm::angleAxis(-HALF_PI, glm::vec3(0.0f, 1.0f, 0.0f));
glm::fquat rotationXY_Z = glm::angleAxis(rotationAngleXY, glm::vec3(0.0f, 0.0f, 1.0f));
glm::fquat rotationXY = rotationXY_Z * rotationXY_Y;

float rotationAngleYZ = glm::orientedAngle(glm::vec3(0.0f, 0.0f, 1.0f), glm::normalize(glm::vec3(0, y, z)), glm::vec3(1.0f, 0.0f, 0.0f));
std::cout << "Rotation Angle YZ = " << rotationAngleYZ << std::endl;
glm::fquat rotationYZ = glm::angleAxis(rotationAngleYZ + HALF_PI, glm::vec3(1.0f, 0.0f, 0.0f));

glm::fquat rotation = glm::normalize(rotationXZ) * glm::normalize(rotationXY) * glm::normalize(rotationYZ);

Answer 1

Your code is using Euler angles (axis-aligned rotations). The wiggles and jumps are because Euler angles are a bad parameterization of the space of 3D rotations. Instead, here are two alternate approaches.

Constructing a Rotation Matrix via a Frame

Assuming the bird is pointing down the x-axis and up in it's own local coordinate system.

Let p = [x y z] be the position of the bird. Let v be its velocity vector. Let

f = v/|v|
up = p/|p|
s = cross(f, up)

Now construct the matrix with rows f, up, s. Specifically:

[  f[0]   f[1]  f[2] ]
[ up[0]  up[1] up[2] ]
[  s[0]   s[1]  s[2] ]

Then generate a quaternion via GLM's quat_cast function.

Avoid gluLookAt because it uses the deprecated fixed function matrix stack.

Building up via Rotations (Quaternions)

Let R0 be the rotation from i to f. (The angle is acos(dot(i,f)) and axis is cross(i,f))

Let R1 be the rotation from R0*j to up. (Using matrix multiplication notation since it is easier in this context)

Let R2 be the rotation from R1*R0*k to s.

The final rotation should be R2*R1*R0. Check that this rotation is equal to the matrix above.

Answer 2

I don't have ready code for you, but how about this idea? Assuming you already have a formula for the x,y,z location of the bird as a function of t (Seiffert's spherical spiral). Then:

eye    = fn(t)
center = fn(t + dt) // where will the bird be in the next time-step
up     = normalize(eye - sphereCenter)

Now, gluLookAt(eye, center, up) will provide a matrix and you should be able to use that to orientate your bird.

This reference might also help: https://gamedev.stackexchange.com/questions/41940/how-does-glulookat-work.

Hope this helps,

--Roger