B-spline (1)

1.什么是B-spline

Bezier Curve的缺点主要是：

• nonlocality，移动一个控制点会影响整条曲线
• 高阶曲线，曲线离控制点很远

对此，可以把N条Bezier Curve组合起来。如下图所示，是两条3次Bezier Curve在D点连接而得。

此时，该曲线是一条$C^0$曲线，共有7个控制点。若想得到一条$C^1$曲线，则需要两条Bezier Curve在D点的导数相同。通过Bezier Curve的导数表达式，可知需要C,D,E点共线并且$|CD|=|DE|$。此时D点失去了控制作用，共用6个控制点，如下图所示。

进一步地，若想得到一条$C^2$曲线，则需要D点的二阶导数相同。那么此时需要又一个控制点失去作用，例如E点，如下图所示。

可见，N条3阶Bezier Curve拼接成一条$C^0$曲线，需要$4N-(N-1)=3N+1$个控制点；一条$C^1$曲线需要$3N+1-(N-1)=2N+2$个控制点；一条$C^2$曲线需要$2N+2-(N-1)=N+3$个控制点。

进而，N条k阶的Bezier Curve拼接成一条$C^{k-1}$连续的B-spline曲线需要N+k个控制点。

但值得注意的是，在上图中，A,B,C,F,G并不是B-spline的控制点，控制点应该是下图的S0-S4。

B-spline具有局部性，例如移动S0，并不会对第二段Bezier Curve产生影响。同样地，移动S4，也不会对第一段Bezier Curve产生影响。B-spline的性质可总结为：

• The curve starts and ends at the first and last control vertices. 
• The curve is tangent to the first and last edges of its control polygon. 
• The B-spline curve is contained in the convex hull of its control polygon. 如下图所示。
• If the B-spline curve is of degree n, it is contained in the union of the convex hulls of every n+1 consecutive vertices. 如下图所示。

• The shape of the B-spline curve is roughly the same as the shape of its control polygon.
• The B-spline curve twists less than its control polygon. 如下图所示。

• Moving a vertex of the control polygon affects only a part of the curve close to the vertex moved.
• The B-spline curve can be translated, rotated or reflected, by translating, rotating, or reflecting the vertices of its control polygon.

2.B-spline表达式

由N条k阶的Bezier Curve组成k阶B-spline曲线有N+k个控制点。B-spline的表达式可以写为：

$$ BS(t)=\sum_{i=0}^{N+k-1}P_iB_{i,k}(t)\space\space\space t\in[0,N]\space\space\space (1) $$

当N=3,k=3时，对应的6个基函数如下图所示。

可见，基函数$B_{0,3}$的有效区间，即值不为零的区间是[0,1]，基函数$B_{1,3}$的有效区间是[0,2]。

The knots of a B-spline curve are the endpoints of the supports of the basis functions that define the curve. 这6个基函数的knots为： $$ \begin{array}{l} supp B_{0,3}=[0,1]\\ supp B_{1,3}=[0,2] \\ supp B_{2,3}=[0,3]\\ supp B_{3,3}=[0,3]\\ supp B_{4,3}=[1,3]\\ supp B_{5,3}=[2,3] \end{array} $$

The knot vector comprises the knots of the basis functions that define a S-spline curve, listed from smallest to largest. 上方B-spline对应的knot vector为$V=\left \{0,0,0,0,1,2,3,3,3,3\right \}$。

N条k阶Bezier Curve组成的k阶B-spline曲线对应的knot vector为：

3.Cox-de boor formula or De Boor's Algorithm

要想根据控制点生成B-spline，需要计算出相应的基函数。De Boor's Algorithm就是一种用knot vector递归计算出基函数的算法。以上方的3阶B-spline为例，他的思想就是通过$B_{0,0}$和$B_{1,0}$得到$B_{0,1}$，以及通过$B_{1,0}$和$B_{2,0}$得到$B_{1,1}$，进而通过$B_{0,1}$和$B_{1,1}$得到$B_{0,2}$，如此递归地得到$B_{0,3}-B_{5,3}$。

De Boor's Algorithm计算公式如下： $$ B_{i,k}(t)=\frac{t-V[i]}{V[i+k]-V[i]}B_{i,k-1}(t)+\frac{V[i+k+1]-t}{V[i+k+1]-V[i+1]}B_{i+1,k-1}(t)\space\space\space(2) $$ 并且 $$ B_{i,0}(t)=\begin{cases} 1,\space\space\space t\in[V[i],V[i+1]) \\ 0,\space\space\space otherwise \end{cases}\space\space\space(3) $$ 以及$B_{i,0}=0\space\space if\space V[i]=V[i+1]$，以及分母为0，则为0。

python代码如下：

import numpy as np
import matplotlib.pyplot as plt

def generateBspline(controlPoints,k,num=100):
    n=len(controlPoints)
    # knot vector
    v0=np.linspace(0,0,k+1)
    v1=np.linspace(1,n-k-1,n-k-1)
    v2=np.linspace(n-k,n-k,k+1)
    knot=np.concatenate((v0,v1,v2))
    
    t=np.linspace(0,n-k,num)

    # order 0
    BF=np.zeros((knot.size-1,num))
    for i in range(knot.size-1):
        if knot[i]!=knot[i+1]:
            for j in range(num):
                time=t[j]
                if i==knot.size-k-2:
                    if knot[i]<=time<=knot[i+1]:BF[i][j]=1
                else:
                    if knot[i]<=time<knot[i+1]:BF[i][j]=1
    
    # recurison
    for order in range(1,k+1):
        tmp=np.zeros((knot.size-1-order,num))
        for i in range(knot.size-1-order):
            for j in range(num):
                time=t[j]
                if knot[i+order]!=knot[i]:
                    tmp[i][j]+=(time-knot[i])*BF[i][j]/(knot[i+order]-knot[i])
                if knot[i+order+1]!=knot[i+1]:
                    tmp[i][j]+=(knot[i+order+1]-time)*BF[i+1][j]/(knot[i+order+1]-knot[i+1])
        BF=tmp
    
    points=np.zeros((num,2))
    for i in range(num):
        for p in range(n):
            point=controlPoints[p]
            points[i]+=BF[p][i]*point
    return points
    
def drawBspline(points,controlPoints):
    fig,ax=plt.subplots()
    ax.plot(points[:,0],points[:,1])
    ax.plot(controlPoints[:,0],controlPoints[:,1],'o--')
    plt.show()

# control points S0-S4
controlPoints=np.array([[0, 0], [2, 3], [6, 3], [9, 0],[10,3]])
points=generateBspline(controlPoints,3)
drawBspline(points,controlPoints)