I am using the Triangle module to generate constrained Delaunay triangulations (available at http://www.lfd.uci.edu/~gohlke/pythonlibs/). In some instances, the function triangle.triangulate crashes, and windows says 'Python.exe has stopped responding'. I have tried using try/except structures as in the following example, which crashes inconsistently. I assume that a part of the triangulate process is random (see 'secondary question', below) but I am not entirely sure.
The inconsistency is rather worrying. I am on Windows.
from shapely.geometry import Polygon, MultiPolygon, LineString
import numpy as np
import triangle
import matplotlib.pyplot as plt
import matplotlib.tri as tri
import random
def meshXOR(polyRef, shape, otherVerts, otherSegs, otherHole):
verts = []
verts3 = []
segs = []
outerLength = len(polyRef[shape[0]])
for i in range(outerLength-1):#-1 because xorpoly duplicates first point into last spot
verts.append((polyRef[shape[0]][i][0],polyRef[shape[0]][i][1])) #append the point to the verts array which will be fed into the delaunay triangulator
if i == outerLength - 2:
segs.append([i,0])
else:
segs.append([i,i+1])
h = []
for cInd in shape[1]:
shift = len(verts)
innerLength = len(polyRef[cInd])
for i in range(innerLength-1):
verts.append((polyRef[cInd][i][0],polyRef[cInd][i][1]))
if i == innerLength - 2:
segs.append([i+shift,shift])
else:
segs.append([i+shift,i+1+shift])
h += list(Polygon(polyRef[cInd]).representative_point().coords)
print 'verts are: ', verts
#output: verts are: [(0.0, 5.0), (0.0, 10.0), (10.0, 10.0), (10.0, 0.0), (0.0, 0.0), (0.0, 5.0), (7.0, 3.0), (7.0, 7.0)]
print 'segs are: ', segs
#output: segs are: [[0, 1], [1, 2], [2, 3], [3, 4], [4, 0], [5, 6], [6, 7], [7, 5]]
print 'holes are: ', h
#output: holes are: [(5.25, 6.0)]
print 'verts: ', verts == otherVerts
print 'segs: ', segs == otherSegs
print 'hole: ', h == otherHole
return verts, segs, h
pA = Polygon([[0.0,0.0],[10.0,0.0],[10.0,10.0],[0.0,10.0]])
pB = Polygon([[0.0,5.0],[7.0,3.0],[7.0,7.0]])
xorPoly = pA.symmetric_difference(pB)
if xorPoly.geom_type == 'Polygon': xorPoly = MultiPolygon([xorPoly])
otherVerts = [(0.0, 5.0), (0.0, 10.0), (10.0, 10.0), (10.0, 0.0), (0.0, 0.0), (0.0, 5.0), (7.0, 3.0), (7.0, 7.0)]
otherSegs = [[0, 1], [1, 2], [2, 3], [3, 4], [4, 0], [5, 6], [6, 7], [7, 5]]
otherHole = [(5.25,6.0)]
xorPolys = []
shapes = []
for poly in xorPoly:
shapes.append([len(xorPolys), [], len(shapes)])
xorPolys.append(list(poly.exterior.coords))
for ip in poly.interiors:
shapes[-1][1].append(len(xorPolys))
xorPolys.append(list(ip.coords))
try:
verts, segs, holes = meshXOR(xorPolys, shapes[0], otherVerts, otherSegs, otherHole) # i even tried placing it here
except:
print 'failed'
if len(holes)>0:
A = dict(vertices = np.asarray(verts), segments = np.asarray(segs), holes = holes)
else:
A = dict(vertices = np.asarray(verts), segments = np.asarray(segs))
print 'about to tri'
B = triangle.triangulate(A, opts = 'pi') #this is the step that try/except doesn't work on
print 'completed tri'
try:
B_t = B["triangles"].tolist()
except:
print 'no trianlges'
if B_t != []:
cols = []
import random
for tri in B_t:
cols.append(random.random())
plt.figure()
plt.gca().set_aspect('equal')
xy = np.asarray(verts)
plt.tripcolor(xy[:,0], xy[:,1], B_t, facecolors=np.array(cols))
#for tri in B_t:
#print 'tri is: ', [verts[t] for t in tri]
plt.show()
else:
print 'no triangles'
Is there a way to do something like a 'Try/Except' structure to catch this error? Alternatively, am I doing something wrong with the triangle module?
Based on the behavior you describe, it seems like the external code is getting caught in a deadlock or an infinite loop. There is no way for you to catch this on the user side, short of spawning an extra thread that kills the current thread if the process takes too long.
However, even that doesn't really solve the problem that you are facing. I looked at the triangle library (the python version is just a wrapper of the version I am linking to), and I found on the website the following. Read through it and see if it explains why your issue is happening:
Triangle doesn't terminate, or just crashes:
Bad things can happen when triangles get so small that the distance between their vertices isn't much larger than the precision of your machine's arithmetic. If you've compiled Triangle for single-precision arithmetic, you might do better by recompiling it for double-precision. Then again, you might just have to settle for more lenient constraints on the minimum angle and the maximum area than you had planned.
You can minimize precision problems by ensuring that the origin lies inside your vertex set, or even inside the densest part of your mesh. If you're triangulating an object whose x coordinates all fall between 6247133 and 6247134, you're not leaving much floating-point precision for Triangle to work with.
Precision problems can occur covertly if the input PSLG contains two segments that meet (or intersect) at an extremely angle, or if such an angle is introduced by the -c switch. If you don't realize that a tiny angle is being formed, you might never discover why Triangle is crashing. To check for this possibility, use the -S switch (with an appropriate limit on the number of Steiner points, found by trial-and-error) to stop Triangle early, and view the output .poly file with Show Me. Look carefully for regions where dense clusters of vertices are forming and for small angles between segments. Zoom in closely, as such segments might look like a single segment from a distance.
If some of the input values are too large, Triangle may suffer a floating exception due to overflow when attempting to perform an orientation or incircle test. (Read the section on exact arithmetic.) Again, I recommend compiling Triangle for double (rather than single) precision arithmetic.
Unexpected problems can arise if you use quality meshing (-q, -a, or -u) with an input that is not segment-bounded - that is, if your input is a vertex set, or you're using the -c switch. If the convex hull of your input vertices has collinear vertices on its boundary, an input vertex that you think lies on the convex hull might actually lie just inside the convex hull. If so, an extremely thin triangle is formed by the vertex and the convex hull edge beside it. When Triangle tries to refine the mesh to enforce angle and area constraints, extremely tiny triangles may be formed, or Triangle may fail because of insufficient floating-point precision