B-spline算法(B样条曲线)

个人博客：https://www.vectormoon.net/
Bezier曲线算法：https://www.vectormoon.net/2020/09/25/Bezier/

算法产生背景

Bezier曲线有以下几个不足点，所以导致出现了B-spline算法：

• 一旦确定特征多边形，就确定了曲线的阶次
• Bezier曲线拼接复杂（需要满足几何连续性，参数连续性等）
• Bezier曲线不能作局部修改（只能整体修改）

B-spline算法是整条曲线用一段一段的曲线连接而成，采用分段连续多段式生成

B-spline曲线定义

B-spline曲线定义为：

$$P(u)=\sum _{i=0}^n{P}_i{B}_{i,k}(u)u\in [u_{k-1},u_{n+1}]$$

其中$P_i$是特征多边形的顶点；$B_{i,k}$称为k阶（k-1次）基函数，B-spline算法阶数是次数加1，这是和Bezier算法的一个不同之处；定义域的解释之后会给出，先给出基函数算法。

B-spline基函数的求出算法应用最广泛的是deBoor-cox递推算法：

$$B_{i,k}(u)=\frac{u-u_i}{u_{i+k-1}-u_i}\ast B_{i,k-1}u+\frac{u_{i+k}-u}{u_{i+k}-u_{i+1}}\ast B_{i+1,k+1}(u)$$

$$B_{i,1}(u)=\{\begin{array}{ll}1& u_i<u<u_{i+1};\\ 0& Otherwise.\end{array}$$

规定0/0=0

B-spline曲线的定义域为$u\in [u_{k-1},u_{n+1}]$。设U为所有节点矢量的集合，显见节点表个数为$n+k+1$个。举例说明，当$n=4,k=4$时，有 U = { u 0 , u 1 , u 2 , u 3 , u 4 , u 5 , u 6 , u 7 , u 8 } U=\{u_0,u_1,u_2,u_3,u_4,u_5,u_6,u_7,u_8\} U={u0​,u1​,u2​,u3​,u4​,u5​,u6​,u7​,u8​}。第一项为$P_0{B}_{0,4}(u)$，由deBoor-cox算法可知，其涉及到$u_0$到$u_4$五个点；第二项同理，设计到$u_1$到$u_5$五个点；其余同理；所以可以画出区间对应坐标轴：

区间合法所需要的条件为：区间内必须有足够基函数与顶点对应，也即区间中基函数覆盖较多的区间才是一个合法区间。所以上例中对应的合法区间为$u\in [u_3,u_5]$也就是$u\in [u_{k-1},u_{n+1}]$

B-spline的类型划分

有均匀B样条曲线，准均匀B样条曲线，分段Bezier曲线，非均匀B样条曲线这几种方法，这里采用均匀B样条曲线的方法绘制曲线。

均匀B样条曲线

定义：当节点沿参数轴均匀等距分布，即$u_{i+1}-u_i=C>0$时，为均匀B样条函数，比如： { 0 , 1 , 2 , 3 , 4 , 5 , 6 } , { 0 , 0.2 , 0.4 , 0.6 , 0.8 , 1 } \{0,1,2,3,4,5,6\},\{0,0.2,0.4,0.6,0.8,1\} {0,1,2,3,4,5,6},{0,0.2,0.4,0.6,0.8,1}

B-spline基函数($B_{i,k}$)的递推公式计算需要用到这里的思想。这里还使用上面当$n=4,k=4$的例子，从上面的定义可知，我们可以把节点集合 U = { u 0 , u 1 , u 2 , u 3 , u 4 , u 5 , u 6 , u 7 , u 8 } U=\{u_0,u_1,u_2,u_3,u_4,u_5,u_6,u_7,u_8\} U={u0​,u1​,u2​,u3​,u4​,u5​,u6​,u7​,u8​}写成这样 U = { 0 , 1 , 2 , 3 , 4 , 5 , 6 , 7 , 8 } U=\{0,1,2,3,4,5,6,7,8\} U={0,1,2,3,4,5,6,7,8}也即$u_0=0,u_1=1,...,u_8=8$，当然U也可以写成其他集合，因为计算的时候是一个比例，所以对计算结果并没有影响，综上就可以简单的算出对应基函数的值。

代码实现

def B_spline(p_list):
	"""
	:param p_list: (list of list of int:[[x0, y0], [x1, y1], ...])point set of p
	result: (list of list of int:[[x0, y0], [x1, y1], ...])point on curve
	绘制三次(四阶)均匀B样条曲线
	"""
	result = []
	n = len(p_list)
	k = 4
	u = k-1
	while (u < n+1):
		x, y = 0, 0
		#calc P(u)
		for i in range(0, n):
			B_ik = deBoor_Cox(u, k, i)
			x += B_ik * p_list[i][0]
			y += B_ik * p_list[i][1]
		result.append((int(x+0.5), int(y+0.5)))
		u += 1/20927 #2020/09/27
	return result

def deBoor_Cox(u, k, i):
	if k==1:
		if i <= u and u <= i+1:
			return 1
		else:
			return 0
	else:
		coef_1, coef_2 = 0, 0
		if (u-i == 0) and (i+k-1-i == 0):
			coef_1 = 0
		else:
			coef_1 = (u-i) / (i+k-1-i)
		if (i+k-u == 0) and (i+k-i-1 == 0)
			coef_2 = 0
		else:
			coef_2 = (i+k-u) / (i+k-i-1)
	return coef_1 * deBoor_Cox(u, k-1, i) + coef_2 * deBoox_Cox(u, k-1, i+1)