So I've been learning about quaternions recently and decided to make my own implementation. I tried to make it simple but I still can't pinpoint my error. x/y/z axis rotation works fine on it's own and y/z rotation work as well, but the second I add x axis to any of the others I get a strange stretching output. I'll attach the important code for the rotations below:(Be warned I'm quite new to cpp).

Here is how I describe a quaternion (as I understand since they are unit quaternions imaginary numbers aren't required):

struct Quaternion {
    float w, x, y, z;
};

The multiplication rules of quaternions:

Quaternion operator* (Quaternion n, Quaternion p) {
    Quaternion o;

    // implements quaternion multiplication rules:
    o.w = n.w * p.w - n.x * p.x - n.y * p.y - n.z * p.z;
    o.x = n.w * p.x   n.x * p.w   n.y * p.z - n.z * p.y;
    o.y = n.w * p.y - n.x * p.z   n.y * p.w   n.z * p.x;
    o.z = n.w * p.z   n.x * p.y - n.y * p.x   n.z * p.w;

    return o;
}

Generating the rotation quaternion to multiply the total rotation by:

Quaternion rotate(float w, float x, float y, float z) {
    Quaternion n;

    n.w = cosf(w/2);
    n.x = x * sinf(w/2);
    n.y = y * sinf(w/2);
    n.z = z * sinf(w/2); 

    return n;
}

And finally, the matrix calculations which turn the quaternion into an x/y/z position:

inline vector<float> quaternion_matrix(Quaternion total, vector<float> vec) {
    float x = vec[0], y = vec[1], z = vec[2];

    // implementation of 3x3 quaternion rotation matrix:
    vec[0] = (1 - 2 * pow(total.y, 2) - 2 * pow(total.z, 2))*x   (2 * total.x * total.y - 2 * total.w * total.z)*y   (2 * total.x * total.z   2 * total.w * total.y)*z;
    vec[1] = (2 * total.x * total.y   2 * total.w * total.z)*x   (1 - 2 * pow(total.x, 2) - 2 * pow(total.z, 2))*y   (2 * total.y * total.z   2 * total.w * total.x)*z;
    vec[2] = (2 * total.x * total.z - 2 * total.w * total.y)*x   (2 * total.y * total.z - 2 * total.w * total.x)*y   (1 - 2 * pow(total.x, 2) - 2 * pow(total.y, 2))*z;

    return vec;
}

That's pretty much it (I also have a normalize function to deal with floating point errors), I initialize all objects quaternion to: w = 1, x = 0, y = 0, z = 0. I rotate a quaternion using an expression like this:

obj.rotation = rotate(angle, x-axis, y-axis, z-axis) * obj.rotation

where obj.rotation is the objects total quaternion rotation value.

I appreciate any help I can get on this issue, if anyone knows what's wrong or has also experienced this issue before. Thanks

EDIT: multiplying total by these quaternions output the expected rotation:

rotate(angle,1,0,0)
rotate(angle,0,1,0)
rotate(angle,0,0,1)
rotate(angle,0,1,1)

However, any rotations such as these make the model stretch oddly:

rotate(angle,1,1,0)
rotate(angle,1,0,1)

EDIT2: here is the normalize function I use to normalize the quaternions:

Quaternion normalize(Quaternion n, double tolerance) {
    // adds all squares of quaternion values, if normalized, total will be 1:
    double total = pow(n.w, 2)   pow(n.x, 2)   pow(n.y, 2)   pow(n.z, 2);
    if (total > 1   tolerance || total < 1 - tolerance) {
        // normalizes value of quaternion if it exceeds a certain tolerance value:
        n.w /= (float) sqrt(total);
        n.x /= (float) sqrt(total);
        n.y /= (float) sqrt(total);
        n.z /= (float) sqrt(total);
    }
    return n;
}

CodePudding user response：

As pointed out in comments, you are not initializing your quaternions correctly.

The following quaternions are not rotations:

rotate(angle,0,1,1)
rotate(angle,1,1,0)
rotate(angle,1,0,1)

The reason is the axis is not normalized e.g., the vector (0,1,1) is not normalized. Also make sure your angles are in radians.

CodePudding user response：

To implement two rotations in sequence you need the quaternion product of the two elementary rotations. Each elementary rotation is specified by an axis and an angle. But in your code you did not make sure you have a unit vector (direction vector) for the axis.

Do the following modification

Quaternion rotate(float w, float x, float y, float z) {
    Quaternion n;
    float f = 1/sqrtf(x*x y*y z*z)

    n.w = cosf(w/2);
    n.x = f * x * sinf(w/2);
    n.y = f * y * sinf(w/2);
    n.z = f * z * sinf(w/2); 

    return n;
}

and then use it as follows

Quaternion n = rotate(angle1,1,0,0) * rotate(angle2,0,1,0)

for the combined rotation of angle1 about the x-axis, and angle2 about the y-axis.