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

1、学习小结：前面学习了 Bezier曲线，B样条基函数和B样条曲线的一些基础知识。掌握关键问题是一条B样条曲线间的多段曲线的光滑连接。因为现在是用多段Bezier曲线来描绘一条B样条曲线，所以问题变为两段Bezier曲线间光滑连接。两段Bezier曲线段(3次)B1和B2光滑连接的条件：(1) .要求B1和B2有共同的连接点，即&连续。(2) .要求B1和B2在连接点处有成比例的一阶导数，即G连续。由端点处的一阶导数B1(1)=3(F3 P* B, 2(0) =3(QQ。)为实现 G 连续，则有：B2(0)=B1 即：Qi-QoP2这也表明, P2,P5(Qo),Qi三点共线。如下图表示

2、了一条3次B样条曲线的所有控制多边形:出R)P(P与拦109/P用p°p图5.3次B样条曲线和所有控制多边形P4图5中，P0至P6为原始3次B样条曲线控制多边形顶点，P。至是计算后最终形成b样条曲线控制多、'T-r ZL-图6.双二次(2x2)B样条曲面6.B样条曲线曲面和 NURBS曲线曲面的C语言实现算法源程序#ifndef _mynu rbs_h#ifndef MYNURBS H#include nglgl.hH#include "math.h"*/确定参数U所在的节点区间下标/n=m-p-1*B样条基函数计算部分/m为节点矢量U的最大下标/P为B样

3、条函数次数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)tif(u<Umid) high=mid; else low=mid;mid=(int)(low+high)/2;Wu>=Umid&&uvUmid+1)找至U u所在的节点区间的下标break; /退出二分搜索 return mid; /返回区间下标/

4、计算所有非零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>=i-k;di-)(if(u>=Udi&&u<Udi+1)tmpNdi0=1;elsetmpNdiO=O;dp+=1；for(j=1 ;j<dp;j+)(tul=Udi+j-Udi;tur=Ud

5、i+j+1-Udi+1;if(tul!=O)left=(u-Udi)/tul;elseleft=O;if(tur!=O)right=(Udi+j+1 -u)/tur;elseright=O;tmpNdi0=left*tmpNdi0-1 +right*tmpNdi+1)Ni-k=tmpNi-kp;/计算基函数的1阶导数并保存在NP中i为参数u所在的节点区间下标/P为B样条函数次数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

6、,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; if(tul!=O) left=(u-Udi)/tul,dl=1/tul;elseleft=0,dl=0;if(tur!=O)right=(Udi+j+1-u)/tur,dr=1/tur;else right=0,dr

7、=0;tmpNdi0=(left*tmpNdi0-1+right*tmpNdi+10-1)；saved=p*(dl*tmpNdi0-1-dr*tmpNdi+10-1)/(p+p-1)；) NPi-k=saved; )Bezier 曲线曲面部*/计算参数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>O&&i<p) BCi=p

8、*(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 tmp3) 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+)(tmpO0+=BCi*pxi;tmpO1+=BCi*pyi;t

9、mpO2+=BCi*pzi;) )/获取P次有理Bezier曲线上的lines个点的值pw,intvoid NBezierPoint(int p,float px,float py,float pz,float lines,float tmp4) (float x,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

10、; w+=BCi*pwi;)tmpjO=x/w;tmpj1=y/w;tmpj2=z/w;tmpj3=w;/ 绘制p次的Bezier曲线void Bezier(int p,float px,float py,float pz,int lines)(float pt1003;intj;BezierPoint(p,px,py,pz,lines,pt); for(j=1 ;j<=lines;j+) (glBegin(GL_LINES);gIVertex3f(ptr10,ptH1,ptH2);glVertex3f(ptj0,ptj1,p52);glEnd(); /绘制P次的有理Bezier曲线vo

11、id NBezier(int p,float px,float py,float pz,float w,int lines)(float pt1004;intj;NBezierPoint(p,px,py,pz,w,lines,pt); for(j=1 ;j<=lines;j+) (glBegin(GL_GNES);glVertex3fglEnd();/计算双p次Bezier曲面上所有的点并保存在 Pt中/u和v分别为曲面(u,v)方向上的网格数void BezierFacePoint(int p,int u,int v,float px4,float pz4,float pt161161

12、3) py4,float (float urx11 161 ,ury11 161 ,urz11 161 ;float tx11,ty11,tz11,tmp1613;int for(j=0;j<p+1 ;j+)(for(i=0;i<p+1 ;i+) txi=pxiD； tyi=pyij; tzi=pz 皿;)BezierPoint(p,tx,ty,tz,v,tmp); for(k=0;k<=v;k+) (urxOk=tmpk0;uryOk=tmpk1;urzOk=tmpk2;)for(i=0;i<=v;i+)(for(k=0;k<p+1 ;k+)(txk=urxki

13、;tyk=uryki;tzk=urzki;)BezierPoint(p,tx,ty,tz,u,tmp);for(j=0;j<=u;j+)(pti0O=tmp0O;ptiO1=tmpO1；ptiO2=tmpO2;)/计算双P次有理Bezier曲面上所有的点并保存在Pt中/u和v分别为曲面(u,v)方向上的网格数void NuBezierFacePoint(int p,int u,int v,float px4,float pz4,float w4,floatpy4,floatpt1611613)float urx11161,ury11161,urz11161,urw11161;float

14、tx11 ,ty11 ,tz11 ,tw11 ,tmp161 4;intfor(j=0;j<p+1 ;j+)(for(i=0;i<p+1 ;i+)(txi=pxiD;tyi=pyi0;tzi=pzi0；twi=wij;NBezierPoint(p,tx5ty,tz,tw5v,tmp);for(k=0;k<=v;k+)urxOk=tmpk0; uryOk=tmpk1; urzjk=tmpk2; urwOk=tmpk3;) for(i=0;i<=v;i+)for(k=0;k<p+1 ;k+)(txk=urxki;tyk=uryki;tzk=urzki;twk=urwk

15、i; )NBezierPoint(p,tx,ty,tz,tw,u,tmp); for(j=0;j<=u;j+) (pti0O=tmp0O;ptij1=tmpO1； ptiO2=tmpO2;)/*B样条曲线曲面部*计算样条曲线的1阶导矢(u所对应的所有点)保存在De山中/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 ;iv=n;i+) (DerBasisFunc(i,p,Ui,U,N);tmp=0;for(j=i;j>=i-

16、p;j-) tmp+=Nj*Pj;Deri-p=tmp;计算曲线上的点(U所对应的所有点)保存在Poi中n=mp-1/p为曲线的次数void BSplinePoint(int n,int p,float U,float P,float Poi) (float N100,tmp;int i,j;for(i=p+1 ;iv=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 B3Spline

17、ControlPoint(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+ad

18、d; 计算2次样条曲线上的所有控制多边形保存在CP口中/m为节点矢量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的最大下标v

19、oid BSpline3L(int m,float ll,float 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

20、 ,pcyi-1 ,pczi-1 );glVertex3f(pcxi,pcyi,pczi);glEnd();glColor3f(1.0f,0.0f,0.0f);7tmcp=m-7;for(i=0;i<=tmcp;i+)for(j=i*3;j<i*3+4;j+)(drxO-i*3=pcxO；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,pcy

21、100,pcz100,drx3,dry3,drz3;int i,j,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;drz0-i*2=pczj;)Bezier(2,drx,dry,drz,20);)/计算双三次(3x3)B样条曲面所有控制多边形顶点，并保存在 pt中/

22、mu,mv分别为节点矢量U,V的最大下标值void BS3FaceControlPoint(int mu,float U,int mv,float V,float py,float pz,floatpx,floatpt1001003)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

23、-dp=pyj;tmzj-dp=pzj;)B3SplineControlPoint(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+)tmxO=tmpxji; tmyO=tmpyOi； tmzj=tmpzji;)B3SplineControlPoint(mv,V,tmx,uvxi);B3SplineControlPoint(mv,V,tmy,uvyi);B3SplineCon

24、trolPoint(mv,VJtmz,uvzi); for(k=0;k<3*mv-17;k+) (ptikO=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,tmpy5

25、0100,tmpz50100;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;

26、j+)tmxO=tmpxji;tmyj=tmpyji；tmzj=tmpzOi；)B2SplineControlPoint(mv,V,tmx,uvxi);B2SplineControlPoint(mv,V,tmy,uvyi);B2SplineControlPoint(mv,V,tmzJuvzi); for(k=0;k<2*mv-7;k+)ptikO=uvxik;ptik1=uvyik;ptik2=uvzik;/设置网格数void SetGridCount(int dtjnt tu,int tmk) int i,tm;tm=tu%dt;for(i=0;i<dt-1 ;i+) tmki=

27、(tu-tm)/dt;tmkdt-1 =tmk0+tm;/计算双三次(3x3次)或双二次(2x2次)B样条曲面上所有的点并保存在bs口皿中/nu,mv分别为节点矢量U,V的最大下标/uk.vk分别为B样条曲面(u,v)方向上的网格数/P为曲面的次数void BSplineFace(int p,int nu,float U,int ukjnt 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,td1611613;flo

28、at 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+) tmx0-i*pl-k*p=tp

29、l0O;tmy0-i*pl-k*p=tpl01； tmz0-i*pl-k*p=tpl02;)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+) bssvhuO=tdsv-i*vdkOhu-k*udkOO;bssvhu1 =tdsv-i*vdk0hu-k*udk01;bssvhu2=tdsv-i*vdk0hu-k*udk02;)/*5*,*,-* Nurbs样条曲线曲面部*计算Nurbs曲线上的点(u所对应的

30、所有点)保存在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 ;iv=n;i+) (BasisFunction(i,p,Ui,U,N); tmp=0;tmw=0;for(j=i;j>=i-p;j-)(tmp+=NO*PO*WO； tmw+=NO*WO； )Poii-p=tmp/tmw; 计算Nurbs曲线的1阶导矢(u所对应的所有点)保存在Der口中n=m-p-1 /p为曲线的次数void Nurb

31、sDer(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 ;iv=n;i+)(DerBasisFunc(i,p,Ui,U,N); tmp=0;tw=0;for(j=i;j>=i-p;j-)(tmp+=N0*P0*W0； tw+=Nj*Wj； iDeri-p=(tmp-tw*CPi-p)/NWi-p; )/计算3次Nurbs曲线上的所有控制多边形保

32、存在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;jfor(i=3;i<cp-1 ;i+=3) (k=(int)i/3;add=Derk/p;CPi=Poik;CPi-1

33、=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=2;i<cp-1;i+=2) (k=(int)i/2;CPi=Poik;CPi+1=Pt

34、m;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(m,U,pz,W,pcz);B3SplineControlPoint(m,U,W,p

35、cw); 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,pcz100,drx3,dry3,drz3; float pcw100,

36、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;dryg-i*2=pcy0;drzj-i*2=pczj;drwj-i*2=pcwj;NBezier(2,drx,dry,drz,drw,20); 中px5flo

37、at计算双三次(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 uvx100100,uvy100100,uvz100100,uvw100100;for(i=0

38、;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,tmx5tmw,tmpxi);Nurbs3ControlPoint(mu5U,tmy,tmw5tmpyi);Nurbs3ControlPoint(mu5U,tmz,tmw5tmpzi);B3SplineControlPoint(mu,U,tmw,tmpwi); for(i=0;i<3*mu-17;i+) (for(j=0;j<mv-3;

39、j+)tmxO=tmpxji;tmyj=tmpyji；tmzO=tmpzOi； 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+)ptikO=uvxik;ptik1=uvyik;ptik2=uvzik;ptik3=uvwik;/计算双二次(2x2) Nurbs样条曲面所有控

40、制多边形顶点，并保存在 pt 中/mu,mv分别为节点矢量U,V的最大下标值void Nurbs2FControlPoint(int mu,float U,int mv,float V,float py,float pz,floatpx,floatW,float pt1001004)int i,j,k,dp;float tmx50 ,tmy 50 ,tmz50 ,tm w50;float tmpx50100,tmpy50100,tmpz50100,tmpw50100; float uvx100100,uvy100100,uvz100100,uvw100100;for(i=0;i<mv-2

41、;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(mu5U4rnx,trnw,tmpxi);Nurbs2ControlRoint(mu,LI,tmy,tmw5tmpyi);Nurbs2ControlPoint(muJU,tmz5tmw,tmpzi);B2SplineControlPoint(mu,U,tmw,tmpwi);for(i=0;i<2*mu-7;i+)(for(j=0;j<mv-2;j+)tmxj=tmp

42、xji； tmyj=tmpyji；tmzO=tmpzOi；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+) ptikO=uvxik;ptik1=uvyik;ptik2=uvzik;ptik3=uvwik;/计算双三次(3x3次)或双二次(2x2次)Nurbs样条曲面上所有的

43、点并保存在/nu,mv分别为节点矢量U,V的最大下标 uk,vk分别为B样条曲面(u,v)方向上的网格数/P为曲面的次数w,floatvoid NurbsFace(int p,int nu,float U,int ukjnt mv,float V,int vk, float px,floatpy,floatpz,floatbs1611613) (int udk20,vdk20,iJ,k,l,hu,sv,du,dv;float tp1001004,td161 1613;float tmx4 4, t my 4 4, t mz4 4, tm w4 4;du=nu-2*p;dv=mv-2*p;Set

44、GridCount(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+) (tmx0-i*pl-k*p=tpl0O;tmy0-i*pl-k*p=tpl01；tmz0-i*pl-k*p=

45、tpl02;显曲面部tmwO-i*pl-k*p=tplD3;)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+)bssvhuO=tdsv-i*vdkOhu-k*udkOO;bssvhu1 =tdsv-i*vdk0hu-k*udk01;bssvhu2=tdsv-i*vdk0hu-k*udk02;* *人帘,mH计算多边形的外法线返回值tmN void getN(float x3,float y3,fl

46、oat 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;float x3,y3,z3,tmn3;for(i=0;i<=

47、v;i+)for(j=0;j<=u;j+)(if(fill!=O)(xO=bsi0O; x1=bsi+10O; x2=bsi+10+1O;yO=bsi01； y1=bsi+101； y2=bsi+10+11；zO=bsi02;z1=bsi+102;z2=bsi+10+12;getN(x,y,z,tmn);glEnable(GL_NORMALIZE);glBegin(GL_QUADS); glNormal3f(tmn0,tmn1,tmn2);if(jvu)(glVertex3f(bsi0O.bsi01,bsi02);glVertex3f(bsi0+1 O,bsi0+1 1 ,bsi0+1

48、2);if(i<v)(glVertex3f(bsi+1 j+1 O,bsi+10+11 ,bsi+1 j+1 glVertex3f(bsi+1D0,bsi+1j1,bsi+1j2);|glEnd();glDisable(GL_NORMALIZE);|elseglBegin(GL_LINES);if(j<u)(glVertex3f(bsi0O,bsi01,bsi02);glVertex3f(bsi0+1 O,bsi0+1 1 ,bsi0+12);|if(i<v)(glVertex3f(bsij0,bsiD1,bsij2);glVertex3f(bsi+1 j0,bsi+1 01

49、 ,bsi+1 02);)glEnd();#endif#endif建立库文件myNurbs.h ”保存在in elude目录下，在文件开始处写上include"myNurbs.h"即可用各函数的功能。例如：#include "windows.h"include "glut.hHinclude "myNurbs.h1'float x=17,17,-17,-17,-17,17,17,10,10, -10,-10,-10,10,10,15,15,-15,-15,-15,15,15,30,30,-30,-30,-30,30,30,20,20,-20,-20,-20,20,20,10,10, -10,-10,-10,10,10,0,0,0,0,0,0,0）；float y=100,100,100,100,