Quick primer on Bowyer-Watson:
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
This was implemented into python 3.7.2 and pygame 1.7~ with relative ease.
After this, I found this comment on another post stating that most Bowyer-Watson algorithms lack calculations for vertices at infinity and implementations to overcome this:
- Check if any of triangles vertices lies at infinity. In other words: check if triangle is sharing some vertices with bounding triangle.
- If it's sharing all three vertices: trivial.
- If it's sharing zero vertices: classical approach - check if distance from point to circumcenter is shorter than circumradius.
- If it's sharing one vertex: check if point lies to the left/right of line defined by the other two vertices. one vertex in infinity
- If it's sharing two vertices: check if point lies to the left/right of line defined by these two vertices but shifted to the third point. In other words: you take only the slope vector from line between these shared vertices and shift it so that the line passes through the third point. two vertices in infinity
Therefore, I created a new Point() class instead of the built-in pygame.Vector2(), with an added isInfinite boolean element.
For the Triangle() class, it takes in 3 Point() and sorts then into 2 lists: finiteVertices and infiniteVertices based on the isInfinite element in each Point and executes the procedure stated above depending on how many elements are in each list.
This would be great if it didn't turn my triangulation into spagetti.
My code is:
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()
In my opinion, the functions that might be causing this problem would be Triangle.IsPointInCircumcircle() and SideOfLineOfPoint(), though it is also equally likely that the original algorithm isn't meant to support infinite-vertices calculation to begin with.
The code works if the whole infinite-vertices calculation is scrapped and normal circumcircle detection is used, though it would be a step backwards from my goal.
I hope someone could find any mistake in my code that would fix this or even just point me in the right direction to begin debugging this mess.
Thanks in advance.
