Convex Hulls

In this tutorial, we show how we can compute convex hulls. The main method for this is to compute them directly from the triangulation, noting that for an unconstrained triangulation the boundary is the convex hull. We also provide a method for computing the convex hull from the point set directly, using the monotone chain algorithm.

Let us first demonstrate how we compute convex hulls from a triangulation. This is only possible from triangulations without a constrained boundary

using DelaunayTriangulation
using CairoMakie
using StableRNGs

We construct a triangulation of a point set.

rng = StableRNG(123)
points = randn(rng, 2, 250)
tri = triangulate(points; rng)
Delaunay Triangulation.
   Number of vertices: 250
   Number of triangles: 489
   Number of edges: 738
   Has boundary nodes: false
   Has ghost triangles: true
   Curve-bounded: false
   Weighted: false
   Constrained: false

To get the convex hull, just use get_convex_hull:

get_convex_hull(tri)
Convex hull.
   Vertices:
10-element Vector{Int64}:
 151
  96
  43
 159
  25
 147
  65
 206
 199
 151

You can also obtain the vertices directly from get_convex_hull_vertices:

get_convex_hull_vertices(tri)
10-element Vector{Int64}:
 151
  96
  43
 159
  25
 147
  65
 206
 199
 151

The vertices refer to indices of points in points, and they are given in counter-clockwise order. To obtain the convex hull directly from the points without constructing the triangulation, use convex_hull:

ch = convex_hull(points)
ch_points = [get_point(tri, i) for i in DelaunayTriangulation.get_vertices(ch)]
fig, ax, sc = lines(ch_points, color = :red, linewidth = 4)
scatter!(ax, points)
fig
Example block output

Just the code

An uncommented version of this example is given below. You can view the source code for this file here.

using DelaunayTriangulation
using CairoMakie
using StableRNGs

rng = StableRNG(123)
points = randn(rng, 2, 250)
tri = triangulate(points; rng)

get_convex_hull(tri)

get_convex_hull_vertices(tri)

ch = convex_hull(points)
ch_points = [get_point(tri, i) for i in DelaunayTriangulation.get_vertices(ch)]
fig, ax, sc = lines(ch_points, color = :red, linewidth = 4)
scatter!(ax, points)
fig

This page was generated using Literate.jl.