辽宁师范大学数学学院，辽宁 大连

收稿日期：2021年11月13日；录用日期：2021年12月9日；发布日期：2021年12月16日

摘要

计算机技术的成熟，使得高维球面上的样条曲线拟合技术在关键帧动画和航空航天飞行器运行轨迹规划等领域都得到了广泛的应用。实际中常需用连续性好，形状可控性高的球面样条曲线来插值高维球面中的数据点。因此，本文的研究目的是提供一种可在任意维球面上构造几何连续的球面插值样条曲线的新方法。首先，将欧氏空间中Beta样条曲线推广到高维球面上，定义了球面Beta样条曲线，并对其连续性进行了讨论；然后，明确了辅助控制顶点与插值点之间的关系；最后，利用已知的插值点来计算辅助控制顶点，从而构造出几何连续的球面Beta插值样条曲线。该方法属于整体构造方法，松弛了连续要求，引入了可以控制样条曲线形状的参数，从而使得构造的样条具有灵活控制性。

关键词

球面Bézier曲线，球面Beta样条，插值，几何连续

Research on the Spherical Beta Spline Curse

Lei Yang, Qian Wang^*, Jiahui Chen, Nan Sun, Yao He

School of Mathematics, Liaoning Normal University, Dalian Liaoning

Received: Nov. 13^th, 2021; accepted: Dec. 9^th, 2021; published: Dec. 16^th, 2021

ABSTRACT

The maturity of computer technology has made the spline curve fitting technology on the high-dimensional spheres widely used in the fields of key frame animation and aerospace vehicle trajectory planning. In practice, it is often necessary to use spherical splines with well continuity and high shape controllability to interpolate points on high-dimensional spheres. Therefore, the purpose of this paper is to provide a new method for constructing geometrically continuous spherical interpolation splines on any dimensional spheres. First, the Beta spline in Euclidean space is extended to the high-dimensional spheres. The spherical Beta spline curve is defined, and its continuity is discussed. Second, the relationship between the auxiliary control vertex and the interpolation point is clarified. Third, the known interpolation points are used to calculate the auxiliary control vertices, all then a geometrically continuous spherical Beta interpolation spline is constructed. This method has overall properties. It relaxes the continuous requirements and introduces parameters that can control the shape of the spline curve, which makes the spherical splines more flexible and controllable.

Keywords:Spherical Bézier Curve, Spherical Beta Spline, Interpolation, Geometric Continuity

Copyright © 2021 by author(s) and Hans Publishers Inc.

This work is licensed under the Creative Commons Attribution International License (CC BY 4.0).

http://creativecommons.org/licenses/by/4.0/

1. 前言

随着计算机技术的快速发展，球面样条曲线在众多研究领域中都发挥其重要作用。研究球面样条曲线的构造方法对关键帧动画，航空航天飞行器运行轨迹规划等问题都有重要的理论支撑和实际意义。这一类问题用欧氏空间中的样条曲线处理效果不佳，因此球面样条曲线的发展至关重要。接下来我们就对如何构造一段光滑的球面样条曲线进行讨论。

受欧氏空间中经典样条曲线构造方法 [1] 的启发，人们提出了多种构造球面插值样条曲线的方法 [2] - [16]，下面介绍一下与本文相关的工作：

1985年，SHOEMAKE [7] 提出了圆弧插值法，将欧氏空间中的De Casteljau算法向单位球面 $S^3$ 进行了推广，从而定义出球面Bézier曲线，并研究了球面Bézier的相关性质，为后续对球面样条曲线的构造研究提供了新工具。NOAKES [8] 将上述算法推广到了m维单位球面 $S^m$ 上，并构造出了 $C^1$ 球面Bézier样条曲线。POPIEL和NOAKES [9] 对球面Bézier曲线在端点处的性质进一步地进行了研究，并利用共轭导矢的概念构造出 $C^2$ 球面Bézier样条。以此为基础可以将Bézier曲线或B样条曲线的相关算法推广到任意维单位球面中。王宁，唐月红在文 [11] 中引入球面NURBS曲线，并给出了等距三次球面B样条曲线保持连续性的插值方法。该球面NURBS曲线不具有分裂性质，即曲线无法分割为多条独立表示的曲线，因此降低了对所构造的球面样条曲线形状的可控性。Nielson在文 [12] 中构造了 $\nu$ -四元数样条曲线，通过迭代法解决了控制顶点的计算，引入的张力参数使得样条曲线具有更加灵活的形状调控能力。罗钟铉，王倩在文 [13] 中解析地给出了球面Bézier曲线段在连接点处 $G^2$ 光滑拼接时，其控制顶点的计算公式，从而提出一种球面几何连续插值的新方法。罗钟铉，王倩等人又在文 [14] 中引入广义有理Bézier曲线，并完成了高维球面上的 $G^2$ Hermite插值问题，提出了一种新的刚体运动方向平滑插值方法。同年，XING等 [15] 将一类广义B样条曲线类推广到单位球面 $S^3$，该类样条曲线具有单调递增且连续可微的基函数，通过修改基函数可改变曲线形状。杨军，赵桔焓 [16] 将 $R^3$ 中一类形状可调Bézier曲线推广到单位球面 $S^3$ 上，利用广义Bernstein基函数及四元数构造出含参数的 $C^1$ 插值球面Bézier曲线。但XING和杨军等人的方法只适用于单位球面 $S^3$，不能推广到高维球面上，具有一定局限性。王倩，潘乐等 [17] [18] 仅利用插值点构造出具有局部性质的 $C^2$ 连续的球面插值样条曲线，扰动个别插值点不会对样条曲线的全局产生影响。

1981年，B. A. Barsky [19] 提出了几何连续的Beta样条，Beta样条曲线具有B样条曲线的若干性质。同时，由于该类样条引入了两组形状参数，可在不改变控制顶点的前提下，通过改变形状参数的值来调整曲线的形状，使样条曲线具有较大的灵活性。基于此，本文的研究目的是利用球面Bézier曲线构造出一种与欧氏空间中Beta样条曲线的几何性质相似的球面Beta样条曲线，由此提供一种在任意维球面中几何连续插值的新方法。

本文结构安排如下，第一节介绍了研究背景，研究意义及内容安排；第二节给出了预备知识；第三节介绍了可在任意维球面上构造几何连续的球面插值样条曲线的新方法；第四节将本文方法应用到二维球面样条曲线的构造设计中，并观察形状参数对样条曲线形状的影响；第五节对全文进行了总结。

2. 预备知识

欧氏空间中的De Casteljau算法是递推定义Bézier曲线的一种算法。该算法可以用几何作图来实现曲线上任一点的计算，也是分析Bézier曲线性质的重要工具。SHOEMAKE将欧氏空间中两点间直线的插值方法推广到球面 $S^3$ 上，提出了圆弧插值 [7]。若 $q_0$ 和 $q_1$ 是球面上的两个点，那么圆弧插值定义为：

$G\left[q_0,q_1\right]\left(t\right)=\frac{\sin \left(\left(1-t\right)\theta \left(q_0,q_1\right)\right)}{\sin \left(\theta \left(q_0,q_1\right)\right)}{q}_0+\frac{\sin \left(t\theta \left(q_0,q_1\right)\right)}{\sin \left(\theta \left(q_0,q_1\right)\right)}{q}_1,t\in \left[0,1\right].$ (1)

其中 $\theta \left(.,.\right)$ 表示两个向量之间的夹角， $G\left[q_0,q_1\right]\left(t\right)$ 可以理解为连接 $q_0$ 和 $q_1$ 的球面直线。利用该插值方法，SHOEMAKE将欧氏空间中的De Casteljau算法向低维球面进行了推广。

进一步地，POPIEL和NOAKES将De Casteljau算法推广到了球面 $S^m$ 上，并给出了球面Bézier曲线的递推定义 [7]。

给定球面 $S^m$ 上的互异点 $q_0,q_1,\cdots ,q_n$，对于 $t\in \left[0,1\right]$，$k=1,2,3,\cdots ,n$，$i=0,1,\cdots ,n-k$，令

$\begin{array}{c}{q}_i^k\left(t\right)=G\left[q_i^{k-1}\left(t\right),q_{i+1}^{k-1}\left(t\right)\right]\left(t\right)\\ =\frac{\sin \left(\left(1-t\right){\theta }_i^{k-1}\right)}{\sin \left({\theta }_i^{k-1}\right)}{q}_i^{k-1}\left(t\right)+\frac{\sin \left(t{\theta }_i^{k-1}\right)}{\sin \left({\theta }_i^{k-1}\right)}{q}_{i+1}^{k-1}\left(t\right).\end{array}$ (2)

其中

$q_i^0=q_i,$

${\theta }_i^{k-1}=\theta \left(q_i^{k-1}\left(t\right),q_{i+1}^{k-1}\left(t\right)\right).$

则称曲线

$Q\left[q_0,q_1,\cdots ,q_n\right]\left(t\right)=q_0^n\left(t\right),$ (3)

是以 $q_0,q_1,\cdots ,q_n$ 为控制顶点的n次球面Bézier曲线。

根据上述定义，POPIEL和NOAKES给出了连接球面 $S^m$ 上两个点 $q_0$ 和 $q_1$ 的一次球面Bézier曲线 $G\left[q_0,q_1\right]$ 的一些性质：

引理1 [9] 对于任意 $t\in \left[0,1\right]$，有

(1) $G\left[q_0,q_1\right]=\theta \left(q_0,q_1\right)$ ；

(2) $\stackrel{\dot{}}{G}\left[q_0,q_1\right]\left(t\right)=\frac{-\theta \left(q_0,q_1\right)\cos \left(\left(1-t\right)\theta \left(q_0,q_1\right)\right)}{\sin \left(\theta \left(q_0,q_1\right)\right)}{q}_0+\frac{\theta \left(q_0,q_1\right)\cos \left(t\theta \left(q_0,q_1\right)\right)}{\sin \left(\theta \left(q_0,q_1\right)\right)}{q}_1$ ；

(3) $\Vert \stackrel{\dot{}}{G}\left[q_0,q_1\right]\left(t\right)\Vert =\theta \left(q_0,q_1\right)$。

此外，POPIEL和NOAKES又给出了n次球面Bézier曲线在端点处的一些性质：

引理2 [9] 一条n次球面Bézier曲线 $t\mapsto q_0^n\left(t\right)$ 在首末端点处有：

(1) $q_0^n\left(0\right)=q_0$ ；

(2) $q_0^n\left(1\right)=q_n$ ；

(3) ${\stackrel{\dot{}}{q}}_0^n\left(0\right)=n\stackrel{\dot{}}{G}\left[q_0,q_1\right]\left(0\right)$ ；

(4) ${\stackrel{\dot{}}{q}}_0^n\left(1\right)=n\stackrel{\dot{}}{G}\left[q_{n-1},q_n\right]\left(1\right)$。

3. 球面Beta样条曲线

3.1. 球面Beta样条曲线的定义

设 $\beta \left(t\right),t\in \left[0,n\right]$ 是一条由n段三次球面Bézier曲线构成的三次样条曲线，根据式(3)，该三次样条曲线的第i段 $\left(i=1,\cdots ,n\right)$ 可以表示为：

${\beta }_i\left(t\right)=Q\left[q_{i-1,0},q_{i-1,1},q_{i-1,2},q_{i-1,3}\right]\left(t-i+1\right),t\in \left[i-1,i\right].$

其控制顶点为 $q_{i-1,k}\left(k=0,1,2,3\right)$，且 $q_{i,0}=q_{i-1,3}$。

如图1，给定一组顶点 $d_0,d_1,\cdots ,d_n$ 和两组形状参数 ${\beta }_{1，0},{\beta }_{1,1},\cdots ,{\beta }_{1,n},{\beta }_{2,0},{\beta }_{2,1},\cdots ,{\beta }_{2,n}$，当 $\beta \left(t\right)$ 的控制顶点 $q_{i-1,k}$ 按如下关系式确定时：

$q_{i-1,1}=G\left[d_{i-1},d_i\right]\left(\frac{{\gamma }_{i-1}}{{\gamma }_{i-1}+1+{\beta }_{1,i}^2{\gamma }_i}\right),i=1,2,\cdots ,n,$ (4)

$q_{i-1,2}=G\left[d_{i-1},d_i\right]\left(\frac{{\gamma }_{i-1}+1}{{\gamma }_{i-1}+1+{\beta }_{1,i}^2{\gamma }_i}\right),i=1,2,\cdots ,n,$ (5)

$q_{i-1,0}=\{\begin{array}{ll}{d}_0,\hfill & i=1,\hfill \\ G\left[q_{i-2,2},q_{i-1,1}\right]\left(\frac{1}{1+{\beta }_{1,i\text{-}1}}\right),\hfill & i=2,\cdots ,n,\hfill \\ d_n,\hfill & i=n+1.\hfill \end{array}$ (6)

其中

${\gamma }_i=\frac{2\left({\beta }_{1,i}+1\right)}{{\beta }_{2,i}+2{\beta }_{1,i}\left({\beta }_{1,i}+1\right)}$ (7)

图1. 三次球面Beta样条曲线的递推定义

则称该样条曲线为三次球面Beta样条曲线，并称给定的顶点 $d_i$ 为辅助控制顶点。实际上，三次球面Beta样条曲线是欧氏空间中Beta样条曲线的球面推广。另外，由引理1(1)可知 $q_{i-1,2},q_{i,0},q_{i,1}$ 之间有如下关系：

$\frac{\theta \left(q_{i-1,2},q_{i,0}\right)}{\theta \left(q_{i,0},q_{i,1}\right)}=\frac{1}{{\beta }_{1,i-1}},i=1,2,\cdots ,n-1.$ (8)

3.2. 球面Beta样条曲线的性质

由球面Bézier曲线的性质可知，当 $i=1,2,\cdots ,n$ 时，球面Beta样条曲线 $\beta \left(t\right)$ 在区间 $\left(i-1,i\right)$ 是无限次可微的，即 $G^{\infty }$ 连续的，且在 $t=i$ 处是 $G^0$ 连续的。下面将进一步讨论样条曲线在 $t=i$ 处的性质：

定理1当 $i=1,2,\cdots ,n$ 时，球面Beta样条曲线 $\beta \left(t\right)$ 在 $t=i$ 处是 $G^1$ 连续的，且满足 ${\stackrel{\dot{}}{\beta }}_{i+1}\left(i\right)={\beta }_{1,i-1}{\stackrel{\dot{}}{\beta }}_i\left(i\right)$。

证明：首先由引理1(3)，引理2(3)和(4)可知，曲线段 ${\beta }_i\left(t\right)$ 与 ${\beta }_{i+1}\left(t\right)$ 在连接点 $q_{i,0}$ 处的单位切方向分别为 $\frac{\stackrel{\dot{}}{G}\left[q_{i-1,2},q_{i,0}\right]\left(1\right)}{\theta \left(q_{i-1,2},q_{i,0}\right)}$ 和 $\frac{\stackrel{\dot{}}{G}\left[q_{i,0},q_{i,1}\right]\left(0\right)}{\theta \left(q_{i,0},q_{i,1}\right)}$。按引理1(2)可将上述式子展开得：

$\frac{\stackrel{\dot{}}{G}\left[q_{i-1,2},q_{i,0}\right]\left(1\right)}{\theta \left(q_{i-1,2},q_{i,0}\right)}=\frac{-1}{\sin \left(\theta \left(q_{i-1,2},q_{i,0}\right)\right)}{q}_{i-1,2}+\frac{\cos \left(\theta \left(q_{i-1,2},q_{i,0}\right)\right)}{\sin \left(\theta \left(q_{i-1,2},q_{i,0}\right)\right)}{q}_{i,0},$ (9)

$\frac{\stackrel{\dot{}}{G}\left[q_{i,0},q_{i,1}\right]\left(0\right)}{\theta \left(q_{i,0},q_{i,1}\right)}=\frac{-\cos \left(\theta \left(q_{i,0},q_{i,1}\right)\right)}{\sin \left(\theta \left(q_{i,0},q_{i,1}\right)\right)}{q}_{i,0}+\frac{1}{\sin \left(\theta \left(q_{i,0},q_{i,1}\right)\right)}{q}_{i,1}.$ (10)

由式(6)可知：

$q_{i,0}=\frac{\sin \left(\left(\frac{{\beta }_{1,i-1}}{1+{\beta }_{1,i-1}}\right)\theta \left(q_{i-1,2},q_{i,1}\right)\right)}{\sin \left(\theta \left(q_{i-1,2},q_{i,1}\right)\right)}{q}_{i-1,2}+\frac{\sin \left(\left(\frac{1}{1+{\beta }_{1,i-1}}\right)\theta \left(q_{i-1,2},q_{i,1}\right)\right)}{\sin \left(\theta \left(q_{i-1,2},q_{i,1}\right)\right)}{q}_{i,1}.$

将式(8)带入上式，通过整理即可得：

$q_{i-1,2}=\frac{\sin \left(\theta \left(q_{i-1,2},q_{i,1}\right)\right)}{\sin \left(\theta \left(q_{i,0},q_{i,1}\right)\right)}{q}_{i,0}-\frac{\sin \left(\theta \left(q_{i-1,2},q_{i,0}\right)\right)}{\sin \left(\theta \left(q_{i,0},q_{i,1}\right)\right)}{q}_{i,1},$ (11)

将式(8)和(11)带入式(9)可得：

$\frac{\stackrel{\dot{}}{G}\left[q_{i-1,2},q_{i,0}\right]\left(1\right)}{\theta \left(q_{i-1,2},q_{i,0}\right)}=\frac{1}{\sin \left(\theta \left(q_{i,0},q_{i,1}\right)\right)}{q}_{i,1}+\left(\frac{-\sin \left(\theta \left(q_{i-1,2},q_{i,1}\right)\right)+\sin \left(\theta \left(q_{i,0},q_{i,1}\right)\right)\cos \left(\theta \left(q_{i-1,2},q_{i,0}\right)\right)}{\sin \left(\theta \left(q_{i-1,2},q_{i,0}\right)\right)\sin \left(\theta \left(q_{i,0},q_{i,1}\right)\right)}\right)q_{i,0},$

利用式(9)以及积化和差公式可将上式化简为：

$\frac{\stackrel{\dot{}}{G}\left[q_{i-1,2},q_{i,0}\right]\left(1\right)}{\theta \left(q_{i-1,2},q_{i,0}\right)}=\frac{1}{\sin \left(\theta \left(q_{i,0},q_{i,1}\right)\right)}{q}_{i,1}+\left(\frac{\sin \left(\left(\frac{1}{1+{\beta }_{1,i-1}}\right)\theta \left(q_{i,0},q_{i,1}\right)\right)-\sin \left(\left(1-\frac{1}{{\beta }_{1,i-1}}\right)\theta \left(q_{i,0},q_{i,1}\right)\right)}{\cos \left(\left(\frac{1}{1+{\beta }_{1,i-1}}\right)\theta \left(q_{i,0},q_{i,1}\right)\right)-\cos \left(\left(1-\frac{1}{{\beta }_{1,i-1}}\right)\theta \left(q_{i,0},q_{i,1}\right)\right)}\right)q_{i,0},$

再利用和差化积公式整理可得：

$\frac{\stackrel{\dot{}}{G}\left[q_{i-1,2},q_{i,0}\right]\left(1\right)}{\theta \left(q_{i-1,2},q_{i,0}\right)}=\frac{-\cos \left(\theta \left(q_{i,0},q_{i,1}\right)\right)}{\sin \left(\theta \left(q_{i,0},q_{i,1}\right)\right)}{q}_{i,0}+\frac{1}{\sin \left(\theta \left(q_{i,0},q_{i,1}\right)\right)}{q}_{i,1}.$ (12)

从而由式(10)和(12)可知 ${\beta }_i\left(t\right)$ 与 ${\beta }_{i+1}\left(t\right)$ 在连接点 $q_{i,0}$ 处的单位切方向相同，且 $\stackrel{\dot{}}{G}\left[q_{i,0},q_{i,1}\right]\left(0\right)={\beta }_{1,i-1}\stackrel{\dot{}}{G}\left[q_{i-1,2},q_{i,0}\right]\left(1\right)$，由此即可证明定理。

4. 插值样条曲线的构造

由上述讨论可知，球面Beta样条曲线的控制顶点是由辅助控制顶点和两组形状参数确定。但在实际应用中，往往已知的是一组点列 $q_0,q_1,\cdots ,q_n$，想要使用球面Beta样条曲线来光滑插值这些点 [20]。那么根据球面Beta样条曲线的定义，首先需要先计算出辅助控制顶点 $d_i,i=0,1,\cdots ,n$，然后才能计算出各控制顶点 $q_{i-1,k},i=1,2,\cdots ,n;k=0,1,2,3$，最终构造出所需曲线。

构造一条通过上述插值点 $q_0,q_1,\cdots ,q_n$ 的球面Beta样条曲线的方法如下。

4.1. 辅助控制顶点与插值点的关系

令 $q_{i,0}=q_i,i=0,1,\cdots ,n$，由引理2可知 $\beta \left(i\right)=q_i$，即样条曲线 $\beta \left(t\right)$ 插值点列 $q_0,q_1,\cdots ,q_n$。下面考虑辅助控制顶点与插值点之间的关系。

定理2球面Beta样条曲线的插值点 $q_0,q_1,\cdots ,q_n$ 与辅助控制顶点 $d_0,d_1,\cdots ,d_n$ 之间存在如下关系：

$\{\begin{array}{l}{d}_0=q_0,\\ \left(\frac{\sin \left(\left(1-{\omega }_{i-1}\right){\phi }_{i-1}\right)\sin \left(\left(1-{\delta }_i\right){\psi }_i\right)}{\sin \left({\phi }_{i-1}\right)\sin \left({\psi }_i\right)}\right)d_{i-1}\\ +\left(\frac{\sin \left({\omega }_{i-1}{\phi }_{i-1}\right)\sin \left(\left(1-{\delta }_i\right){\psi }_i\right)}{\sin \left({\phi }_{i-1}\right)\sin \left({\psi }_i\right)}+\frac{\sin \left(\left(1-{\eta }_i\right){\phi }_i\right)\sin \left({\delta }_i{\psi }_i\right)}{\sin \left({\phi }_i\right)\sin \left({\psi }_i\right)}\right)d_i\\ +\left(\frac{\sin \left({\eta }_i{\phi }_i\right)\sin \left({\delta }_i{\psi }_i\right)}{\sin \left({\phi }_i\right)\sin \left({\psi }_i\right)}\right)d_{i+1}=q_i,i=1,2,\cdots ,n-1,\\ d_n=q_n.\end{array}$ (13)

其中

$\begin{array}{l}{\phi }_i=\theta \left(d_i,d_{i+1}\right),{\psi }_i=\theta \left(q_{i-1,2},q_{i,1}\right),{\delta }_i=\frac{1}{1+{\beta }_{1,i}},\\ {\omega }_i=\frac{{\gamma }_i+1}{{\gamma }_i+1+{\beta }_{1,i+1}^2{\gamma }_{i+1}},{\eta }_i=\frac{{\gamma }_i}{{\gamma }_i+1+{\beta }_{1,i+1}^2{\gamma }_{i+1}}.\end{array}$ (14)

证明：由三次球面Beta样条的定义可知， $d_0=q_0,d_n=q_n$，利用式(1)将式(4)~(6)分别展开，即可得到：

$q_{i-1,2}=\frac{\sin \left(\left(1-{\omega }_{i-1}\right){\phi }_{i-1}\right)}{\sin \left({\phi }_{i-1}\right)}{d}_{i-1}+\frac{\sin \left({\omega }_{i-1}{\phi }_{i-1}\right)}{\sin \left({\phi }_{i-1}\right)}{d}_i,$ (15)

$q_{i,1}=\frac{\sin \left(\left(1-{\eta }_i\right){\phi }_i\right)}{\sin \left({\phi }_i\right)}{d}_i+\frac{\sin \left({\eta }_i{\phi }_i\right)}{\sin \left({\phi }_i\right)}{d}_{i+1},$ (16)

$q_{i,0}=\frac{\sin \left(\left(1-{\delta }_i\right){\psi }_i\right)}{\sin \left({\psi }_i\right)}{q}_{i-1,2}+\frac{\sin \left({\delta }_i{\psi }_i\right)}{\sin \left({\psi }_i\right)}{q}_{i,1},$ (17)

将式(15) (16)代入式(17)即可得控制顶点 $q_{i,0}$ 与辅助控制顶点 $d_i$ 之间的关系：

$\begin{array}{c}{q}_{i,0}=\left(\frac{\sin \left(\left(1-{\omega }_{i-1}\right){\phi }_{i-1}\right)\sin \left(\left(1-{\delta }_i\right){\psi }_i\right)}{\sin \left({\phi }_{i-1}\right)\sin \left({\psi }_i\right)}\right)d_{i-1}\\ +\left(\frac{\sin \left({\omega }_{i-1}{\phi }_{i-1}\right)\sin \left(\left(1-{\delta }_i\right){\psi }_i\right)}{\sin \left({\phi }_{i-1}\right)\sin \left({\psi }_i\right)}+\frac{\sin \left(\left(1-{\eta }_i\right){\phi }_i\right)\sin \left({\delta }_i{\psi }_i\right)}{\sin \left({\phi }_i\right)\sin \left({\psi }_i\right)}\right)d_i\\ +\left(\frac{\sin \left({\eta }_i{\phi }_i\right)\sin \left({\delta }_i{\psi }_i\right)}{\sin \left({\phi }_i\right)\sin \left({\psi }_i\right)}\right)d_{i+1}.\end{array}$ (18)

由此即可证明定理。

4.2. 辅助控制顶点的计算

由上述讨论可知，当样条曲线插值于给定点列 $q_0,q_1,\cdots ,q_n$ 时，其辅助控制顶点 $d_0,d_1,\cdots ,d_n$ 是由 $n+1$ 个方程构成的 $n+1$ 元非线性方程组(13)的解。我们采用迭代法来近似计算非线性方程组(13)的解，具体迭代方法如下：

输入：插值点 $q_0,q_1,\cdots ,q_n$ ；

形状参数 ${\beta }_{1，0},{\beta }_{1,1},\cdots ,{\beta }_{1,n},{\beta }_{2,0},{\beta }_{2,1},\cdots ,{\beta }_{2,n}$。

输出：辅助控制顶点 $d_0,d_1,\cdots ,d_n$。

Step1 选择初始辅助控制顶点 $d_i^0=q_i,i=0,1,\cdots ,n$ ；

Step2令 $k=1$ ；

Step3根据式(14)计算 ${\phi }_i^{k-1},{\psi }_i^{k-1},{\delta }_i^{k-1}$ ；

Step4令 $d_0^k=q_0$，$d_n^k=q_n$，对于 $i=1,2,\cdots ,n-1$，根据定理2计算辅助控制顶点 $d_i^k$ ：

$\begin{array}{c}{d}_i^k=\frac{\sin \left({\psi }_i^{k-1}\right)}{{\lambda }_{i-1}^{k-1}}{q}_i-\frac{\sin \left(\left(1-{\omega }_{i-1}^{k-1}\right){\phi }_{i-1}^{k-1}\right)\sin \left(\left(1-{\delta }_i^{k-1}\right){\psi }_i^{k-1}\right)}{{\lambda }_{i-1}^{k-1}\sin \left({\phi }_{i-1}^{k-1}\right)}{d}_{i-1}^{k-1}\\ -\frac{\sin \left({\eta }_i^{k-1}{\phi }_i^{k-1}\right)\sin \left({\delta }_i^{k-1}{\psi }_i^{k-1}\right)}{{\lambda }_{i-1}^{k-1}\sin \left({\phi }_i^{k-1}\right)}{d}_{i+1}^{k-1}.\end{array}$ (19)

其中

${\lambda }_{i-1}^{k-1}=\frac{\sin \left({\omega }_{i-1}^{k-1}{\phi }_{i-1}^{k-1}\right)\sin \left(\left(1-{\delta }_i^{k-1}\right){\psi }_i^{k-1}\right)}{\sin \left({\phi }_{i-1}^{k-1}\right)}+\frac{\sin \left(\left(1-{\eta }_i^{k-1}\right){\phi }_i^{k-1}\right)\sin \left({\delta }_i^{k-1}{\psi }_i^{k-1}\right)}{\sin \left({\phi }_i^{k-1}\right)}.$

Step5如果 ${\displaystyle \underset{i=0}{\overset{n}{\sum }}{\Vert d_i^k-d_i^{k-1}\Vert }^2}<{10}^{-5}$，则算法终止；否则，令 $k:=k+1$，并转向Step3。

利用上述算法可以近似计算出球面Beta样条曲线的辅助顶点 $d_0,d_1,\cdots ,d_n$，进而计算出样条曲线的控制顶点 $q_{i-1,k},i=1,2,\cdots ,n;k=0,1,2,3$，从而构造出插值给定点列 $q_0,q_1,\cdots ,q_n$ 的球面Beta样条曲线。

5. 实验

在本节中，我们将通过几个数值例子来说明本文方法的优势。

例1随机给定球面 $S^2$ 上的四个插值点：

$q_0=\left(0.6296,-0.1171,-0.7680\right),$

$q_1=\begin{array}{l}\left(-0.8142,-0.2856,0.5055\right)\hfill \end{array},$

$q_2=\begin{array}{l}\left(-0.4037,0.6855,0.6059\right)\hfill \end{array},$

$q_3=\begin{array}{l}\left(-0.7206,-0.5462,-0.4270\right)\hfill \end{array}.$

令形状参数 ${\beta }_{1,j}=1$，${\beta }_{2,j}=0$ $\left(j=1,2,3,4\right)$，利用本文方法，即可构造出一条参数连续的球面插值样条曲线，如图2所示。

图2. 基于给定数据的球面插值样条曲线

例2随机给定球面 $S^2$ 上的四个插值点：

$q_0=\left(\text{1,0,0}\right),$

$q_1=\begin{array}{l}\left(0.6852,-0.5424,0.4861\right)\hfill \end{array},$

$q_2=\begin{array}{l}\left(\text{0}\text{.7071},0,\text{0}\text{.7071}\right)\hfill \end{array},$

$q_3=\begin{array}{l}\left(\text{0,1,0}\right)\hfill \end{array}.$

令形状参数 ${\beta }_{1,j}=1$，${\beta }_{2,j}=0$ $\left(j=1,2,3,4\right)$，类似的可构造出一条参数连续的球面插值样条曲线，如图3黑色曲线所示。进一步地，改变形状参数 ${\beta }_{1,2}$ 的值，从图3中可以看出样条曲线的形状发生变化。随着 ${\beta }_{1,2}$ 值的减小，第一段曲线曲率增大，第二段曲线曲率减小。由此可知，改变形状参数的值就可以对样条曲线段产生影响，从而实现对样条曲线形状的控制。

图3. 形状参数 ${\beta }_{1,2}$ 的改变对样条曲线形状的影响。黑、蓝、红色曲线分别对应 ${\beta }_{1,2}$ 的取值为1、0.5、0.1

6. 结论

给定若干个球面插值点数据，本文提出了一种构造几何连续的球面插值样条曲线的新方法。该方法是欧氏空间中Beta样条曲线的推广，属于整体构造方法，松弛了连续要求，引入了可以控制样条曲线形状的参数，从而使得构造的样条具有更加灵活的形状控制性。在后续的研究中，我们将会对该球面样条曲线的高阶几何连续性，以及定量分析形状参数对曲线形状的影响等问题进一步展开工作。

文章引用

杨 雷,王 倩,陈佳惠,孙 楠,何 耀. 球面Beta样条曲线相关问题研究
Research on the Spherical Beta Spline Curse[J]. 应用数学进展, 2021, 10(12): 4209-4217. https://doi.org/10.12677/AAM.2021.1012447

参考文献