Delaunay 三角分割


Voronoi 與 Delaunay 三角分割是一體的兩面,在〈玩轉 p5.js〉中就談過兩者的關係。

在〈玩玩 spatial 模組〉可以看到,Scipy 的 spatial 也提供了 Delaunay,可以建立 Delaunay 三角分割,給予隨機點後,建立的 Delaunay 物件,會有個 simplices 特性,其中包含了三角形的索引:

import numpy as np
from scipy import spatial
import cadquery as cq

def polygon(points):
    return cq.Workplane().polyline(points).close()

n = 50
points = np.random.rand(n, 2)
delaunay  = spatial.Delaunay(points)

r = cq.Workplane()
for tri in delaunay.simplices:
    r.add(polygon(points[tri]))

索引可用來取得點座標,上例會繪製出以下的結果:

Delaunay 三角分割

類似地,Delaunay 也可以進行三維以上的分割,例如:

import numpy as np
from scipy import spatial
import cadquery as cq

from random import randint
from scipy.spatial import ConvexHull
from cadquery import Vector, Edge, Wire, Solid, Shell, Face

def polyhedron(points, faces):
    def _edges(vectors, face_indices):
        leng_vertices = len(face_indices)   
        return (
            Edge.makeLine(
                vectors[face_indices[i]], 
                vectors[face_indices[(i + 1) % leng_vertices]]
            ) 
            for i in range(leng_vertices)
        )

    vectors = [Vector(*p) for p in points]

    return Solid.makeSolid(
        Shell.makeShell(
            Face.makeFromWires(
                Wire.assembleEdges(
                    _edges(vectors, face_indices)
                )
            )
            for face_indices in faces
        )
    )

def hull3D(points):
    hull = ConvexHull(points)
    vertices = [points[i] for i in hull.vertices]
    v_i_lookup = {v: i for i, v in enumerate(vertices)}
    faces = [
        [v_i_lookup[points[i]] for i in face]
        for face in hull.simplices
    ]

    return polyhedron(vertices, faces)

def random_color():
    return (randint(0, 255), randint(0, 255), randint(0, 255))

# 3D Delaunay

n = 10
points = np.random.rand(n, 3)
delaunay = spatial.Delaunay(points)

r = cq.Workplane()
for tetrahedron in delaunay.simplices:
    p = hull3D([tuple(p) for p in points[tetrahedron]])
    r.add(p)

二維是三角分割,那三維呢?若三角形多一個平面外的頂點,會是四面體,四面體的外接球,絕不會包含另一個四面體外接球的球心。

若只檢視上例繪製出來的模型框線,會看到以下的成果:

Delaunay 三角分割