非均匀有理B样条.docx
- 文档编号:27054476
- 上传时间:2023-06-26
- 格式:DOCX
- 页数:36
- 大小:103.16KB
非均匀有理B样条.docx
《非均匀有理B样条.docx》由会员分享,可在线阅读,更多相关《非均匀有理B样条.docx(36页珍藏版)》请在冰豆网上搜索。
非均匀有理B样条
非均匀有理B样条
(学习记录和上机练习)
非均匀有理B样条,通常简称为NURBS(Non-UniformRationalB-Splines)。
NURBS是非有理B样条、有理以及非有理Bezier曲线曲面的推广。
一、Bezier曲线
1、Bezier曲线
一条n次Bezier曲线可以表示为
0≤
≤1
(1)
其中,基函数(也称为混合函数)
是著名的n次Bernstein多项式,其定义为
(2)
(1)式中的几何系数
称为控制点。
//计算参数u的n次基函数值并存在B[]中
voidBernstein(intn,doubleu,floatB[])
{
doubleu1;
intj,k;
B[0]=1.0;
u1=1.0-u;
for(j=1;j<=n;j++)
{
floatsaved=0.0;
for(k=0;k { floattemp=B[k]; B[k]=saved+u1*temp; saved=u*temp; } B[j]=saved; } } 2、有理Bezier曲线 一条n次有理Bezier曲线的定义为 0≤ ≤1(3) 其中 是标量,称为权因子。 图1、3次Bezier曲线图2、3次有理Bezier曲线 图3、n=3Bernstein图形图4、n=9Bernstein图形 图5、Bezier双三次曲面片图6、有理Bezier双三次曲面片 二、B样条基函数 1、B样条基函数的定义和性质 令 是一个实数序列,即 ≤ 。 其中, 称为节点, 称为节点矢量,用 表示第i个p次(p+1阶)B样条基函数,其定义为: (4) 由此可知: (1) 是一个阶梯函数,它在半开区间 外都为零; (2)当p>0时, 是两个p-1次基函数的线性组合; (3)计算一组基函数时需要事先指定节点矢量 和次数p; (4)上式中可能出现0/0,我们规定0/0=0; (5) 是定义在整个实数轴上的分段多项式函数,它在区间 上有意义; (6)半开区间 称为第i个节点区间,它的长度可以为零,因为相邻节点可以是相同的; (7)计算p次基函数的过程生成一个三角阵列。 性质: 1局部支撑性,如果 ,则 =0。 2非负性,对于所有的i,p和u,有 ≥0。 3规范性,对于任意的节点区间 ,当 时 。 4可微性,在节点区间内部, 是无限次可微的(在每个节点区间内部,它是一个多项式),在节点处 是p-k次连续可微的,其中k是节点的重复度。 2、B样条基函数的导数 1基函数的求导公式为: (5) 对基函数 求导得到一般的求导公式: (6) 为了完整性,另一个计算B样条基函数各阶导数的公式(参考[Butt76]): (7) //------------------------------------------------------------------------- //计算所有非零B样条基函数并返回其值 //i为参数u所在的节点区间下标 voidBasisFunction(inti,intp,floatu,floatU[],floatN[]) { intj,di,dp,k; floattul,tur,left,right; floattmpN[50][50]; for(k=0;k<=p;k++) { dp=0; for(di=i+p-k;di>=i-k;di--) { if(u>=U[di]&&u tmpN[di][0]=1; else tmpN[di][0]=0; dp+=1; for(j=1;j { tul=U[di+j]-U[di]; tur=U[di+j+1]-U[di+1]; if(tul! =0) left=(u-U[di])/tul; else left=0; if(tur! =0) right=(U[di+j+1]-u)/tur; else right=0; tmpN[di][j]=left*tmpN[di][j-1]+right*tmpN[di+1][j-1]; } } N[i-k]=tmpN[i-k][p]; } } //----------------------------------------------------------------------- //计算基函数的1阶导数并保存在NP[]中 //i为参数u所在的节点区间下标 //p为B样条函数次数P>2 voidDerBasisFunc(inti,intp,floatu,floatU[],floatNP[]) { intj,di,dp,k; floattul,tur,left,right,saved,dl,dr; floattmpN[50][50]; for(k=0;k<=p;k++) { dp=0; for(di=i+p-k;di>=i-k;di--) { if(u>=U[di]&&u tmpN[di][0]=1; else tmpN[di][0]=0; dp+=1; for(j=1;j { tul=U[di+j]-U[di]; tur=U[di+j+1]-U[di+1]; if(tul! =0) left=(u-U[di])/tul,dl=1/tul; else left=0,dl=0; if(tur! =0) right=(U[di+j+1]-u)/tur,dr=1/tur; else right=0,dr=0; tmpN[di][j]=left*tmpN[di][j-1]+right*tmpN[di+1][j-1]; saved=p*(dl*tmpN[di][j-1]-dr*tmpN[di+1][j-1])/(p+1); } } NP[i-k]=saved; } } 图7、三次B样条曲线 三、B样条曲线曲面 1B样条曲线的定义 P次B样条曲线为: ≤ ≤ (8) 这里 是控制点, 是定义在非周期节点矢量上的p次B样条基函数。 除非另外声明,我们假定a=0,b=1,并且a和b的个数为(p+1)个。 2性质: (1)如果n=p,U={0,…,0,1,…,1},那么C(u)是Bezier曲线(如下图所示)。 P1P3 P0P2 图8、定义在U[]={0,0,0,0,1,1,1,1};上的三次B样条曲线是三次Bezier曲线 (2)C(u)是分段多项式曲线(因为 是分段多项式函数);次数p,控制点个数n+1和节点个数m+1满足关系: m=n+p+1。 (3)端点插值性: C(0)=P0,C (1)=Pn。 (4)仿射不变性: 对B样条曲线进行仿射变换,所得曲线不变。 (5)强凸包性: 曲线C(u)包含在它的控制多边形的凸包内。 (6)局部修改性: 移动 只改变C(u)在区间 上的形状。 这是因为对于如果 ,则 =0。 P6 P_3 P3 P0 图9、移动P3(到P_3)局部修改曲线 (7)控制多边形是对B样条曲线的一个分段线性逼近。 (8)基函数 只与对应的区间有意义。 (9)C(u)在节点区间内部是无限次可微的。 (10)变差减少性: 任何一个平面与曲线的交点个数不多于它和控制多边形的交点个数。 (11)利用重控制顶点是可能的。 3B样条曲线的导矢 (9) //计算样条曲线的1阶导矢(u所对应的所有点)保存在Der[]中 //n=m-p-1 //p为曲线的次数 voidBSplineDer(intn,intp,floatU[],floatP[],floatDer[]) { floatN[100],tmp; inti,j; for(i=p+1;i<=n;i++) { DerBasisFunc(i,p,U[i],U,N); tmp=0; for(j=i;j>=i-p;j--) tmp+=N[j]*P[j]; Der[i-p]=tmp; } } //计算曲线上的点(u所对应的所有点)保存在Poi[]中 //n=m-p-1 //p为曲线的次数 voidBSplinePoint(intn,intp,floatU[],floatP[],floatPoi[]) { floatN[100],tmp; inti,j; for(i=p+1;i<=n;i++) { BasisFunction(i,p,U[i],U,N); tmp=0; for(j=i;j>=i-p;j--) tmp+=N[j]*P[j]; Poi[i-p]=tmp; } } 4B样条曲面的定义 B样条曲面由两个方向的控制点网格,两个节点矢量和单变量B样条基函数的乘积来定义,其方程为 (10) 节点矢量为 和 U中含有r+1个节点,V中含有s+1个节点(r=n+p+1,s=m+q+1)。 图10、双三次B样条曲面 5上机练习 下面用双二次B样条曲面生成一个简单飞机模型为例说明曲面生成的过程。 图11、双二次B样条曲面生成的飞机模型实物图 图12、双二次B样条曲面生成的飞机模型网格图 程序源如下: #include"glut.h" #include"math.h" //飞机机身头部数据 floathx[]={-360,-360,-360,-360,-360,-360,-360, -350,-350,-350,-350,-350,-350,-350, -300,-300,-300,-300,-300,-300,-300, -250,-250,-250,-250,-250,-250,-250, -200,-200,-200,-200,-200,-200,-200, -50,-50,-50,-50,-50,-50,-50, 0,0,0,0,0,0,0}; floathy[]={0,0,0,0,0,0,0, 0,20,20,0,-20,-20,0, 0,40,40,0,-40,-40,0, 0,60,60,0,-60,-60,0, 0,120,120,0,-60,-60,0, 0,110,110,0,-60,-60,0, 0,100,100,0,-60,-60,0}; floathz[]={0,0,0,0,0,0,0, -20,-20,20,20,20,-20,-20, -40,-40,40,40,40,-40,-40, -60,-60,60,60,60,-60,-60, -60,-60,60,60,60,-60,-60, -60,-60,60,60,60,-60,-60, -60,-60,60,60,60,-60,-60}; floatHU[]={0,0,0,0.2,0.4,0.6,0.8,1,1,1}; //飞机机身尾部数据 floatlx[]={0,0,0,0,0,0,0, 150,150,150,150,150,150,150, 200,200,200,200,200,200,200, 250,250,250,250,250,250,250, 300,300,300,300,300,300,300, 350,350,350,350,350,350,350, 400,400,400,400,400,400,400}; floatly[]={0,100,100,0,-60,-60,0, 0,90,90,0,-60,-60,0, 0,80,80,0,-60,-60,0, 0,70,70,0,-60,-60,0, 0,60,60,0,-60,-60,0, 0,60,60,0,-60,-60,0, 0,40,40,0,-40,-40,0}; floatlz[]={-60,-60,60,60,60,-60,-60, -60,-60,60,60,60,-60,-60, -60,-60,60,60,60,-60,-60, -60,-60,60,60,60,-60,-60, -60,-60,60,60,60,-60,-60, -60,-60,60,60,60,-60,-60, -40,-40,40,40,40,-40,-40}; floatLU[]={0,0,0,0.2,0.4,0.6,0.8,1,1,1}; //前机翼数据 floatrx[]={25,-50,-50,100,100,25, 40,-30,-30,110,110,40, 55,-10,-10,120,120,55, 60,10,10,120,120,60, 80,30,30,130,130,80, 95,95,95,95,95,95}; floatry[]={10,10,-10,-10,10,10, 10,10,-10,-10,10,10, 10,10,-10,-10,10,10, 10,10,-10,-10,10,10, 10,10,-10,-10,10,10, 0,0,0,0,0,0}; floatrz[]={-60,-60,-60,-60,-60,-60, -120,-120,-120,-120,-120,-120, -180,-180,-180,-180,-180,-180, -240,-240,-240,-240,-240,-240, -300,-300,-300,-300,-300,-300, -310,-310,-310,-310,-310,-310}; floatrz1[]={60,60,60,60,60,60, 120,120,120,120,120,120, 180,180,180,180,180,180, 240,240,240,240,240,240, 300,300,300,300,300,300, 310,310,310,310,310,310}; floatRU[]={0,0,0,0.25,0.50,0.75,1,1,1}; //后机翼数据 floatex[]={275,200,200,350,350,275, 290,220,220,360,360,290, 305,240,240,370,370,305, 320,260,260,380,380,320, 335,280,280,390,390,335, 350,350,350,350,350,350}; floatey[]={10,10,-10,-10,10,10, 10,10,-10,-10,10,10, 10,10,-10,-10,10,10, 10,10,-10,-10,10,10, 10,10,-10,-10,10,10, 0,0,0,0,0,0}; floatez[]={-60,-60,-60,-60,-60,-60, -90,-90,-90,-90,-90,-90, -120,-120,-120,-120,-120,-120, -150,-150,-150,-150,-150,-150, -180,-180,-180,-180,-180,-180, -190,-190,-190,-190,-190,-190}; floatez1[]={60,60,60,60,60,60, 90,90,90,90,90,90, 120,120,120,120,120,120, 150,150,150,150,150,150, 180,180,180,180,180,180, 190,190,190,190,190,190}; floatEU[]={0,0,0,0.25,0.50,0.75,1,1,1}; //后上机翼数据 floattx[]={265,180,180,350,350,265, 280,200,200,360,360,280, 295,220,220,370,370,295, 310,240,240,380,380,310, 325,260,260,390,390,325, 340,340,340,340,340,340}; floatty[]={60,60,60,50,50,60, 85,85,85,85,85,85, 110,110,110,110,110,110, 135,135,135,135,135,135, 160,160,160,160,160,160, 170,170,170,170,170,170}; floattz[]={-10,-10,10,10,-10,-10, -10,-10,10,10,-10,-10, -10,-10,10,10,-10,-10, -10,-10,10,10,-10,-10, -10,-10,10,10,-10,-10, 0,0,0,0,0,0}; floatTU[]={0,0,0,0.25,0.50,0.75,1,1,1}; floatbst[161][161][3]; floatbst1[161][161][3]; floatbst2[161][161][3]; floatbst3[161][161][3]; floatbst4[161][161][3]; floatbst5[161][161][3]; floatbst6[161][161][3]; //计算所有非零B样条基函数并返回其值 //i为参数u所在的节点区间下标 voidBasisFunction(inti,intp,floatu,floatU[],floatN[]) { intj,di,dp,k; floattul,tur,left,right; floattmpN[50][50]; for(k=0;k<=p;k++) { dp=0; for(di=i+p-k;di>=i-k;di--) { if(u>=U[di]&&u tmpN[di][0]=1; else tmpN[di][0]=0; dp+=1; for(j=1;j { tul=U[di+j]-U[di]; tur=U[di+j+1]-U[di+1]; if(tul! =0) left=(u-U[di])/tul; else left=0; if(tur! =0) right=(U[di+j+1]-u)/tur; else right=0; tmpN[di][j]=left*tmpN[di][j-1]+right*tmpN[di+1][j-1]; } } N[i-k]=tmpN[i-k][p]; } } //计算参数u的p次基函数值并存在BC[]中 voidBernsteinFunc(intn,doubleu,floatB[]) { doubleu1; intj,k; B[0]=1.0; u1=1.0-u; for(j=1;j<=n;j++) { floatsaved=0.0; for(k=0;k { floattemp=B[k]; B[k]=saved+u1*temp; saved=u*temp; } B[j]=saved; } } //获取p次Bezier曲线上的lines个点的值 voidBezierPoint(intp,floatpx[],floatpy[],floatpz[],intlines,floattmp[][3]) { floatBC[20]; inti,j; for(j=0;j<=lines;j++) { doublet=j/(float)lines; BernsteinFunc(p,t,BC); tmp[j][0]=tmp[j][1]=tmp[j][2]=0; for(i=0;i { tmp[j][0]+=BC[i]*px[i]; tmp[j][1]+=BC[i]*py[i]; tmp[j][2]+=BC[i]*pz[i]; } } } //--------------------------------------------------------------------------------- //计算双p次Bezier曲面上所有的点并保存在Pt[][][]中 //u和v分别为曲面(u,v)方向上的网格数 voidBezierFacePoint(intp,intu,intv,floatpx[][4],floatpy[][4],floatpz[][4],floatpt[161][161][3]) { floaturx[11][161],ury[11][161],urz[11][161]; floattx[11],ty[11],tz[11],tmp[161][3]; inti,j,k; for(j=0;j { for(i=0;i { tx[i]=px[i][j]; ty[i]=py[i][j]; tz[i]=pz[i][j]; } BezierPoint(p,tx,ty,tz,v,tmp); for(k=0;k<=v;k++) { urx[j][k]=tmp[k][0]; ury[j][k]=tmp[k][1]; urz[j][k]=tmp[k][2]; } } for(i=0;i<=v;i++) { for(k=0;k { tx[k]=urx[k][i]; ty[k]=ury[k][i]; tz[k]=urz[k][i]; } BezierPoint(p,tx,ty,tz,u,tmp); for(j=0;j<=u;j++) { pt
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 均匀 有理