API reference

InfinitePolygon

This is the type we use to represent dual cells / Voronoi cells. It's a possibly infinite polygon stored as an array of points with optional head and tail rays.

Methods

• segments to get the line segments for the polygon.

PointsFromMatrix(A [,i1=1,i2=2])

This implicitly extracts 2d point tuples from a matrix using row indices i1 and i2 for the coordinates. The columns of the matrix because individual points.

PointsFromMatrix(A) == vec(reinterpret(Tuple{Int,Int},A[1:2,:]))

This function can be used to transparently map a matrix into a Delaunator set of points. (There is no copying involved).

Example

A = reshape(1:20, 4, 5)
rval = triangulate(PointsFromMatrix(A))) # uses A[1:2,:]

contains(p::InfinitePolygon, pt)

Test if the infinite polygon contains a point. The point type pt must be able to be iterated to a pair of numbers

eachindex(t::AbstractDelaunatorData)

Return the indicies of each point in the dataset.

isfinite(poly)
isinfinite(poly)

Test if the polygon is infinite or finite.

Walk the corners by their region codes to find points on the bounding box to add as we move from outside the bbox back inside...

  # This is the order of codes in counter-clockwise order. 
  # Top-Left        Top         Top-Right
  # 1001     <-     1000    <-  1010
  #   |                           ^
  #   v                           |
  # 0001                         0010
  # Left                         Right 
  #   |                           ^
  #   v                           |  
  # 0101     ->    0100     ->  0110
  # Bottom Left   Bottom       Bottom right 
  # 
  # bottom right means x >= xmax (0010) and y <= ymin (0100)

Create a function to simplify type management.

Calculate upper and lowerbounds on the coordinates in points.

Get the largest distance along x or y dimension.

Create a wrapper to simplify type management.

bt, cdata = basictriangulation(IntType=Int32,] [FloatType=Float64,] [points];[sizehint=length(points),] [tol])

Allocate a basic triangulation structure and associated compute data in order to implement the Deluantor method.

This method does not actaully compute a triangulation, but only allocates the data. See triangulate or update! for the computational methods.

The triangulation data structure

These data structures are explained at https://mapbox.github.io/delaunator/ (but here, all the indices have been modified a little).

• triangles: a length T array of 3 tuples, where each tuple is a triangle
• halfedges: The halfedge index for the edge in the triangle array. The halfedges for triangle t are 3(t-1)+1, 3(t-1)+2, 3(t-2)+3. Each halfedge index gives the entry of the other halfedge.
• points a copy of the input set of points

Compute intersections between the bbox and a point and ray combo. Or on the line between two points. For a line intersection, use tmin=0, tmax=1. For a ray intersection, use tmin=0, tmax=Inf.

This will return the origin point (pt) if there are no intersections with the bbox. This will return a duplicated point if there is only a single intersection.

These return a coordinate that is guaranteed to be on the bbox, unless the return value is pt.

minpt, maxpt = bounds(t::AbstractDelaunatorData)

Return the coordinate bounds on the points. All of the points (as of the computation of the algorithm) lie within minpt <= pt <= maxpt.

cellarea(p)

Compute the area of a possibly infinite polygon.

circumcenters!(array, bt; [collinearthresh=0])

Write information on the circumcenters into array. Array must have length(bt.triangles) allocated. collinearthresh is an option that

clippedpoly(p::InfinitePolygon, bbox; [closed=true])
clippedpoly!(pts, p::InfinitePolygon, bbox)

returns an empty array if the poly is entirely outside the bounding box. Otherwise, return a set of points that represent edges of the infinite polygon clipped to the bounding box. The set of points will be closed (where the first point is equal to the last) if the closed=true option is set. The set of points may be closed even if this isn't set. Using closed=true results in simpler behavior. This is not an option on the mutating version, see below for how to get this functionality.

The mutating version will update the pts array by using

- `push!(pts, <newpt>)` 
- `last(pts)`
- `isempty(pts)`

It will return the input type pts.

Example code

# generate polygon regions for all of the dualcells 
# clipped to a 5% expansion of the point bounding box
# as a list of NaN separated paths, with all 
# polygons closed. 
using GeometryBasics
t = triangulate(rand(Point2f, 15))
ppts = Point2f[]
for i in eachindex(t)
    ind = lastindex(ppts)
    clippedpoly!(ppts, dualcell(t, i), margin_bbox(t, 0.05))
    # check if the polygon was closed... 
    if lastindex(ppts) > ind # then we added point
        if ppts[ind+1] != ppts[end] # check if we are closed
            push!(ppts, ppts[ind+1]) # close the polygon
        end 
    end 
    push!(ppts, (NaN,NaN)) # add the NaN separator 
end 

See also dualcell.

clippedpoly(p::InfinitePolygon, bbox; [closed=true])
clippedpoly!(pts, p::InfinitePolygon, bbox)

returns an empty array if the poly is entirely outside the bounding box. Otherwise, return a set of points that represent edges of the infinite polygon clipped to the bounding box. The set of points will be closed (where the first point is equal to the last) if the closed=true option is set. The set of points may be closed even if this isn't set. Using closed=true results in simpler behavior. This is not an option on the mutating version, see below for how to get this functionality.

The mutating version will update the pts array by using

- `push!(pts, <newpt>)` 
- `last(pts)`
- `isempty(pts)`

It will return the input type pts.

Example code

# generate polygon regions for all of the dualcells 
# clipped to a 5% expansion of the point bounding box
# as a list of NaN separated paths, with all 
# polygons closed. 
using GeometryBasics
t = triangulate(rand(Point2f, 15))
ppts = Point2f[]
for i in eachindex(t)
    ind = lastindex(ppts)
    clippedpoly!(ppts, dualcell(t, i), margin_bbox(t, 0.05))
    # check if the polygon was closed... 
    if lastindex(ppts) > ind # then we added point
        if ppts[ind+1] != ppts[end] # check if we are closed
            push!(ppts, ppts[ind+1]) # close the polygon
        end 
    end 
    push!(ppts, (NaN,NaN)) # add the NaN separator 
end 

See also dualcell.

dualcell(t, i)
dualcell(t, centers, i)

Return the finite or infinite polygon description of the dual cell to a given point index in the triangulation. The dual cell is the Voronoi cell to a Delaunay triangulation.

The default usage uses the circumcenters of the Delaunay triangles as the coordinates of the Voroni vertices. However, you can override this by giving another collection of centers. The number of centers must be equal to the number of triangles.

edgelines(t::Triangulation)

Return a generator that can be used with Makie's linesegments function to display the edges of the triangulation. Each edge is only drawn once.

Example

using GLMakie
t = triangulate(rand(StableRNG(1), Point2f, 10))
linesegments(collect(edgelines(rval)))

Sets h = e + f.  See the long version of my paper for details. 

If round-to-even is used (as with IEEE 754), maintains the strongly
nonoverlapping property.  (That is, if e is strongly nonoverlapping, h
will be also.)  Does NOT maintain the nonoverlapping or nonadjacent
properties.

hullpoly(t)

Return the coordinates of the convex hull suitable for plotting as a polygon.

Example

t = triangulate(rand(StableRNG(1), Point2f, 10))
f = scatter(t.points)
poly!(f.axis,collect(hullpoly(t)),color=:transparent, strokewidth=1)
f

hullvertices!(hull, bt, cdata)

After an update! call, we can use the data structures to get a simple list of vertices on the convex hull. This routine will run push!(hull, v) for each vertex on the convex hull. The order will be in order around the hull. This returns the hull variable.

This is part of the advanced interface. See triangulate which returns a more complete data structure including the hull information automatically.

Sample calls

julia> using Delaunator, StableRNGs, GeometryBasics
julia> bt, cdata = basictriangulation(pts)
julia> hull = hullvertices!(Vector{Int}, bt, cdata )

Store the first time a point occurs in halfedges into the index.

inhull(t::Triangulation, i::Integer)

Return the index of the vertex in the list of hull vertices (t.hull) if it's in the hull, or zero if the vertex is not in the hull.

isfinite(poly)
isinfinite(poly)

Test if the polygon is infinite or finite.

bbox = margin_bbox(t, [margin])
bbox = margin_bbox(t, xmargin, ymargin)

Returns a bounding box (bbox) for the triangulation after applying a margin. The default value of margin is 0.05.

Return a positive value if the points pa, pb, and pc occur in counterclockwise order a negative value if they occur in clockwise order and zero if they are collinear. The result is also a rough approximation of twice the signed area of the triangle defined by the three points.

points(t::AbstractDelaunatorData)

Return the point coordinates that were given as input to the algorithm. Note that changing these does not dynamically change the triangulation.

raystart, rayend = rays(t)

This returns arrays that are indexed by the index of a point on the convex hull. So to get the infinite rays associated with the nearest point cell use:

function rays_for_point(t)
    rs, re = rays(t) # only need to compute this one 
    hullindex = inhull(t, i)
    if hullindex > 0 
        return rs[hullindex],re[hullindex]
    else
        return 
    end 
end 

segments(p::InfinitePolygon; [dist = eps()])

Return an iterator over linesegments involved in the polygon. If the polygon is infinite, then this will not be closed and the first two points will be along the incoming ray and the last two points will be along the outgoing ray (each will be an arbitrary length along this ray)

If the polygon is finite, then it will be closed.

triangles(t::AbstractDelaunatorData)

Return the point indices for each triangle in the triangulation.

triangles(t, i)

Given a Triangulation and a point index i, return an iterator over the indices of triangles that include point i. These are returned in counter-clockwise order.

triangulate([Int32,] [FloatType=Float64,] points; [tol=eps(FloatType),])

Computes a triangulation of a set of points using the Delaunator algorithm. This is designed for quick graphics applications and speed rather than exact computational geometry.

Inputs

• points is any type that has integer indexing and length supported. In addition, p = points[i] should be a type where p[1] and p[2] are the x, y coordinates of p. Or you need to define the functions Delaunator.getX(p), Delaunator.getY(p) for your own type p.

  If you wish to use a a matrix to give the point information, PointsFromMatrix

• tol is used to determine when points are sufficiently close not to include

Return value

A triangulation, with methods to explore edges, hull points, dual cells.

See also basictriangulation

truncatedcircumcenter

For a nearly collinear triangle, then the circumcenter can be off at a point near infinity. Since the goal of this library is not computational geometry, a pragmatic choice is to truncate these wildly divergent near infinite circumcenters.

function from d3-delaunay / Voronoi.js

bt, cdata = update!(points, bt, cdata)

Reuse the memory and arrays to update the triangulation. Note that the resulting return values may be new as this routine can allocate new arrays if needed to handle the updated set of points. This implements the key step of the Delaunator method.