从凹多边形生成三角形网格的高效算法

Question

我正在从事一个项目，该项目涉及从潜在的凹多边形生成三角形网格。与最佳解决方案相比，我对高性能解决方案更感兴趣。我遇到的大部分内容都涉及使用 OpenGL（不是一个选项）或专注于最佳或近似最佳地分解为凸多边形，使用 O(n³) 或 O(n⁴)。剪耳似乎是最直接的解决方案，应该是 O(n²):

for /* ever */ {
    for each(vertex) {
        candidate := Triangle{previousVertex, currentVertex, nextVertex}
        if !candidate.IsInterior() || candidate.ContainsAnyOtherVertex() {
            continue
        }

        result.Push(candidate)
        vertex = vertex.Except(currentVertex)
    }

    if len(vertex) == 3 {
        result.Push(Triangle(vertex))
        break
    }
}

我不明白这怎么可能是 O(n²)。如我所见，这将必须涉及三个嵌套循环，每个循环都与 N 成比例：上面的两个显式循环和 candidate.ContainsAnyOtherVertex()。我错过了什么吗？有没有办法检测候选三角形是否包含任何其他剩余顶点而不涉及迭代所有剩余顶点？

或者，我应该使用其他算法吗？我对分解成凸多边形的一个很好，但是 A) 我不在乎解决方案是否是最优的，并且 B) 我对阅读 CGAL 的源代码或学术文献不感兴趣。

Answer 1

这样剪耳是O(n²):

首先，对于每个角度 < 180 度的顶点，计算其角度内的顶点数。 (O(n²))
您可以剪裁 + 移除任何此类顶点数为 0 的顶点。您将执行此操作 O(n) 次。当你移除一个 vetex 时：
- 首先将它从它所在的任何三角形的计数中删除 (O(n))；和
- 计算您通过裁剪形成的至多 2 个新三角形中的其他顶点（也为 O(n)）。

Answer 2

对于遇到此问题的任何其他人，这是我正在使用的内容的浓缩版：

package mesh

import "sort"

type Vec []float32

type Face []Vec

type FaceTri struct {
    face  Face
    index [3]int
}

func triangulate(face Face, offset int) [][3]int {
    normal := face.Normal()

    index := make([]int, len(face))
    convex := map[int]bool{}
    reflex := map[int]bool{}
    ear := map[int]bool{}

    // Mark convex and reflex
    for i := range index {
        index[i] = i
        tri := face.tri(i, i+1, i+2)

        // Skip colinear vertices
        if tri.Area().Len() == 0 {
            continue
        }

        if tri.Area().Dot(normal) > 0 {
            convex[tri.index[1]] = true
        } else {
            reflex[tri.index[1]] = true
        }
    }

    // Mark ears
    for i := range convex {
        if isEar(face.tri(i-1, i, i+1), reflex) {
            ear[i] = true
        }
    }

    var tris [][3]int
    for len(reflex) > 0 || len(index) > 3 {
        ni := len(index)

        // Next ear
        var j int
        for j = range ear {
            break
        }

        // Find J in the list of unclipped vertices and get 2 neighbors to each side
        x := sort.SearchInts(index, j)
        h, i, k, l := (ni+x-2)%ni, (ni+x-1)%ni, (x+1)%ni, (x+2)%ni
        h, i, k, l = index[h], index[i], index[k], index[l]

        // Clip (I,J,K)
        index = append(index[:x], index[x+1:]...)
        tris = append(tris, [3]int{offset + i, offset + j, offset + k})
        delete(ear, j)
        delete(convex, j)

        // Update prior vertex
        update(face.tri(h, i, k), normal, reflex, convex, ear)

        // Update later vertex
        if h != l {
            update(face.tri(i, k, l), normal, reflex, convex, ear)
        }
    }

    tris = append(tris, [3]int{
        offset + index[0],
        offset + index[1],
        offset + index[2],
    })
    return tris
}

func update(tri *FaceTri, faceNormal Vec, reflex, convex, ear map[int]bool) {
    idx := tri.index[1]
    wasConvex, wasEar := convex[idx], ear[idx]

    isConvex := wasConvex || tri.Area().Dot(faceNormal) > 0
    if !wasConvex && isConvex {
        convex[idx] = true
        delete(reflex, idx)
    }

    if !wasEar && isConvex && isEar(tri, reflex) {
        ear[idx] = true
    }
}

func isEar(tri *FaceTri, reflex map[int]bool) bool {
    // It is sufficient to only check reflex vertices - a convex vertex
    // cannot be contained without a reflex vertex also being contained
    for j := range reflex {
        if tri.HasIndex(j) {
            continue
        }

        if tri.ContainsPoint(tri.face[j]) {
            return false
        }
    }

    return true
}

// Euclidean length of V
func (v Vec) Len() float32

// Dot product of A and B
func (a Vec) Dot(b Vec) float32

// Calculates the normal vector to the face - for concave polygons, this
// is the summation of the normals of each vertex (normalized)
func (Face) Normal() Vec

// Tri returns a FaceTri for I, J, K, modulo len(f)
func (Face) tri(i, j, k int) *FaceTri

// Area returns the cross product of the vector from v1 to v0 and the vector
// from v1 to v2
func (*FaceTri) Area() Vec

// Returns true if v is in f.index
func (f *FaceTri) HasIndex(v int) bool

// Returns true if P is contained within the described triangle
func (*FaceTri) ContainsPoint(p Vec) bool

这是基于 https://www.geometrictools.com/Documentation/TriangulationByEarClipping.pdf 的描述。我认为这与马特描述的本质上是一样的。

为了清楚和简洁起见，我省略了一些方法的内容，以及我以降低可读性的方式重用结果的情况。在我看来，我遗漏的唯一有趣的部分是 ContainsPoint。在 2D 中，确定三角形 ABC 是否包含点 P 很简单。在 3D 中，情况并非如此，因为 P 不一定与 ABC 共面：

V = vector from B to A
U = vector from B to C
Q = vector from B to P

M = MatrixFromColumns(V, U, V x U)
X = (inverse of M) * Q

Q === (V * X[0]) + (U * X[1]) + ((V x U) * X[2])

If X[2] is zero, then Q has a component out of the plane of the triangle

P is in ABC if and only if X[0] and X[1] are positive and X[0] + X[1] <= 1

从凹多边形生成三角形网格的高效算法

Efficient algorithm for generating a triangle mesh from a concave polygon

algorithm

graphics

concave