This notebook is part of https://github.com/AudioSceneDescriptionFormat/splines, see also http://splines.readthedocs.io/.

# Derivation of Non-Uniform Catmull--Rom Splines

Recursive algorithm developed by
<cite data-cite="barry1988recursive">Barry and Goldman (1988)</cite>,
according to
<cite data-cite="yuksel2011parameterization">Yuksel et al. (2011)</cite>, figure 3.


```python
import sympy as sp
sp.init_printing()
```


```python
from utility import NamedExpression, NamedMatrix
```


```python
x_1, x0, x1, x2 = sp.symbols('xbm_-1 xbm:3')
```


```python
t, t_1, t0, t1, t2 = sp.symbols('t t_-1 t:3')
```


```python
p_10 = NamedExpression('pbm_-1,0', x_1 * (t0 - t) / (t0 - t_1) + x0 * (t - t_1) / (t0 - t_1))
p_10
```


```python
p01 = NamedExpression('pbm_0,1', x0 * (t1 - t) / (t1 - t0) + x1 * (t - t0) / (t1 - t0))
p01
```


```python
p12 = NamedExpression('pbm_1,2', x1 * (t2 - t) / (t2 - t1) + x2 * (t - t1) / (t2 - t1))
p12
```


```python
p_101 = NamedExpression('pbm_-1,0,1', p_10.name * (t1 - t) / (t1 - t_1) + p01.name * (t - t_1) / (t1 - t_1))
p_101
```


```python
p012 = NamedExpression('pbm_0,1,2', p01.name * (t2 - t) / (t2 - t0) + p12.name * (t - t0) / (t2 - t0))
p012
```


```python
p = NamedExpression('pbm', p_101.name * (t1 - t) / (t1 - t0) + p012.name * (t - t0) / (t1 - t0))
p
```


```python
p = p.subs([p_101, p012]).subs([p_10, p01, p12])
p
```


```python
p_normalized = p.expr.subs(t, t * (t1 - t0) + t0)
```


```python
M_CR = NamedMatrix(
    r'{M_\text{CR}}',
    sp.Matrix([[c.expand().coeff(x).factor() for x in (x_1, x0, x1, x2)]
               for c in p_normalized.as_poly(t).all_coeffs()]))
```


```python
deltas = [
    (t_1, -sp.Symbol('Delta_-1')),
    (t0, 0),
    (t1, sp.Symbol('Delta0')),
    (t2, sp.Symbol('Delta0') + sp.Symbol('Delta1'))
]
```


```python
M_CR.simplify().subs(deltas).factor()
```


```python
uniform = [
    (sp.Symbol('Delta_-1'), 1),
    (sp.Symbol('Delta0') , 1),
    (sp.Symbol('Delta1') , 1),
]
```


```python
M_CR.subs(deltas).subs(uniform).pull_out(sp.S.Half).expr
```


```python
velocity = p.expr.diff(t)
```


```python
velocity.subs(t, t0).subs(deltas).factor()
```


```python
velocity.subs(t, t1).subs(deltas).factor()
```

in general:

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

You might encounter another way to write the equation for $\boldsymbol{\dot{x}}_0$
(e.g. at https://stackoverflow.com/a/23980479/):


```python
(x0 - x_1) / (t0 - t_1) - (x1 - x_1) / (t1 - t_1) + (x1 - x0) / (t1 - t0)
```

... but this is equivalent to the equation shown above:


```python
_.subs(deltas).factor()
```

Yet another way to skin this cat -- sometimes referred to as Bessel--Overhauser -- is to define the velocity of the left and right chords:


```python
v_left = (x0 - x_1) / (t0 - t_1)
v_right = (x1 - x0) / (t1 - t0)
```

... and then combine them in this way:


```python
((t1 - t0) * v_left + (t0 - t_1) * v_right) / (t1 - t_1)
```

Again, that's the same as we had above:


```python
_.subs(deltas).factor()
```

c9ca7563c421f0e3d40ffc577b12b7d0150afb4f
ipynb
Jupyter Notebook
doc/catmull-rom-non-uniform.ipynb
mgeier/splines
f54b09479d98bf13f00a183fd9d664b5783e3864
MIT
doc/catmull-rom-non-uniform.ipynb
mgeier/splines
f54b09479d98bf13f00a183fd9d664b5783e3864
MIT
doc/catmull-rom-non-uniform.ipynb
mgeier/splines
f54b09479d98bf13f00a183fd9d664b5783e3864
MIT
Qwen/Qwen-72B
1. YES
2. YES