This interpolation method uses a cubic curve between each pair of control points.
It's much smoother than linear interpolation.
Desirable properties:
- curve should interpolate the control points;
- curve should be $C^1$ (i.e. continuous in its first derivative).
How to make each cubic piece?
Define $Q_i(u) =$ curve between control points $q_i$ and $q_{i+1}$.
Note: On $Q_i(u)$, $t$ varies in $[t_i, t_{i+1}]$, but $u$ varies in $[0,1]$. So $$u = {t - t_i \over t_{i+1} - t_i}$$
There are four constraints on $Q_i(u)$:
First, $Q_i(u)$ must interpolate the control points:
- $Q_i(0) = q_i$
- $Q_i(1) = q_{i+1}$
Second, $Q_i(u)$ must be $C^1$ at the control points:
- ${d \over du} Q_{i-1}(1) = {d \over du} Q_i(0)$
- ${d \over du} Q_i(1) = {d \over du} Q_{i_{i+1}}(0)$
From now on:
We will write $Q'(u)$ for ${d \over du} Q(u)$.
Note on degrees of freedom
Given control points $q_1 \ldots q_n$, we have $n-1$ segments, $Q_1 \ldots Q_{n-1}$ between them.
Each segment, $Q_i$, is a cubic: $$Q_i(u) = a_i\ u^3 + b_i\ u^2 + c_i\ u + d_i$$
for as-yet unknown coefficients $a_i, b_i, c_i, d_i$.
Therefore, the whole curve has $4 (n-1)$ unknowns, or "degrees of freedom".
How many constraints do we have?
There are $2(n-1)$ constraints for $$Q_i(0) = q_i$$ $$Q_i(1) = q_{i+1}$$
Note that each $C^1$ constraint is shared by two segments. If we assume a "virtual segment" before $q_1$ and after $q_n$ for a total of $n+1$ segments, then we have $n$ continuity constraints, one at each control point: $$Q'_{i-1}(1) = Q'_i(0)$$ $$Q'_i(1) = Q'_{i+1}(0)$$
(The "virtual segments" are used to control the behaviour of the curve at its ends.)
So we have $4(n-1)$ degrees of freedom, but only $2(n-1)+n$ constraints, and $n-2$ degrees of freedom left after taking into account the constraints: $$\begin{array}{rl} 4(n-1) & \textrm{DOF} \\ 2(n-1) + n & \textrm{constraints} \\ \hline n-2 & \textrm{DOF left over} \end{array}$$
We could take advantage of this by instead using quadratic curves in each section, instead of cubic curves, and dropping the two $C^1$ constraints on the first and last points. Then we would have $$\begin{array}{rl} 3(n-1) & \textrm{DOF in quadratic curves} \\ 2(n-1) + (n-2) & \textrm{constraints} \\ \hline 1 & \textrm{DOF left over} \end{array}$$
But the quadratic curves have the big disadvantage that a change at one part of the curve will propagate (through the constraints) to affect all parts of the curve.
The cubic curves, on the other hand, do not have this problem. They have the property of "local control" in which a change in one control point does not affect the whole curve.
$\mathbf{C^1}$ continuity
The key idea with Catmull-Rom curves is that the derivatives at end of one curve, $Q'_{i-1}(1)$, and at the beginning of the next curve, $Q'_i(0)$ are set to a particular value.
This is a stronger constraint than saying $Q'_{i-1}(1) = Q'_i(0)$.
It's saying $$Q'_{i-1}(1) = Q'_i(0) = v$$
for some value, $v$.
With Catmull-Rom curves, that derivative is chosen to be the first finite difference of the surrounding control points: $$\displaystyle Q'_{i-1}(1) = Q'_i(0) = {q_{i+1} - q_{i-1} \over t_{i+1} - t_{i-1}}$$
This is the slope between the surrounding points, as shown below:
So the contstraints on $Q_i(u)$ are:
- $Q_i(0) = q_i$
- $Q_i(1) = q_{i+1}$
- $\displaystyle Q'_i(0) = {q_{i+1} - q_{i-1} \over t_{i+1} - t_{i-1}}$
- $\displaystyle Q'_i(1) = {q_{i+2} - q_i \over t_{i+2} - t_i}$
This ensures that $Q_i(u)$ is only affected by the control points $q_{i-1}, q_i, q_{i+1}, q_{i+2}$. So distant parts of the overall curve have no effect on the shape of $Q_i(u)$. This property is called local control, as only local control points affect the shape.
Change-of-basis matrix
Since each $Q_i(u)$ is $$Q_i(u) = a_i\ u^3 + b_i\ u^2 + c_i\ u + d_i$$
we have four constraints (above) and four unknowns ($a_i, b_i, c_i, d_i$).
So we can solve for these unknowns and, in doing so, determine the Catmull-Rom change-of-basis matrix.
We need four equations relating the $q$'s and the coefficients. Let's drop the $i$ subscript in $Q$ to get $$Q(u) = a u^3 + b u^2 + c u + d$$ $$Q'(u) = 3 a u^2 + 2 b u + c$$
Then we get the following equations:
Interpolation equations
Continuity equations
To simply things, we will assume that the spacing between each pair of adjacent keyframes is one unit: $$\forall i\ \ t_{i+1} - t_i = 1$$
Then the continuity equations are simplified to $$q_{i-1} = a+b-c+d$$ $$q_{i+2} = 6a+4b+2c+d$$
The four equations and four unknowns can now be written in matrix form: $$\left[\begin{array}{c} q_{i-1} \\ q_i \\ q_{i+1} \\ q_{i+2} \end{array}\right] = \underbrace{\left[\begin{array}{rrrr} 1 & 1 & -1 & 1 \\ 0 & 0 & 0 & 1 \\ 1 & 1 & 1 & 1 \\ 6 & 4 & 2 & 1 \end{array}\right]}_{M^{-1}} \left[\begin{array}{c} a \\ b \\ c \\ d \end{array}\right]$$
The matrix, $M^{-1}$, above is the inverse of the change-of-basis matrix. Inverting it, we get $$\left[\begin{array}{c} a \\ b \\ c \\ d \end{array}\right] = \underbrace{\frac12 \left[\begin{array}{rrrr} -1 & 3 & -3 & 1 \\ 2 & -5 & 4 & -1 \\ -1 & 0 & 1 & 0 \\ 0 & 2 & 0 & 0 \end{array}\right]}_{M}\ \ \left[\begin{array}{c} q_{i-1} \\ q_i \\ q_{i+1} \\ q_{i+2} \end{array}\right]$$
Note that the Catmull-Rom change-of-basis matrix, $M$, is constant. (This requires the assumption of even spacing in the keyframes.)
Summary
The original problem was to interpolate a control value, $q(t)$, between the control points at keyframes.
If $t \in [t_i,t_{i+1}]$, we interpolate with the cubic curve $Q_i(u) = a_i\ u^3 + b_i\ u^2 + c_i\ u + d_i$ for $u = {t - t_i \over t_{i+1} - t_i}$.
The coefficients of $Q_i(u)$ are computed from the control points as $$\left[\begin{array}{c} a_i \\ b_i \\ c_i \\ d_i \end{array}\right] = M\ \left[\begin{array}{c} q_{i-1} \\ q_i \\ q_{i+1} \\ q_{i+2} \end{array}\right]$$
This is done separately for each curve: $Q_1(u), Q_2(u), \ldots, Q_{n-1}(u)$.
Note that the use of $q_{i-1}$ and $q_{i+2}$ above requires that we have a point $q_0$ before the first point, $q_1$, and that we have a point $q_{n+1}$ after the last point, $q_n$.
These "phantom points" are usually chosen with $q_0 = q_1$ and $q_n = q_{n+1}$.