基于PSO的B样条曲线光顺重构算法

刘武飞 张旭

关键词：B样条;曲线光顺;曲率变化;粒子群优化算法

摘要：针对工程实践中存在的曲线重构技术很难同时考虑曲线误差和曲线光顺性的问题，提出了一种基于粒子群优化（PSO）的B样条曲线光顺重构算法.该算法利用PSO算法同时调整影响曲率坏点、坏区，以及最坏点处的主、副等多个控制顶点，找出控制点位置的最优解，优先对曲线上曲率符号不一致的坏点或坏区进行光顺，以避免曲线上出现多余拐点，而后对曲率变化剧烈的区域进行光顺，迭代更新生成最优曲线.实验结果表明，该算法有效地提升了光顺效率，得到了更好的光顺效果，且能够满足任意给定的误差精度，验证了其应用于工程实践的可行性.

Abstract：In order to solve the problem that the curve reconstruction technology existing in engineering practice was difficult to consider both curve error and curve smoothness， a fairing reconstruction algorithm of B-spline curve was proposed based on particle swarm optimization （PSO）.The algorithm used PSO to adjust the main and secondary control vertices which affect the curvature bad points， bad region andthe main control and subcontrol points of the worst points at the same time， so as to find the optimal solution for the position of the control points.The bad points or bad region with inconsistent curvature symbols on the curve were faired first to avoid the occurrence of redundant inflection points on the curve.Then the region with sharp curvature changes were faired and the optimal curve was generated iteratively.The experimental results showed that the algorithm effectively improved the fairing efficiency， obtained better fairing effect， and could satisfy any given error accuracy， which proved its feasibility in engineering practice.

0 引言

曲線重构工作经常出现在计算机辅助设计（CAD）、计算机辅助制造（CAM）等领域，工程技术人员一般通过插值的方法进行曲线重构得到高精度的曲线.但在数据获取阶段，人工或设备因素会导致原始数据点掺杂了一些误差和噪声，使得通过插值方法构造的高精度曲线难有较好的光顺性，无法满足目前航空、造船、模具等工业领域对产品外观、产品功能等方面越来越高的光顺要求，因此国内外业内技术人员对曲线光顺方法展开了深入研究.

龙小平[1]根据型值点信息筛选坏型值点，遵循最优化原则，通过每次调整1个控制点来调整局部的控制顶点，得到局部能量最小的光顺曲线.局部光顺算法虽然计算速度较快，但很难保证光顺后曲线的整体效果.因此，罗卫兰等[2]通过增加权重规定控制顶点的扰动量，利用最小化原则求出新的控制顶点，进而对曲线进行光顺;张莉等[3]建立了广义B样条曲线的对偶基，用最佳逼近的方法求出新的曲线控制顶点，以实现曲线光顺.全局光顺算法不需要选点，但当数据点过于复杂时，运算速度较慢，且优化方程可能无解.为了提高光顺效率和重构质量，研究人员引进了小波分析，A.Ceruti等[4]基于二进小波光顺，把多分辨分析应用在A级曲线曲面上，对坏控制顶点进行调整，得到了较好的光顺效果.A.Z.Wang等[5-6]提出了光顺性指标的概念，发现经典的二进小波光顺算法的使用会受控制顶点数目的约束.为满足任意控制顶点数目曲线曲面的光顺要求，潘洋宇等[7]通过增加控制顶点的数目使其得以实现，给出了基于第2代小波变换的多分辨率光顺算法.R.Pan等[8]通过节点插入构造出双正交非均匀B样条小波.纪小刚等[9]提出了一种新的能在任何二进尺度下进行连续小波变换的光顺算法.该类算法生成的曲线光顺性较好，但随着光顺的进行，控制顶点的减少会导致曲线产生很大的误差.

随着人工智能算法的兴起，基于人工智能的曲线光顺研究日益受到重视.E.Ulker等[10]采用人工免疫系统来优化节点矢量，对曲线进行光顺，使曲线有较好的光顺性;X.Y.Zhao等[11]基于改进的k-means算法，建立节点的高斯混合模型，再根据其高斯概率求解新节点位置光顺曲线，但该方法只适用于闭曲线;A.Galvez等[12]提出粒子群光顺算法，该算法对于具有奇点、尖点的曲线也能生成较好的结果;H.M.Kang等[13]采用稀疏优化的模型求解节点向量，这会产生较少数量的节点，使曲线拥有较好的光顺性;胡良臣等[14]使用带有法向约束的粒子群优化（PSO）算法迭代求解曲线的节点向量.这类算法由于参数的设定不同，使得算法有可能收敛于局部最优解，导致曲线的光顺效果不理想.但由于人工智能算法对复杂、非线性、多目标问题可以获得全局最优解，且不同的智能算法还有各自独到的优点，所以基于人工智能算法的曲线光顺重构，依旧是业界当前的研究热点.

在实际建模过程中，设计人员根据点云数据来生成曲线，再逐一调整单个控制点，以改变曲线与点云之间的偏差，使得生成的曲线与点云贴合;画出曲线的曲率梳，再逐一调整单个控制点，以控制曲线曲率的变化波动，使得曲率变化平缓，曲线光顺.但曲线诸多控制点的位置变动又会导致曲线误差发生改变，需再进行曲线误差的调整，确保曲线误差在允许范围内，如此循环往复，需要不停地反复操作，直到在给定的误差范围内，得到一条曲率变化平缓、光顺性较好的曲线.但当面对数据量较大且包含信息复杂的点云数据时，这种人工操作需要很长时间，会耗费大量的人力物力，且不能保证得到最好的结果.

针对现有算法对曲线曲面光顺效果不够理想和实际建模过程中手动调整曲线曲面光顺会耗费大量人力物力的问题，本文拟以目前普遍使用的3次B样条曲线为研究对象，提出在给定误差约束下、基于PSO的B样条曲线光顺重构算法，以期提高B样条曲线的光顺效果，有效控制光顺后生成曲线的误差，以满足工程需求.

1 曲线光顺性影响因素分析

1.1 曲线的表示

根据式①和B样条曲线控制多边形的性质，通过直接修改B样条曲线控制顶点位置进行曲线光顺，可使光顺过程更加直接.3次B样条曲线上每个点的值由4个控制顶点和4个B样条基函数决定，其中最大的B样条基函数值所对应的控制点，称之为主控制顶点.将与主控制顶点相邻的两个控制顶点称为副控制顶点，当主控制顶点为4个控制顶点中的第一个或最后一个时，只有1个副控制顶点.

基于通用的光顺准则[15-16]，本文通过观察曲线曲率变化来评价曲线的光顺性.如果曲线光顺性比较好，那么曲线的曲率梳应是连续的，且曲率梳由尽可能少的单调曲线组成[17-18].3次B样条曲线的曲率可由下式求得：

1.2 曲率坏点调整

越简单的曲线光顺性越好，复杂的曲线可以拆分成简单的曲线进行光顺.

光顺性良好的曲线上每点的曲率符号应该保持一致，与相邻点曲率符号不一致的点可通过下式找出：

这些点称为曲率坏点，与其下标对应的曲线上的投影点称为曲线坏点.将曲线坏点提取出来，根据下标由小到大重新排列组成序列{Qi1，Qi2，…，Qis}，其中s是坏点的总个数，给出与之对应的主控制顶点的序列对主控制顶点序列中具有相同下标的控制顶点进行处理，使其只出现1次，进而构造新的主控制顶点序列本文通过调整主控制顶点序列中出现的控制顶点来处理曲线坏点.

1.3 曲率坏区调整

对提取的曲线坏点进行曲率符号调整后，曲线的曲率符号可能仍不一致.根据式②得到的曲线曲率符号一般分为两种情况：一种是曲线曲率经过一次变号，小部分曲线上点的曲率符号与其他曲线上点的曲率符号相反，呈现“正负”或“负正”的趋势;另一種是曲线的曲率经过两次变号，呈现“正负正”或“负正负”的趋势.曲率符号多变的曲线可以

简单

划分为这两种情况进行处理.

统计曲率符号为正的投影点个数，若数量多于曲线符号为负的投影点数量，就将曲率符号为正的区域称为曲率坏区;反之，就将曲率符号为负的区域称为曲率坏区.

与曲率坏区对应曲线上的点构成的区域是曲线坏区.提取出给定数据点在曲线坏区上的所有投影点，根据其下标，由小到大地重新排列组成序列{Qj1，Qj2，…，Qjl}，其中l是曲线坏区包含投影点的总个数，推出与其下标对应的主控制顶点的序列{Pi1，Pi2，…，Pil}.对具有相同下标的控制顶点进行同样处理，得到新的主控制顶点序列{Pi1，…，Pim}，m≤l，对其中的控制顶点位置进行优化调整以使曲线各点曲率符号相同.

1.4 曲率变化调整

光顺准则中要求曲率变化较均匀，设曲线上第i个内部点左右曲率变化率差为从上式可知，Δki小，表示该点左右曲率变化率相差不大，也表示该处曲线曲率变化平缓;当Δki=0时，表示曲线上该点和相邻两点（第i-1，i，i+1个数据点处）的曲率变化均匀，曲线光顺.以Δki最大值作为光顺性指标会使曲率高的区域光顺次数较多，曲率低的区域光顺次数较少.为弥补这种不足，选用Δdi作为光顺性指标：

其中，分母是曲线上该点和相邻两点曲率绝对值中的最大值.计算所有内部点的Δdi，将Δdi值最大的点作为待光顺的最坏点.

曲线上某点处的B样条基函数值的大小，决定了对应控制顶点对该点的影响程度.一般是选择移动影响程度最大的控制顶点来光顺曲线，即通过调整最坏点处主控制顶点的位置来光顺曲线，使得该点处曲率变化比之前平缓.然而，在实际光顺过程中，往往会出现一些特殊的情况，譬如最坏点处的B样条基函数的最大值与次大值之间相差不大，这就表示与之对应的控制顶点对最坏点的影响程度相当，在这种情况下，采取只移动主控制顶点的光顺策略可能得不到理想的光顺结果，光顺效率也不高;又因为B样条基函数中的最小值都远小于最大值，所以与之对应的控制顶点对最坏点的影响程度远小于主控制顶点，导致调整全部控制顶点的光顺策略也可能得不到理想的光顺结果.因此，本文同时调整主控制顶点和副控制顶点位置对最坏点进行光顺，这样同时调整两个或者三个控制顶点的位置可以大大提升曲线光顺的效率，缩短光顺时间.

2 光顺重构算法设计

针对上文分析得出的导致重构曲线光顺性差的3个原因（曲率坏点;曲率坏区;曲率变化剧烈），分别通过式③和④，以及曲率符号变化情况，找出曲线光顺性差的点或区域，在光顺过程中不断调整影响该点或区域的多个控制顶点位置，找出控制点位置最优解，完成曲线光顺.由于PSO算法采用实数编码，相对于其他智能算法来说，PSO算法更易于操作和实现，且算法中的参数可以根据经验值来设定，不需要过多的人为干预.因此，本文基于PSO进行光顺算法设计.

2.1 PSO算法

PSO算法[19]，简单来说，就是模拟一群鸟在一片区域内找寻唯一一块食物的过程，给出了每只鸟到食物的距离，搜索与食物距离最近的那只鸟附近的区域，即为这群鸟找到食物的最佳办法.该算法将优化问题可能的解都表示成搜索空间中不计重量和体积、具有一定飞行速度的微粒，相当于区域中飞行的鸟，优化问题的最优解就是区域中唯一的那块食物，由目标函数计算得到的适应度值就是鸟到食物的距离.

PSO算法流程[20]如下：

步骤1 初始化种群，设定初始位置和速度.

步骤2 计算每个微粒的适应度，并初始化pbest和gbest.

步骤3 按照式⑤和式⑥更新微粒的速度和位置.

步骤4 计算更新后各微粒的适应度，更新pbest和gbest.

步骤5 如果满足终止条件（最大迭代次数或满意的适应度），则结束，输出当前最优解;否则，返回步骤3.

计算各点离曲线的最近距离，再取其最大值作为曲线的误差，考虑到光顺进程会增大曲线的误差，可以将给定的误差精度作为PSO算法的约束条件，限定种群中微粒的飞行空间，能确保光顺后输出的曲线满足误差精度要求，而误差精度可以根据设计的需求进行设定.

2.2 光顺性量化标准

通过曲率图只能得到曲线光顺性的定性评价，为了更好地进行数据分析，还需要对曲线光顺性进行量化评价.通常是将能量函数的值作为评判光顺性的标准，但是这种方法不一定能对曲线整体的光顺性给出准确的评价，例如半径相同但圆心角不同的两条弧，其光顺程度相当，但圆心角大的弧能量更大.因此，本文基于光顺性指标Δdi量化曲线的光顺性，设光顺性度量为.

其中，ΔD表示曲线上所有内部数据点Δdi之和，将ΔD作为PSO算法的适应度函数进行种群的更新.随着曲线光顺次数的增多，ΔD一般会减少直至不再变化，中间偶尔会出现增大再减少的情况.为避免过早收敛，现在一般的做法是连续记录4次曲线光顺后的ΔD值，当其中最小值与最大值的比值大于某个固定值时停止光顺进程.在本文算法中，根据参考文献[2]和经验，设定固定值为0.995.

2.3 光顺重构算法流程

针对前述影响曲线光顺性的3个原因，分别给出对应的曲率光顺方法：曲率坏点光顺算法、曲率坏区光顺算法和曲率变化光顺算法.考虑实现这3种光顺方法的先后顺序对曲线完成光顺进程效率和质量的影响，本文光顺重构算法的流程如图1所示.

光顺重构算法具体步骤如下.

步骤1 对于给定数据点，首先利用最小二乘法重构得到一条3次B样条曲线，及其给出的一系列控制顶点Pi，计算曲线上各点的曲率k.

步骤2 通过式③找出曲线上的所有坏点，构造相应的主控制顶点序列{Pj1，Pj2，…，Pjn}.

步骤3 将主控制顶点序列{Pj1，Pj2，…，Pjn}作为PSO算法的初始种群，通过迭代重新计算曲线的误差、曲率和曲线的光顺性度量ΔD.当满足给定的误差且曲线无坏点时，输出曲线新的控制顶点序列和曲线曲率k，进入下一步，否则重复步骤3.

步骤4 确定曲线上的坏区，得到主控制顶点序列{Pi1，Pi2，…，Pim}

步骤5 将主控制顶点序列{Pi1，

作为PSO算法的初始种群，通过迭代重新计算曲线的误差、曲率和曲线的光顺性度量ΔD.当满足给定的误差和曲线无坏区（曲线曲率符号相同），输出曲线新的控制顶点序列曲线曲率k，否则重复步骤5.

步骤6 通过式④确定曲线上的最坏点，得到需要调整的控制顶点序列 {Pi-1，Pi，Pi+1}.

步骤7 将控制顶点序列{Pi-1，Pi，Pi+1}作为PSO算法的初始种群，通过迭代重新计算曲线的误差、曲率和曲线的光顺性度量ΔD.

步骤8 若满足PSO算法终止条件，且不再满足给定的误差精度或达到最大迭代次数时，输出新的控制顶点，记录ΔD值，ΔD值不足4个时，返回步骤6，不少于4个时，进入下一步;否则返回步骤7.

步骤9 检查光顺进程的终止条件，当连续的4个ΔD值中最小值与最大值的比值大于0.995时，停止光顺进程，输出控制顶点构造曲线;否则返回步骤6.

3 实验结果与分析

3.1 相同误差精度光顺结果与分析

为验证本文光顺重构算法的性能，在给定相同误差精度（E=0.15 mm）的条件下，基于Matlab（2016b）软件，使用传统选点光顺算法与本文算法，对给定的数据点分别进行10次光顺重构，对比这两种算法输出曲线的曲率梳.对某模型底部轮廓线的部分数据点（431个）进行拟合，得到的原始曲线及其原始曲率如图2所示.

从图2b）可以看出，原始曲线存在两处曲率波动较剧烈的区域，需要进行光顺处理.另外，在原始曲率的后半部分還存在曲率符号不一致的区域（图中圈出部分），这表示原始曲线上存在多余拐点，需要先进行曲率符号调整.本文提出的曲率坏点、坏区光顺算法光顺前后的曲率对比见图3.由图3可以看出，光顺之前原始曲率上的两处坏区于光顺之后都得到了调整，现在整条曲线的曲率符号保持一致，即曲线上不存在拐点，且曲线的曲率波动也较之前有所减少.

在同等误差精度（E=0.15 mm）的约束下，传统选点光顺算法和本文光顺重构算法的比较结果见图4.由图4可以看出，在原始曲率剧烈波动的区域，本文算法的变化波动较传统选点光顺算法要小得多，且本文算法生成的曲线要更光顺一些.将图4中曲线尾部的第二处曲率坏区（图中圈出部分）放大得到图5.由图5可以看出，经本文算法光顺后的曲线不存在曲率坏区，而采用传统选点光顺法光顺后的曲线上仍存在多余的拐点，即存有曲率坏区，这进一步说明本文算法的光顺效果要优于传统算法.

曲线内部数据点Δd值的变化情况可以反映曲率变化情况的好坏：Δd值变化越大，表示曲线的曲率变化越剧烈，曲线越不光顺.选择原始曲率剧烈波动的第30—55个数据点，本文算法、原始曲率和传统算法的内部数据点Δd值见图6.由图6可以看出，本文算法和传统算法光顺后的曲线曲率变化率都较原曲线大幅下降，且本文算法的曲率变化更平缓.表1为传统算法与本文算法的光顺性评价数据.由表1可知，本文算法的光顺性度量ΔD值最小，表明曲率变化最平缓，虽然传统算法与本文算法基于能量标准的曲率积分和十分相近，但是本文算法完成整个曲线光顺过程所需次数要少于传统算法.综上说明，在同等误差精度约束下，本文算法比传统算法光顺效果更好，即生成的曲线更光顺，且能显著提高曲线的光顺效率.

3.2 不同误差精度光顺结果与分析

在不同误差精度的条件下，验证本文光顺重构算法的光顺效果.图7是某叶片横截面轮廓原始曲线及其原始曲率.在误差精度分别为E=0.20 mm，E=0.15 mm，E=0.10 mm的条件下，采用本文算法对原始曲线进行18次光顺，结果如图8所示.本文算法在不同误差精度下的光顺性评价见表2.

由图8可见，在允许的误差范围内，经本文算法光顺后的曲线都较原曲线更光顺，光顺后曲线的曲率变化也都较原曲线的曲率更平缓，且给定的误差精度值越大，光顺效果越好，输出的曲线越光顺.这是因为给定的误差精度作为PSO算法中种群更新的约束条件，限制了种群中粒子的飞行空间，从而有效控制了算法输出曲线的误差.由表2可知，经本文算法光顺后的曲线误差均处于允许的误差范围内，说明本文算法实现了对光顺后生成的曲线误差的有效控制;较原始曲线，本文算法的光顺性度量ΔD值和基于能量标准的曲率积分和更优.综上可知，本文算法在不同误差精度条件下都能获得更好的光顺效果，且可有效控制光顺后生成曲线的曲线误差.

4 结语

本文提出了一种基于PSO的B样条曲线光顺重构算法.该算法提出了利用PSO算法同时调整多个控制点的策略，首先对曲线上曲率符号不一致的点和区域进行了预处理，使得曲线上点的曲率符号保持一致，以避免出现冗余的拐点.而后对曲率变化剧烈区域进行光顺，迭代更新生成最优曲线.实验结果表明，本文算法拥有较好的光顺效果，提升了光顺效率，缩短了光顺进程，能很好地满足工程上对曲线光顺性的要求，且算法对生成曲线误差的控制在工程上也将极大地降低建模工作所需耗费的人力物力.本文算法针对的是平面B样条曲线的光顺重构，而将其推广到空间曲线或曲面的光顺重构，是未来工作的一部分.

参考文献：

[1] 龙小平.局部能量最优法与曲线曲面的光顺[J].计算机辅助设计与图形学学报， 2002， 14（12）：1109.

[2] 罗卫兰，杨勋年，郑建民.B样条曲线的约束光顺算法[J].浙江大学学报（理学版），2004，31（1）：51.

[3] 张莉，葛先玉，檀结庆.广义B样条曲线的节点去除与光顺算法[J].计算机辅助设计与图形学报，2016，28（4）：540.

[4] CERUTI A， LIVERANI A， CALIGIANA G.Fairing with neighbourhood LOD filtering to upgrade interactively B-Spline into Class-A curve[J].IJIDeM： International Journal on Interactive Design and Manufacturing， 2014， 8（2）：67.

[5] WANG A Z， ZHAO G， LI Y D.Fairness degree based fairness criterion and fairing algorithm[J].Applied Mathematics and Computation，2015，253：184.

[6] 王爱增，赵罡，穆国旺.基于数字化光顺性指标的NURBS曲线自适应光顺[J].计算机学报，2011，34（8）：1548.

[7] 潘洋宇，姜福祥.任意控制点曲线小波光顺方法研究[J].机械设计，2009，26（11）：12.

[8] PAN R， YAO Z.Biorthogonal nonuniform B-spline wavelets based on a discrete norm[J].Com-puter Aided Geometric Design， 2009， 26（4）：480.

[9] 纪小刚，杨艳，薛杰.基于多分辨技术的任意控制顶点曲面光顺[J].机械工程学报，2015，51（11）：159.

[10]ULKER E， ARSLAN A.Automatic knot adjustment using an artificial immune system for B-spline curve approximation[J].Information Sciences：An International Journal， 2009， 179（10）：1483.

[11]ZHAO X Y，YANG B，ZHANG C M，et al.Adaptive knot placement using a GMM-based continuous optimization algorithm in B-spline curve approximation[J].Computer-Aided Design， 2011， 43（6）：598.

[12]GALVEZ A， IGLESIAS A.Efficient particle swarm optimization approach for data fitting with free knot B-splines[J].Computer-Aided Design， 2011， 43（12）：1683.

[13]KANG H M， CHEN F L，LI Y S， et al.Knot calculation for spline fitting via sparse optimization[J].Computer-Aided Design， 2015， 58：179.

[14]胡良臣，壽华好.PSO求解带法向约束的B样条曲线逼近问题[J].计算机辅助设计与图形学学报，2016，28（9）：1443.

[15]苏步青，刘鼎元.计算几何[M].上海：上海科学技术出版社，1981：297.

[16]王士玮，刘利刚，张举勇，等.基于稀疏模型的曲线光顺算法[J].计算机辅助设计与图形学学报，2016，28（12）：2043.

[17]章虎冬.基于局部能量的三次B样条曲线自动光顺算法[J].西安航空学院学报，2016，34（1）：79.

[18]王爱增，何川，赵罡，等.基于几何方法的曲率单调Bézier曲线的一个充分必要准则[J].计算机辅助设计与图形学学报，2019，31（9）：1617.

[19]郑晓月.用快速收敛粒子群优化算法解决函数优化问题[J].轻工学报，2016，31（3）：89.

[20]李巧燕，全海燕.基于改进粒子群的独立分量分析算法研究[J].轻工学报，2016，31（2）：103.