Example codes of Bezier and B-spline surfaces (Python)

The following codes in Python are only tested with Grasshopper-Python and not with pure Python provided with Rhinoceros.

Bezier curve

Input: P (List Access)

Input parameter P is a list of sequential points.
2.png

import rhinoscriptsyntax as rs
import Rhino.Geometry as rg
import math
dim=len(P)-1
def u(x,nu):
    return x**nu
def v(x,nu):
    return (1-x)**(dim-nu)
def nv(nu):
    return math.factorial(dim)/(math.factorial(dim-nu)*math.factorial(nu))
 
tt=list()
for i in rs.frange(0,1.0,0.01):
    f=list()
    tt.append(f)
    for j in range(0,dim+1):
        f.append(nv(j)*u(i,j)*v(i,j))
a=list()
print tt
for s in tt:
    PP=rg.Point3d(0,0,0)
    for d in range(0,dim+1):
        PP=PP+s[d]*P[d]
    a.append(PP)

1.png

Bezier Surface

Input: P (List Access)
Input: nU
Input: nV
Input parameter P is a list of sequential points.
nU*nV must be equal to the total number of the points.

3.png 
import rhinoscriptsyntax as rs
import Rhino.Geometry as rg
import math
dimU=nU-1
dimV=nV-1
def u(x,nn):
    return x**nn
def v(x,nn,dim):
    return (1-x)**(dim-nn)
def nv(nn,dim):
    return math.factorial(dim)/(math.factorial(dim-nn)*math.factorial(nn))
 
def main():
    tt=list()
    for x in rs.frange(0,1.0,0.1):
        for y in rs.frange(0,1.0,0.1):
            f=list()
            tt.append(f)
            for j in range(0,nV):
                for k in range(0,nU):
                    s1=nv(k,dimU)*u(x,k)*v(x,k,dimU)
                    s2=nv(j,dimV)*u(y,j)*v(y,j,dimV)
                    f.append(s1*s2)
        out=list()
    print tt[120]
 
    for s in tt:
        PP=rg.Point3d(0,0,0)
        for d in range(0,nU*nV):
            PP=PP+s[d]*P[d]
        print PP
        out.append(PP)
    return out
if len(P)==(nU)*(nV):
    a=main()

4.png

B-spline curve

Input: P (List access)
P is a point list.
B-spline curve is 99% same as Nurbs curve.
5.png

import rhinoscriptsyntax as rs
import Rhino.Geometry as rg
import math
ddim=dim-1
N=len(P)
nSeg=N-1
knot=range(N-ddim+1)
startK=knot[0]
endK=knot[knot.Count-1]
for i in range(ddim):
    knot.insert(0,startK)
for i in range(ddim):
    knot.insert(knot.Count,endK)
def binominal(m,n):
    return math.factorial(m)/math.factorial(m-n)/math.factorial(n)
def func(i,j):
    k=dim
    T=binominal(k-1,i)/math.factorial(k-1)
    value=0
    for l in range(j,k):
        value+=binominal(k,l-j)*(k-(l+1))**i*(-1)**(l-j)
    return value*T
def fN(_i,_k,_dim):
    if _dim==1:
        F=list()
        for i in range(dim):
            F.append(0)
        if _k==_i:
            F[dim-1]=1
        return F
    S1=fN(_i,_k,_dim-1)
    S2=fN(_i,_k+1,_dim-1)
    E1=knot[_k+_dim-2]-knot[_k-1]
    E2=knot[_k+_dim-1]-knot[_k]
    D1=list()
    D2=list()
    D1.append(0)
    D1.append(0)
    D2.append(0)
    D2.append(0)
    if  E1>0:
        D1[0]=1/E1
        D1[1]=-knot[_k-1]/E1
    if E2>0:
        D2[0]=-1/E2
        D2[1]=knot[_k+_dim-1]/E2
    F=list()
    for i in range(dim):
        F.append(0)
    for i in range(1,dim):
        F[i-1]=F[i-1]+S1[i]*D1[0];
        F[i]=F[i]+S1[i]*D1[1];
        F[i-1]=F[i-1]+S2[i]*D2[0];
        F[i]=F[i]+S2[i]*D2[1];
    return F
def fM(i):
    M=rg.Matrix(dim,dim)
    M.Zero()
    for k in range(i,dim+i):
        D=fN(i+ddim,k,dim)
        for n in range(dim):
            M[n,k-i]=D[n]
 
    S=rg.Matrix(dim,dim)
    S.Zero()
    for n in range(1,dim+1):
        for k in range(1+n,dim+2):
            if n==dim:
                for t in range(n-1):
                   S[t,n-1]=0
                S[n-1,n-1]=1
            else:
                S[k-2,n-1]=binominal(dim-n,dim+1-k)*((i-1)**(k-1-n))
    G=S*M
    return G
def draw():
    out=list()
    segment=list()
    for i in range(1,N-ddim+1):
        M=fM(i)            
        points=list()
        for n in range(dim):
            points.append(P[i-1+n])
        segment.append((M,points))
    T=rg.Matrix(1,dim)
    for seg in segment:
        M=seg[0]
        points=seg[1]
        for t in rs.frange(0,1.0,0.2):
            for i in range(dim):
                T[0,i]=t**(ddim-i)
            F=T*M
            PP=rg.Point3d(0,0,0)
            for i in range(dim):
                PP=PP+F[0,i]*points[i]
            out.append(PP)
    return out
if N>dim-1:
    a=draw()

6.png
An equivalent rhinoscriptsyntax is 
import Rhino.Geometry as rg
import rhinoscriptsyntax as rs
a=rg.NurbsCurve.CreateControlPointCurve(P,int(dim)-1)

7.png
Note that the same result is obtained.

B-spline surface

Input: P (List access)
Input: nU
Input: nV
Input: uDim
Input: vDim
nU*nV must be equal to the total number of points.
8.png

import rhinoscriptsyntax as rs
import Rhino.Geometry as rg
import math
def split_seq(seq, size):
    return [seq[i:i+size] for i in range(0, len(seq), size)]
 
def binominal(m,n):
    return math.factorial(m)/math.factorial(m-n)/math.factorial(n)
def func(i,j):
    k=dim
    T=binominal(k-1,i)/math.factorial(k-1)
    value=0
    for l in range(j,k):
        value+=binominal(k,l-j)*(k-(l+1))**i*(-1)**(l-j)
    return value*T
def fN(_i,_k,_dim,dim,knot):
    if _dim==1:
        F=list()
        for i in range(dim):
            F.append(0)
        if _k==_i:
            F[dim-1]=1
        return F
    S1=fN(_i,_k,_dim-1,dim,knot)
    S2=fN(_i,_k+1,_dim-1,dim,knot)
    E1=knot[_k+_dim-2]-knot[_k-1]
    E2=knot[_k+_dim-1]-knot[_k]
    D1=list()
    D2=list()
    D1.append(0)
    D1.append(0)
    D2.append(0)
    D2.append(0)
    if  E1>0:
        D1[0]=1/E1
        D1[1]=-knot[_k-1]/E1
    if E2>0:
        D2[0]=-1/E2
        D2[1]=knot[_k+_dim-1]/E2
    F=list()
    for i in range(dim):
        F.append(0)
    for i in range(1,dim):
        F[i-1]=F[i-1]+S1[i]*D1[0];
        F[i]=F[i]+S1[i]*D1[1];
        F[i-1]=F[i-1]+S2[i]*D2[0];
        F[i]=F[i]+S2[i]*D2[1];
    return F
def fM(i,dim,ddim,knot):
    M=rg.Matrix(dim,dim)
    M.Zero()
    for k in range(i,dim+i):
        D=fN(i+ddim,k,dim,dim,knot)
        for n in range(dim):
            M[n,k-i]=D[n]
 
    S=rg.Matrix(dim,dim)
    S.Zero()
    for n in range(1,dim+1):
        for k in range(1+n,dim+2):
            if n==dim:
                for t in range(n-1):
                   S[t,n-1]=0
                S[n-1,n-1]=1
            else:
                S[k-2,n-1]=binominal(dim-n,dim+1-k)*((i-1)**(k-1-n))
    G=S*M
    return G
def draw(P):
    out=list()
    segment=list()
    for j in range(1,nV-vDdim+1):
        vM=fM(j,vDim,vDdim,vKnot)
        for i in range(1,nU-uDdim+1):
            uM=fM(i,uDim,uDdim,uKnot)
            points=list()
            for n in range(vDim):
                uPoints=list()
                points.append(uPoints)
                for m in range(uDim):
                    uPoints.append(P[j-1+n][i-1+m])
            segment.append((uM,vM,points))
    uT=rg.Matrix(1,uDim)
    vT=rg.Matrix(1,vDim)
    for seg in segment:
        (uM,vM,points)=seg
        for ut in rs.frange(0,1.0,0.1):
            for i in range(uDim):
                uT[0,i]=ut**(uDdim-i)
            for vt in rs.frange(0,1.0,0.1):
                for i in range(vDim):
                    vT[0,i]=vt**(vDdim-i)
                uF=uT*uM
                vF=vT*vM
                PP=rg.Point3d(0,0,0)
                for i in range(uDim):
                    for j in range(vDim):
                        PP=PP+uF[0,i]*vF[0,j]*points[j][i]
                out.append(PP)
    return out
uDdim=uDim-1
vDdim=vDim-1
N=len(P)
uNSeg=nU-1
vNSeg=nV-1
uKnot=range(nU-uDdim+1)
vKnot=range(nV-vDdim+1)
uStartK=uKnot[0]
uEndK=uKnot[uKnot.Count-1]
vStartK=vKnot[0]
vEndK=vKnot[vKnot.Count-1]
for i in range(uDdim):
    uKnot.insert(0,uStartK)
for i in range(uDdim):
    uKnot.insert(uKnot.Count,uEndK)
for i in range(vDdim):
    vKnot.insert(0,vStartK)
for i in range(vDdim):
    vKnot.insert(vKnot.Count,vEndK)
 
print vKnot
if nU*nV==N and nU>uDim-1 and nV>vDim-1:
    P=split_seq(P,nU)
    a=draw(P)

9.png

An equivalent rhinoscriptsyntax is

import Rhino.Geometry as rg
a=rg.NurbsSurface.CreateFromPoints(P,nV,nU,vDim-1,uDim-1)

10.png
Note that the same result is obtained.
page revision: 24, last edited: 04 Nov 2013 04:57
Edit Tags History Files Print Site tools + Options