GLM 4.5 Air

Quaternion Interpolation: Complete SLERP Guide

Learn how to properly interpolate quaternions using SLERP for accurate 3D rotations. Complete guide with code examples and best practices for constant angular velocity interpolation.

#3d-rotations#computer-graphics#game-development#interpolation#programming-tutorial#quaternion-math#slerp
Question

How to calculate an intermediate quaternion between two quaternions based on a ratio value?

I’m working with quaternions that represent rotations, and I need to find an intermediate quaternion between two given quaternions (q1 and q2) based on a ratio value between 0 and 1.

What I’ve Tried

Approach 1: Linear Interpolation

I first tried a naive linear interpolation approach:

cpp
q3 = {
    q1[0] + (q2[0] - q1[0])*ratio,
    q1[1] + (q2[1] - q1[1])*ratio,
    q1[2] + (q2[2] - q1[2])*ratio,
    q1[3] + (q2[3] - q1[3])*ratio,
};

Approach 2: SLERP

Then I tried using a SLERP (Spherical Linear Interpolation) function:

cpp
vec4 quaternionSlerp(vec4 q1, vec4 q2, float t)
{
    // Based on http://www.euclideanspace.com/maths/algebra/realNormedAlgebra/quaternions/slerp/index.htm
    float cosHalfTheta = q1[3]*q2[3] + q1[0]*q2[0] + q1[1]*q2[1] + q1[2]*q2[2];

    // if q1=q2 or q1=-q2 then theta = 0 and we can return q1
    if (std::fabs(cosHalfTheta) >= 1.0F) {
        return q1;
    }

    float halfTheta = std::acos(cosHalfTheta);
    float sinHalfTheta = std::sqrt(1.0F - cosHalfTheta*cosHalfTheta);

    // If theta = 180 degrees then result is not fully defined
    // We could rotate around any axis normal to q1 or q2
    if (std::fabs(sinHalfTheta) < 1e-4F) {
        return {
            q1[0] * 0.5F + q2[0] * 0.5F,
            q1[1] * 0.5F + q2[1] * 0.5F,
            q1[2] * 0.5F + q2[2] * 0.5F,
            q1[3] * 0.5F + q2[3] * 0.5F
        };
    }
    else {
        // Default quaternion calculation
        float ratioA = std::sin((1.0F - t) * halfTheta) / sinHalfTheta;
        float ratioB = std::sin(t * halfTheta) / sinHalfTheta;

        return {
            q1[0] * ratioA + q2[0] * ratioB,
            q1[1] * ratioA + q2[1] * ratioB,
            q1[2] * ratioA + q2[2] * ratioB,
            q1[3] * ratioA + q2[3] * ratioB
        };
    }
}

Approach 3: Quaternion Multiplication

After some research, I learned that quaternion addition is actually multiplication, and subtraction involves multiplying by the conjugate. I tried this approach:

cpp
vec4 multiplyQuaternions(vec4 q1, vec4 q2)
{
    float w1 = q1[3], x1 = q1[0], y1 = q1[1], z1 = q1[2];
    float w2 = q2[3], x2 = q2[0], y2 = q2[1], z2 = q2[2];

    return {
        w1 * x2 + x1 * w2 + y1 * z2 - z1 * y2,  // X
        w1 * y2 - x1 * z2 + y1 * w2 + z1 * x2,  // Y
        w1 * z2 + x1 * y2 - y1 * x2 + z1 * w2,  // Z
        w1 * w2 - x1 * x2 - y1 * y2 - z1 * z2   // W
    };
}

vec4 intermediateQuaternion(vec4 q1, vec4 q2, float ratio)
{
    float cosHalfTheta = q1[3]*q2[3] + q1[0]*q2[0] + q1[1]*q2[1] + q1[2]*q2[2];
    float halfTheta = std::acos(cosHalfTheta);
    float sinHalfTheta = std::sqrt(1.0F - cosHalfTheta*cosHalfTheta);

    ratio = std::sin(ratio * halfTheta) / sinHalfTheta;

    vec4 qDifference = multiplyQuaternions(q1, { -q2[0], -q2[1], -q2[2], q2[3] });
    return multiplyQuaternions(q1, { qDifference[0]*ratio, qDifference[1]*ratio, qDifference[2]*ratio, qDifference[3]*ratio });
}

None of these approaches have produced the desired result. I need to find a proper way to calculate an intermediate quaternion between two quaternions based on a ratio value.

Context: I’m using vec4 as an alias for std::array<float, 4> where index 0 is x, 1 is y, 2 is z, and 3 is w. The quaternions represent rotations. I’m not rendering anything, just need the mathematical interpolation.

GLM 4.5 Air

How to Calculate an Intermediate Quaternion Between Two Quaternions Based on a Ratio Value

The proper way to calculate an intermediate quaternion between two quaternions is to use Spherical Linear Interpolation (SLERP), which maintains constant angular velocity and follows the shortest path on the quaternion unit sphere. Your SLERP implementation is actually the correct approach, but there might be some issues with how you’re using it or with your quaternion representation.

Contents

Understanding Quaternion Interpolation

Quaternions represent rotations in 3D space and exist on the surface of a 4D unit sphere. When interpolating between two rotations, we need to find the shortest path on this sphere, which is a great circle arc.

Unlike simple linear interpolation (LERP), which would take a straight-line path through 4D space and result in non-uniform rotation speed, SLERP follows the surface of the sphere and maintains constant angular velocity.

Key Insight: Quaternions q and -q represent the same rotation, so when interpolating, we should always choose the shortest path between them. If the dot product of two quaternions is negative, we should negate one of them before interpolation.

Why SLERP is the Correct Approach

Your SLERP implementation follows the standard mathematical formula for spherical interpolation. Let’s break down why it works:

  1. Dot Product Calculation: cosHalfTheta = q1[3]*q2[3] + q1[0]*q2[0] + q1[1]*q2[1] + q1[2]*q2[2]
    This calculates the cosine of the angle between the two quaternions on the 4D unit sphere.

  2. Handling Special Cases: When quaternions are nearly identical or opposite, the interpolation needs special handling to avoid numerical instability.

  3. Spherical Interpolation Formula: The core SLERP formula uses sine and cosine functions to ensure the interpolation follows the great circle arc on the quaternion sphere.

Your SLERP implementation is mathematically correct. If it’s not producing the expected results, the issue might be with:

  • The quaternion format/order
  • Not handling the q/-q case properly
  • Numerical precision issues
  • How you’re applying the ratio value

Correct SLERP Implementation

Here’s a robust SLERP implementation that handles all edge cases:

cpp
#include <cmath>
#include <algorithm>

vec4 quaternionSlerp(vec4 q1, vec4 q2, float t)
{
    // Ensure quaternions are normalized
    float norm1 = std::sqrt(q1[0]*q1[0] + q1[1]*q1[1] + q1[2]*q1[2] + q1[3]*q1[3]);
    float norm2 = std::sqrt(q2[0]*q2[0] + q2[1]*q2[1] + q2[2]*q2[2] + q2[3]*q2[3]);
    
    q1[0] /= norm1; q1[1] /= norm1; q1[2] /= norm1; q1[3] /= norm1;
    q2[0] /= norm2; q2[1] /= norm2; q2[2] /= norm2; q2[3] /= norm2;

    // Calculate cosine of angle between quaternions
    float cosHalfTheta = q1[3]*q2[3] + q1[0]*q2[0] + q1[1]*q2[1] + q1[2]*q2[2];
    
    // If negative, use -q2 to take the shorter path
    if (cosHalfTheta < 0.0f) {
        q2[0] = -q2[0];
        q2[1] = -q2[1];
        q2[2] = -q2[2];
        q2[3] = -q2[3];
        cosHalfTheta = -cosHalfTheta;
    }
    
    // If quaternions are very close, use linear interpolation
    if (cosHalfTheta > 0.9995f) {
        float t_inv = 1.0f - t;
        return {
            q1[0] * t_inv + q2[0] * t,
            q1[1] * t_inv + q2[1] * t,
            q1[2] * t_inv + q2[2] * t,
            q1[3] * t_inv + q2[3] * t
        };
    }
    
    // Calculate interpolation parameters
    float halfTheta = std::acos(cosHalfTheta);
    float sinHalfTheta = std::sqrt(1.0f - cosHalfTheta*cosHalfTheta);
    float ratioA = std::sin((1.0f - t) * halfTheta) / sinHalfTheta;
    float ratioB = std::sin(t * halfTheta) / sinHalfTheta;
    
    // Perform interpolation
    return {
        q1[0] * ratioA + q2[0] * ratioB,
        q1[1] * ratioA + q2[1] * ratioB,
        q1[2] * ratioA + q2[2] * ratioB,
        q1[3] * ratioA + q2[3] * ratioB
    };
}

Common Issues and Solutions

  1. Quaternion Format: Your code uses [x, y, z, w] format. Ensure this is consistent throughout your application.

  2. Not Handling q/-q Case: The most common issue is not checking if the dot product is negative. When it is, you should use -q2 for the shortest path.

  3. Non-normalized Quaternions: If your input quaternions aren’t normalized, the interpolation won’t work correctly.

  4. Numerical Precision: For very small angles, the standard SLERP formula can suffer from numerical instability. The implementation above includes a fallback to linear interpolation when quaternions are very close.

  5. Ratio Value: Ensure your ratio (t) is between 0.0 and 1.0, inclusive.

Alternative: LERP with Normalization

While not as accurate as SLERP for large angle differences, a normalized LERP can be useful for performance-critical applications:

cpp
vec4 quaternionNlerp(vec4 q1, vec4 q2, float t)
{
    // Ensure quaternions are normalized
    float norm1 = std::sqrt(q1[0]*q1[0] + q1[1]*q1[1] + q1[2]*q1[2] + q1[3]*q1[3]);
    float norm2 = std::sqrt(q2[0]*q2[0] + q2[1]*q2[1] + q2[2]*q2[2] + q2[3]*q2[3]);
    
    q1[0] /= norm1; q1[1] /= norm1; q1[2] /= norm1; q1[3] /= norm1;
    q2[0] /= norm2; q2[1] /= norm2; q2[2] /= norm2; q2[3] /= norm2;
    
    // Take the shortest path
    float cosHalfTheta = q1[3]*q2[3] + q1[0]*q2[0] + q1[1]*q2[1] + q1[2]*q2[2];
    if (cosHalfTheta < 0.0f) {
        q2[0] = -q2[0];
        q2[1] = -q2[1];
        q2[2] = -q2[2];
        q2[3] = -q2[3];
    }
    
    // Linear interpolation
    vec4 result = {
        q1[0] + (q2[0] - q1[0]) * t,
        q1[1] + (q2[1] - q1[1]) * t,
        q1[2] + (q2[2] - q1[2]) * t,
        q1[3] + (q2[3] - q1[3]) * t
    };
    
    // Normalize the result
    float norm = std::sqrt(result[0]*result[0] + result[1]*result[1] + 
                          result[2]*result[2] + result[3]*result[3]);
    result[0] /= norm; result[1] /= norm; result[2] /= norm; result[3] /= norm;
    
    return result;
}

Edge Cases and Special Considerations

  1. Identical Quaternions: When q1 and q2 represent the same rotation (or opposite rotations), the interpolation should return one of them consistently.

  2. 180-degree Difference: When quaternions are exactly opposite (180-degree difference), the interpolation axis isn’t uniquely defined. The implementation above handles this by using linear interpolation in this special case.

  3. Very Small Ratios: For very small ratio values, ensure numerical stability by checking against a small epsilon value.

  4. Quaternion Normalization: Always ensure your input quaternions are normalized before interpolation.

Performance Considerations

Method Accuracy Performance Best Use Case
SLERP Highest Moderate Animation, camera movements
Normalized LERP Good Highest Many small interpolations, performance-critical code
Linear Interpolation Poor Highest Quick approximations, debugging

SLERP is computationally more expensive than LERP but provides the most accurate interpolation. For most applications, especially those involving animation or camera movements, SLERP is the preferred method.

Conclusion and Best Practices

To calculate an intermediate quaternion between two quaternions based on a ratio value:

  1. Use SLERP for the most accurate interpolation with constant angular velocity
  2. Handle the q/-q case by checking if the dot product is negative
  3. Normalize your quaternions before interpolation
  4. Consider edge cases like identical quaternions or 180-degree differences
  5. For performance-critical applications, consider normalized LERP as a faster alternative

Your SLERP implementation was on the right track, but the complete solution needs to handle the quaternion negation case and ensure proper normalization. With these considerations, you should get the correct intermediate quaternion for any ratio value between 0 and 1.

What's the difference between SLERP and LERP for quaternion interpolation?How to handle edge cases when interpolating quaternions?When should I use normalized LERP instead of SLERP for performance?How to interpolate multiple quaternions in a sequence for smooth animation?What are the common mistakes when implementing quaternion interpolation?How to visualize quaternion interpolation paths in 3D space?
Ask GLM 4.5 Air