B样条曲线曲面和NURBS曲线曲1

B样条曲线曲面和NURBS曲线曲面(学习笔记和上机练习)非均匀有理B样条，通常简称为NURBS(Non-UniformRationalB-Splines)。NURBS是非有理B样条、有理以及非有理Bézier曲线曲面的推广。一、要对B样条曲线曲面和NURBS曲线曲面有所了解应先了解B样条基函数B样条基函数的定义和性质令muuuU,...,,10是一个实数列，即iu≤1iu，i=0,1,…,m-1。其中，iu称为节点，U称为节点矢量，用)(,uNpi表示第i个p次（p+1阶）B样条基函数，其定义为,0,1)(0,uNi若iu≤u＜1iu值为1，其他值为0)()()(1,11111,,uNuuuuuNuuuuuNpiipipipiipiipi(2)由(2)式可知：（1）)(0,uNi是一个阶梯函数，它在半开区间),[1iiuuu外都为零；（2）当p＞0时，)(,uNpi是两个p-1次基函数的线性组合；（3）计算一组基函数时需要事先指定节点矢量U和次数p；（4）（2）式中可能出现0/0，我们规定0/0=0；（5）)(,uNpi是定义在整个实数轴上的分段多项式函数，只在区间][,0muu上有意义；（6）半开区间),[1iiuu称为第i个节点区间，长度可以为零，因相邻节点可以相同；B样条基函数的一些重要性质：1如果),[1piiuuu，则)(,uNpi=0。2对于所有的pi,和u，有)(,uNpi≥0.3对于任意的节点区间),[1iiuu，当),[1iiuuu时，ipijpjuN1)(,。4在节点区间内部，)(,uNpi是无限次可微的。在节点处)(,uNpi是kp次连续可微的，其中k是节点的重复度。5除0p的情况外，)(,uNpi严格的达到一次最大值。B样条基函数的导数基函数的求导公式为)()(1,1111,,uNuupuNuupNpiipipiipipi（4）根据求导法则gfgffg)(，对基函数)()()(1,11111,,uNuuuuuNuuuuuNpiipipipiipiipi求导得到一般的求导公式：)(11)1(1,1)1(1,)(,ipikpiipikpikpiuuNuuNpN（5）另一个计算B样条基函数各阶导数的公式（参考[Butt76]）：)()(1,1111)(1,)(,kpiipipikpiipiikpiNuuuuNuuuukppN，,1,0k…,1p(6)B样条基函数的一些重要计算a–确定节点(u)在节点矢量(U)的区间下标(据),[1iiuuu区间)//确定参数u所在的节点区间下标//n=m-p-1//m为节点矢量U[]的最大下标//p为B样条函数次数intBuspan(intn,intp,floatu,floatU[]){intlow,high,mid;if(u==U[n+1])//特殊情况returnn;//进行二分搜索low=p;high=n+1;mid=(int)(low+high)/2;while(uU[mid]||uU[mid]){if(uU[mid])high=mid;elselow=mid;mid=(int)(low+high)/2;if(u=U[mid]&&uU[mid+1])//找到u所在的节点区间的下标break;//退出二分搜索}returnmid;//返回区间下标}b–计算节点参数u的所有非零B样条基函数值(据定义(2)式)//计算所有非零B样条基函数并返回其值,floatN[]//i为参数u所在的节点区间下标voidAllBNipu(inti,intp,floatu,floatU[],floatN[]){floattm1,tm2,TN[100][100];floatleft[100][100],right[100][100];for(intti=0;ti=i+p;ti++)for(intj=0;j=p;j++){tm1=U[ti+j]-U[ti];tm2=U[ti+j+1]-U[ti+1];if(tm1!=0)left[ti][j]=(u-U[ti])/tm1;elseleft[ti][j]=0.0;if(tm2!=0)right[ti][j]=(U[ti+j+1]-u)/tm2;elseright[ti][j]=0.0;if(u=U[ti]&&uU[ti+1])TN[ti][0]=1.0;elseTN[ti][0]=0.0;}N[0]=TN[i][0];for(intt=i+p-1;t=i-p;t--)for(intk=1;kp;k++)TN[t][k]=left[t][k]*TN[t][k-1]+right[t][k]*TN[t+1][k-1];for(intr=i;r=i-p;r--)N[r]=left[r][p]*TN[r][p-1]+right[r][p]*TN[r+1][p-1];}c–计算B样条基函数的导数(据(6)式)//计算基函数的1阶导数并保存在NP[]中//i为参数u所在的节点区间下标//p为B样条函数次数P2voidBNipip(inti,intp,floatu,floatU[],floatNP[]){floattm1,tm2,TN[100][100];floatleft[100][100],right[100][100];for(intti=0;ti=i+p+1;ti++)for(intj=0;j=p;j++){tm1=U[ti+j]-U[ti];tm2=U[ti+j+1]-U[ti+1];if(tm1!=0)left[ti][j]=(u-U[ti])/tm1;elseleft[ti][j]=0.0;if(tm2!=0)right[ti][j]=(U[ti+j+1]-u)/tm2;elseright[ti][j]=0.0;if(u=U[ti]&&uU[ti+1])TN[ti][0]=1.0;elseTN[ti][0]=0.0;}for(intt=i+p-1;t=i-p;t--)for(intk=1;kp;k++)TN[t][k]=left[t][k]*TN[t][k-1]+right[t][k]*TN[t+1][k-1];for(intl=i;l=i-p;l--)NP[l]=p*((u-U[l])*TN[l][p-1]/(U[l+p]-U[l])-(U[l+p+1]-u)*TN[l+1][p-1]/(U[l+p+1]-U[l+1]))/(p-1);}二、B样条曲线曲面B样条曲线的定义p次B样条曲线的定义为niipiPuNuC0,)()(，a≤u≤b(1)其中iP是控制点，)(,uNpi是定义在非周期（并且非均匀）节点矢量bbuuaaUpmp,...,,...,,...11，(a,…a)和(b,…b)为(p+1)个节点，(包含m+1个节点)上的p次B样条基函数。一般情况下，a=0,b=1。如下定义一条3次B样条曲线：U[11]={0.0,0.0,0.0,0.0,0.25,0.5,0.75,1.0,1.0,1.0,1.0};下图为定义在上节点矢量的一条3次B样条曲线(P0,P1,…,P6为控制多边形)。P1P2P60PP3P5P4图1、用三次Bezier曲线形式表示分段四次多项式3次B样条曲线B样条曲线的性质（1）如果U[]={0.0,…,0.0,1.0,…,1.0};那么C(u)是Bezier曲线。如下图2所示P1P2P0P3图2、定义在U[8]={0.0,0.0,0.0,0.0,1.0,1.0,1.0,1.0};上的3次B样条曲线（2）C(u)是分段多项式曲线，因)(,uNpi是分段多项式函数。次数p、控制点个数n+1和节点个数m+1应满足以下条件m=n+p+1（3）端点插值性：C(0)=P0，C(1)=Pn。（4）仿射不变性：对B样条曲线进行仿射变换，得到的仍为B样条曲线。（5）强凸包性：曲线C(u)包含在它的控制多边形的凸包内（见图1和图2）。（6）局部修改性：移动Pi只改变C(u)在区间),[1piiuu上的形状。原因为),[1piiuuu，则)(,uNpi=0。如下图所示定义在U[10]={0.0,0.0,0.0,0.0,0.34,0.7,1.0,1.0,1.0,1.0};上的三次曲线P1P2P41P5P0P3P4图3、移动P4到P41（7）控制多边形是对B样条曲线的一个分段性逼近。（8）沿着曲线C(u)从u=0移动到u=1，基函数)(,uNpi只与对应的u起作用。（9）C(u)的连续性和可微性（10）变差减少性（11）重控制顶点的利用B样条曲线的导矢niikpikPuNuC0)(,)()()(（7）其中)()(uCk表示)(uC的k阶导矢。B样条曲线在端点处的一阶导矢：)()0(011PPupCp，)(1)1(11nnpmPPupC（8）B样条曲线在端点处的二阶导矢：2221121101)2()()1()0(pppppppuPuuPuuuPuppC，222112111)2(1)1)(1()2(11)1()1(pmnpmpmnpmpmpmnpmuPuuPuuuPuppC（9）对于Bezier曲线上面（8）、（9）式可简化为下面对应的表达式：)()0(01PPnC)()1(1nnPPnC)2)(1()0(012PPPnnC)2)(1()1(21nnnPPPnnCB样条曲线的一些计算//计算B样条曲线上的点并返回其值//n=m-p-1//m为节点矢量U[]的最大下标//p为B样条函数次数floatBCurveP(intn,intp,floatU[],floatP[],floatu){intspan;floatN[100],Cp;span=Buspan(n,p,u,U);AllBNipu(span,p,u,U,N);Cp=0.0;for(inti=span-1;i=span-p;i--)Cp+=N[i]*P[i];returnCp;}//计算曲线上的1阶导矢并返回1阶的导矢//n为参数u所在的区间下标//p为次数floatBCurveCN(intn,intp,floatu,floatU[],floatP[]){floatCN,NP[100];BNipip(n,p,u,U,NP);CN=0.0;for(inti=n-1;i=n-p;i--)CN+=NP[i]*P[i];returnCN;}B样条曲面的定义B样条曲面由两个方向的控制点网格、两个节点矢量和单变量B样条基函数的乘积来定义，方程如下nimjjiqjpiPvNuNvuS00,,,)()(),(（10）节点矢量为}1,...,1,,...,0,...,0{1pruU和}1,...,1,,...,0,...,0{1qsvVU中含有r+1个节点，V中有s+1个节点，据B样条曲线性质（2）有r=n+p+1,s=m+q+1B样条曲面的性质（1）曲面插值于控制网格的四个角点：S(0,0)=P0,0,S(1,0)=Pn,0,S(0,1)=P0,m,S(1,1)=Pn,m得到1)1()1()1()0()0()1()0()0(,,,,0,0,,0,0qmpnqmpqpnqpNNNNNNNN（2）仿射不变性（3）局部修改性（4）凸包性B样条曲面的偏导矢jilqjnimjkpilklkPvNuNvuSvu,)(,00)(,)()(),(（11）角点偏导矢公式（u,v）=(0,0):)()0,0(0,00,11PPupSpu)()0,0(0,01,01PPvqSqv)()0,0(0,00,1