[Bf-committers] Surface and Volume modeling with new developed Bezier-Surface math
Roland Adorni
rol at solnet.ch
Sun Sep 29 18:46:14 CEST 2013
Hello,
I am nowhere with the web page that I wanted to create. However I
started to dig a bit into blender scripting and so I created a script
which makes the thing more visible.
The script creates tangential planes and handles on the vertices of
the selected mesh. However it does not create bezier surfaces or volume
yet.
The name is: Bezier Volume Operator
Sry if the code may look funny. I never programmed phyton before nor I
know what blender offers (math functions and data structures).
regards
rad
--------------------
import bpy
import math
class XYZPoint:
def __init__(self,x,y,z):
self.x = x
self.y = y
self.z = z
description = "XYZPoint stores the x,y and z coord of a point or a
vector"
author = "rad"
class HNFPlane:
def __init__(self,a,b,c,d):
self.a = a
self.b = b
self.c = c
self.d = d
description = "HNFPlane a*x + b*y + c* z + d = 0, whereas
a*a+b*b+c*c =1"
author = "rad"
def calcdistance(self,x,y,z):
dist = self.a * x + self.b * y + self.c * z + self.d
return dist
def flipnormal(self):
self.a = -self.a
self.b = -self.b
self.c = -self.c
self.d = -self.d
def flipnormaltowardpoint(self, x,y,z):
dist = self.calcdistance(x,y,z)
if dist < 0:
self.flipnormal()
def moveplanetopoint(self, x,y,z):
dist = self.calcdistance(x,y,z)
newd = self.d - dist * (self.a*self.a + self.b*self.b +
self.c*self.c)
self.d = newd
def CalculateXYZGlobalAngles(self):
vBackTurnAngleY = - math.atan2(self.a, self.c)
vNewZ = math.sqrt(self.a * self.a + self.c * self.c);
vBackTurnAngleX = math.atan2(self.b, vNewZ)
GlobalAngleX = - vBackTurnAngleX;
GlobalAngleY = - vBackTurnAngleY;
return (GlobalAngleX, GlobalAngleY, 0)
class BezierVolumeOperator(bpy.types.Operator):
bl_idname = "object.beziervolume_operator"
bl_label = "Bezier Volume Operator"
def execute(self, context):
print("Start")
HandleLengthFactor = 0.333 #Factor of edge lengh
SourceData = context.active_object.data
VertIndexOfInterest = 0;
for VertOfInterest in SourceData.vertices:
TangentialPlaneVerts =
ExtractTangentialVerticesArray(SourceData, VertIndexOfInterest)
TangentialHNFPlane = CreateTangentialPlane(VertOfInterest,
TangentialPlaneVerts)
HandelsVectors =
CalculateHandleVectors(TangentialHNFPlane, VertOfInterest,
TangentialPlaneVerts, HandleLengthFactor)
drawplane(context, VertOfInterest, TangentialHNFPlane,
HandelsVectors)
#print(TangentialPlaneVerts)
VertIndexOfInterest += 1
print("Stop")
return {'FINISHED'}
def drawplane(context, destver, iHNFPlane, iHandelsVectors):
vGlobalPlaneAngles = iHNFPlane.CalculateXYZGlobalAngles();
#rotate the handles to global plane
globalizedHhandels = []
for handle in iHandelsVectors:
globalhandle = TurnToGlobal(handle, vGlobalPlaneAngles)
globalizedHhandels.append(globalhandle)
maxhandlelength = 0
for gbh in globalizedHhandels :
length = math.sqrt(gbh[0] * gbh[0] + gbh[1] * gbh[1] +
gbh[2] * gbh[2])
if length > maxhandlelength:
maxhandlelength = length
coords = []
radius = 0.33 * maxhandlelength
vAngles = (0,45,90,135,180,225,270,315)
for vAngle in vAngles:
vX = radius * math.cos(3.1415926 * vAngle/180)
vY = radius * math.sin(3.1415926 * vAngle/180)
coords.append((vX,vY,0))
#coords = [(-1,-1,0), (1,-1,0), (1,1,0), (-1,1,0)]
faces = [(0,1,2,3,4,5,6,7)]
vertexcount = 8
vtrianglebasefactor = 0.1
for globhandle in globalizedHhandels:
point1 = ( vtrianglebasefactor*globhandle[1], -
vtrianglebasefactor*globhandle[0], 0)
point2 = (-vtrianglebasefactor*globhandle[1],
vtrianglebasefactor*globhandle[0], 0)
coords.append(point1)
coords.append(globhandle)
coords.append(point2)
faces.append((vertexcount,vertexcount+1,vertexcount+2))
vertexcount += 3
#print(coords)
#print(faces)
me = bpy.data.meshes.new("TangPlaneMesh")
ob = bpy.data.objects.new("TangPlaneObject", me)
ob.rotation_euler = vGlobalPlaneAngles
ob.location = (destver.co[0],destver.co[1],destver.co[2])
bpy.context.scene.objects.link(ob)
ob.parent = context.active_object
me.from_pydata(coords,[],faces)
me.update(calc_edges=True)
return {'FINISHED'}
def TurnToGlobal(handle, vGlobalPlaneAngles):
vX = handle[0]
vY = handle[1]
vZ = handle[2]
# Z first
vAlfaSin = math.sin(-vGlobalPlaneAngles[2]);
vAlfaCos = math.cos(-vGlobalPlaneAngles[2]);
vNewX = vAlfaCos * vX - vAlfaSin * vY;
vNewY = vAlfaSin * vX + vAlfaCos * vY;
vX = vNewX
vY = vNewY
# Y first
vAlfaSin = math.sin(-vGlobalPlaneAngles[1]);
vAlfaCos = math.cos(-vGlobalPlaneAngles[1]);
vNewZ = vAlfaCos * vZ - vAlfaSin * vX;
vNewX = vAlfaSin * vZ + vAlfaCos * vX;
vZ = vNewZ
vX = vNewX
# X after
vAlfaSin = math.sin(-vGlobalPlaneAngles[0]);
vAlfaCos = math.cos(-vGlobalPlaneAngles[0]);
vNewY = vAlfaCos * vY - vAlfaSin * vZ;
vNewZ = vAlfaSin * vY + vAlfaCos * vZ;
vY = vNewY
vZ = vNewZ
return (vX,vY,vZ)
def ExtractTangentialVerticesArray(SourceData, indexofinterest):
verticesindex = []
arrvertices = []
for ed in SourceData.edges:
found = False
for vert in ed.vertices:
if indexofinterest == vert:
found = True
if (found):
for vert in ed.vertices:
if indexofinterest != vert:
verticesindex.append(vert)
for vert in verticesindex:
arrvertices.append(SourceData.vertices[vert])
return arrvertices
def CreateTangentialPlane(destver, planevertarray):
cEpsilon = 0.00000001
arraysize = 0.0
for vert in planevertarray:
arraysize += 1.0
if arraysize < cEpsilon:
return HNFPlane(0,0,1,0)
vDN = arraysize
vSx = 0.0
vSy = 0.0
vSx2 = 0.0
vSy2 = 0.0
vSz2 = 0.0
vSxy = 0.0
vSz = 0.0
vSxz = 0.0
vSyz = 0.0
for vert in planevertarray:
x = vert.co[0]
y = vert.co[1]
z = vert.co[2]
vSx += x
vSy += y
vSz += z
vSx2 += x*x
vSy2 += y*y
vSz2 += z*z
vSxy += x*y
vSxz += x*z
vSyz += y*z
#[a,b,c]
vNennerABC = (vSx2 * (vSyz * vSyz - vSy2 * vSz2) + vSxy * vSxy *
vSz2 - 2.0 * vSxy * vSxz * vSyz + vSxz * vSxz * vSy2)
#[b,c,d]
vNennerBCD = (vSy2 * (vSz * vSz - vSz2 * vDN) + vSyz * vSyz * vDN +
vSy * vSy * vSz2 - 2.0 * vSy * vSyz * vSz)
#[a,c,d]
vNennerACD = (vSx2 * (vSz * vSz - vSz2 * vDN) + vSxz * vSxz * vDN +
vSx * vSx * vSz2 - 2.0 * vSx * vSxz * vSz)
#[a,b,d]
vNennerABD = (vSx2 * (vSy * vSy - vSy2 * vDN) + vSxy * vSxy * vDN +
vSx * vSx * vSy2 - 2.0 * vSx * vSxy * vSy)
vMaxNenner = abs(vNennerABC)
if abs(vNennerBCD) > vMaxNenner:
vMaxNenner = abs(vNennerBCD)
if abs(vNennerACD) > vMaxNenner:
vMaxNenner = abs(vNennerACD)
if abs(vNennerABD) > vMaxNenner:
vMaxNenner = abs(vNennerABD)
if vMaxNenner< cEpsilon:
return HNFPlane(0,0,1,0)
vA = 0.0
vB = 0.0
vC = 0.0
vD = 0.0
if abs(vNennerABC) >= vMaxNenner:
vFactor = 1.0 / vNennerABC
vD = 1.0
vA = vFactor * (vD * vSx * (vSy2 * vSz2 - vSyz * vSyz) + vD *
vSxy * (vSyz * vSz - vSy * vSz2) + vD * vSxz * (vSy * vSyz - vSy2 * vSz))
vB = vFactor * -(vD * vSx2 * (vSyz * vSz - vSy * vSz2) + vD *
vSx * vSxy * vSz2 + vSxz * (-vD * vSxy * vSz - vD * vSx * vSyz) + vD *
vSxz * vSxz * vSy)
vC = vFactor * (vD * vSx2 * (vSy2 * vSz - vSy * vSyz) - vD *
vSxy * vSxy * vSz + vD * vSx * vSxy * vSyz + vSxz * (vD * vSxy * vSy -
vD * vSx * vSy2))
elif abs(vNennerBCD) >= vMaxNenner:
vFactor = 1.0 / vNennerBCD
vA = 1.0
vB = vFactor * (vA * vSxy * (vSz2 * vDN - vSz * vSz) + vA * vSyz
* (vSx * vSz - vSxz * vDN) + vA * vSy * (vSxz * vSz - vSx * vSz2))
vC = vFactor * (vA * vSy2 * (vSxz * vDN - vSx * vSz) + vSyz *
(vA * vSx * vSy - vA * vSxy * vDN) + vA * vSxy * vSy * vSz - vA * vSxz *
vSy * vSy)
vD = vFactor * -(vA * vSy2 * (vSxz * vSz - vSx * vSz2) + vA *
vSxy * vSy * vSz2 + vSyz * (-vA * vSxy * vSz - vA * vSxz * vSy) + vA *
vSx * vSyz * vSyz)
elif abs(vNennerACD) >= vMaxNenner:
vFactor = 1.0 / vNennerACD
vB = 1.0
vA = vFactor * (vB * vSxy * (vSz2 * vDN - vSz * vSz) + vB * vSxz
* (vSy * vSz - vSyz * vDN) + vB * vSx * (vSyz * vSz - vSy * vSz2))
vC = vFactor * -(vB * vSx2 * (vSy * vSz - vSyz * vDN) + vSxz *
(vB * vSxy * vDN - vB * vSx * vSy) - vB * vSx * vSxy * vSz + vB * vSx *
vSx * vSyz)
vD = vFactor * -(vB * vSx2 * (vSyz * vSz - vSy * vSz2) + vB * vSx
* vSxy * vSz2 + vSxz * (-vB * vSxy * vSz - vB * vSx * vSyz) + vB * vSxz
* vSxz * vSy)
elif abs(vNennerABD) >= vMaxNenner:
vFactor = 1.0 / vNennerABD
vC = 1.0;
vA = vFactor * (vC * vSxy * (vSy * vSz - vSyz * vDN) + vC * vSxz
* (vSy2 * vDN - vSy * vSy) + vC * vSx * (vSy * vSyz - vSy2 * vSz))
vB = vFactor * -(vC * vSx2 * (vSy * vSz - vSyz * vDN) + vSxz *
(vC * vSxy * vDN - vC * vSx * vSy) - vC * vSx * vSxy * vSz + vC * vSx *
vSx * vSyz)
vD = vFactor * (vC * vSx2 * (vSy2 * vSz - vSy * vSyz) - vC * vSxy
* vSxy * vSz + vC * vSx * vSxy * vSyz + vSxz * (vC * vSxy * vSy - vC *
vSx * vSy2))
vNormalize = math.sqrt(vA * vA + vB * vB + vC * vC);
if vNormalize < cEpsilon:
return HNFPlane(0,0,1,0)
vNormalizeFactor = 1.0 / vNormalize
vA = vNormalizeFactor * vA
vB = vNormalizeFactor * vB
vC = vNormalizeFactor * vC
vD = vNormalizeFactor * vD
vHNFPlane = HNFPlane(vA,vB,vC,vD)
vHNFPlane.flipnormaltowardpoint(destver.co[0], destver.co[1],
destver.co[2])
vHNFPlane.moveplanetopoint(destver.co[0], destver.co[1],
destver.co[2])
#print(vHNFPlane.a)
#print(vHNFPlane.b)
#print(vHNFPlane.c)
#print(vHNFPlane.d)
return vHNFPlane
def CalculateHandleVectors(iTangentialHNFPlane, iVertOfInterest,
iTangentialPlaneVerts, iHandleLengthFactor):
Handles = []
for vert in iTangentialPlaneVerts:
dist = iTangentialHNFPlane.calcdistance(vert.co[0], vert.co[1],
vert.co[2])
planehitpoint = (vert.co[0] - dist*iTangentialHNFPlane.a,
vert.co[1] - dist*iTangentialHNFPlane.b, vert.co[2] -
dist*iTangentialHNFPlane.c)
handlevector = (planehitpoint[0] - iVertOfInterest.co[0],
planehitpoint[1] - iVertOfInterest.co[1], planehitpoint[2] -
iVertOfInterest.co[2])
reducedhandle = (iHandleLengthFactor*handlevector[0],
iHandleLengthFactor*handlevector[1], iHandleLengthFactor*handlevector[2])
Handles.append(reducedhandle)
return Handles
def register():
bpy.utils.register_class(BezierVolumeOperator)
if __name__ == "__main__":
register()
More information about the Bf-committers
mailing list