delaunay {RCDT} | R Documentation |
2D Delaunay triangulation
Description
Performs a (constrained) Delaunay triangulation of a set of 2d points.
Usage
delaunay(points, edges = NULL, elevation = FALSE)
Arguments
points |
a numeric matrix with two or three columns (three colums for an elevated Delaunay triangulation) |
edges |
the edges for the constrained Delaunay triangulation,
an integer matrix with two columns; |
elevation |
Boolean, whether to perform an elevated Delaunay triangulation (also known as 2.5D Delaunay triangulation) |
Value
A list. There are three possibilities. #'
-
If the dimension is 2 and
edges=NULL
, the returned value is a list with three fields:vertices
,mesh
andedges
. Thevertices
field contains the given vertices. Themesh
field is an object of classmesh3d
, ready for plotting with the rgl package. Theedges
field provides the indices of the vertices of the edges, given by the first two columns of a three-columns integer matrix. The third column, namedborder
, only contains some zeros and some ones; a border (exterior) edge is labelled by a1
. -
If the dimension is 2 and
edges
is notNULL
, the returned value is a list with four fields:vertices
,mesh
,edges
, andconstraints
. Thevertices
field contains the vertices of the triangulation. They coincide with the given vertices if the constraint edges do not intersect; otherwise there are the intersections in addition to the given vertices. Themesh
andedges
fields are similar to the previous case, the unconstrained Delaunay triangulation. Theconstraints
field is an integer matrix with two columns, it represents the constraint edges. They are not the same as the ones provided by the user if these ones intersect. If they do not intersect, then in general these are the same, but not always, in some rare corner cases. -
If
elevation=TRUE
, the returned value is a list with five fields:vertices
,mesh
,edges
,volume
, andsurface
. Thevertices
field contains the given vertices. Themesh
field is an object of classmesh3d
, ready for plotting with the rgl package. Theedges
field is similar to the previous cases. Thevolume
field provides a number, the sum of the volumes under the Delaunay triangles, that is to say the total volume under the triangulated surface. Finally, thesurface
field provides the sum of the areas of all triangles, thereby approximating the area of the triangulated surface.
Note
The triangulation can depend on the order of the points; this is shown in the examples.
Examples
library(RCDT)
# random points in a square ####
set.seed(314)
library(uniformly)
square <- rbind(
c(-1, 1), c(1, 1), c(1, -1), c(-1, -1)
)
pts_in_square <- runif_in_cube(10L, d = 2L)
pts <- rbind(square, pts_in_square)
del <- delaunay(pts)
opar <- par(mar = c(0, 0, 0, 0))
plotDelaunay(
del, type = "n", xlab = NA, ylab = NA, asp = 1,
fillcolor = "random", luminosity = "light", lty_edges = "dashed"
)
par(opar)
# the order of the points matters ####
# the Delaunay triangulation is not unique in general;
# it can depend on the order of the points
points <- cbind(
c(1, 2, 1, 3, 2, 1, 4, 3, 2, 1, 4, 3, 2, 4, 3, 4),
c(1, 1, 2, 1, 2, 3, 1, 2, 3, 4, 2, 3, 4, 3, 4, 4)
)
del <- delaunay(points)
opar <- par(mar = c(0, 0, 0, 0))
plotDelaunay(
del, type = "p", pch = 19, xlab = NA, ylab = NA, axes = FALSE,
asp = 1, lwd_edges = 2, lwd_borders = 3
)
par(opar)
# now we randomize the order of the points
set.seed(666L)
points2 <- points[sample.int(nrow(points)), ]
del2 <- delaunay(points2)
opar <- par(mar = c(0, 0, 0, 0))
plotDelaunay(
del2, type = "p", pch = 19, xlab = NA, ylab = NA, axes = FALSE,
asp = 1, lwd_edges = 2, lwd_borders = 3
)
par(opar)
# a constrained Delaunay triangulation: outer and inner dodecagons ####
# points
nsides <- 12L
angles <- seq(0, 2*pi, length.out = nsides+1L)[-1L]
points <- cbind(cos(angles), sin(angles))
points <- rbind(points, points/1.5)
# constraint edges
indices <- 1L:nsides
edges_outer <- cbind(
indices, c(indices[-1L], indices[1L])
)
edges_inner <- edges_outer + nsides
edges <- rbind(edges_outer, edges_inner)
# constrained Delaunay triangulation
del <- delaunay(points, edges)
# plot
opar <- par(mar = c(0, 0, 0, 0))
plotDelaunay(
del, type = "n", fillcolor = "yellow", lwd_borders = 2, asp = 1,
axes = FALSE, xlab = NA, ylab = NA
)
par(opar)
# another constrained Delaunay triangulation: a face ####
V <- read.table(
system.file("extdata", "face_vertices.txt", package = "RCDT")
)
E <- read.table(
system.file("extdata", "face_edges.txt", package = "RCDT")
)
del <- delaunay(
points = as.matrix(V)[, c(2L, 3L)], edges = as.matrix(E)[, c(2L, 3L)]
)
opar <- par(mar = c(0, 0, 0, 0))
plotDelaunay(
del, type="n", col_edges = NULL, fillcolor = "salmon",
col_borders = "black", col_constraints = "purple",
lwd_borders = 3, lwd_constraints = 3,
asp = 1, axes = FALSE, xlab = NA, ylab = NA
)
par(opar)