B样条曲线曲面和NURBS曲线曲面C语言算法源程序

1、学习小结：前面学习了Bezier曲线,B样条基函数和 B样条曲线的一些根底知识.掌握关键问题是一条B样条曲线间的多段曲线的光滑连接.由于现在是用多段Bezier曲线来描绘一条B样条曲线,所以问题变为两段 Bezier曲线间光滑连接.两段 Bezier曲线段(3次) B1和B2光滑连接的条件：(1) .要求B1和B2有共同的连接点,即 片连续.(2) .要求B1和B2在连接点处有成比例的一阶导数,即 G连续.由端点处的一阶导数1B1=3(P3 P2), B'2(0) =3(Q1 Qo),为实现 G连续,那么有：B20)=B1 即：Q1Qo=RP2这也说明,P2,巳(Qo),Q1三点共线.

2、如下列图表示了一条 3次B样条曲线的所有限制多边形：13P p/4 (P 11)P5 -/ / PP却 P0PP/先/ P 7P 二? 8图5.3次B样条曲线和所有限制多边形图5中,P0至P6为原始3次B样条曲线限制多边形顶点, 由线限制多边形顶点.厂占府(P12)109P4P0至P12是计算后最终形成B丹(P1)P=2"图6.双二6.B样条曲线曲面和 NURBSffi线曲面的次(2x2) B样条曲面C语言实现算法源程序#ifndef _mynurbs_h#ifndef _MYNURBS_H#include "glgl.h"#include "math.

3、h"*-*-*-*-*-*-*-*-*-*-*-*-*-*_*B 样条基函数计算音 B 分*-*-*-*-*-*-*-*-*-*-*-*-*/确定参数u所在的节点区间下标/n=m-p-1/m为节点矢量U口的最大下标/p为B样条函数次数int FindSource(int n,int p,float u,float U)int low,high,mid;if(u=Un+1)特殊情况return n;/进行二分搜索low=p;high=n+1;mid=(int)(low+high)/2;while(u<Umid|u>Umid)if(u<Umid)high=mid;els

4、elow=mid;mid=(int)(low+high)/2;if(u>=Umid&&u<Umid+1)找到u所在的节点区间的下标break; /退出二分搜索return mid; / 返回区间下标计算所有非零B样条基函数并返回其值/i为参数u所在的节点区间下标void BasisFunction(int i,int p,float u,float U,float N) int j,di,dp,k;float tul,tur,left,right;float tmpN5050;for(k=0;k<=p;k+) dp=0;for(di=i+p-k;di>=

5、i-k;di-) if(u>=Udi&&u<Udi+1)tmpNdi0=1;elsetmpNdi0=0;dp+=1;for(j=1;j<dp;j+)(tul=Udi+j-Udi;tur=Udi+j+1-Udi+1;if(tul!=0)left=(u-Udi)/tul;elseleft=0;if(tur!=0)right=(Udi+j+1卜u)/tur;elseright=0;tmpNdij=left*tmpNdij-1+right*tmpNdi+1j-1;Ni-k=tmpNi-kp;/计算基函数的1阶导数并保存在NP中/i为参数u所在的节点区间下标/p为B样条

6、函数次数P>2void DerBasisFunc(int i,int p,float u,float U,float NP)(int j,di,dp,k;float tul,tur,left,right,saved,dl,dr;float tmpN5050;for(k=0;k<=p;k+)(dp=0;for(di=i+p-k;di>=i-k;di-)(if(u>=Udi&&u<Udi+1)tmpNdi0=1;elsetmpNdi0=0;dp+=1;for(j=1;j<dp;j+)(tul=Udi+j-Udi;tur=Udi+j+1-Udi+1

7、;if(tul!=0)left=(u-Udi)/tul,dl=1/tul; else left=0,dl=0; if(tur!=0)right=(Udi+j+1卜u)/tur,dr=1/tur; elseright=0,dr=0;tmpNdij=(left*tmpNdij-1+right*tmpNdi+1j-1);saved=p*(dl*tmpNdij-1卜dr*tmpNdi+1j-1)/(p+p-1); NPi-k=saved; /*-*-*-*-*-*-*-*-*-*-*-*-*-*Bezier 曲 线 曲 面 部 分*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*/计算

8、参数u的p次基函数值并存在BC口中void BernsteinFunc(int p,double t,float BC) for(int i=0;i<=p;i+)if(i=0)BC0=(float)pow(1-t,p);if(i=p)BCp=(float)pow(t,p);if(i>0&&i<p)BCi=p*(float)pow(t,i)*(float)pow(1-t,p-i);/获取p次Bezier曲线上的lines 个点的值void BezierPoint(int p,float px,float py,float pz,int lines,float t

9、mp3)float BC20;int i,j;for(j=0;j<=lines;j+)double t=j/(float)lines;BernsteinFunc(p,t,BC);tmpj0=tmpj1=tmpj2=0;for(i=0;i<p+1;i+)(tmpj0+=BCi*pxi;tmpj1+=BCi*pyi;tmpj2+=BCi*pzi; ) ) )/获取p次有理Bezier曲线上的lines个点的值pw,intvoid NBezierPoint(int p,float px,float py,float pz,float lines,float tmp4) (float x,

10、y,z,w,BC20; int i,j;for(j=0;j<=lines;j+) (double t=j/(float)lines;BernsteinFunc(p,t,BC);x=y=z=w=0;for(i=0;i<p+1;i+) ( x+=BCi*pxi*pwi;y+=BCi*pyi*pwi;z+=BCi*pzi*pwi;w+=BCi*pwi; ) tmpj0=x/w; tmpj1=y/w; tmpj2=z/w; tmpj3=w; ) )/ /绘制p次的Bezier曲线void Bezier(int p,float px,float py口,float pz,int lines)

11、 (float pt1003;int j;BezierPoint(p,px,py,pz,lines,pt);for(j=1;j<=lines;j+)(glBegin(GL_LINES);glVertex3f(ptj-10,ptj-11,ptj-12);glVertex3f(ptj0,ptj1,ptj2);glEnd(); ) ) /-/绘制p次的有理Bezier曲线void NBezier(int p,float px,float py,float pz,float w,int lines) float pt1004; int j;NBezierPoint(p,px,py,pz,w,li

12、nes,pt);for(j=1;j<=lines;j+) glBegin(GL_LINES);glVertex3f(ptj-10,ptj-11,ptj-12);glVertex3f(ptj0,ptj1,ptj2);glEnd(); ) ) /计算双p次Bezier曲面上所有的点并保存在Pt 中/u和v分别为曲面(u,v)方向上的网格数py4,floatvoid BezierFacePoint(int p,int u,int v,float px4,float pz4,float pt1611613) float urx11161,ury11161,urz11161;float tx11,

13、ty11,tz11,tmp1613;int i,j,k;for(j=0;j<p+1;j+) for(i=0;i<p+1;i+) txi=pxij; tyi=pyij; tzi=pzij; )BezierPoint(p,tx,ty,tz,v,tmp);for(k=0;k<=v;k+) urxjk=tmpk0;uryjk=tmpk1;urzjk=tmpk2; ) )for(i=0;i<=v;i+) ( for(k=0;k<p+1;k+) (txk=urxki;tyk=uryki;tzk=urzki; )BezierPoint(p,tx,ty,tz,u,tmp);for

14、(j=0;j<=u;j+) (ptij0=tmpj0;ptij1=tmpj1;ptij2=tmpj2;) ) ) /计算双p次有理Bezier曲面上所有的点并保存在Pt 中/u和v分别为曲面(u,v)方向上的网格数py4,floatvoid NuBezierFacePoint(int p,int u,int v,float px4,float pz4,float w4,float pt1611613) (float urx11161,ury11161,urz11161,urw11161;float tx11,ty11,tz11,tw11,tmp1614;int i,j,k;for(j=0

15、;j<p+1;j+) (for(i=0;i<p+1;i+) ( txi=pxij; tyi=pyij; tzi=pzij; twi=wij; )NBezierPoint(p,tx,ty,tz,tw,v,tmp);for(k=0;k<=v;k+)urxjk=tmpk0;uryjk=tmpk1;urzjk=tmpk2;urwjk=tmpk3;)for(i=0;i<=v;i+)(for(k=0;k<p+1;k+)(txk=urxki;tyk=uryki;tzk=urzki; twk=urwki;)NBezierPoint(p,tx,ty,tz,tw,u,tmp); fo

16、r(j=0;j<=u;j+) (ptij0=tmpj0;ptij1=tmpj1;ptij2=tmpj2;)曲 线 曲 面 部 分Der中-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*/计算样条曲线的1阶导矢(u所对应的所有点)保存在/n=m-p-1/p为曲线的次数void BSplineDer(int n,int p,float U,float P,float Der)(float N100,tmp;int i,j;for(i=p+1;i<=n;i+)(DerBasisFunc(i,p,Ui,U,N);tmp=0;for(j=i;j>=i-p;j-)tmp+

17、=Nj*Pj;Deri-p=tmp;/计算曲线上的点(u所对应的所有点)保存在 Poi中 /n=m-p-1/p为曲线的次数void BSplinePoint(int n,int p,float U,float P,float Poi)(float N100,tmp;int i,j;for(i=p+1;i<=n;i+)(BasisFunction(i,p,Ui,U,N);tmp=0;for(j=i;j>=i-p;j-)tmp+=Nj*Pj;Poii-p=tmp;/计算3次样条曲线上的所有限制多边形保存在CP口中/m为节点矢量U的最大下标void B3SplineControlPoin

18、t(int m,float U,float P,float CP)(int n,k,i,cp,p;float Poi100,Der100,add;p=3;n=m-p-1;BSplinePoint(n,p,U,P,Poi);BSplineDer(n,p,U,P,Der);cp=(n-p)*3+p;for(i=0;i<2;i+)(CPi=Pi;CPcp-i=Pn-i;for(i=3;i<cp-1;i+=3)(k=(int)i/3;add=Derk/p;CPi=Poik;CPi-1=CPi卜add;CPi+1=CPi+add;/计算2次样条曲线上的所有限制多边形保存在CP口中/m为节点矢

19、量U的最大下标void B2SplineControlPoint(int m,float U,float P,float CP) (int n,k,tm,i,cp,p;float Poi100;p=2;n=m-p-1;BSplinePoint(n,p,U,P,Poi);cp=(n-p)*2+p;for(i=0;i<2;i+)CPi=Pi;CPcp=Pn;tm=2;for(i=2;i<cp-1;i+=2)(k=(int)i/2;CPi=Poik;CPi+1=Ptm;tm+;/绘制3次B样条曲线/m为节点矢量U的最大下标void BSpline3L(int m,float U,floa

20、t px,float py,float pz)(float pcx100,pcy100,pcz100,drx4,dry4,drz4;int i,j,tmcp;B3SplineControlPoint(m,U,px,pcx);B3SplineControlPoint(m,U,py,pcy);B3SplineControlPoint(m,U,pz,pcz);/*glColor3f(0.0f,0.0f,0.0f);for(i=1;i<3*m-17;i+)(glBegin(GL_LINES);glVertex3f(pcxi-1,pcyi-1,pczi-1);glVertex3f(pcxi,pcy

21、i,pczi);glEnd();glColor3f(1.0f,0.0f,0.0f);*/tmcp=m-7;for(i=0;i<=tmcp;i+)for(j=i*3;j<i*3+4;j+)(drxj-i*3=pcxj;dryj-i*3=pcyj;drzj-i*3=pczj;)Bezier(3,drx,dry,drz,20);)/绘制2次B样条曲线/m为节点矢量U口的最大下标void BSpline2L(int m,float U,float px,float py,float pz)(float pcx100,pcy100,pcz100,drx3,dry3,drz3;int i,j,

22、tmcp;B2SplineControlPoint(m,U,px,pcx);B2SplineControlPoint(m,U,py,pcy);B2SplineControlPoint(m,U,pz,pcz);tmcp=m-5;for(i=0;i<=tmcp;i+)(for(j=i*2;j<i*2+3;j+)(drxj-i*2=pcxj;dryj-i*2=pcyj;drzj-i*2=pczj;)Bezier(2,drx,dry,drz,20);)计算双三次(3x3) B样条曲面所有限制多边形顶点,并保存在pt 中/mu,mv分别为节点矢量 U,V的最大下标值px,floatvoid

23、BS3FaceControlPoint(int mu,float U,int mv,float V,float py,float pz口,float pt1001003)(int i,j,k,dp;float tmx50,tmy50,tmz50;float tmpx50100,tmpy50100,tmpz50100;float uvx100100,uvy100100,uvz100100;for(i=0;i<mv-3;i+) (dp=i*(mu-3);for(j=dp;j<mu-3+dp;j+) ( tmxj-dp=pxj; tmyj-dp=pyj; tmzj-dp=pzj;)B3S

24、plineControlPoint(mu,U,tmx,tmpxi);B3SplineControlPoint(mu,U,tmy,tmpyi);B3SplineControlPoint(mu,U,tmz,tmpzi);)for(i=0;i<3*mu-17;i+) (for(j=0;j<mv-3;j+) ( tmxj=tmpxji; tmyj=tmpyji; tmzj=tmpzji;)B3SplineControlPoint(mv,V,tmx,uvxi);B3SplineControlPoint(mv,V,tmy,uvyi);B3SplineControlPoint(mv,V,tmz,

25、uvzi); for(k=0;k<3*mv-17;k+) ( ptik0=uvxik;ptik1=uvyik;ptik2=uvzik; ) ) )计算双二次(2x2) B样条曲面所有限制多边形顶点,并保存在pt 中/mu,mv分别为节点矢量 U,V的最大下标值px,floatvoid BS2FaceControlPoint(int mu,float U,int mv,float V,float py,float pz口,float pt1001003) (int i,j,k,dp;float tmx50,tmy50,tmz50;float tmpx50100,tmpy50100,tmpz

26、50100;float uvx100100,uvy100100,uvz100100;for(i=0;i<mv-2;i+) ( dp=i*(mu-2); for(j=dp;j<mu-2+dp;j+)tmxj-dp=pxj;tmyj-dp=pyj;tmzj-dp=pzj;)B2SplineControlPoint(mu,U,tmx,tmpxi);B2SplineControlPoint(mu,U,tmy,tmpyi);B2SplineControlPoint(mu,U,tmz,tmpzi);)for(i=0;i<2*mu-7;i+)(for(j=0;j<mv-2;j+)(t

27、mxj=tmpxji;tmyj=tmpyji;tmzj=tmpzji;)B2SplineControlPoint(mv,V,tmx,uvxi);B2SplineControlPoint(mv,V,tmy,uvyi);B2SplineControlPoint(mv,V,tmz,uvzi); for(k=0;k<2*mv-7;k+) (ptik0=uvxik;ptik1=uvyik;ptik2=uvzik;)/-/设置网格数void SetGridCount(int dt,int tu,int tmk) (int i,tm;tm=tu%dt;for(i=0;i<dt-1;i+) tmk

28、i=(tu-tm)/dt;tmkdt-1=tmk0+tm;)/-bs中/计算双三次(3x3次)或双二次(2x2次)B样条曲面上所有的点并保存在 /nu,mv分别为节点矢量 U,V的最大下标/uk,vk 分别为B样条曲面(u,v)方向上的网格数/p为曲面的次数void BSplineFace(int p,int nu,float U,int uk,int mv,float V,int vk,float px,float py,float pz口,float bs1611613)int udk20,vdk20,i,j,k,l,hu,sv,du,dv;float tp1001003,td161161

29、3;float tmx44,tmy44,tmz44;du=nu-2*p;dv=mv-2*p;SetGridCount(du,uk,udk);SetGridCount(dv,vk,vdk);if(p=3)BS3FaceControlPoint(nu,U,mv,V,px,py,pz,tp);if(p=2)BS2FaceControlPoint(nu,U,mv,V,px,py,pz,tp);for(i=0;i<dv;i+)for(k=0;k<du;k+)for(j=i*p;j<p+1+i*p;j+)for(l=k*p;l<p+1+k*p;l+)tmxj-i*pl-k*p=tp

30、lj0;tmyj-i*pl-k*p=tplj1;tmzj-i*pl-k*p=tplj2;BezierFacePoint(p,udkk,vdki,tmx,tmy,tmz,td);for(sv=i*vdk0;sv<=vdki+i*vdk0;sv+)for(hu=k*udk0;hu<=udkk+k*udk0;hu+) bssvhu0=tdsv-i*vdk0hu-k*udk00; bssvhu1=tdsv-i*vdk0hu-k*udk01; bssvhu2=tdsv-i*vdk0hu-k*udk02;/-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*Nurbs 样 条 曲

31、线 曲 面 部 分-*-*-*-*-*-*-*-*-*-*-*-*-*-*-/计算Nurbs曲线上的点(u所对应的所有点)保存在 Poi口中/n=m-p-1/p为曲线的次数void NurbsPoint(int n,int p,float U,float P,float W,float Poi)(float N100,tmp,tmw;int i,j;for(i=p+1;i<=n;i+)(BasisFunction(i,p,Ui,U,N);tmp=0;tmw=0;for(j=i;j>=i-p;j-) (tmp+=Nj*Pj*Wj;tmw+=Nj*Wj;Poii-p=tmp/tmw;/

32、计算Nurbs曲线的1阶导矢(u所对应的所有点)保存在Der口中/n=m-p-1/p为曲线的次数void NurbsDer(int n,int p,float U,float P,float W,float Der)(float N100,CP100,NW100,tmp,tw;int i,j;NurbsPoint(n,p,U,P,W,CP);BSplinePoint(n,p,U,W,NW);for(i=p+1;i<=n;i+) (DerBasisFunc(i,p,Ui,U,N);tmp=0;tw=0;for(j=i;j>=i-p;j-) (tmp+=Nj*Pj*Wj;tw+=Nj*

33、Wj;Deri-p=(tmp-tw*CPi-p)/NWi-p;/计算3次Nurbs曲线上的所有限制多边形保存在CP口中/m为节点矢量U的最大下标void Nurbs3ControlPoint(int m,float U,float P,float W,float CP)(int n,k,i,cp,p;float Poi100,Der100,add;p=3;n=m-p-1;NurbsPoint(n,p,U,P,W,Poi);NurbsDer(n,p,U,P,W,Der);cp=(n-p)*3+p;for(i=0;i<2;i+)CPi=Pi;CPcp-i=Pn-i;for(i=3;i<

34、cp-1;i+=3)k=(int)i/3;add=Derk/p;CPi=Poik;CPi-1=CPi卜add;CPi+1=CPi+add;/计算2次Nurbs曲线上的所有限制多边形保存在CP口中/m为节点矢量U的最大下标void Nurbs2ControlPoint(int m,float U,float P,float W,float CP) int n,k,tm,i,cp,p;float Poi100;p=2;n=m-p-1;NurbsPoint(n,p,U,P,W,Poi);cp=(n-p)*2+p;for(i=0;i<2;i+)CPi=Pi;CPcp=Pn;tm=2;for(i=

35、2;i<cp-1;i+=2)k=(int)i/2;CPi=Poik;CPi+1=Ptm;tm+;/绘制3次Nurbs样条曲线/m为节点矢量U的最大下标void Nurbs3L(int m,float U,float px,float py,float pz,float W) (float pcx100,pcy100,pcz100,drx4,dry4,drz4;float pcw100,drw4;int i,j,tmcp;Nurbs3ControlPoint(m,U,px,W,pcx);Nurbs3ControlPoint(m,U,py,W,pcy);Nurbs3ControlPoint(

36、m,U,pz,W,pcz);B3SplineControlPoint(m,U,W,pcw);tmcp=m-7;for(i=0;i<=tmcp;i+)(for(j=i*3;j<i*3+4;j+)(drxj-i*3=pcxj;dryj-i*3=pcyj;drzj-i*3=pczj;drwj-i*3=pcwj;NBezier(3,drx,dry,drz,drw,20);/绘制2次Nurbs样条曲线/m为节点矢量U的最大下标void Nurbs2L(int m,float U,float px,float py,float pz,float W) (float pcx100,pcy100,

37、pcz100,drx3,dry3,drz3;float pcw100,drw3;int i,j,tmcp;Nurbs2ControlPoint(m,U,px,W,pcx);Nurbs2ControlPoint(m,U,py,W,pcy);Nurbs2ControlPoint(m,U,pz,W,pcz);B2SplineControlPoint(m,U,W,pcw);tmcp=m-5;for(i=0;i<=tmcp;i+)(for(j=i*2;j<i*2+3;j+)(drxj-i*2=pcxj;dryj-i*2=pcyj;drzj-i*2=pczj;drwj-i*2=pcwj; )N

38、Bezier(2,drx,dry,drz,drw,20);)中px口,float计算双三次(3x3) Nurbs样条曲面所有限制多边形顶点,并保存在pt/mu,mv分别为节点矢量 U,V口的最大下标值void Nurbs3FControlPoint(int mu,float U,int mv,float V,float py,float pz口,float W,float pt1001004)int i,j,k,dp;float tmx50,tmy50,tmz50,tmw50;float tmpx50100,tmpy50100,tmpz50100,tmpw50100;float uvx1001

39、00,uvy100100,uvz100100,uvw100100;for(i=0;i<mv-3;i+) dp=i*(mu-3); for(j=dp;j<mu-3+dp;j+) tmxj-dp=pxj;tmyj-dp=pyj;tmzj-dp=pzj;tmwj-dp=Wj; )Nurbs3ControlPoint(mu,U,tmx,tmw,tmpxi);Nurbs3ControlPoint(mu,U,tmy,tmw,tmpyi);Nurbs3ControlPoint(mu,U,tmz,tmw,tmpzi);B3SplineControlPoint(mu,U,tmw,tmpwi);)fo

40、r(i=0;i<3*mu-17;i+)for(j=0;j<mv-3;j+) tmxj=tmpxji; tmyj=tmpyji; tmzj=tmpzji; tmwj=tmpwji;)Nurbs3ControlPoint(mv,V,tmx,tmw,uvxi);Nurbs3ControlPoint(mv,V,tmy,tmw,uvyi);Nurbs3ControlPoint(mv,V,tmz,tmw,uvzi);B3SplineControlPoint(mv,V,tmw,uvwi);for(k=0;k<3*mv-17;k+)(ptik0=uvxik;ptik1=uvyik;ptik2

41、=uvzik;ptik3=uvwik;中px,float计算双二次(2x2) Nurbs样条曲面所有限制多边形顶点,并保存在pt/mu,mv分别为节点矢量 U,V口的最大下标值void Nurbs2FControlPoint(int mu,float U,int mv,float V,float py,float pz口,float W,float pt1001004)(int i,j,k,dp;float tmx50,tmy50,tmz50,tmw50;float tmpx50100,tmpy50100,tmpz50100,tmpw50100;float uvx100100,uvy10010

42、0,uvz100100,uvw100100;for(i=0;i<mv-2;i+)(dp=i*(mu-2);for(j=dp;j<mu-2+dp;j+)(tmxj-dp=pxj;tmyj-dp=pyj;tmzj-dp=pzj;tmwj-dp=Wj;Nurbs2ControlPoint(mu,U,tmx,tmw,tmpxi);Nurbs2ControlPoint(mu,U,tmy,tmw,tmpyi);Nurbs2ControlPoint(mu,U,tmz,tmw,tmpzi);B2SplineControlPoint(mu,U,tmw,tmpwi);for(i=0;i<2*mu

43、-7;i+)(for(j=0;j<mv-2;j+)(tmxj=tmpxji;tmyj=tmpyji;tmzj=tmpzji;tmwj=tmpwji;Nurbs2ControlPoint(mv,V,tmx,tmw,uvxi);Nurbs2ControlPoint(mv,V,tmy,tmw,uvyi);Nurbs2ControlPoint(mv,V,tmz,tmw,uvzi);B2SplineControlPoint(mv,V,tmw,uvwi);for(k=0;k<2*mv-7;k+) ptik0=uvxik;ptik1=uvyik;ptik2=uvzik;ptik3=uvwik;

44、/-/计算双三次(3x3次)或双二次(2x2次)Nurbs样条曲面上所有的点并保存在bs 中/nu,mv分别为节点矢量U,V口的最大下标 /uk,vk 分别为B样条曲面(u,v)方向上的网格数 /p为曲面的次数void NurbsFace(int p,int nu,float U,int uk,int mv,float V,int vk,float px,float py,float pz,float w,float bs1611613) int udk20,vdk20,i,j,k,l,hu,sv,du,dv;float tp1001004,td1611613;float tmx44,tmy4

45、4,tmz44,tmw44;du=nu-2*p;dv=mv-2*p;SetGridCount(du,uk,udk);SetGridCount(dv,vk,vdk);if(p=3)Nurbs3FControlPoint(nu,U,mv,V,px,py,pz,w,tp);if(p=2)Nurbs2FControlPoint(nu,U,mv,V,px,py,pz,w,tp);for(i=0;i<dv;i+) for(k=0;k<du;k+) for(j=i*p;j<p+1+i*p;j+) for(l=k*p;l<p+1+k*p;l+) tmxj-i*pl-k*p=tplj0;

46、tmyj-i*pl-k*p=tplj1;tmzj-i*pl-k*p=tplj2;)NuBezierFacePoint(p,udkk,vdki,tmx,tmy,tmz,tmw,td);for(sv=i*vdk0;sv<=vdki+i*vdk0;sv+)for(hu=k*udk0;hu<=udkk+k*udk0;hu+)bssvhu0=tdsv-i*vdk0hu-k*udk00;bssvhu1=tdsv-i*vdk0hu-k*udk01;bssvhu2=tdsv-i*vdk0hu-k*udk02;) /*-*-*-*-*-*-*-*-*-*-*-*-*-*-*/计算多边形的外法线返回值

47、tmN口void getN(float x3,float y3,float z3,float tmN3) float p1,p2,p3,q1,q2,q3;float nx,ny,nz;p1=x1-x0;p2=y1-y0;p3=z1-z0;q1=x2-x1;q2=y2-y1;q3=z2-z1;nx=p2*q3-q2*p3;ny=q1*p3-p1*q3;nz=p1*q2-p2*q1;tmN0=nx;tmN1=ny;tmN2=nz;)/显示B样条曲面/fill 取值为0或1void ShowSurface(int u,int v,float bs1611613,int fill) int i,j;f

48、loat x3,y3,z3,tmn3;for(i=0;i<=v;i+)for(j=0;j<=u;j+)(if(fill!=0)(x0=bsij0;x1=bsi+1j0;x2=bsi+1j+10;y0=bsij1;y1=bsi+1j1;y2=bsi+1j+11;z0=bsij2;z1=bsi+1j2;z2=bsi+1j+12;getN(x,y,z,tmn);glEnable(GL_NORMALIZE);glBegin(GL_QUADS);glNorma13f(tmn0,tmn1,tmn2);if(j<u)(glVertex3f(bsij0,bsij1,bsij2);glVertex3f(bsij+10,bsij+11,bsij+12);if(i<v)(glVertex3f(bsi+1j+10,bsi+1j+11,bs