Clipped Voronoi Tessellations
Clipping to a Rectangle
In the previous tutorial, we demonstrated how we can clip to the convex hull of the point set. However, it is often useful to clip to a rectangle, for example if you want to clip to a region of interest in a simulation. Here we obtain the coordinates just by looping over all the polygons, but you can also use this by providing a clip_polygon
into voronoi
, as described in this tutorial.
Let us now demonstrate. First, we construct a tessellation of some example point set.
using DelaunayTriangulation
using CairoMakie
A = (-3.0, 7.0)
B = (1.0, 6.0)
C = (-1.0, 3.0)
D = (-2.0, 4.0)
E = (3.0, -2.0)
F = (5.0, 5.0)
G = (-4.0, -3.0)
H = (3.0, 8.0)
points = [A, B, C, D, E, F, G, H]
tri = triangulate(points)
vorn = voronoi(tri)
Voronoi Tessellation.
Number of generators: 8
Number of polygon vertices: 9
Number of polygons: 8
Weighted: false
Let us show the tessellation, and the rectangle we want to clip the tessellation to.
fig, ax, sc = voronoiplot(vorn)
a, b, c, d = -2.0, 3.0, 0.0, 7.0
lines!(ax, [(a, c), (b, c), (b, d), (a, d), (a, c)], color = :black, linewidth = 4)
fig
To apply this clipping, we need to provide a bounding box of the form (xmin, xmax, ymin, ymax)
. Here, we will use
bounding_box = (a, b, c, d)
(-2.0, 3.0, 0.0, 7.0)
You can obtain some reasonable defaults for this bounding box using DelaunayTriangulation.polygon_bounds(vorn). The coordinates for each polygon clipped to this box can be obtained as follows.
clipped_coords = Vector{Vector{NTuple{2, Float64}}}(undef, num_polygons(vorn))
for i in each_polygon_index(vorn)
clipped_coords[i] = get_polygon_coordinates(vorn, i, bounding_box)
end
clipped_coords
8-element Vector{Vector{Tuple{Float64, Float64}}}:
[(-2.000000000000001, 7.0), (-2.0, 5.666666666666667), (-1.1363636363636365, 5.954545454545455), (-0.8750000000000009, 7.0), (-2.000000000000001, 7.0)]
[(-0.30000000000000004, 4.7), (2.357142857142857, 2.9285714285714284), (3.0, 5.499999999999999), (3.0, 6.0), (2.0, 7.0), (-0.8750000000000009, 7.0), (-1.1363636363636365, 5.954545454545455), (-0.30000000000000004, 4.7)]
[(0.375, 0.0), (2.710526315789474, 1.868421052631579), (2.357142857142857, 2.9285714285714284), (-0.30000000000000004, 4.7), (-2.0, 3.0), (-2.0, 0.0), (0.375, 0.0)]
[(-2.0, 3.0), (-0.30000000000000004, 4.7), (-1.1363636363636365, 5.954545454545455), (-2.0, 5.666666666666667), (-2.0, 3.0)]
[(3.0, 0.0), (3.0, 1.7857142857142858), (2.710526315789474, 1.868421052631579), (0.375, 0.0), (3.0, 0.0)]
[(3.0, 1.7857142857142858), (3.0, 5.499999999999999), (2.357142857142857, 2.9285714285714284), (2.710526315789474, 1.868421052631579), (3.0, 1.7857142857142858)]
[]
[(3.0, 7.0), (2.0, 7.0), (3.0, 6.0), (3.0, 7.0)]
Now let's plot these.
fig, ax, sc = poly(clipped_coords, color = :white, strokewidth = 4)
fig
As we can see, the polygons have been clipped to the rectangle. Note that if you just want this for plotting, you can also call voronoiplot
with the bounding_box
keyword argument.
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
A = (-3.0, 7.0)
B = (1.0, 6.0)
C = (-1.0, 3.0)
D = (-2.0, 4.0)
E = (3.0, -2.0)
F = (5.0, 5.0)
G = (-4.0, -3.0)
H = (3.0, 8.0)
points = [A, B, C, D, E, F, G, H]
tri = triangulate(points)
vorn = voronoi(tri)
fig, ax, sc = voronoiplot(vorn)
a, b, c, d = -2.0, 3.0, 0.0, 7.0
lines!(ax, [(a, c), (b, c), (b, d), (a, d), (a, c)], color = :black, linewidth = 4)
fig
bounding_box = (a, b, c, d)
clipped_coords = Vector{Vector{NTuple{2, Float64}}}(undef, num_polygons(vorn))
for i in each_polygon_index(vorn)
clipped_coords[i] = get_polygon_coordinates(vorn, i, bounding_box)
end
clipped_coords
fig, ax, sc = poly(clipped_coords, color = :white, strokewidth = 4)
fig
This page was generated using Literate.jl.