鲍耶·沃森快速入门：
function BowyerWatson (pointList)
// pointList is a set of coordinates defining the points to be triangulated
triangulation := empty triangle mesh data structure
add super-triangle to triangulation // must be large enough to completely contain all the points in pointList
for each point in pointList do // add all the points one at a time to the triangulation
badTriangles := empty set
for each triangle in triangulation do // first find all the triangles that are no longer valid due to the insertion
if point is inside circumcircle of triangle
add triangle to badTriangles
polygon := empty set
for each triangle in badTriangles do // find the boundary of the polygonal hole
for each edge in triangle do
if edge is not shared by any other triangles in badTriangles
add edge to polygon
for each triangle in badTriangles do // remove them from the data structure
remove triangle from triangulation
for each edge in polygon do // re-triangulate the polygonal hole
newTri := form a triangle from edge to point
add newTri to triangulation
for each triangle in triangulation // done inserting points, now clean up
if triangle contains a vertex from original super-triangle
remove triangle from triangulation
return triangulation
相对容易地实现了python 3.7.2和pygame 1.7〜。
[此后，我在另一篇文章中发现this comment，其中指出大多数Bowyer-Watson算法都缺乏对无穷远点的顶点的计算以及克服这一问题的实现方法：
- 检查三角形的顶点是否位于无穷大处。换一种说法：检查三角形是否与边界三角形共享某些顶点。
- [如果共享所有三个顶点：琐碎。
- [如果共享零个顶点：经典方法-检查是否点到外接心的距离短于外接半径。
- 如果共享一个顶点，请检查点是否位于顶点的左侧/右侧由其他两个顶点定义的线。 one vertex in infinity
- 如果共享两个顶点，请检查点是否位于左侧/右侧：这两个顶点定义的线的线，但移至第三个点。换句话说：您仅从直线获取斜率矢量在这些共享顶点之间移动它，以便线通过通过第三点。 two vertices in infinity
因此，我创建了一个新的Point（）类，而不是内置的pygame.Vector2（），并添加了isInfinite布尔元素。
对于Triangle（）类，它接收3个Point（），然后将其分类为2个列表：finiteVertices和infiniteVertices基于每个Point中的isInfinite元素，并根据每个列表中有多少个元素来执行上述过程。
如果不把我的三角剖分变成意大利细面条，那就太好了。
我的代码是：
import pygame
import pygame.gfxdraw
import math
import random
pygame.init()
def Sign(value):
if value < 0:
return -1
if value > 0:
return 1
if value == 0:
return 0
def SideOfLineOfPoint(x1,y1,x2,y2,posX,posY):
d = (posX-x1)*(y2-y1) - (posY-y1)*(x2-x1)
return Sign(d)
def LineIsEqual(line1,line2): # Detect congruence of line, no matter which end is defined first
if (line1[0] == line2[0] and line1[1] == line2[1]) or (line1[0] == line2[1] and line1[1] == line2[0]):
return True
return False
class Point:
def __init__(self,x,y,isInfinite):
self.x = x
self.y = y
self.isInfinite = isInfinite
def distanceTo(self,other):
return math.sqrt( (self.x-other.x)**2 + (self.y-other.y)**2 )
class Triangle:
def __init__(self,a,b,c):
self.vertices = [a,b,c] # a,b,c are vertices defining the triangle
self.edges = [[a,b],
[b,c],
[c,a]] # Manually defining all edges of triangle ([])
self.CalculateCircumcenter()
self.infiniteVertices = []
self.finiteVertices = []
for vertex in self.vertices:
if vertex.isInfinite:
self.infiniteVertices.append(vertex)
else:
self.finiteVertices.append(vertex)
def CalculateCircumcenter(self): # Copied from Delaunator page
a = [self.vertices[0].x , self.vertices[0].y]
b = [self.vertices[1].x , self.vertices[1].y]
c = [self.vertices[2].x , self.vertices[2].y]
ad = a[0] * a[0] + a[1] * a[1]
bd = b[0] * b[0] + b[1] * b[1]
cd = c[0] * c[0] + c[1] * c[1]
D = 2 * (a[0] * (b[1] - c[1]) + b[0] * (c[1] - a[1]) + c[0] * (a[1] - b[1]))
self.circumcenter = Point(1 / D * (ad * (b[1] - c[1]) + bd * (c[1] - a[1]) + cd * (a[1] - b[1])),
1 / D * (ad * (c[0] - b[0]) + bd * (a[0] - c[0]) + cd * (b[0] - a[0])),
False)
def IsPointInCircumcircle(self,point):
if len(self.infiniteVertices) == 3:
return True # Any point is within the circumcircle if all therr vertices are infinite
elif len(self.infiniteVertices) == 2: # If two infinite vertices: check if point lies to the left/right of line defined by these two vertices but shifted to the third point.
x1 = self.finiteVertices[0].x
y1 = self.finiteVertices[0].y
x2 = self.infiniteVertices[0].x - self.infiniteVertices[1].x + x1
y2 = self.infiniteVertices[0].y - self.infiniteVertices[1].y + y1
sideOfLineOfVertex = SideOfLineOfPoint(x1,y1,x2,y2,point.x,point.y)
sideOfLineOfPoint = SideOfLineOfPoint(x1,y1,x2,y2,self.infiniteVertices[0].x,self.infiniteVertices[0].x)
if sideOfLineOfVertex == sideOfLineOfPoint:
return False
else:
return True
elif len(self.infiniteVertices) == 1: # If one infinite vertex: check if point lies to the left/right of line defined by the other two vertices.
x1 = self.finiteVertices[0].x
y1 = self.finiteVertices[0].y
x2 = self.finiteVertices[1].x
y2 = self.finiteVertices[1].y
sideOfLineOfVertex = SideOfLineOfPoint(x1,y1,x2,y2,point.x,point.y)
sideOfLineOfPoint = SideOfLineOfPoint(x1,y1,x2,y2,self.infiniteVertices[0].x,self.infiniteVertices[0].x)
if sideOfLineOfVertex == sideOfLineOfPoint:
return False
else:
return True
elif len(self.infiniteVertices) == 0: # For triangle with finite vertices
if self.vertices[0].distanceTo(self.circumcenter) > point.distanceTo(self.circumcenter):
return True # If point is closer to circumcenter than any vertices, point is in circumcircle
else:
return False
def HasVertex(self,point):
if point in self.vertices:
return True
return False
def Show(self,screen,colour):
for edge in self.edges:
pygame.draw.aaline(screen,colour,(edge[0].x,edge[0].y),(edge[1].x,edge[1].y))
class DelaunayTriangulation:
def __init__(self,points,width,height):
self.triangulation = [] # Create empty list
self.superTriangleA = Point(-100,-100,True)
self.superTriangleB = Point(2*width+100,-100,True)
self.superTriangleC = Point(-100,2*height+100,True)
superTriangle = Triangle(self.superTriangleA,self.superTriangleB,self.superTriangleC)
self.triangulation.append(superTriangle) # Create super-triangle
for point in points: # For every single point to be triangulated
self.addPoint(point)
def addPoint(self,point):
invalidTriangles = [] # Invalid triangle list
for triangle in self.triangulation: # For every existing triangle
if triangle.IsPointInCircumcircle(point): # If new point is in the circumcenter of triangle
invalidTriangles.append(triangle) # Triangle is invalid and added to invalid triangle list
polygon = [] # List for remaining edges after removal of invalid triangles
for triangle in invalidTriangles: # For every invalid triangle
for triangleEdge in triangle.edges: # For each invalid triangle's edges
isShared = False # Assume no edges are shared
for other in invalidTriangles: # For every other invalid triangle
if triangle == other: # If both are the same triangle
continue
for otherEdge in other.edges: # For every edge in other triangle
if LineIsEqual(triangleEdge,otherEdge):
isShared = True # Congruent edges are shared
if isShared == False: # Only append edges not shared by invalid triangles to polygon hole
polygon.append(triangleEdge)
for triangle in invalidTriangles: # Remove invalid triangles
self.triangulation.remove(triangle)
for edge in polygon:
newTriangle = Triangle(edge[0],edge[1],point) # Retriangulate based on edges of polygon hole and point
self.triangulation.append(newTriangle)
def Show(self,screen,colour):
self.shownTriangulation = self.triangulation
superTriangles = [] # List for triangles that are part of super-triangle
for triangle in self.triangulation:
if (triangle.HasVertex(self.superTriangleA) or triangle.HasVertex(self.superTriangleB) or triangle.HasVertex(self.superTriangleC)) and (triangle in self.triangulation):
superTriangles.append(triangle) # Add triangles that have any super-triangle vertex
for triangle in superTriangles:
self.triangulation.remove(triangle) # Remove super-triangles
for triangle in self.shownTriangulation:
triangle.Show(screen,colour)
background = 20,40,100
red = 255,0,0
white = 255,255,255
width = int(500)
height = int(500)
amount = int(5)
screen = pygame.display.set_mode((width,height))
screen.fill(background)
points = []
for i in range(amount):
x = random.randint(1,width-1)
y = random.randint(1,height-1)
points.append(Point(x,y,False))
delaunay = DelaunayTriangulation(points,width,height)
delaunay.Show(screen,white)
pygame.display.update()
我认为，可能引起此问题的函数是Triangle.IsPointInCircumcircle（）和SideOfLineOfPoint（），尽管同样有可能原始算法并不打算一开始就支持无限顶点计算。
如果取消整个无限顶点计算并使用正常的外接圆检测，该代码将起作用，尽管这比我的目标要落后一步。
我希望有人能在我的代码中找到任何错误来解决此问题，甚至只是将我指向正确的方向以开始调试此问题。
提前感谢。
[通过避免昂贵的math.sqrt操作可能获得最大的性能改进。
而不是比较两点之间的Euclidean distances，而是比较距离的平方：
class Point:
# [...]
def distanceToSquare(self,other):
dx = self.x-other.x
dy = self.y-other.y
return dx*dx + dy*dy
class Triangle:
# [...]
def IsPointInCircumcircle(self,point):
return (self.vertices[0].distanceToSquare(self.circumcenter) >
point.distanceToSquare(self.circumcenter))
此外，您的代码中有简单的错别字。 x分量infiniteVertices[0]两次传递给方法SideOfLineOfPoint，但是y分量丢失（在两种情况下）：
SideOfLineOfPoint(...,self.infiniteVertices[0].x,self.infiniteVertices[0].x)
此外，如果点在同一侧，则方法IsPointInCircumcircle必须返回True。您执行相反操作：
if sideOfLineOfVertex == sideOfLineOfPoint: return False else: return True