Newer
Older
# -*- coding:utf-8 -*-
# ##### BEGIN GPL LICENSE BLOCK #####
#
# This program is free software; you can redistribute it and/or
# modify it under the terms of the GNU General Public License
# as published by the Free Software Foundation; either version 2
# of the License, or (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program; if not, write to the Free Software Foundation,
# Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
#
# ##### END GPL LICENSE BLOCK #####
bl_info = {
"name": "Delaunay Voronoi",
"description": "Points cloud Delaunay triangulation in 2.5D "
"(suitable for terrain modelling) or Voronoi diagram in 2D",
"author": "Domlysz, Oscurart",
"version": (1, 3),
"blender": (2, 7, 0),
"location": "3D View > Toolshelf > Create > Delaunay Voronoi",
"warning": "",
"wiki_url": "https://github.com/domlysz/BlenderGIS/wiki",
"category": "Add Mesh"
}
import bpy
from .DelaunayVoronoi import (
computeVoronoiDiagram,
computeDelaunayTriangulation,
)
from bpy.types import (
Operator,
Panel,
)
from bpy.props import EnumProperty
try:
from scipy.spatial import Delaunay
import bmesh
import numpy as np
HAS_SCIPY = True
except:
HAS_SCIPY = False
pass
# Globals
# set to True to enable debug_prints
DEBUG = False
def debug_prints(text=""):
global DEBUG
if DEBUG and text:
print(text)
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
class Point:
def __init__(self, x, y, z):
self.x, self.y, self.z = x, y, z
def unique(L):
"""Return a list of unhashable elements in s, but without duplicates.
[[1, 2], [2, 3], [1, 2]] >>> [[1, 2], [2, 3]]"""
# For unhashable objects, you can sort the sequence and
# then scan from the end of the list, deleting duplicates as you go
nDupli = 0
nZcolinear = 0
# sort() brings the equal elements together; then duplicates
# are easy to weed out in a single pass
L.sort()
last = L[-1]
for i in range(len(L) - 2, -1, -1):
if last[:2] == L[i][:2]: # XY coordinates compararison
if last[2] == L[i][2]: # Z coordinates compararison
nDupli += 1 # duplicates vertices
else: # Z colinear
nZcolinear += 1
del L[i]
else:
last = L[i]
# list data type is mutable, input list will automatically update
# and doesn't need to be returned
return (nDupli, nZcolinear)
def checkEqual(lst):
return lst[1:] == lst[:-1]
class ToolsPanelDelaunay(Panel):
bl_category = "Create"
bl_label = "Delaunay Voronoi"
bl_space_type = "VIEW_3D"
bl_context = "objectmode"
bl_region_type = "TOOLS"
bl_options = {"DEFAULT_CLOSED"}
def draw(self, context):
layout = self.layout
box = layout.box()
col = box.column(align=True)
col.label("Point Cloud:")
col.operator("delaunay.triangulation")
col.operator("voronoi.tesselation")
class OBJECT_OT_TriangulateButton(Operator):
bl_idname = "delaunay.triangulation"
bl_label = "Triangulation"
bl_description = ("Terrain points cloud Delaunay triangulation in 2.5D\n"
"Needs an existing Active Mesh Object")
bl_options = {"REGISTER", "UNDO"}
@classmethod
def poll(cls, context):
obj = context.active_object
return (obj is not None and obj.type == "MESH")
def execute(self, context):
# move the check into the poll
obj = context.active_object
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
if HAS_SCIPY:
# Use scipy when present (~18 x faster)
bpy.ops.object.mode_set(mode='EDIT')
bm = bmesh.from_edit_mesh(obj.data)
points_3D = [list(v.co) for v in bm.verts]
points_2D = np.array([[v[0], v[1]] for v in points_3D])
print("Triangulate " + str(len(points_3D)) + " points...")
# Triangulate
tri = Delaunay(points_2D)
faces = tri.simplices.tolist()
# Create new mesh structure
print("Create mesh...")
bpy.ops.object.mode_set(mode='OBJECT')
mesh = bpy.data.meshes.new("TIN")
mesh.from_pydata(points_3D, [], faces)
mesh.update(calc_edges=True)
my = bpy.data.objects.new("TIN", mesh)
context.scene.objects.link(my)
my.matrix_world = obj.matrix_world.copy()
obj.select = False
my.select = True
context.scene.objects.active = my
self.report({'INFO'}, "Mesh created (" + str(len(faces)) + " triangles)")
print("Total :%s faces %s verts" % (len(faces), len(points_3D)))
return {'FINISHED'}
# Get points coodinates
r = obj.rotation_euler
s = obj.scale
mesh = obj.data
vertsPts = [vertex.co for vertex in mesh.vertices]
# Remove duplicate
verts = [[vert.x, vert.y, vert.z] for vert in vertsPts]
nDupli, nZcolinear = unique(verts)
nVerts = len(verts)
debug_prints(text=str(nDupli) + " duplicate points ignored")
debug_prints(str(nZcolinear) + " z colinear points excluded")
self.report({"WARNING"},
"Not enough points to continue. Operation Cancelled")
return {"CANCELLED"}
# Check colinear
xValues = [pt[0] for pt in verts]
yValues = [pt[1] for pt in verts]
if checkEqual(xValues) or checkEqual(yValues):
self.report({'ERROR'}, "Points are colinear")
return {'FINISHED'}
# Triangulate
debug_prints(text="Triangulate " + str(nVerts) + " points...")
vertsPts = [Point(vert[0], vert[1], vert[2]) for vert in verts]
triangles = computeDelaunayTriangulation(vertsPts)
# reverse point order --> if all triangles are specified anticlockwise then all faces up
triangles = [tuple(reversed(tri)) for tri in triangles]
debug_prints(text=str(len(triangles)) + " triangles")
# Create new mesh structure
debug_prints(text="Create mesh...")
tinMesh = bpy.data.meshes.new("TIN") # create a new mesh
tinMesh.from_pydata(verts, [], triangles) # Fill the mesh with triangles
tinMesh.update(calc_edges=True) # Update mesh with new data
# Create an object with that mesh
tinObj = bpy.data.objects.new("TIN", tinMesh)
# Place object
tinObj.location = obj.location.copy()
tinObj.rotation_euler = r
tinObj.scale = s
# Update scene
bpy.context.scene.objects.link(tinObj) # Link object to scene
bpy.context.scene.objects.active = tinObj
tinObj.select = True
obj.select = False
self.report({"INFO"},
"Mesh created (" + str(len(triangles)) + " triangles)")
return {'FINISHED'}
class OBJECT_OT_VoronoiButton(Operator):
bl_idname = "voronoi.tesselation"
bl_label = "Diagram"
bl_description = ("Points cloud Voronoi diagram in 2D\n"
"Needs an existing Active Mesh Object")
bl_options = {"REGISTER", "UNDO"}
meshType = EnumProperty(
items=[('Edges', "Edges", "Edges Only - do not fill Faces"),
('Faces', "Faces", "Fill Faces in the new Object")],
name="Mesh type",
description="Type of geometry to generate"
)
@classmethod
def poll(cls, context):
obj = context.active_object
return (obj is not None and obj.type == "MESH")
def execute(self, context):
# move the check into the poll
obj = context.active_object
# Get points coodinates
r = obj.rotation_euler
s = obj.scale
mesh = obj.data
vertsPts = [vertex.co for vertex in mesh.vertices]
# Remove duplicate
verts = [[vert.x, vert.y, vert.z] for vert in vertsPts]
nDupli, nZcolinear = unique(verts)
nVerts = len(verts)
debug_prints(text=str(nDupli) + " duplicates points ignored")
debug_prints(text=str(nZcolinear) + " z colinear points excluded")
self.report({"WARNING"},
"Not enough points to continue. Operation Cancelled")
return {"CANCELLED"}
# Check colinear
xValues = [pt[0] for pt in verts]
yValues = [pt[1] for pt in verts]
if checkEqual(xValues) or checkEqual(yValues):
self.report({"WARNING"},
"Points are colinear. Operation Cancelled")
return {"CANCELLED"}
debug_prints(text="Tesselation... (" + str(nVerts) + " points)")
xbuff, ybuff = 5, 5
zPosition = 0
vertsPts = [Point(vert[0], vert[1], vert[2]) for vert in verts]
if self.meshType == "Edges":
pts, edgesIdx = computeVoronoiDiagram(
vertsPts, xbuff, ybuff,
polygonsOutput=False, formatOutput=True
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
)
else:
pts, polyIdx = computeVoronoiDiagram(
vertsPts, xbuff, ybuff, polygonsOutput=True,
formatOutput=True, closePoly=False
)
pts = [[pt[0], pt[1], zPosition] for pt in pts]
# Create new mesh structure
voronoiDiagram = bpy.data.meshes.new("VoronoiDiagram") # create a new mesh
if self.meshType == "Edges":
# Fill the mesh with triangles
voronoiDiagram.from_pydata(pts, edgesIdx, [])
else:
# Fill the mesh with triangles
voronoiDiagram.from_pydata(pts, [], list(polyIdx.values()))
voronoiDiagram.update(calc_edges=True) # Update mesh with new data
# create an object with that mesh
voronoiObj = bpy.data.objects.new("VoronoiDiagram", voronoiDiagram)
# place object
voronoiObj.location = obj.location.copy()
voronoiObj.rotation_euler = r
voronoiObj.scale = s
# update scene
bpy.context.scene.objects.link(voronoiObj) # Link object to scene
bpy.context.scene.objects.active = voronoiObj
voronoiObj.select = True
obj.select = False
# Report
if self.meshType == "Edges":
self.report({"INFO"}, "Mesh created (" + str(len(edgesIdx)) + " edges)")
self.report({"INFO"}, "Mesh created (" + str(len(polyIdx)) + " polygons)")
return {'FINISHED'}
# Register
def register():
bpy.utils.register_class(OBJECT_OT_VoronoiButton)
bpy.utils.register_class(OBJECT_OT_TriangulateButton)
bpy.utils.register_class(ToolsPanelDelaunay)
def unregister():
bpy.utils.unregister_class(OBJECT_OT_VoronoiButton)
bpy.utils.unregister_class(OBJECT_OT_TriangulateButton)
bpy.utils.unregister_class(ToolsPanelDelaunay)
if __name__ == "__main__":
register()