Nurbs样条线算法推导和python实现

1、百度文库让每个人平等地捉升口我 Nurbs样条线算法推导及python实现 注：在现在的建模软件中Nurbs曲线还引入一个跨度的概念，在这里并不体现，这里是有几 个控制点就是控制点数阶，最后好像到36阶曲线就发散了耶 三个控制点： a二(ax, ay) b=(bx, by) C=(cx, cy) OWtWl e为ab线上一点则： ae = t * ab e=(ex, ey)=(t*(bx-ax)+axz t*(by-ay)+ay) 同理，同一 t值，f为be线上一点为： f=(fx, fy)=(t*(cx-bx)+bx/ t*(cy-by)+by) g点为动线ef上的动点为： g=(t*(fx

2、-ex)+ex, t*(fy-ey)+ey)=(gx, gy) gx= tA2*(cx-2*bx+ax)+2*t*(bx-ax)+ax gy= tA2* (cy-2 * by+ay)+2*t*( by-ay)+ay 动点g的运动方程入上图所示 0： (xo，yo); l：(xvyi); 2： (x3#y3); x = t2 (x2 - 2 xx + x0) + t 2 (xx - x0) + x0; y = t?(y? 2 + y。)+1 2 ( y) + y0； 四个控制点 OWtWl a=(ax, ay) (xo，yo) b=(bx, by) (xpyi) C=(cx, cy) (X22)

3、 d=(dx, dy) 仪33) e=(ex/ ey)=(t*(bx-ax)+ax, t*(by-ay)+ay) t(Xi_Xo) + x0/. f = (fx, fy)=(t*(cx-bx)+bx, t*(cy-by)+by) t(x2 - Xi) +xlf g= (gx, gy)=(t*(dx-cx)+cx, t*(dy-cy)+cy) t(x3 -x2) +x2/. h=(hx, hy)=(t*(fx-ex)+ex, t*(fy-ey)+ey) t2 (x2 - 2 xx + x0) + 2 t (xx - x0) + x0,. I = (ix, iy )=(t*(gx-fx)+fx,

4、t*(gy-fy)+fy) t2 (x3 2 “2 + xj + 2 t (“2 _ Xx) + Xp. j = (jx, jy )=(t*(ix-hx)+hxz t*(iy-hy)+hy) t3 (x3 3 X2 + 3 X x0) + 3 t2 (x2 2 X + x0) + 3 t (xx x0) + x0,. t3 (x3 - 3 - (x2 一 xj X。) + 3 (X2 2 Xi + X。) + 3 t 一 x0) + x0,. N个控制点： def cv(Nz t): n= Point(t*( Ni+l.x- Nix)+Nix t*( Ni+l.y- Ni.y)+Ni.y fo

5、r i in range(len(N)-l) 讦 len(n)=l: return n0 else: return cv(n, t) OWtWl N=pOz pl, p2.pn print cv(N, t) python代码实现：(按分段数计算样条线上的点) def cv(N, seg): step=l/seg sP=NO eP=N-l P=sP oldN=N for i in range(step, 1, step): t=i*step count= len(N) while countl: n= Point(t*( N(i+l.x- Ni.x)+Ni.x z t*( Ni+l.y- Ni.

6、y)+Ni.y) for i in range(count-1) N 二n cou N=oldN p.append(nO) p.appe nd(eP) return p 用多项式的方法优化运算: 各项常数确泄：Cx, Cy （xoo） （xpYi） （X22） （*33） （X44） U t （X _ x0） + x0,. t- （x2 - Xi） + xlf. t-（x3 -X2） +X2, t （x4 -x3） +x3/. 因为y和x的运算 过程相同的故下 面只求X 17 t2 （x2 - xj - （xx - x0） + 2 t （X - X。）+ X。 t2 （x3 一 x2） 一（X

7、2 一 Xi） + 2 t （x2 - xj + xx t2 （X4 _ X3） _（X3 _ X2） + 2 t （X3 _ X2） + X2 t3 （X3 - x2） 一（X2 - Xi） 一 （x2 一 Xi） - （Xi 一 x0） + 3 t2 （x2 一 Xi） 一 （xx 一 x0） + 3 t （X1 -Xo） +x0 t3 （x4 一 x3） 一（X3 一 x2） 一 （x3 一 x2） -（x2 一 Xi） + 3 t2 （x3 一 x2） 一 （x2 一 xj） + 3 t（X2-Xi）+Xi （x2 _ xi） _ （xi _ X。） + 4 2 （x3 _ x2）-

8、（x2 _ xx） _ （x2 _ xQ _ （xx _ X。） + 6 t2 （x2 一 xx） 一 （xx 一 x0） + 4 t （X1 一 X。） + x0; Co = Xq 4-C1=4-(x1-C0) 6-C2 = 6-(x2-x1)-C1) 4 C3 = 4 (x3 一 x2) 一 (x2 一 xx) 一 C2) C4 = ( (X4 _ x3) _ (x3 _ X2) _ (x3 _ x2) _ (x2 _ xj) _ C3) U Dl = xl o2 =(X2-DJ D3 = (x3 -X2)- D2) D4 = (X4 _ X3) -(X3 - X2) D3) U E2 =

9、x2 E3 = (x3 -e2) E4 = (X4 -x3)- E3) U F3 = X3 F4 = (x4 一 F3) U G4 = x4 U F4 = (G4-F3)G4 = F3+F4 E4 =(F4 - E3) u F4 = E3 + E4 E3 = (F3-E2)F3 = E2 + E3 d4 = (e4 _ d3) e4 = d3 + d4 D3 = (E3 D?) E3 = D2 + D3 Dz = (E2 - Di) E2 = Di + D2 C4 = ( D4 C3) D4 = C3 + c4 C3 = ( D3 - C2) u D3 = C2 + C3 C2 = ( D?

10、C D2 = C + C2 C = ( D Cq) u D = Co + Ci G4 F3M4 Eq, E3, E4 DD2* D3/ D4 C(),Ci，C2, Cf C4 0 0 1 2 3 4 X4 X3 F4 *2, E3, E4 X1/ Dq* D3, D4 各项系数增量：scale Izl1 1,2,1 2 13,3,1 3 146414 各项系数增M scale: scales = scalen.LH1 + scalen_14： (scalen,0 = l,scalen/n = 1) scale, = scaleni def getScales(n): 各项常系数的缩放值 14n

11、=l 1,2,1 n=2 1,3,34n=3 1,464n=4 返回第n行数据 n为控制点数目 III if *2: return None n=n2 scales=lzl for i in range(n): scales=slj+slj+l for j in range(len(sl)-l) scales.insertfO,!) scales.append(l) sl=scales return scales 各项常系数C Cq3 = Cn,i： (Cn?o = X。或y 或 Z() G = scale, - Cn4 def getCs(scales,points): n=0 n=l “各

12、项常系数 x_4 x_3,F_4 n=2 n=3 n=5 x_2，E_3,E_4 x_lQ_2,D_3,D_4 x_0,C_l,C_2,C_3,C_4 points为控制点 （三角最左边对应各控制点同一轴向坐标， 顺序自下向上对应控制点输入顺序） 返回各轴轴常系数-第n行数据（杨辉三角底行） III points=list(points) points.reverse() axises=zip(*points) Cs= for axis in axises: Cl=axisO for i in range(ljen(axis): C= C.appe nd(axisi) for j in ran

13、ge(len(Cl): C.append(Clj-CO) C1=C C=tuple(map(lambda c,s:c*s, C, scales) Cs.append(C) return Cs 各分段各项U值： def getTs(n，seg): hi n控制点数目，seg分段数 按分段数获取曲线上所有分段点对应 0=t=l范围内t*i的所有数据 III step=l/seg Ts= for i in range(seg+l): tl= s=i*step for j in range(n): tl.append(s*j) Ts.append(tl) Ts.append(l*n) return T

14、s 样条线多项式: n 1 X =tiq： （n为控制点数目） i=0 def cv(Tsz Cs): hi 返回样条线上所有分段点坐标 III n=len(CsO) nAxis=len(Cs) ps= for i in range(len(Ts): sumR=0 for i in range(nAxis) for j in range(n): for k in range(nAxis): sumRk+=Tsij*Cskj ps.appe nd(sumR) return ps 样条线长度：未完待续 判断一点是否在样条线内：未完待续 Python 代码： #nurbs # coding: utf

15、-8 -*-# from copy import deepcopy def getScales(self, n): 各项常系数的缩放值 1,1 n=l 1,2,1 n=2 1,33,1 n=3 14641n=4 返回第n行数据 n为控制点数目 III if *2: return None n=n2 scales=lzl for i in range(n): scales=sl(j+slj+l for j in range(len(sl)-l) scales.insert(O,l) scales.append(l) sl=scales return scales def getCs(selt s

16、cales, points): “各项常系数 x_4 n=0 x_3,F_4 n=l x_2,E_3,E_4 n=2 xJLQ_2,D_3,D_4 n=3 x_0,C_l,C_2,C_3,C_4 n=5 points为控制点 (三角最左边对应各控制点同一轴向坐标. 顺序自下向上对应控制点输入顺序) 返回各轴轴常系数-第n行数据(杨辉三角底行) III if scales is None: return None # points=list(points)#points.reverse()#axises=zip(* points) axises=zip(*reversed(poi nts) Cs

17、= for axis in axises: Cl=axisO for i in range(ljen(axis): c= C.appe nd(axisi) for j in range(len(Cl): C.append(Clj-CO) C1=C #print (t) C=tuple(map(lambda c,s:c*s, C, scales) Cs.appe nd(C) return Cs def getTsfself, n, seg): in n控制点数目，seg分段数 按分段数获取曲线上所有分段点对应 0=t=l范围内t*i的所有数据 III if nl: return None ste

18、p=l/seg Ts= for i in range(seg+l): s=i*step #tl=#for j in range(n):#tl.append(s*j) tl=s*j for j in range(n) Ts.append(tl) Ts.append(l*n) return Ts def getCvps(selfz Ts, Cs): in 返回样条线上所有分段点坐标 III if Cs is None: return None n=len(Cs0) n Axis=len(Cs) ps= for i in range(len(Ts)-l): sumR=O for i in range

19、(nAxis) for j in range(n): for k in range(nAxis): sumRk+=lsilU*CskU ps.appe nd(sumR) return ps class Nurbs(object): def _init_(self, points=(), segs=100): self._seg=segs self._points=poi nts self._n=len(self_points) self._scales=getScales(self/ self._n) self._Cs=getCs(self/ self._scales, self._point

20、s) self._Ts=getTs(self, self._n, self._seg) self._cvps=getCvps(self, self._Ts, self.Cs) #get or set ctrPoints property def ctrPoints(self): return deepcopy(self._points) ctrPoints.setter def ctrPoints(selt points): if points!=self._points: self._poin ts=deepcopy(poi nts) n=len(self._points) if n!=se

21、lf._n: self._n=n self._scales=getScales(self, self._n) self._Ts=getTs(self, self._n, self._seg) self._Cs=getCs(selfz self._scales, self._points) self._cvps=getCvps(selfz self._Ts, self.Cs) #get or set ctrPoints (property def seg(self): return self._seg seg.setter def seg(selfz segs): if segs!=self._

22、seg: self._seg=segs self._Ts=getTs(self, self._n, self._seg) self._cvps=getCvps(self/ self._Ts, self.Cs) (property def scales(self): return deepcopy(self._scales) (property def Cs(self): return deepcopy(self._Cs) property def Ts(self): return deepcopyself._Ts) (property def cvps(self): return deepco

23、py(self._cvps) #get ctrPoints amount def_len_(self): return self._n #tkinter交互程序 #-*- coding:utf-8 -*-# import tkinterz time, nurbs class Cv(tkinter.Canvas): def getObjs(self, tag, p, r): objs=self.find_overlapping(p(0-r; pl-c p0+r, pl+r) return tuple(filter(lambda obj:tag in self.gettags(obj), objs

24、) def createPoint(self, pos, radius, color二Non已 lineColor=None): return self.create_oval(pos0-radius, posllj-radius, pos0+radius，posl+radius, fill二color, outlinG=lineColor; tags=,ctrpl) def modifyPointfselt id, pos=None, color=None, lineColor二None): if pos: xO, yO, xl, yl=self.coords(id) wmz hm=(xl-

25、x0)/2z (yl-y0)/2 self.coords(idz pos0卜wm,poslhm，pos0+wm, posl+hm) if color: self.itemconfig(idz fill=color) if lineColor: self.itemconfigtid, outline=lineColor) def unpackCoordsfself, I): return j for i in I for j in i def sortCtrPsfself, die, I): return dici for i in I def _init_(self, root, size=(

26、640/ 480): super()._init_(root, bgwhite： width=sizeO, height=sizel) self.lineDict=#id:nb self.curLineld=N on e#cvld self.dotLi neld=None#ctrpCvld self.curCtrPsDict=#ctrpld:pos self.curCtrPslnd=#ctrpld.J self.curCtrPld=None#ctrpld self. need2update 二 False self.editState=, new self.r un ning=l root.b

27、ind(,/ self.eQuit) self.bind(, self.eLeftClick) self.bindfBl-Motion, self.eLeftDrag) self.bind(l,/ self.eLeftRelease) self.bind(： self.eMidClick) self.bind(, self.eRightClick) self.colorNormal=,darkorange self.colorPick=gree rV self.colorE nd=,blue self.colorDashLi ne 二orange self.colorCurve=darkblu

28、e, self.stateText=self.create_text(100/100, text=self.editState) self.pack(side=,left,) self.segScale=tkinter.Scale(root, command=self.setSeg, width=10, length=480, sliderlength=15, resolution=l)#滑块 self.segScale.pack(side=,right1) def setSeg(self, seg):#设It分段数 if self.curLineld: self.lineDictself.c

29、urLineld.seg=max(l/int(seg) self.need2update=True #print (event) def eQuit(self, event): self.running=0 def eLeftClick(selt event): pos=event.x/ event.y if self.editState=new,: getCtrps=self.getObjs(,ctrp,/ pos, 3) if len(getCtrps): self.curCtrPld=getCtrpsO self.modifyPoint(self.curCtrPld/ lineColor

30、=self.colorPick) else: self.curCtrPld=self.createPoint(pos, 3, lineColor=self.colorPick) self.curCtrPsDictself.curCtrPld=pos self.curCtrPslnd.append(self.curCtrPld) if self.curLineld is None and len(self.curCtrPslnd)l: ctrps=self.sortCtrPs(self.curCtrPsDict, self.curCtrPsInd) curve 二n urbs.Nurbs(ctr

31、ps) self.curLineld二selfcreate_line(curvecvps, tags二curve： fill=self.colorCurve) self.lineDictself.curLineld=curve self.dotLineld=self.createJine(ctrps, dash=l/ fill=self.colorDashLine) elif self.editState=,select,: getCurve=self.getObjs(,curve/ pos, 3) if len(getCurve): self.curLineld=getCurveO #sel

32、f.curCtrPsDict.clear() #self.curCtrPslnd.clear() ctrps=self.lineDictself.curLineld.ctrPoints for pos in ctrps: point=self.createPoint(pos/ 3, lineColor=self.colorNormal) self.curCtrPsDictpoint=pos self.curCtrPsI nd.append(point) self.modifyPointfself.curCtrPsIndl-l, lineColor=self.colorEnd) self.dot

33、Lineld=self.createine(ctrps, dash=lz fill=self.colorDashLine) self. need2update 二 True self.segScale.set(self.lineDictself.curLineld.seg)# 获取分段数 self. editState二new self.itemconfig(self.stateTextz text=self.editState) def eLeftDrag(selfz event): pos=event.x/ event.y if not self.curCtrPId is None: se

34、lf.curCtrPsDictself.curCtrPld=pos self.modifyPoint(self.curCtrPld, pos) self.need2update=True def eLeftReleasefself, event): pos=event.x, event.y if not self.curCtrPId is None: print (self.curCtrPId) self.curCtrPsDictself.curCtrPld=pos self.modifyPoint(self.curCtrPld, pos) if self.curCtrPId 二二 self.curCtrPsInd 卜 1: if len(self.curCtrPslnd)l: self.modifyPoint(self.curCtrPslnd-2z lineColor=self.colorNormal) self.modifyPoint(self.curCtrPld/ lineColor=self.colorEnd) else: self.modifyPoint(self.curCtrPldz lineColor=self.colorNormal) self.curCtrPld=None self.need2update=True def eMidClick(self,