这篇论文 [1] 比较基础,在很多与 B 样条有关的论文中都能找到对它的引用。
B 样条曲线是参数方程,使用参数 \(t\) 的多项式表达连续曲线。
多项式的计算是递归计算的,k 阶 B 样条的多项式的参数是 k 阶的,需要从 k-1 阶参数计算。可以看回形针视频 [2] 有一个简单的认识。
注意一下 spline 的 degree 与 order。对于多项式而言 degree 与 order 是相等的,是同义词。而对于 spline,order = degree + 1。[3]
B-splines of order \(p+1\) are connected piece-wise polynomial functions of degree \(p\).
所以在 [2] 中提到的 “4 个控制点 3 次 B 样条”的 order 为 4,degree 为 3。B 样条上的任意一点受到 order 个控制点的影响,可以看做是 3 次多项式。
作者在这篇论文提出了一种 B 样条曲线使用 Matrix 表达的方式,这种方式在计算上有优势。
回想到斐波那契数列,也是递归计算,其中有一种计算方式也是用矩阵表达[4]。甚至可以从这个矩阵的表达中得出斐波那契数列的通项公式。
以下正文中的公式引用使用两种形式:本文的公式,如(1);原文 [1] 的公式,如(1*)。
以下推导过程应该会很啰嗦,我的目的是我自己能够看得懂,通过此文能够毫不费力地了解到论文中所有没有明说的细节。
1. Toeplitz Matrix 表达多项式
Toeplitz Matrix 的定义在 [1] 中有,2.1 的第一句话。
The Toeplitz matrix is one whose elements on any line parallel to the main diagonal are all equal.
如果 Toeplitz Matrix 是一个 lower triangular matrix 的话,它就能表达一个多项式。我尝试用二次型的方式去解释这个矩阵,但是失败了。反正 [2] 2.1 的第二个矩阵 \(\mathbf{T}\) 可以表示一个多项式。
\[\begin{align}
f(x) = a_0 + a_1 x + a_x x^2 + \dots + a_{n-1} x^{n-1} (a_{n-1} \neq 0)
\end{align}\]
你问我这个 \(\mathbf{T}\) 是多少维的方阵,我也不知道,按需构建。这不重要。重点在 [1] 2.2 对它的应用中。
假设有两个多项式。
\[\begin{align}
g(x) &= c_0 + c_1 x + c_x x^2 + \dots + c_{m-1} x^{m-1} (c_{m-1} \neq 0) \\
q(x) &= d_0 + d_1 x + d_x x^2 + \dots + d_{n-1} x^{n-1} (d_{n-1} \neq 0)
\end{align}\]
\(f(x) = g(x) q(x)\),\(f(x)\) 的最高次应当是 \(m + n - 2\)。
所以 [1] 2.2 中的矩阵具体的维数应当如下。
\[\begin{align}
f(x) &= g(x) q(x) \\
&= \mathbf{X} \begin{bmatrix}
c_0 & \phantom{} & \phantom{} & \phantom{} & \phantom{} & \phantom{} & 0 \\
c_1 & c_0 \\
\vdots & \ddots & \ddots \\
c_{m - 1} & \cdots & \ddots \\
\phantom{} & \ddots & \cdots & \ddots & \ddots \\
\phantom{} & \phantom{} & c_{m-1} & \cdots & \ddots & c_0 \\
0 & \phantom{} & \phantom{} & c_{m-1} & \cdots & c_1 & c_0
\end{bmatrix}_{(m + n - 1) \times (m + n - 1)} \begin{bmatrix}
d_0 \\ d_1 \\ \vdots \\ d_{n - 1} \\ 0 \\ \vdots \\ 0
\end{bmatrix}_{(m + n - 1) \times 1} \\
&= \mathbf{X} \begin{bmatrix}
c_0 & \phantom{} & \phantom{} & 0 \\
c_1 & c_0 \\
\vdots & \ddots & \ddots \\
\vdots & \ddots & \ddots & c_0\\
\vdots & \ddots & \ddots & \vdots\\
c_{m - 1} & c_{m-2} & \cdots & c_{m-n} \\
\phantom{} & \ddots & \cdots & \vdots \\
\phantom{} & \phantom{} & c_{m-1} & c_{m-2} \\
0 & \phantom{} & \phantom{} & c_{m-1}
\end{bmatrix}_{(m + n - 1) \times n} \begin{bmatrix}
d_0 \\ d_1 \\ \vdots \\ d_{n - 1}
\end{bmatrix}_{n \times 1}
\end{align}\]
[5] 给了一个具体的例子,可以看一下。
2. Toeplitz Matrix 表达 B 样条的 basis functions
我对 B 样条不熟,各种乱七八糟的多项式,不清晰。希望看完这个能够清晰一些。
(1*) 给出了 B 样条各阶 basis function 的递归公式。
作者说 By means of basis translation from B-spline to power basis, ...
,所以 basis function 就从以参数 \(t\) 作为变量,转换到以 \(u\) 作为变量。作者说参考的两个参考文献都是书,而且其中的 Practical Guide to Splines 还只能找到扫描版。不容易搜索关键字 power basis
。所以我另外找了一篇参考文献 [6],也是 De Boor 写的。看 [6] 的 8. Condition
,是精度的问题,用 power basis 转换变量的取值空间能够更容易控制精度。注意 Condition
这个词,联想到了 Matrix Condition Number
([7] 2.2.2 与 4.3.2)。Matrix 用于表示一个线性系统,Condition Number 表示 forward error / backward error
,forward error 是 x
的误差,而 backward error 是 f(x)
的误差。希望 condition number 尽可能小,即 f(x)
是观察值,用观察值去估计实际关心的状态值 x
,希望观察值小的抖动,只会造成状态值小的抖动,这样就容许了观察值的观测误差,对测量条件没有这么苛刻。[6] 的 8. Condition
也是这个意思。
\[\begin{align}
B_{j,k}(t) &= {t - t_j \over t_{j + k - 1} - t_j} B_{j,k-1}(t) + {t_{j + k} - t \over t_{j + k} - t_{j + 1}} B_{j + 1, k -1}(t) \\
&= ({t_i - t_j \over t_{j + k - 1} - t_j} + {t - t_i \over t_{i + 1} - t_i}{t_{i + 1} - t_i \over t_{j + k - 1} - t_j}) B_{j,k-1}(t) \\
&\phantom{=} + ({t_{j + k} - t_i \over t_{j + k} - t_{j + 1}} - {t - t_i \over t_{i + 1} - t_i}{t_{i + 1} - t_i \over t_{j + k} - t_{j + 1}}) B_{j + 1, k -1}(t) \\
B_{j,k}(u) &= (d_{0,j} + ud_{1,j}) B_{j,k-1}(u) + (h_{0,j} + uh_{1,j}) B_{j + 1, k -1}(u) \\
u &= {t - t_i \over t_{i+1} - t_i}, u\in [0, 1)
\end{align}\]
上式中出现的 index \(i\) 是什么意思?公式 (3*) 上方的文字做了解释,\(i\) 表示当前关注的 \(t\in[t_i,t_{i+1})\),这个区间内的 B 样条与周围 \(k\) 个 basis functions 以及它们对应的 control points 有关,\(j\) 是这 \(k\) 个部分的 index。
作者给出 \(B_{j,k-1}(u)\) 的定义。
\[\begin{align}
B_{j,k-1}(u) &= \begin{bmatrix}
1 & u & \dots & u^{k-2}
\end{bmatrix} \begin{bmatrix}
N_{0,j}^{k-1} \\ N_{1,j}^{k-1} \\ \vdots \\ N_{k-2,j}^{k-1}
\end{bmatrix} \\
u &= {t - t_i \over t_{i+1} - t_i}, u\in [0, 1)
\end{align}\]
用 Toeplitz Matrix 表示这个迭代过程,\(u\) 从 \(k-2\) 此上升到 \(k-1\) 次。
\[\begin{align}
B_{j,k}(u) &= (d_{0,j} + ud_{1,j}) B_{j,k-1}(u) + (h_{0,j} + uh_{1,j}) B_{j + 1, k -1}(u) \\
&= \begin{bmatrix}
1 & u & \dots & u^{k-1}
\end{bmatrix} \left\{ \begin{bmatrix}
N_{0,j}^{k-1} & 0 \\
N_{1,j}^{k-1} & N_{0,j}^{k-1} \\
\vdots & N_{1,j}^{k-1} \\
N_{k-2,j}^{k-1} & \vdots \\
0 & N_{k-2,j}^{k-1}
\end{bmatrix} \begin{bmatrix}
d_{0,j} \\ d_{1,j}
\end{bmatrix} + \begin{bmatrix}
N_{0,j+1}^{k-1} & 0 \\
N_{1,j+1}^{k-1} & N_{0,j+1}^{k-1} \\
\vdots & N_{1,j+1}^{k-1} \\
N_{k-2,j+1}^{k-1} & \vdots \\
0 & N_{k-2,j+1}^{k-1}
\end{bmatrix} \begin{bmatrix}
h_{0,j} \\ h_{1,j}
\end{bmatrix} \right\}
\end{align}\]
如公式 (4*),在 \(t\in[t_i,t_{i+1})\) 这个区间内,点的坐标可以写作。
\[\begin{align}
\mathbf{c}_{i-k+1} &= \begin{bmatrix}
B_{i-k+1,k}(u) & B_{i-k+2,k}(u) & \cdots & B_{i,k}(u)
\end{bmatrix} \times \begin{bmatrix}
\mathbf{V}_{i-k+1} \\ \mathbf{V}_{i-k+2} \\ \vdots \\ \mathbf{V}_{i}
\end{bmatrix}
\end{align}\]
这 \(k\) 个 basis function 都是 \(k-1\) 次的。扩展开。
\[\begin{align}
&\phantom{=}\begin{bmatrix}
B_{i-k+1,k}(u) & B_{i-k+2,k}(u) & \cdots & B_{i,k}(u)
\end{bmatrix} \\ &= \begin{bmatrix}
1 & u & \dots & u^{k-1}
\end{bmatrix} \mathbf{M}^k(i) \\
\mathbf{M}^k(i) &= \begin{bmatrix}
N_{0,i-k+1}^{k} & N_{0,i-k+2}^{k} & \dots & N_{0,i}^{k} \\ N_{1,i-k+1}^{k} & N_{1,i-k+2}^{k} & \dots & N_{1,i}^{k} \\ \vdots & \vdots & & \vdots \\ N_{k-1,i-k+1}^{k} & N_{k-1,i-k+2}^{k} & \dots & N_{k-1,i}^{k}
\end{bmatrix}
\end{align}\]
所以,B 样条公式的关注重点在 \(\mathbf{M}^k(i)\),\(\mathbf{M}^k(i)\) 是变化的部分,其他部分相对固定。
3. \(\mathbf{M}^k(i)\) 的迭代计算
\(\mathbf{M}^k(i)\) 第 1 列。
\[\begin{align}
\begin{bmatrix}
N_{0,i-k+1}^{k} \\ N_{1,i-k+1}^{k} \\ \vdots \\ \vdots \\ N_{k-1,i-k+1}^{k}
\end{bmatrix} &= \begin{bmatrix}
N_{0,i-k+1}^{k-1} & 0 \\ N_{1,i-k+1}^{k-1} & N_{0,i-k+1}^{k-1} \\ \vdots & N_{1,i-k+1}^{k-1} \\ N_{k-2,i-k+1}^{k-1} & \vdots \\ 0 & N_{k-2,i-k+1}^{k-1}
\end{bmatrix} \begin{bmatrix}
d_{0,i-k+1} \\ d_{1,i-k+1}
\end{bmatrix} + \begin{bmatrix}
N_{0,i-k+2}^{k-1} & 0 \\ N_{1,i-k+2}^{k-1} & N_{0,i-k+2}^{k-1} \\ \vdots & N_{1,i-k+2}^{k-1} \\ N_{k-2,i-k+2}^{k-1} & \vdots \\ 0 & N_{k-2,i-k+2}^{k-1}
\end{bmatrix} \begin{bmatrix}
h_{0,i-k+1} \\ h_{1,i-k+1}
\end{bmatrix}
\end{align}\]
\(\mathbf{M}^k(i)\) 可以每一列分别计算再求和,记得补 \(\mathbf{0}\) 列。
\[\begin{align}
\begin{bmatrix}
N_{0,i-k+1}^{k} & \overbrace{0 \cdots 0}^{k-1} \\ N_{1,i-k+1}^{k} & 0 \cdots 0 \\ \vdots & \vdots \\ \vdots & \vdots \\ N_{k-1,i-k+1}^{k} & 0 \cdots 0
\end{bmatrix} &= \begin{bmatrix}
N_{0,i-k+1}^{k-1} & 0 \\ N_{1,i-k+1}^{k-1} & N_{0,i-k+1}^{k-1} \\ \vdots & N_{1,i-k+1}^{k-1} \\ N_{k-2,i-k+1}^{k-1} & \vdots \\ 0 & N_{k-2,i-k+1}^{k-1}
\end{bmatrix} \begin{bmatrix}
d_{0,i-k+1} & \overbrace{0 \cdots 0}^{k-1} \\ d_{1,i-k+1} & 0 \cdots 0
\end{bmatrix} + \begin{bmatrix}
N_{0,i-k+2}^{k-1} & 0 \\ N_{1,i-k+2}^{k-1} & N_{0,i-k+2}^{k-1} \\ \vdots & N_{1,i-k+2}^{k-1} \\ N_{k-2,i-k+2}^{k-1} & \vdots \\ 0 & N_{k-2,i-k+2}^{k-1}
\end{bmatrix} \begin{bmatrix}
h_{0,i-k+1} & \overbrace{0 \cdots 0}^{k-1} \\ h_{1,i-k+1} & 0 \cdots 0
\end{bmatrix}
\end{align}\]
下面使用了矩阵分块。
\[\begin{align}
\mathbf{M}^k(i) &= \begin{bmatrix}
N_{0,i-k+1}^{k} & N_{0,i-k+2}^{k} & \dots & N_{0,i}^{k} \\ N_{1,i-k+1}^{k} & N_{1,i-k+2}^{k} & \dots & N_{1,i}^{k} \\ \vdots & \vdots & & \vdots \\ N_{k-1,i-k+1}^{k} & N_{k-1,i-k+2}^{k} & \dots & N_{k-1,i}^{k}
\end{bmatrix} \\
&= \begin{bmatrix}
N_{0,i-k+1}^{k} & 0 & \dots & 0 \\ N_{1,i-k+1}^{k} & 0 & \dots & 0 \\ \vdots & \vdots & & \vdots \\ N_{k-1,i-k+1}^{k} & 0 & \dots & 0
\end{bmatrix} + \begin{bmatrix}
0 & N_{0,i-k+2}^{k} & \dots & 0 \\ 0 & N_{1,i-k+2}^{k} & \dots & 0 \\ \vdots & \vdots & & \vdots \\ 0 & N_{k-1,i-k+2}^{k} & \dots & 0
\end{bmatrix} + \cdots + \begin{bmatrix}
0 & 0 & \dots & N_{0,i}^{k} \\ 0 & 0 & \dots & N_{1,i}^{k} \\ \vdots & \vdots & & \vdots \\ 0 & 0 & \dots & N_{k-1,i}^{k}
\end{bmatrix} \\
&= \begin{bmatrix}
N_{0,i-k+1}^{k-1} & 0 \\ N_{1,i-k+1}^{k-1} & N_{0,i-k+1}^{k-1} \\ \vdots & N_{1,i-k+1}^{k-1} \\ N_{k-2,i-k+1}^{k-1} & \vdots \\ 0 & N_{k-2,i-k+1}^{k-1}
\end{bmatrix} \begin{bmatrix}
d_{0,i-k+1} & \overbrace{0 \cdots 0}^{k-1} \\ d_{1,i-k+1} & 0 \cdots 0
\end{bmatrix} + \begin{bmatrix}
N_{0,i-k+2}^{k-1} & 0 \\ N_{1,i-k+2}^{k-1} & N_{0,i-k+2}^{k-1} \\ \vdots & N_{1,i-k+2}^{k-1} \\ N_{k-2,i-k+2}^{k-1} & \vdots \\ 0 & N_{k-2,i-k+2}^{k-1}
\end{bmatrix} \begin{bmatrix}
0 & d_{0,i-k+2} & \overbrace{0 \cdots 0}^{k-2} \\ 0 & d_{1,i-k+2} & 0 \cdots 0
\end{bmatrix} \\
&+ \cdots + \begin{bmatrix}
N_{0,i}^{k-1} & 0 \\ N_{1,i}^{k-1} & N_{0,i}^{k-1} \\ \vdots & N_{1,i}^{k-1} \\ N_{k-2,i}^{k-1} & \vdots \\ 0 & N_{k-2,i}^{k-1}
\end{bmatrix} \begin{bmatrix}
\overbrace{0 \cdots 0}^{k-1} & d_{0,i} \\ 0 \cdots 0 & d_{1,i}
\end{bmatrix} \\
&+ \begin{bmatrix}
N_{0,i-k+2}^{k-1} & 0 \\ N_{1,i-k+2}^{k-1} & N_{0,i-k+2}^{k-1} \\ \vdots & N_{1,i-k+2}^{k-1} \\ N_{k-2,i-k+2}^{k-1} & \vdots \\ 0 & N_{k-2,i-k+2}^{k-1}
\end{bmatrix} \begin{bmatrix}
h_{0,i-k+1} & \overbrace{0 \cdots 0}^{k-1} \\ h_{1,i-k+1} & 0 \cdots 0
\end{bmatrix} + \begin{bmatrix}
N_{0,i-k+3}^{k-1} & 0 \\ N_{1,i-k+3}^{k-1} & N_{0,i-k+3}^{k-1} \\ \vdots & N_{1,i-k+3}^{k-1} \\ N_{k-2,i-k+3}^{k-1} & \vdots \\ 0 & N_{k-2,i-k+3}^{k-1}
\end{bmatrix} \begin{bmatrix}
0 & h_{0,i-k+2} & \overbrace{0 \cdots 0}^{k-2} \\ 0 & h_{1,i-k+2} & 0 \cdots 0
\end{bmatrix} \\
&+ \cdots + \begin{bmatrix}
N_{0,i+1}^{k-1} & 0 \\ N_{1,i+1}^{k-1} & N_{0,i+1}^{k-1} \\ \vdots & N_{1,i+1}^{k-1} \\ N_{k-2,i+1}^{k-1} & \vdots \\ 0 & N_{k-2,i+1}^{k-1}
\end{bmatrix} \begin{bmatrix}
\overbrace{0 \cdots 0}^{k-1} & h_{0,i} \\ 0 \cdots 0 & h_{1,i}
\end{bmatrix} \\
&= \begin{bmatrix}
N_{0,i-k+1}^{k-1} & 0 \\ N_{1,i-k+1}^{k-1} & N_{0,i-k+1}^{k-1} \\ \vdots & N_{1,i-k+1}^{k-1} \\ N_{k-2,i-k+1}^{k-1} & \vdots \\ 0 & N_{k-2,i-k+1}^{k-1}
\end{bmatrix} \begin{bmatrix}
d_{0,i-k+1} & \overbrace{0 \cdots 0}^{k-1} \\ d_{1,i-k+1} & 0 \cdots 0
\end{bmatrix} \\
&+ \begin{bmatrix}
N_{0,i-k+2}^{k-1} & 0 \\ N_{1,i-k+2}^{k-1} & N_{0,i-k+2}^{k-1} \\ \vdots & N_{1,i-k+2}^{k-1} \\ N_{k-2,i-k+2}^{k-1} & \vdots \\ 0 & N_{k-2,i-k+2}^{k-1}
\end{bmatrix} \begin{bmatrix}
h_{0,i-k+1} & d_{0,i-k+2} & \overbrace{0 \cdots 0}^{k-2} \\ h_{1,i-k+1} & d_{1,i-k+2} & 0 \cdots 0
\end{bmatrix} \\
&+ \begin{bmatrix}
N_{0,i-k+3}^{k-1} & 0 \\ N_{1,i-k+3}^{k-1} & N_{0,i-k+3}^{k-1} \\ \vdots & N_{1,i-k+3}^{k-1} \\ N_{k-2,i-k+3}^{k-1} & \vdots \\ 0 & N_{k-2,i-k+3}^{k-1}
\end{bmatrix} \begin{bmatrix}
0 & h_{0,i-k+2} & d_{0,i-k+3} & \overbrace{0 \cdots 0}^{k-3} \\ 0 & h_{1,i-k+2} & d_{1,i-k+3} & 0 \cdots 0
\end{bmatrix} \\
&+ \cdots \\
&+ \begin{bmatrix}
N_{0,i}^{k-1} & 0 \\ N_{1,i}^{k-1} & N_{0,i}^{k-1} \\ \vdots & N_{1,i}^{k-1} \\ N_{k-2,i}^{k-1} & \vdots \\ 0 & N_{k-2,i}^{k-1}
\end{bmatrix} \begin{bmatrix}
\overbrace{0 \cdots 0}^{k-2} & h_{0,i-1} & d_{0,i} \\ 0 \cdots 0 & h_{1,i-1} & d_{1,i}
\end{bmatrix} \\
&+ \begin{bmatrix}
N_{0,i+1}^{k-1} & 0 \\ N_{1,i+1}^{k-1} & N_{0,i+1}^{k-1} \\ \vdots & N_{1,i+1}^{k-1} \\ N_{k-2,i+1}^{k-1} & \vdots \\ 0 & N_{k-2,i+1}^{k-1}
\end{bmatrix} \begin{bmatrix}
\overbrace{0 \cdots 0}^{k-1} & h_{0,i} \\ 0 \cdots 0 & h_{1,i}
\end{bmatrix} \\
\end{align}\]
\[\begin{align}
&=^{(**)} \begin{bmatrix}
N_{0,i-k+2}^{k-1} \\ N_{1,i-k+2}^{k-1} \\ \vdots \\ N_{k-2,i-k+2}^{k-1} \\ 0
\end{bmatrix} \begin{bmatrix}
h_{0,i-k+1} & d_{0,i-k+2} & \overbrace{0 \cdots 0}^{k-2}
\end{bmatrix} + \begin{bmatrix}
0 \\ N_{0,i-k+2}^{k-1} \\ N_{1,i-k+2}^{k-1} \\ \vdots \\ N_{k-2,i-k+2}^{k-1}
\end{bmatrix} \begin{bmatrix}
h_{1,i-k+1} & d_{1,i-k+2} & \overbrace{0 \cdots 0}^{k-2}
\end{bmatrix} \\
&+ \begin{bmatrix}
N_{0,i-k+3}^{k-1} \\ N_{1,i-k+3}^{k-1} \\ \vdots \\ N_{k-2,i-k+3}^{k-1} \\ 0
\end{bmatrix} \begin{bmatrix}
0 & h_{0,i-k+2} & d_{0,i-k+3} & \overbrace{0 \cdots 0}^{k-3}
\end{bmatrix} + \begin{bmatrix}
0 \\ N_{0,i-k+3}^{k-1} \\ N_{1,i-k+3}^{k-1} \\ \vdots \\ N_{k-2,i-k+3}^{k-1}
\end{bmatrix} \begin{bmatrix}
0 & h_{1,i-k+2} & d_{1,i-k+3} & \overbrace{0 \cdots 0}^{k-3}
\end{bmatrix} \\
&+ \cdots \\
&+ \begin{bmatrix}
N_{0,i}^{k-1} \\ N_{1,i}^{k-1} \\ \vdots \\ N_{k-2,i}^{k-1} \\ 0
\end{bmatrix} \begin{bmatrix}
\overbrace{0 \cdots 0}^{k-2} & h_{0,i-1} & d_{0,i}
\end{bmatrix} + \begin{bmatrix}
0 \\ N_{0,i}^{k-1} \\ N_{1,i}^{k-1} \\ \vdots \\ N_{k-2,i}^{k-1}
\end{bmatrix} \begin{bmatrix}
\overbrace{0 \cdots 0}^{k-2} & h_{1,i-1} & d_{1,i}
\end{bmatrix} \\
&= \begin{bmatrix}
N_{0,i-k+2}^{k-1} & N_{0,i-k+3}^{k-1} & \cdots & N_{0,i}^{k-1} \\ N_{1,i-k+2}^{k-1} & N_{1,i-k+3}^{k-1} & \cdots & N_{1,i}^{k-1} \\ \vdots & \vdots & & \vdots \\ N_{k-2,i-k+2}^{k-1} & N_{k-2,i-k+3}^{k-1} & \cdots & N_{k-2,i}^{k-1} \\ 0 & 0 & \cdots & 0
\end{bmatrix} \begin{bmatrix}
h_{0,i-k+1} & d_{0,i-k+2} & & & 0 \\
& h_{0,i-k+2} & d_{0,i-k+3} & & \\
& & \ddots & \ddots & \\
0 & & & h_{0,i-1} & d_{0,i} \\
\end{bmatrix} \\
&+ \begin{bmatrix}
0 & 0 & \cdots & 0 \\ N_{0,i-k+2}^{k-1} & N_{0,i-k+3}^{k-1} & \cdots & N_{0,i}^{k-1} \\ N_{1,i-k+2}^{k-1} & N_{1,i-k+3}^{k-1} & \cdots & N_{1,i}^{k-1} \\ \vdots & \vdots & & \vdots \\ N_{k-2,i-k+2}^{k-1} & N_{k-2,i-k+3}^{k-1} & \cdots & N_{k-2,i}^{k-1}
\end{bmatrix} \begin{bmatrix}
h_{1,i-k+1} & d_{1,i-k+2} & & & 0 \\
& h_{1,i-k+2} & d_{1,i-k+3} & & \\
& & \ddots & \ddots & \\
0 & & & h_{1,i-1} & d_{1,i} \\
\end{bmatrix} \\
&= \begin{bmatrix}
\mathbf{M}^{k-1}(i) \\ \mathbf{0}^T
\end{bmatrix} \begin{bmatrix}
h_{0,i-k+1} & d_{0,i-k+2} & & & 0 \\
& h_{0,i-k+2} & d_{0,i-k+3} & & \\
& & \ddots & \ddots & \\
0 & & & h_{0,i-1} & d_{0,i} \\
\end{bmatrix} \\
&+ \begin{bmatrix}
\mathbf{0}^T \\ \mathbf{M}^{k-1}(i)
\end{bmatrix} \begin{bmatrix}
h_{1,i-k+1} & d_{1,i-k+2} & & & 0 \\
& h_{1,i-k+2} & d_{1,i-k+3} & & \\
& & \ddots & \ddots & \\
0 & & & h_{1,i-1} & d_{1,i} \\
\end{bmatrix} \\
&= \begin{bmatrix}
\mathbf{M}^{k-1}(i) \\ \mathbf{0}^T
\end{bmatrix} \begin{bmatrix}
1-d_{0,i-k+2} & d_{0,i-k+2} & & & 0 \\
& 1-d_{0,i-k+3} & d_{0,i-k+3} & & \\
& & \ddots & \ddots & \\
0 & & & 1-d_{0,i} & d_{0,i} \\
\end{bmatrix} \\
&+ \begin{bmatrix}
\mathbf{0}^T \\ \mathbf{M}^{k-1}(i)
\end{bmatrix} \begin{bmatrix}
-d_{1,i-k+2} & d_{1,i-k+2} & & & 0 \\
& -d_{1,i-k+3} & d_{1,i-k+3} & & \\
& & \ddots & \ddots & \\
0 & & & -d_{1,i} & d_{1,i} \\
\end{bmatrix}
\end{align}\]
(**): 为什么第一项与最后一项都为 0?这个等号前方的式子第一项与最后一项的“\(N\) 的部分”分别是 \(B_{i-k+1,k-1}(u), B_{i+1,k-1}(u)\) 的系数部分,是 \(k-1\) 阶的系数。结合公式 (3*),知道 \(k\) 阶的区间涉及范围是 \([i-k+1, i]\)。所以 \(k-1\) 阶的区间涉及范围是 \([i-k+2, i]\),没有包含 \(i-k+1,i+1\),所以在 \(k-1\) 阶时这两项为 0。
\(\mathbf{M}^k(i)\) 的初始值,一阶只含有一个 \(N^1_{0,i}\),\(\mathbf{M}^1(i)=\begin{bmatrix}1\end{bmatrix}\),即常值函数。
Uniform B-splines 的每一个区间的长度相等,假设为 1。可以将 \(\mathbf{M}^k(i)\) 的 \(d_{0,*}, d_{1,*}\) 写成分数,并且与 \(i\) 无关系,所以不管有多少结点,都只用算一次 \(\mathbf{M}^k(i)\),也就可以直接写作 \(\mathbf{M}^k\)。得到的结果如公式 (7*)。
\[\begin{align}
\mathbf{M}^k &= {1 \over k-1} \left\{ \begin{bmatrix}
\mathbf{M}^{k-1} \\ \mathbf{0}^T
\end{bmatrix} \begin{bmatrix}
1 & k-2 & & & 0 \\
& 2 & k-3 & & \\
& & \ddots & \ddots & \\
0 & & & k-1 & 0 \\
\end{bmatrix} + \begin{bmatrix}
\mathbf{0}^T \\ \mathbf{M}^{k-1}
\end{bmatrix} \begin{bmatrix}
-1 & 1 & & & 0 \\
& -1 & 1 & & \\
& & \ddots & \ddots & \\
0 & & & -1 & 1 \\
\end{bmatrix} \right\}
\end{align}\]
通项公式怎么计算,我不知道。
论文之后计算了一些偏导,但是符号系统实在让人看不懂。需要用的时候自然能推导出。
Reference
[1] Qin, Kaihuai. "General matrix representations for B-splines." The Visual Computer 16.3-4 (2000): 177-186.
[2] 【回形针PaperClip】如何设计一个逼真的三维模型.
[3] Wikipedia: De Boor's algorithm.
[4] Wikipedia: Fibonacci number 3.3 Matrix form.
[5] Unknown: Matrixces and Polynomials.
[6] De Boor, Carl. B (asic)-Spline Basics. No. MRC-TSR-2952. WISCONSIN UNIV-MADISON MATHEMATICS RESEARCH CENTER, 1986.
[7] Solomon, Justin. Numerical algorithms: methods for computer vision, machine learning, and graphics. CRC press, 2015.