This page was generated from doc/rotation/catmull-rom-uniform.ipynb. Interactive online version: .

Uniform Catmull–Rom-Like Quaternion Splines

see notebook about De Casteljau’s algorithm (with Slerp)

notebook about Euclidean Catmull–Rom splines

\begin{equation*} \boldsymbol{\dot{x}}_i = \frac{\boldsymbol{x}_{i+1} - \boldsymbol{x}_{i-1}}{2} = \frac{(\boldsymbol{x}_i - \boldsymbol{x}_{i-1}) + (\boldsymbol{x}_{i+1} - \boldsymbol{x}_i) }{2} = \frac{\boldsymbol{x}_i - \boldsymbol{x}_{i-1}}{2} + \frac{\boldsymbol{x}_{i+1} - \boldsymbol{x}_i}{2} \end{equation*}

section about cubic Euclidean Bézier splines: division by 3

\begin{equation*} q_{i,\text{tangent}} \overset{?}{=} \begin{cases} \left(q_{i+1} {q_{i-1}}^{-1}\right)^{\frac{1}{2} \frac{1}{3}} \\ \left(q_\text{out} q_\text{in} \right)^{\frac{1}{2} \frac{1}{3}} \\ \left( {q_\text{out}}^{\frac{1}{2}} {q_\text{in}}^{\frac{1}{2}} \right)^\frac{1}{3} \\ {q_\text{out}}^{\frac{1}{2} \frac{1}{3}} {q_\text{in}}^{\frac{1}{2} \frac{1}{3}} \\ \left( \left( {q_\text{out}}^{\frac{1}{2} \frac{1}{3}} {q_\text{in}}^{-1 \frac{1}{2} \frac{1}{3}} \right)^{\frac{1}{2}} {q_\text{in}}^{\frac{1}{2} \frac{1}{3}} \right)^2 = \left( \left( {q_\text{in}}^{\frac{1}{2} \frac{1}{3}} {q_\text{out}}^{-1 \frac{1}{2} \frac{1}{3}} \right)^{\frac{1}{2}} {q_\text{out}}^{\frac{1}{2} \frac{1}{3}} \right)^2 \\ \exp\left( \frac{1}{2} \frac{1}{3} (\ln(q_\text{in}) + \ln(q_\text{out})) \right) \end{cases} \end{equation*}

where \(q_\text{in} = q_i {q_{i-1}}^{-1}\) and \(q_\text{out} = q_{i+1} {q_i}^{-1}\).

The first four options are quite certainly wrong (although not by much), so we are only looking at the last two here:

def tangent1(q_1, q0, q1):
    q_in = (q0 * q_1.inverse())**(1 / (2 * 3))
    q_out = (q1 * q0.inverse())**(1 / (2 * 3))
    return ((q_out * q_in.inverse())**(1 / 2) * q_in)**2

def tangent2(q_1, q0, q1):
    q_in = q0 * q_1.inverse()
    q_out = q1 * q0.inverse()
    return UnitQuaternion.exp_map((q_in.log_map() + q_out.log_map()) / (2 * 3))

import numpy as np

from splines.quaternion import DeCasteljau, UnitQuaternion, canonicalized

helper.py

from helper import angles2quat, animate_rotations, display_animation

def create_closed_curve(rotations, tangent_func):
    rotations = list(canonicalized(rotations + rotations[:2]))
    control_points = []
    for q_1, q0, q1 in zip(rotations, rotations[1:], rotations[2:]):
        q_tangent = tangent_func(q_1, q0, q1)
        control_points.extend([q_tangent.inverse() * q0, q0, q0, q_tangent * q0])
    control_points = control_points[-2:] + control_points[:-2]
    segments = list(zip(*[iter(control_points)] * 4))
    return DeCasteljau(segments)

rotations = [
    angles2quat(0, 0, 180),
    angles2quat(0, 45, 90),
    angles2quat(90, 45, 0),
    angles2quat(90, 90, -90),
    angles2quat(180, 0, -180),
    angles2quat(-90, -45, 180),
]

s1 = create_closed_curve(rotations, tangent1)
s2 = create_closed_curve(rotations, tangent2)

times = np.linspace(0, len(rotations), 200, endpoint=False)

ani = animate_rotations({
    '1': s1.evaluate(times),
    '2': s2.evaluate(times),
}, figsize=(4, 2))

display_animation(ani, default_mode='loop')

Once Loop Reflect

The results are very similar, but not quite identical:

max(max(map(abs, q.xyzw)) for q in (s1.evaluate(times) - s2.evaluate(times)))

0.0002228419691587824