A (long) while ago I posted about "soon-to-be-released" JTS code for polygon triangulation using Ear Clipping. It turned out it was actually in the category of "never-to-be-released". However, later I worked with a student, Dan Tong, on a coding exercise sponsored by Facebook. We decided to tackle polygon triangulation. By the end of the project he implemented a functional algorithm, including handling holes (via the hole-joining technique described by Eberly). To make it release-ready the code needed performance and structural improvements. This has finally happened and the code is out!
The new codebase has been substantially rewritten to improve modularity and the API. A new data structure to represent triangulations has been introduced, called a Tri It is simpler, more memory efficient, and easier to work with than the QuadEdge data structure previously used in JTS. It provides an excellent basis for developing further algorithms based on triangulations (of which more below)
Most of the reworking involved implementing key performance optimizations, utilizing existing and new JTS spatial index structures:
- the Hole Joining phase now uses a spatial index to optimize the test for whether a join line is internal to a polygon
- the Ear-Clipping phase has to check whether a candidate ear contains any vertices of the polygon shell being clipped. This has been made performant by using a new spatial index structure called the VertexSequencePackedRtree
|Dataset||Size||JTS Time||Earcut4j Time|
|Lake Superior||3,478 vertices / 28 holes||18 ms||13 ms|
|Canada||10,522 vertices||35 ms||26 ms|
|Complex polygon||44,051 vertices / 125 holes||1.1 s||0.64 s|
|North America with Lakes||115,206 vertices / 1,660 holes||13.3 s||5.4 s|
So JTS ear clipping is not quite as fast as Earcut4j. But it's still plenty fast enough for production use!
The triangulation algorithm is exposed as the PolygonTriangulator class (deliberately named after the effect, rather than the algorithm producing it). Here's an example of the output:
Constrained Delaunay Triangulation
The output of PolygonTriangulator is guaranteed to be a valid triangulation, but in the interests of performance it does not attempt to produce a high-quality one. Triangulation quality is measured in terms of minimizing the number of "skinny" triangles (or equivalently, maximizing the sum of the interior angles). Optimal quality is provided by Constrained Delaunay Triangulation (CDT). In addition to being optimal in terms of triangle shape, the CDT has the nice properties of being (essentially) unique, and of accurately representing the structure of the input polygon. It also removes the influence that Hole Joining holes has on the raw triangulation provided by Ear-Clipping.
I originally proposed that a CDT can be achieved by a Delaunay Improvement algorithm based on iterated triangle flipping. It turned out that this technique is not only effective, but also pleasingly performant. This is implemented in the ConstrainedDelaunayTriangulator class. Running it on the example above shows the improvement obtained. Note how there are fewer narrow triangles, and their arrangement more closely represents the structure of the polygon.
This code will appear in the next JTS version (1.18.3 or 1.19), and has already been released in GEOS 3.10.
It seems possible to adapt the machinery for Ear Clipping to provide two other useful constructions:
- Polygon Concave Hull ("Outer Hull") - this requires a triangulation of the area between the polygon and its convex hull. That can be obtained by running Ear Clipping outwards, rather than inwards. A Concave Hull can be constructed by limiting the number of triangles generated based on some criteria such as area, longest edge, or number of vertices. The result is a simplified polygon with fewer vertices which is guaranteed to contain the original polygon.
- Polygon Inner Hull - the triangulation constructed by Ear Clipping can be "eroded" to provide a simpler polygon which is guaranteed to be contained in the original.
- Approximate Skeleton construction (also known as the Medial Axis) by creating edges through the CDT triangles
- Equal-Area Polygon Subdivision, by partitioning the graph induced by the CDT
- Concave Hull of Point Set
- Constrained Delaunay Triangulation of Lines and Points