Option Explicit
‘ 3d Voronoi AKA Project Cell
‘ (c) Gabe Smedresman 2005
‘ All Rights Reserved.
”’ global naming variables
Dim XX: XX = 0
Dim YY: YY = 1
Dim ZZ: ZZ = 2
GenerateVoronoiCells
Sub GenerateVoronoiCells()
Randomize
Dim arrBoundingVolume ‘ array of polysurfaces indicating bounding volume
Dim arrPoints ‘ array of string references to voronoi points
Dim i,j
Dim startTime : startTime = Now
‘ select objects
arrBoundingVolume = Rhino.GetObjects(“Select a Bounding Volume”,16+8,vbTrue,vbTrue)
If IsNull(arrBoundingVolume) Then Exit Sub
HideObjects arrBoundingVolume
arrPoints = Rhino.GetObjects(“Select Cell Points”,1,vbFalse,vbFalse)
ShowObjects arrBoundingVolume
If IsNull(arrPoints) Then Exit Sub
‘ go through each point in the list, and create a cell, name it, color it, And update the progress report.
Dim strCell, strName
Rhino.Print “Beginning cell divisions: ” & UBound(arrPoints)+1 & ” cells total.”
For i = 0 To UBound(arrPoints) ‘ for each point to make a cell around ‘ one point at a time
strCell = GenerateCell(arrPoints(i),arrPoints,arrBoundingVolume)
Dim value : value = Int(Rnd()*255)
Rhino.objectColor strCell, RGB(255,255,value)
strName = “Cell #” & i
Rhino.objectname strCell, strName
Dim strTime : strTime = GetTimeDescription(startTime, (i+1) * 1.0 / (UBound(arrPoints)+1) )
Rhino.Print i+1 & ” of ” & UBound(arrPoints)+1 & ” cells (“& Int((i+1) * 100 / (UBound(arrPoints)+1)) & “%) completed. ” & strTime
Next ‘i
‘ hide source geometry to reveal the cells
HideObjects arrPoints
HideObjects arrBoundingVolume
End Sub
” ———————————————————————-
” given the center, an array of points, and the bounding volume,
” perform the operation, make sure all intermediate geometries
” have been deleted, then return the voronoi cell.
” ———————————————————————-
Function GenerateCell(centerPoint,arrPoints,arrBoundingVolume)
Dim arrBlocks
arrBlocks = CreateBlocks(centerPoint, arrPoints)
Dim strCell
strCell = IntersectBlocks(arrBlocks,arrBoundingVolume)
Dim m: For m = 0 To UBound(arrBlocks)
If IsPolySurface(arrBlocks(m)) Then DeleteObject(arrBlocks(m))
Next
GenerateCell = strCell
End Function
” ———————————————————————-
” using the array of point strings, and one point reference for the center
” create a group of blocks, which, when intersected, will result in the
” voronoi volume.
”
” this is done by, for each point, generating a plane perpendicular to the
” line between the center and the test point, at the midpoint of that line.
” this plane faces the center point and is equidistant from the two points
” across the entire surface.
”
” then extrude the plane towards the center point, theoretically an infinite
” amount but really to the end of the point cloud. meaning, this extrusion
” is closer to the center than it is to the test point.
”
” intersect all of these volumes and you will have the area closest
” to the center.
” ———————————————————————-
Function CreateBlocks(centerPoint, arrPoints)
Dim arrBlocks
ReDim arrBlocks(UBound(arrPoints))
Dim newPlane
Dim numBlocks : numBlocks = 0
Dim i
Dim midpoint, normal
Dim greatestdiagonalspread
greatestdiagonalspread = FindGreatestDiagonalSpread(arrPoints)
Dim centercoords,pointcoords
centercoords = Rhino.PointCoordinates(centerPoint)
‘ —- CREATE BOUNDING PLANES
For i = 0 To UBound(arrPoints) ‘ for each point to consider boundary to
If Not centerPoint = arrPoints(i) Then ‘ as long as the point isn’t the central point
pointcoords = Rhino.PointCoordinates(arrPoints(i))
Dim strPlane
strPlane = CreateBisectingPlane(centercoords, pointcoords, greatestdiagonalspread)
If(Rhino.IsSurface(strPlane)) Then ‘if it’s a valid object, add To the array
‘newPlane has been created
midpoint = VectorMidpoint(centercoords, pointcoords)
normal = VectorUnitize(VectorSubtract(centercoords, pointcoords))
’strPath is now created
Dim strPath
strPath = Rhino.AddLine(midpoint, VectorAdd(midpoint,VectorScale(normal, greatestdiagonalspread)))
Dim strExtrusion
strExtrusion = Rhino.ExtrudeSurface( strPlane, strPath )
Rhino.objectColor strExtrusion, RGB(128,128,128)
arrBlocks(numBlocks) = strExtrusion
numBlocks = numBlocks + 1
If(Rhino.IsObject(strPath)) Then Rhino.DeleteObject(strPath)
If(Rhino.IsObject(strPlane)) Then Rhino.DeleteObject(strPlane)
End If ‘end if it’s a surface
End If ‘end if considering different points
Next ‘i
ReDim Preserve arrBlocks(numBlocks-1)
‘Rhino.Print(numBlocks & ” blocks created.”)
CreateBlocks = arrBlocks
End Function
” ———————————————————————-
” Takes two point coordinate arrays, finds the midpoint
” makes a coordinate system based on the line between these two points
” and uses this coordinate system to make the plane that marks the 3d bisector
” of that line.
” IE this plane is the halfway boundary betweensc the two points
” ———————————————————————-
Function CreateBisectingPlane(arrPtOne, arrPtTwo, reach)
Dim i
Dim center
center = VectorMidpoint(arrPtOne, arrPtTwo)
‘ make new coordinate system p,r,s: p is the line between 1 and 2, r Is To the side, s Is up(ish) from the line
Dim p : p = VectorSubtract(arrPtTwo,arrPtOne) : p = VectorUnitize(p)
Dim up : up = Array(0,0,1)
Dim r : r = VectorCrossProduct(p,up) : If IsVectorZero(r) Then r = Array(0,1,0)
r = VectorUnitize(r) ” points to the right
Dim s : s = VectorCrossProduct(p,r) : s = VectorUnitize(s) ” points To perpendicular To p (forward) And r (side)
‘ now find four points, 1 up 2 left 3 down 4 right (looking from one to two)
Dim arrCorners(3)
arrCorners(0) = VectorAdd(center,VectorScale(s,reach))
arrCorners(1) = VectorAdd(center,VectorScale(r,reach * -1))
arrCorners(2) = VectorAdd(center,VectorScale(s,reach * -1))
arrCorners(3) = VectorAdd(center,VectorScale(r,reach))
Dim strPlane
strPlane = Rhino.AddSrfPt( arrCorners )
CreateBisectingPlane = strPlane
End Function
” ———————————————————————-
” given the array of generated blocks, perform a boolean intersection
” with the bounding volume and each block, one by one.
” delete each source as you iterate. some blocks remain: if
” the intersection fails (no areas intersect), the sources
” will not be deleted.
” ———————————————————————-
Function IntersectBlocks(arrBlocks, arrBoundingVolume)
Dim i
IntersectBlocks = Null
i = 0
Dim results, pendingresults
Dim strBlock
Dim group2
Do
strBlock = arrBlocks(i)
group2 = array(strBlock)
pendingresults = Rhino.BooleanIntersection(arrBoundingVolume,group2,vbFalse)
i = i + 1
Loop While Not IsArray(pendingresults)
results = pendingresults
Dim j
For j = i To UBound(arrBlocks)
strBlock = arrBlocks(j)
group2 = array(strBlock)
pendingresults = Rhino.BooleanIntersection(results,group2)
If IsArray(pendingresults) Then results = pendingresults
Next
IntersectBlocks = results(0)
End Function
” ———————————————————————-
” find the maximum 3d diagonal length between one extreme corner
” of a point cloud and the other
” ———————————————————————-
Function FindGreatestDiagonalSpread(arrPoints)
If(Not IsArray(arrPoints)) Then FindGreatestDiagonalSpread = 0
Dim max : max = Rhino.PointCoordinates(arrPoints(0))
Dim min : min = Rhino.PointCoordinates(arrPoints(0))
Dim i, pt, spread
For i = 0 To UBound(arrPoints)
pt = Rhino.PointCoordinates(arrPoints(i))
If(max(XX) < pt(XX)) Then max(XX) = pt(XX)
If(max(YY) < pt(YY)) Then max(YY) = pt(YY)
If(max(ZZ) pt(XX)) Then min(XX) = pt(XX)
If(min(YY) > pt(YY)) Then min(YY) = pt(YY)
If(min(ZZ) > pt(ZZ)) Then min(ZZ) = pt(ZZ)
Next ‘i
FindGreatestDiagonalSpread = VectorLength(VectorSubtract(max,min))
End Function
” ———————————————————————-
” given the starting time and fraction complete, estimate time to
” completion and then generate a sentence of the format
” “3 minutes 4 seconds elapsed: should be finished at 04:00 PM today.”
” ———————————————————————-
Function GetTimeDescription(startTime, fractioncomplete)
Dim strDescription
strDescription = “”
Dim elapsedseconds : elapsedseconds = DateDiff(“s”,startTime,Now)
Dim m,h,s : s = elapsedseconds
m = Int(s/60) : s = s – m * 60
h = Int(m/60) : m = m – h * 60
If h = 1 Then strDescription = strDescription & ” ” & h & ” hour”
If h > 1 Then strDescription = strDescription & ” ” & h & ” hours”
If m = 1 Then strDescription = strDescription & ” ” & m & ” minute”
If m > 1 Then strDescription = strDescription & ” ” & m & ” minutes”
If s = 1 Then strDescription = strDescription & ” ” & s & ” second”
If s > 1 Then strDescription = strDescription & ” ” & s & ” seconds”
strDescription = strDescription & ” elapsed:”
Dim secondstogo : secondstogo = elapsedseconds / fractioncomplete * (1-fractioncomplete)
Dim ETA : ETA = DateAdd(“s”,secondstogo,Now)
If(fractioncomplete < 1) Then
strDescription = strDescription & " should be finished at " & FormatDateTime(ETA,4)
If(DatePart("d",ETA) = DatePart("d",Now)) Then
strDescription = strDescription & " today."
Else
strDescription = strDescription & " in " & DateDiff("d",Now,ETA) & " day(s)."
End If
Else
strDescription = strDescription & " finished " & Now & "!"
End If
GetTimeDescription = strDescription
End Function
''' ***************************************************************************
''' Rhinoscript Vector Functions
''' ***************************************************************************
''' —————————————————————————
''' Make a vector from two 3D points
''' —————————————————————————
Public Function VectorCreate(p1, p2)
VectorCreate = Array(p2(0) – p1(0), p2(1) – p1(1), p2(2) – p1(2))
End Function
''' —————————————————————————
''' Unitize a 3D vector
''' —————————————————————————
Public Function VectorUnitize(v)
VectorUnitize = Null
Dim dist, x, y, z, x2, y2, z2
x = v(XX) : y = v(YY) : z = v(ZZ)
x2 = x * x : y2 = y * y : z2 = z * z
dist = x2 + y2 + z2
If dist <= 0.0 Then Exit Function
dist = Sqr(dist)
x = x / dist
y = y / dist
z = z / dist
VectorUnitize = Array(x, y, z)
End Function
''' —————————————————————————
''' Return the length of a 3D vector
''' —————————————————————————
Public Function VectorLength(v)
VectorLength = Null
Dim dist, x, y, z, x2, y2, z2
x = v(XX) : y = v(YY) : z = v(ZZ)
x2 = x * x : y2 = y * y : z2 = z * z
dist = x2 + y2 + z2
VectorLength = Sqr(dist)
End Function
''' —————————————————————————
''' Return the dot product of two 3D vectors
''' —————————————————————————
Public Function VectorDotProduct(v1, v2)
VectorDotProduct = v1(XX) * v2(XX) + v1(YY) * v2(YY) + v1(ZZ) * v2(ZZ)
End Function
''' —————————————————————————
''' Return the cross product of two 3D vectors
''' —————————————————————————
Public Function VectorCrossProduct(v1, v2)
VectorCrossProduct = Null
Dim x, y, z
x = v1(YY) * v2(ZZ) – v1(ZZ) * v2(YY)
y = v1(ZZ) * v2(XX) – v1(XX) * v2(ZZ)
z = v1(XX) * v2(YY) – v1(YY) * v2(XX)
VectorCrossProduct = Array(x, y, z)
End Function
''' —————————————————————————
''' Add two 3D vectors
''' —————————————————————————
Public Function VectorAdd(v1, v2)
VectorAdd = Null
VectorAdd = Array(v1(XX) + v2(XX), v1(YY) + v2(YY), v1(ZZ) + v2(ZZ))
End Function
''' —————————————————————————
''' Subtract two 3D vectors
''' —————————————————————————
Public Function VectorSubtract(v1, v2)
VectorSubtract = Null
VectorSubtract = Array(v1(XX) – v2(XX), v1(YY) – v2(YY), v1(ZZ) – v2(ZZ))
End Function
''' —————————————————————————
''' Multiply two 3D vectors
''' —————————————————————————
Public Function VectorMultiply(v1, v2)
VectorMultiply = Null
VectorMultiply = Array(v1(XX) * v2(XX), v1(YY) * v2(YY), v1(ZZ) * v2(ZZ))
End Function
''' —————————————————————————
''' Scale a 3D vectors by a value
''' —————————————————————————
Public Function VectorScale(v, d)
VectorScale = Null
VectorScale = Array(v(XX) * d, v(YY) * d, v(ZZ) * d)
End Function
''' —————————————————————————
''' Compare two 3D vectors for equality
''' —————————————————————————
Public Function VectorCompare(v1, v2)
VectorCompare = vbFalse
If v1(XX) = v2(XX) And v1(YY) = v2(YY) And v1(ZZ) = v2(ZZ) Then
VectorCompare = vbTrue
End If
End Function
''' —————————————————————————
''' Negate a 3D vector
''' —————————————————————————
Public Function VectorNegate(v)
VectorNegate = Null
VectorNegate = Array(-v(XX), -v(YY), -v(ZZ))
End Function
''' —————————————————————————
''' Tiny vector test
''' —————————————————————————
Public Function IsVectorTiny(v)
IsVectorTiny = vbFalse
Dim tol : tol = 1.0e-12 ' ON_ZERO_TOLERANCE
If (Abs(v(XX)) <= tol) And (Abs(v(YY)) <= tol) And (Abs(v(ZZ)) <= tol) Then
IsVectorTiny = vbTrue
End If
End Function
''' —————————————————————————
''' Zero vector test
''' —————————————————————————
Public Function IsVectorZero(v)
IsVectorZero = vbFalse
If (v(XX) = 0.0) And (v(YY) = 0.0) And (v(ZZ) = 0.0) Then IsVectorZero = vbTrue
End Function
''' My more specialized functions
''' —————————————————————————
''' Find midpoint between two vectors
''' —————————————————————————
Public Function VectorMidpoint(v1, v2)
VectorMidpoint = Null
VectorMidpoint = Array((v1(XX) + v2(XX))/2, (v1(YY) + v2(YY))/2,(v1(ZZ) + v2(ZZ))/2)
End Function
''' —————————————————————————
''' Return the length of a 3D vector
''' —————————————————————————
Public Function VectorSquaredDistance(v1,v2)
VectorSquaredDistance = Null
Dim dist, x, y, z, x2, y2, z2
x = v1(XX) – v2(XX)
y = v1(YY) – v2(YY)
z = v1(ZZ) – v2(ZZ)
x2 = x * x
y2 = y * y
z2 = z * z
dist = x2 + y2 + z2
VectorSquaredDistance = dist
End Function
