Title: | Triangulation of Irregularly Spaced Data |
---|---|
Description: | A constrained two-dimensional Delaunay triangulation package providing both triangulation and generation of voronoi mosaics of irregular spaced data. Please note that most of the functions are now also covered in package interp, which is a re-implementation from scratch under a free license based on a different triangulation algorithm. |
Authors: | Albrecht Gebhardt [aut, cre, cph] (R functions), R. J. Renka [aut] (author of underlying Fortran code), Stephen Eglen [aut, cph] (contributions), Sergej Zuyev [aut, cph] (contributions), Denis White [aut, cph] (contributions) |
Maintainer: | Albrecht Gebhardt <[email protected]> |
License: | ACM | file LICENSE |
Version: | 1.3-9.2 |
Built: | 2024-11-20 06:56:36 UTC |
Source: | CRAN |
This subroutine provides for creation of a constrained
Delaunay triangulation which, in some sense, covers an
arbitrary connected region R rather than the convex hull
of the nodes. This is achieved simply by forcing the
presence of certain adjacencies (triangulation arcs) corresponding
to constraint curves. The union of triangles
coincides with the convex hull of the nodes, but triangles
in R can be distinguished from those outside of R. The
only modification required to generalize the definition of
the Delaunay triangulation is replacement of property 5
(refer to tri.mesh
by the following:
5') If a node is contained in the interior of the circumcircle of a triangle, then every interior point of the triangle is separated from the node by a constraint arc.
In order to be explicit, we make the following definitions. A constraint region is the open interior of a simple closed positively oriented polygonal curve defined by an ordered sequence of three or more distinct nodes (constraint nodes) P(1),P(2),...,P(K), such that P(I) is adjacent to P(I+1) for I = 1,...,K with P(K+1) = P(1). Thus, the constraint region is on the left (and may have nonfinite area) as the sequence of constraint nodes is traversed in the specified order. The constraint regions must not contain nodes and must not overlap. The region R is the convex hull of the nodes with constraint regions excluded.
Note that the terms boundary node and boundary arc are reserved for nodes and arcs on the boundary of the convex hull of the nodes.
The algorithm is as follows: given a triangulation which includes one or more sets of constraint nodes, the corresponding adjacencies (constraint arcs) are forced to be present (Fortran subroutine EDGE). Any additional new arcs required are chosen to be locally optimal (satisfy the modified circumcircle property).
add.constraint(tri.obj,cstx,csty,reverse=FALSE)
add.constraint(tri.obj,cstx,csty,reverse=FALSE)
tri.obj |
object of class |
cstx |
vector containing x coordinates of the constraint curve. |
csty |
vector containing y coordinates of the constraint curve. |
reverse |
if |
An new object of class "tri"
.
R. J. Renka (1996). Algorithm 751: TRIPACK: a constrained two-dimensional Delaunay triangulation package. ACM Transactions on Mathematical Software. 22, 1-8.
tri
, print.tri
, plot.tri
, summary.tri
, triangles
, convex.hull
.
# we will use the simple test data from TRIPACK: data(tritest) tritest.tr<-tri.mesh(tritest) opar<-par(mfrow=c(2,2)) plot(tritest.tr) # include all points in a big triangle: tritest.tr<-add.constraint(tritest.tr,c(-0.1,2,-0.1), c(-3,0.5,3),reverse=TRUE) # insert a small cube: tritest.tr <- add.constraint(tritest.tr, c(0.4, 0.4,0.6, 0.6), c(0.6, 0.4,0.4, 0.6), reverse = FALSE) par(opar)
# we will use the simple test data from TRIPACK: data(tritest) tritest.tr<-tri.mesh(tritest) opar<-par(mfrow=c(2,2)) plot(tritest.tr) # include all points in a big triangle: tritest.tr<-add.constraint(tritest.tr,c(-0.1,2,-0.1), c(-3,0.5,3),reverse=TRUE) # insert a small cube: tritest.tr <- add.constraint(tritest.tr, c(0.4, 0.4,0.6, 0.6), c(0.6, 0.4,0.4, 0.6), reverse = FALSE) par(opar)
This function returns some info about the cells of a voronoi mosaic, including the coordinates of the vertices and the cell area.
cells(voronoi.obj)
cells(voronoi.obj)
voronoi.obj |
object of class |
The function calculates the neighbourhood relations between the underlying triangulation and translates it into the neighbourhood relations between the voronoi cells.
retruns a list of lists, one entry for each voronoi cell which contains
cell |
cell index |
center |
cell 'center' |
neighbours |
neighbour cell indices |
nodes |
2 times |
area |
cell area |
outer cells have area=NA
, currently also nodes=NA
which is not really useful – to be done later
A. Gebhardt
data(tritest) tritest.vm <- voronoi.mosaic(tritest$x,tritest$y) tritest.cells <- cells(tritest.vm) # higlight cell 12: plot(tritest.vm) polygon(t(tritest.cells[[12]]$nodes),col="green") # put cell area into cell center: text(tritest.cells[[12]]$center[1], tritest.cells[[12]]$center[2], tritest.cells[[12]]$area)
data(tritest) tritest.vm <- voronoi.mosaic(tritest$x,tritest$y) tritest.cells <- cells(tritest.vm) # higlight cell 12: plot(tritest.vm) polygon(t(tritest.cells[[12]]$nodes),col="green") # put cell area into cell center: text(tritest.cells[[12]]$center[1], tritest.cells[[12]]$center[2], tritest.cells[[12]]$area)
This function plots circles at given locations with given radii.
circles(x, y, r, ...)
circles(x, y, r, ...)
x |
vector of x coordinates |
y |
vector of y coordinates |
r |
vactor of radii |
... |
additional graphic parameters will be passed through |
This function needs a previous plot where it adds the circles.
A. Gebhardt
x<-rnorm(10) y<-rnorm(10) r<-runif(10,0,0.5) plot(x,y, xlim=c(-3,3), ylim=c(-3,3), pch="+") circles(x,y,r)
x<-rnorm(10) y<-rnorm(10) r<-runif(10,0,0.5) plot(x,y, xlim=c(-3,3), ylim=c(-3,3), pch="+") circles(x,y,r)
Sample data for the link{circumcircle}
function.
circtest2
are points sampled from a circle with some jitter
added, i.e. they represent the most complicated case for the
link{circumcircle}
function.
This function returns the circumcircle of a triangle.
circum(x, y)
circum(x, y)
x |
Vector of three elements, giving the x coordinatres of the triangle nodes. |
y |
Vector of three elements, giving the y coordinatres of the triangle nodes. |
This is an interface to the Fortran function CIRCUM found in TRIPACK.
x |
'x' coordinate of center |
y |
'y' coordinate of center |
radius |
circumcircle radius |
signed.area |
signed area of riangle (positive iff nodes are numbered counter clock wise) |
aspect.ratio |
ratio "radius of inscribed circle"/"radius of circumcircle", varies between 0 and 0.5 0 means collinear points, 0.5 equilateral trangle. |
This function is mainly intended to be used by circumcircle
.
Fortran code: R. J. Renka, R code: A. Gebhardt
R. J. Renka (1996). Algorithm 751: TRIPACK: a constrained two-dimensional Delaunay triangulation package. ACM Transactions on Mathematical Software. 22, 1-8.
circum(c(0,1,0),c(0,0,1))
circum(c(0,1,0),c(0,0,1))
This function returns the (smallest) circumcircle of a set of n points
circumcircle(x, y = NULL, num.touch=2, plot = FALSE, debug = FALSE)
circumcircle(x, y = NULL, num.touch=2, plot = FALSE, debug = FALSE)
x |
vector containing x coordinates of the data. If |
y |
vector containing y coordinates of the data. |
num.touch |
How often should the resulting circle touch the convex hull of the given points? default: 2 possible values: 2 or 3 Note: The circumcircle of a triangle is usually defined to touch
at 3 points, this function searches by default the minimum circle,
which may be only touching at 2 points. Set parameter
|
plot |
Logical, produce a simple plot of the result. default: |
debug |
Logical, more plots, only needed for debugging. default: |
This is a (naive implemented) algorithm which determines the smallest circumcircle of n points:
First step: Take the convex hull.
Second step: Determine two points on the convex hull with maximum distance for the diameter of the set.
Third step: Check if the circumcircle of these two points already contains all other points (of the convex hull and hence all other points).
If not or if 3 or more touching points are desired
(num.touch=3
),
search a point with minimum enclosing circumcircle among the
remaining points of the convex hull.
If such a point cannot be found (e.g. for data(circtest2)
),
search the remaining triangle combinations of points from the convex
hull until an enclosing circle with minimum radius is found.
The last search uses an upper and lower bound for the desired miniumum radius:
Any enclosing rectangle and its circumcircle gives an upper bound (the axis-parallel rectangle is used).
Half the diameter of the set from step 1 is a lower bound.
x |
'x' coordinate of circumcircle center |
y |
'y' coordinate of circumcircle center |
radius |
radius of circumcircle |
Albrecht Gebhardt
data(circtest) # smallest circle: circumcircle(circtest,num.touch=2,plot=TRUE) # smallest circle with maximum touching points (3): circumcircle(circtest,num.touch=3,plot=TRUE) # some stress test for this function, data(circtest2) # circtest2 was generated by: # 100 random points almost one a circle: # alpha <- runif(100,0,2*pi) # x <- cos(alpha) # y <- sin(alpha) # circtest2<-list(x=cos(alpha)+runif(100,0,0.1), # y=sin(alpha)+runif(100,0,0.1)) # circumcircle(circtest2,plot=TRUE)
data(circtest) # smallest circle: circumcircle(circtest,num.touch=2,plot=TRUE) # smallest circle with maximum touching points (3): circumcircle(circtest,num.touch=3,plot=TRUE) # some stress test for this function, data(circtest2) # circtest2 was generated by: # 100 random points almost one a circle: # alpha <- runif(100,0,2*pi) # x <- cos(alpha) # y <- sin(alpha) # circtest2<-list(x=cos(alpha)+runif(100,0,0.1), # y=sin(alpha)+runif(100,0,0.1)) # circumcircle(circtest2,plot=TRUE)
Given a triangulation tri.obj
of points in the plane, this
subroutine returns two vectors containing the coordinates
of the nodes on the boundary of the convex hull.
convex.hull(tri.obj, plot.it=FALSE, add=FALSE,...)
convex.hull(tri.obj, plot.it=FALSE, add=FALSE,...)
tri.obj |
object of class |
plot.it |
logical, if |
add |
logical. if |
... |
additional plot arguments |
x |
x coordinates of boundary nodes. |
y |
y coordinates of boundary nodes. |
A. Gebhardt
R. J. Renka (1996). Algorithm 751: TRIPACK: a constrained two-dimensional Delaunay triangulation package. ACM Transactions on Mathematical Software. 22, 1-8.
tri
, print.tri
, plot.tri
, summary.tri
, triangles
, add.constraint
.
# rather simple example from TRIPACK: data(tritest) tr<-tri.mesh(tritest$x,tritest$y) convex.hull(tr,plot.it=TRUE) # random points: rand.tr<-tri.mesh(runif(10),runif(10)) plot(rand.tr) rand.ch<-convex.hull(rand.tr, plot.it=TRUE, add=TRUE, col="red") # use a part of the quakes data set: data(quakes) quakes.part<-quakes[(quakes[,1]<=-17 & quakes[,1]>=-19.0 & quakes[,2]<=182.0 & quakes[,2]>=180.0),] quakes.tri<-tri.mesh(quakes.part$lon, quakes.part$lat, duplicate="remove") plot(quakes.tri) convex.hull(quakes.tri, plot.it=TRUE, add=TRUE, col="red")
# rather simple example from TRIPACK: data(tritest) tr<-tri.mesh(tritest$x,tritest$y) convex.hull(tr,plot.it=TRUE) # random points: rand.tr<-tri.mesh(runif(10),runif(10)) plot(rand.tr) rand.ch<-convex.hull(rand.tr, plot.it=TRUE, add=TRUE, col="red") # use a part of the quakes data set: data(quakes) quakes.part<-quakes[(quakes[,1]<=-17 & quakes[,1]>=-19.0 & quakes[,2]<=182.0 & quakes[,2]>=180.0),] quakes.tri<-tri.mesh(quakes.part$lon, quakes.part$lat, duplicate="remove") plot(quakes.tri) convex.hull(quakes.tri, plot.it=TRUE, add=TRUE, col="red")
Identify points in a plot of "x"
with its
coordinates. The plot of "x"
must be generated with plot.tri
.
## S3 method for class 'tri' identify(x,...)
## S3 method for class 'tri' identify(x,...)
x |
object of class |
... |
additional paramters for |
an integer vector containing the indexes of the identified points.
A. Gebhardt
tri
, print.tri
, plot.tri
, summary.tri
data(tritest) tritest.tr<-tri.mesh(tritest$x,tritest$y) plot(tritest.tr) identify.tri(tritest.tr)
data(tritest) tritest.tr<-tri.mesh(tritest$x,tritest$y) plot(tritest.tr) identify.tri(tritest.tr)
Given a triangulation tri.obj
of points in the plane, this
subroutine returns a logical vector indicating if the points
are contained within the convex hull of
tri.obj
.
in.convex.hull(tri.obj, x, y)
in.convex.hull(tri.obj, x, y)
tri.obj |
object of class |
x |
vector of x-coordinates of points to locate |
y |
vector of y-coordinates of points to locate |
Logical vector.
A. Gebhardt
R. J. Renka (1996). Algorithm 751: TRIPACK: a constrained two-dimensional Delaunay triangulation package. ACM Transactions on Mathematical Software. 22, 1-8.
tri
, print.tri
, plot.tri
,
summary.tri
, triangles
,
add.constraint
, convex.hull
.
# example from TRIPACK: data(tritest) tr<-tri.mesh(tritest$x,tritest$y) in.convex.hull(tr,0.5,0.5) in.convex.hull(tr,c(0.5,-1,1),c(0.5,1,1)) # use a part of the quakes data set: data(quakes) quakes.part<-quakes[(quakes[,1]<=-10.78 & quakes[,1]>=-19.4 & quakes[,2]<=182.29 & quakes[,2]>=165.77),] q.tri<-tri.mesh(quakes.part$lon, quakes.part$lat, duplicate="remove") in.convex.hull(q.tri,quakes$lon[990:1000],quakes$lat[990:1000])
# example from TRIPACK: data(tritest) tr<-tri.mesh(tritest$x,tritest$y) in.convex.hull(tr,0.5,0.5) in.convex.hull(tr,c(0.5,-1,1),c(0.5,1,1)) # use a part of the quakes data set: data(quakes) quakes.part<-quakes[(quakes[,1]<=-10.78 & quakes[,1]>=-19.4 & quakes[,2]<=182.29 & quakes[,2]>=165.77),] q.tri<-tri.mesh(quakes.part$lon, quakes.part$lat, duplicate="remove") in.convex.hull(q.tri,quakes$lon[990:1000],quakes$lat[990:1000])
This function returns a logical vector indicating which elements of the given points P0 are left of the directed edge P1->P2.
left(x0, y0, x1, y1, x2, y2)
left(x0, y0, x1, y1, x2, y2)
x0 |
Numeric vector, 'x' coordinates of points P0 to check |
y0 |
Numeric vector, 'y' coordinates of points P0 to check, same length as 'x'. |
x1 |
'x' coordinate of point P1 |
y1 |
'y' coordinate of point P1 |
x2 |
'x' coordinate of point P2 |
y2 |
'y' coordinate of point P2 |
Logical vector.
This is an interface to the Fortran function VLEFT, wich is modeled after TRIPACKs LEFT function but accepts more than one point P0.
A. Gebhardt
left(c(0,0,1,1),c(0,1,0,1),0,0,1,1)
left(c(0,0,1,1),c(0,1,0,1),0,0,1,1)
Extract a list of neighbours from a triangulation object
neighbours(tri.obj)
neighbours(tri.obj)
tri.obj |
object of class |
nested list of neighbours per point
A. Gebhardt
R. J. Renka (1996). Algorithm 751: TRIPACK: a constrained two-dimensional Delaunay triangulation package. ACM Transactions on Mathematical Software. 22, 1-8.
tri
, print.tri
, plot.tri
, summary.tri
, triangles
data(tritest) tritest.tr<-tri.mesh(tritest$x,tritest$y) tritest.nb<-neighbours(tritest.tr)
data(tritest) tritest.tr<-tri.mesh(tritest$x,tritest$y) tritest.nb<-neighbours(tritest.tr)
Given a triangulation tri.obj
of points in the plane, this
subroutine returns a logical vector indicating if the points
lay on the convex hull of
tri.obj
.
on.convex.hull(tri.obj, x, y)
on.convex.hull(tri.obj, x, y)
tri.obj |
object of class |
x |
vector of x-coordinates of points to locate |
y |
vector of y-coordinates of points to locate |
Logical vector.
A. Gebhardt
R. J. Renka (1996). Algorithm 751: TRIPACK: a constrained two-dimensional Delaunay triangulation package. ACM Transactions on Mathematical Software. 22, 1-8.
tri
, print.tri
, plot.tri
,
summary.tri
, triangles
,
add.constraint
, convex.hull
, in.convex.hull
.
# example from TRIPACK: data(tritest) tr<-tri.mesh(tritest$x,tritest$y) on.convex.hull(tr,0.5,0.5) on.convex.hull(tr,c(0.5,-1,1),c(0.5,1,1)) # use a part of the quakes data set: data(quakes) quakes.part<-quakes[(quakes[,1]<=-10.78 & quakes[,1]>=-19.4 & quakes[,2]<=182.29 & quakes[,2]>=165.77),] q.tri<-tri.mesh(quakes.part$lon, quakes.part$lat, duplicate="remove") on.convex.hull(q.tri,quakes.part$lon[1:20],quakes.part$lat[1:20])
# example from TRIPACK: data(tritest) tr<-tri.mesh(tritest$x,tritest$y) on.convex.hull(tr,0.5,0.5) on.convex.hull(tr,c(0.5,-1,1),c(0.5,1,1)) # use a part of the quakes data set: data(quakes) quakes.part<-quakes[(quakes[,1]<=-10.78 & quakes[,1]>=-19.4 & quakes[,2]<=182.29 & quakes[,2]>=165.77),] q.tri<-tri.mesh(quakes.part$lon, quakes.part$lat, duplicate="remove") on.convex.hull(q.tri,quakes.part$lon[1:20],quakes.part$lat[1:20])
This version of outer
evaluates FUN
only on that part of the grid that is enclosed within
the convex hull of the points (
px
,py
).
This can be useful for spatial estimation if no extrapolation is wanted.
outer.convhull(cx,cy,px,py,FUN,duplicate="remove",...)
outer.convhull(cx,cy,px,py,FUN,duplicate="remove",...)
cx |
x cordinates of grid |
cy |
y cordinates of grid |
px |
vector of x coordinates of points |
py |
vector of y coordinates of points |
FUN |
function to be evaluated over the grid |
duplicate |
indicates what to do with duplicate
|
... |
additional arguments for |
Matrix with values of FUN
(NA
s if outside the
convex hull).
A. Gebhardt
x<-runif(20) y<-runif(20) z<-runif(20) z.lm<-lm(z~x+y) f.pred<-function(x,y) {predict(z.lm,data.frame(x=as.vector(x),y=as.vector(y)))} xg<-seq(0,1,0.05) yg<-seq(0,1,0.05) image(xg,yg,outer.convhull(xg,yg,x,y,f.pred)) points(x,y)
x<-runif(20) y<-runif(20) z<-runif(20) z.lm<-lm(z~x+y) f.pred<-function(x,y) {predict(z.lm,data.frame(x=as.vector(x),y=as.vector(y)))} xg<-seq(0,1,0.05) yg<-seq(0,1,0.05) image(xg,yg,outer.convhull(xg,yg,x,y,f.pred)) points(x,y)
plots the triangulation "x"
## S3 method for class 'tri' plot(x, add=FALSE,xlim=range(x$x),ylim=range(x$y), do.points=TRUE, do.labels = FALSE, isometric=FALSE,...)
## S3 method for class 'tri' plot(x, add=FALSE,xlim=range(x$x),ylim=range(x$y), do.points=TRUE, do.labels = FALSE, isometric=FALSE,...)
x |
object of class |
add |
logical, if |
do.points |
logical, indicates if points should be plotted. |
do.labels |
logical, indicates if points should be labelled |
xlim , ylim
|
x/y ranges for plot |
isometric |
generate an isometric plot (default |
... |
additional plot parameters |
None
A. Gebhardt
R. J. Renka (1996). Algorithm 751: TRIPACK: a constrained two-dimensional Delaunay triangulation package. ACM Transactions on Mathematical Software. 22, 1-8.
# random points plot(tri.mesh(rpois(100,lambda=20),rpois(100,lambda=20),duplicate="remove")) # use a part of the quakes data set: data(quakes) quakes.part<-quakes[(quakes[,1]<=-10.78 & quakes[,1]>=-19.4 & quakes[,2]<=182.29 & quakes[,2]>=165.77),] quakes.tri<-tri.mesh(quakes.part$lon, quakes.part$lat, duplicate="remove") plot(quakes.tri) # use the whole quakes data set # (will not work with standard memory settings, hence commented out) #plot(tri.mesh(quakes$lon, quakes$lat, duplicate="remove"), do.points=F)
# random points plot(tri.mesh(rpois(100,lambda=20),rpois(100,lambda=20),duplicate="remove")) # use a part of the quakes data set: data(quakes) quakes.part<-quakes[(quakes[,1]<=-10.78 & quakes[,1]>=-19.4 & quakes[,2]<=182.29 & quakes[,2]>=165.77),] quakes.tri<-tri.mesh(quakes.part$lon, quakes.part$lat, duplicate="remove") plot(quakes.tri) # use the whole quakes data set # (will not work with standard memory settings, hence commented out) #plot(tri.mesh(quakes$lon, quakes$lat, duplicate="remove"), do.points=F)
Plots the mosaic "x"
.
Dashed lines are used for outer tiles of the mosaic.
## S3 method for class 'voronoi' plot(x,add=FALSE, xlim=c(min(x$tri$x)- 0.1*diff(range(x$tri$x)), max(x$tri$x)+ 0.1*diff(range(x$tri$x))), ylim=c(min(x$tri$y)- 0.1*diff(range(x$tri$y)), max(x$tri$y)+ 0.1*diff(range(x$tri$y))), all=FALSE, do.points=TRUE, main="Voronoi mosaic", sub=deparse(substitute(x)), isometric=FALSE, ...)
## S3 method for class 'voronoi' plot(x,add=FALSE, xlim=c(min(x$tri$x)- 0.1*diff(range(x$tri$x)), max(x$tri$x)+ 0.1*diff(range(x$tri$x))), ylim=c(min(x$tri$y)- 0.1*diff(range(x$tri$y)), max(x$tri$y)+ 0.1*diff(range(x$tri$y))), all=FALSE, do.points=TRUE, main="Voronoi mosaic", sub=deparse(substitute(x)), isometric=FALSE, ...)
x |
object of class |
add |
logical, if |
xlim |
x plot ranges, by default modified to hide dummy points outside of the plot |
ylim |
y plot ranges, by default modified to hide dummy points outside of the plot |
all |
show all (including dummy points in the plot |
do.points |
logical, indicates if points should be plotted. |
main |
plot title |
sub |
plot subtitle |
isometric |
generate an isometric plot (default |
... |
additional plot parameters |
None
A. Gebhardt
R. J. Renka (1996). Algorithm 751: TRIPACK: a constrained two-dimensional Delaunay triangulation package. ACM Transactions on Mathematical Software. 22, 1-8.
voronoi
, print.voronoi
,
summary.voronoi
# plot a random mosaic plot(voronoi.mosaic(runif(100),runif(100),duplicate="remove")) # use isometric=TRUE and all=TRUE to see the complete mosaic # including extreme outlier points: plot(voronoi.mosaic(runif(100),runif(100),duplicate="remove"), all=TRUE, isometric=TRUE) # use a part of the quakes data set: data(quakes) quakes.part<-quakes[(quakes[,1]<=-17 & quakes[,1]>=-19.0 & quakes[,2]<=182.0 & quakes[,2]>=180.0),] quakes.vm<-voronoi.mosaic(quakes.part$lon, quakes.part$lat, duplicate="remove") plot(quakes.vm, isometric=TRUE) # use the whole quakes data set # (will not work with standard memory settings, hence commented out here) #plot(voronoi.mosaic(quakes$lon, quakes$lat, duplicate="remove"), isometric=TRUE)
# plot a random mosaic plot(voronoi.mosaic(runif(100),runif(100),duplicate="remove")) # use isometric=TRUE and all=TRUE to see the complete mosaic # including extreme outlier points: plot(voronoi.mosaic(runif(100),runif(100),duplicate="remove"), all=TRUE, isometric=TRUE) # use a part of the quakes data set: data(quakes) quakes.part<-quakes[(quakes[,1]<=-17 & quakes[,1]>=-19.0 & quakes[,2]<=182.0 & quakes[,2]>=180.0),] quakes.vm<-voronoi.mosaic(quakes.part$lon, quakes.part$lat, duplicate="remove") plot(quakes.vm, isometric=TRUE) # use the whole quakes data set # (will not work with standard memory settings, hence commented out here) #plot(voronoi.mosaic(quakes$lon, quakes$lat, duplicate="remove"), isometric=TRUE)
plots an voronoi.polygons
object
## S3 method for class 'voronoi.polygons' plot(x, which, color=TRUE, ...)
## S3 method for class 'voronoi.polygons' plot(x, which, color=TRUE, ...)
x |
object of class |
which |
index vector selecting which polygons to plot |
color |
logical, determines if plot should be colored, default: |
... |
additional plot arguments |
A. Gebhardt
##---- Should be DIRECTLY executable !! ---- ##-- ==> Define data, use random, ##-- or do help(data=index) for the standard data sets. data(tritest) tritest.vm <- voronoi.mosaic(tritest$x,tritest$y) tritest.vp <- voronoi.polygons(tritest.vm) plot(tritest.vp) plot(tritest.vp,which=c(1,3,5))
##---- Should be DIRECTLY executable !! ---- ##-- ==> Define data, use random, ##-- or do help(data=index) for the standard data sets. data(tritest) tritest.vm <- voronoi.mosaic(tritest$x,tritest$y) tritest.vp <- voronoi.polygons(tritest.vm) plot(tritest.vp) plot(tritest.vp,which=c(1,3,5))
Prints some information about tri.obj
## S3 method for class 'summary.tri' print(x, ...)
## S3 method for class 'summary.tri' print(x, ...)
x |
object of class |
... |
additional paramters for |
None
A. Gebhardt
R. J. Renka (1996). Algorithm 751: TRIPACK: a constrained two-dimensional Delaunay triangulation package. ACM Transactions on Mathematical Software. 22, 1-8.
tri
,tri.mesh
,
print.tri
, plot.tri
,
summary.tri
.
Prints some information about x
## S3 method for class 'summary.voronoi' print(x, ...)
## S3 method for class 'summary.voronoi' print(x, ...)
x |
object of class |
... |
additional paramters for |
None
A. Gebhardt
R. J. Renka (1996). Algorithm 751: TRIPACK: a constrained two-dimensional Delaunay triangulation package. ACM Transactions on Mathematical Software. 22, 1-8.
voronoi
,voronoi.mosaic
,
print.voronoi
, plot.voronoi
,
summary.voronoi
.
prints a adjacency list of "x"
## S3 method for class 'tri' print(x,...)
## S3 method for class 'tri' print(x,...)
x |
object of class |
... |
additional paramters for |
None
A. Gebhardt
R. J. Renka (1996). Algorithm 751: TRIPACK: a constrained two-dimensional Delaunay triangulation package. ACM Transactions on Mathematical Software. 22, 1-8.
prints a summary of "x"
## S3 method for class 'voronoi' print(x,...)
## S3 method for class 'voronoi' print(x,...)
x |
object of class |
... |
additional paramters for |
None
A. Gebhardt
R. J. Renka (1996). Algorithm 751: TRIPACK: a constrained two-dimensional Delaunay triangulation package. ACM Transactions on Mathematical Software. 22, 1-8.
voronoi
,
plot.voronoi
,
summary.voronoi
Returns some information (number of nodes, triangles, arcs, boundary
nodes and constraints) about object
.
## S3 method for class 'tri' summary(object,...)
## S3 method for class 'tri' summary(object,...)
object |
object of class |
... |
additional paramters for |
An objekt of class "summary.tri"
, to be printed by
print.summary.tri
.
It contains the number of nodes (n
), of arcs (na
), of
boundary nodes (nb
), of triangles (nt
) and constraints
(nc
).
A. Gebhardt
R. J. Renka (1996). Algorithm 751: TRIPACK: a constrained two-dimensional Delaunay triangulation package. ACM Transactions on Mathematical Software. 22, 1-8.
tri
, print.tri
, plot.tri
,
print.summary.tri
.
Returns some information about object
## S3 method for class 'voronoi' summary(object,...)
## S3 method for class 'voronoi' summary(object,...)
object |
object of class |
... |
additional parameters for |
Object of class "summary.voronoi"
.
It contains the number of nodes (nn
) and dummy nodes (nd
).
A. Gebhardt
R. J. Renka (1996). Algorithm 751: TRIPACK: a constrained two-dimensional Delaunay triangulation package. ACM Transactions on Mathematical Software. 22, 1-8.
voronoi
,voronoi.mosaic
,
print.voronoi
, plot.voronoi
,
print.summary.voronoi
.
R object that represents the triangulation of a set of 2D points,
generated by tri.mesh
or add.constraint
.
n |
Number of nodes |
x |
x coordinates of the triangulation nodes |
y |
y coordinates of the triangulation nodes |
tlist |
Set of nodal indexes which, along with |
tlptr |
Set of pointers in one-to-one
correspondence with the elements of |
tlend |
Set of pointers to adjacency lists. |
tlnew |
Pointer to the first empty location in |
nc |
number of constraints |
lc |
starting indices of constraints in |
call |
call, which generated this object |
The elements tlist
, tlptr
, tlend
and tlnew
are mainly intended for internal use in the appropriate Fortran
routines.
A. Gebhardt
R. J. Renka (1996). Algorithm 751: TRIPACK: a constrained two-dimensional Delaunay triangulation package. ACM Transactions on Mathematical Software. 22, 1-8.
tri.mesh
, print.tri
, plot.tri
, summary.tri
Return a vector of Delaunay segment lengths for the voronoi object.
The Delaunay triangles connected to sites contained in exceptions
vector are ignored (unless inverse
is TRUE, when only those
Delaunay triangles are accepted).
The exceptions
vector is provided so that sites at the border
of a region can be removed, as these tend to bias the distribution of
Delaunay segment lengths. exceptions
can be created by
voronoi.findrejectsites
.
tri.dellens(voronoi.obj, exceptions = NULL, inverse = FALSE)
tri.dellens(voronoi.obj, exceptions = NULL, inverse = FALSE)
voronoi.obj |
object of class |
exceptions |
a numerical vector |
inverse |
Logical |
A vector of Delaunay segment lengths.
S. J. Eglen
voronoi.findrejectsites
, voronoi.mosaic
,
data(tritest) tritest.vm <- voronoi.mosaic(tritest$x,tritest$y) tritest.vm.rejects <- voronoi.findrejectsites(tritest.vm, 0,1, 0, 1) trilens.all <- tri.dellens(tritest.vm) trilens.acc <- tri.dellens(tritest.vm, tritest.vm.rejects) trilens.rej <- tri.dellens(tritest.vm, tritest.vm.rejects, inverse=TRUE) par(mfrow=c(3,1)) dotchart(trilens.all, main="all Delaunay segment lengths") dotchart(trilens.acc, main="excluding border sites") dotchart(trilens.rej, main="only border sites")
data(tritest) tritest.vm <- voronoi.mosaic(tritest$x,tritest$y) tritest.vm.rejects <- voronoi.findrejectsites(tritest.vm, 0,1, 0, 1) trilens.all <- tri.dellens(tritest.vm) trilens.acc <- tri.dellens(tritest.vm, tritest.vm.rejects) trilens.rej <- tri.dellens(tritest.vm, tritest.vm.rejects, inverse=TRUE) par(mfrow=c(3,1)) dotchart(trilens.all, main="all Delaunay segment lengths") dotchart(trilens.acc, main="excluding border sites") dotchart(trilens.rej, main="only border sites")
This subroutine locates a point P=(x
,y
) relative to a triangulation
created by tri.mesh
. If P is
contained in a triangle, the three vertex indexes are
returned. Otherwise, the indexes of the rightmost and
leftmost visible boundary nodes are returned.
tri.find(tri.obj,x,y)
tri.find(tri.obj,x,y)
tri.obj |
an triangulation object |
x |
x-coordinate of the point |
y |
y-coordinate of the point |
A list with elements i1
,i2
,i3
containing
nodal indexes, in counterclockwise order,
of the vertices of a triangle containing
P=(x
,y
), or, if P is not contained in the convex
hull of the nodes, i1
indexes the
rightmost visible boundary node, i2
indexes
the leftmost visible boundary node,
and i3
= 0. Rightmost and leftmost are
defined from the perspective of P, and a
pair of points are visible from each
other if and only if the line segment
joining them intersects no triangulation
arc. If P and all of the nodes lie on a
common line, then i1
=i2
=i3
= 0 on
output.
A. Gebhardt
R. J. Renka (1996). Algorithm 751: TRIPACK: a constrained two-dimensional Delaunay triangulation package. ACM Transactions on Mathematical Software. 22, 1-8.
tri
, print.tri
, plot.tri
,
summary.tri
, triangles
,
convex.hull
data(tritest) tritest.tr<-tri.mesh(tritest$x,tritest$y) plot(tritest.tr) pnt<-list(x=0.3,y=0.4) triangle.with.pnt<-tri.find(tritest.tr,pnt$x,pnt$y) attach(triangle.with.pnt) lines(tritest$x[c(i1,i2,i3,i1)],tritest$y[c(i1,i2,i3,i1)],col="red") points(pnt$x,pnt$y)
data(tritest) tritest.tr<-tri.mesh(tritest$x,tritest$y) plot(tritest.tr) pnt<-list(x=0.3,y=0.4) triangle.with.pnt<-tri.find(tritest.tr,pnt$x,pnt$y) attach(triangle.with.pnt) lines(tritest$x[c(i1,i2,i3,i1)],tritest$y[c(i1,i2,i3,i1)],col="red") points(pnt$x,pnt$y)
This subroutine creates a Delaunay triangulation of a set of N arbitrarily distributed points in the plane referred to as nodes. The Delaunay triangulation is defined as a set of triangles with the following five properties:
1) The triangle vertices are nodes.
2) No triangle contains a node other than its vertices.
3) The interiors of the triangles are pairwise disjoint.
4) The union of triangles is the convex hull of the set of nodes (the smallest convex set which contains the nodes).
5) The interior of the circumcircle of each triangle contains no node.
The first four properties define a triangulation, and the last property results in a triangulation which is as close as possible to equiangular in a certain sense and which is uniquely defined unless four or more nodes lie on a common circle. This property makes the triangulation well-suited for solving closest point problems and for triangle-based interpolation.
The triangulation can be generalized to a constrained
Delaunay triangulation by a call to add.constraint
.
This allows for user-specified boundaries defining a nonconvex
and/or multiply connected region.
The operation count for constructing the triangulation is close to O(N) if the nodes are presorted on X or Y components. Also, since the algorithm proceeds by adding nodes incrementally, the triangulation may be updated with the addition (or deletion) of a node very efficiently. The adjacency information representing the triangulation is stored as a linked list requiring approximately 13N storage locations.
tri.mesh(x, y = NULL, duplicate = "error", jitter = 10^-12, jitter.iter = 6, jitter.random = FALSE)
tri.mesh(x, y = NULL, duplicate = "error", jitter = 10^-12, jitter.iter = 6, jitter.random = FALSE)
x |
vector containing x coordinates of the data. If |
y |
vector containing y coordinates of the data. |
duplicate |
flag indicating how to handle duplicate elements.
Possible values are: |
jitter |
Jitter of amount of Note that the jitter is not generated randomly unless
|
jitter.iter |
number of iterations to retry with jitter, amount
will be increased in each iteration by |
jitter.random |
logical, see |
An object of class "tri"
There is re-implementation of this function available in
package interp
under the same name with the same arguments.
But the return value is of a different class. So returned objects
from this function can not be used by functions of same name in
package interp
. But code written to use functions from this
package can be reused with the new package unless a constrained
triangulation is wanted. This is the only thing missing in the new
implementation.
R. J. Renka (1996). Algorithm 751: TRIPACK: a constrained two-dimensional Delaunay triangulation package. ACM Transactions on Mathematical Software. 22, 1-8.
tri
, print.tri
, plot.tri
,
summary.tri
, triangles
,
convex.hull
, neighbours
,
add.constraint
.
data(tritest) tritest.tr<-tri.mesh(tritest$x,tritest$y) tritest.tr
data(tritest) tritest.tr<-tri.mesh(tritest$x,tritest$y) tritest.tr
This function extracts a triangulation data structure
from an triangulation object created by tri.mesh
.
The vertices in the returned matrix (let's denote it with
retval
) are ordered
counterclockwise with the first vertex taken
to be the one with smallest index. Thus,
retval[i,"node2"]
and retval[i,"node3"]
are larger
than
retval[i,"node3"]
and index adjacent neighbors of
node retval[i,"node1"]
. The columns trx
and
arcx
, x=1,2,3 index the triangle and arc,
respectively, which are opposite (not shared
by) node nodex
, with trix
= 0 if
arcx
indexes a boundary arc. Vertex
indexes range from 1 to N, triangle indexes
from 0 to NT, and, if included, arc indexes
from 1 to NA = NT+N-1. The triangles are
ordered on first (smallest) vertex indexes,
except that the sets of constraint triangles
(triangles contained in the closure of a constraint
region) follow the non-constraint
triangles.
triangles(tri.obj)
triangles(tri.obj)
tri.obj |
object of class |
A matrix with columns node1
,node2
,node3
,
representing the vertex nodal indexes,
tr1
,tr2
,tr3
, representing neighboring triangle
indexes and arc1
,arc2
,arc3
reresenting arc indexes.
Each row represents one triangle.
A. Gebhardt
R. J. Renka (1996). Algorithm 751: TRIPACK: a constrained two-dimensional Delaunay triangulation package. ACM Transactions on Mathematical Software. 22, 1-8.
tri
, print.tri
, plot.tri
, summary.tri
, triangles
# use a slighlty modified version of data(tritest) data(tritest2) tritest2.tr<-tri.mesh(tritest2$x,tritest2$y) triangles(tritest2.tr)
# use a slighlty modified version of data(tritest) data(tritest2) tritest2.tr<-tri.mesh(tritest2$x,tritest2$y) triangles(tritest2.tr)
Internal tripack functions
These functions are not intended to be called by the user.
A very simply set set of points to test the tripack functions, taken
from the FORTRAN original. tritest2
is a slight modification by
adding runif(,-0.1,0.1)
random numbers to the coordinates.
R. J. Renka (1996). Algorithm 751: TRIPACK: a constrained two-dimensional Delaunay triangulation package. ACM Transactions on Mathematical Software. 22, 1-8.
An voronoi
object is created with voronoi.mosaic
x , y
|
x and y coordinates of nodes of the voronoi mosaic. Each node is a circumcircle center of some triangle from the Delaunay triangulation. |
node |
logical vector, indicating real nodes of the voronoi mosaic. These nodes are the centers of circumcircles of triangles with positive area of the delaunay triangulation. If |
n1 , n2 , n3
|
indices of neighbour nodes. Negative indices indicate dummy points as neighbours. |
tri |
triangulation object, see |
area |
area of triangle i. |
ratio |
aspect ratio (inscribed radius/circumradius) of triangle i. |
radius |
circumradius of triangle i. |
dummy.x , dummy.y
|
x and y coordinates of dummy points. They are used for plotting of unbounded tiles. |
A. Gebhardt
Computes the area of each Voronoi polygon. For some sites at the edge of the region, the Voronoi polygon is not bounded, and so the area of those sites cannot be calculated, and hence will be NA.
voronoi.area(voronoi.obj)
voronoi.area(voronoi.obj)
voronoi.obj |
object of class |
A vector of polygon areas.
S. J. Eglen
data(tritest) tritest.vm <- voronoi.mosaic(tritest$x,tritest$y) tritest.vm.areas <- voronoi.area(tritest.vm) plot(tritest.vm) text(tritest$x, tritest$y, tritest.vm.areas)
data(tritest) tritest.vm <- voronoi.mosaic(tritest$x,tritest$y) tritest.vm.areas <- voronoi.area(tritest.vm) plot(tritest.vm) text(tritest$x, tritest$y, tritest.vm.areas)
Find the sites in the Voronoi tesselation that lie at the edge of the region. A site is at the edge if any of the vertices of its Voronoi polygon lie outside the rectangle with corners (xmin,ymin) and (xmax,ymax).
voronoi.findrejectsites(voronoi.obj, xmin, xmax, ymin, ymax)
voronoi.findrejectsites(voronoi.obj, xmin, xmax, ymin, ymax)
voronoi.obj |
object of class |
xmin |
minimum x-coordinate of sites in the region |
xmax |
maximum x-coordinate of sites in the region |
ymin |
minimum y-coordinate of sites in the region |
ymax |
maximum y-coordinate of sites in the region |
A logical vector of the same length as the number of sites. If the site is a reject, the corresponding element of the vector is set to TRUE.
S. J. Eglen
This function creates a Voronoi mosaic.
It creates first a Delaunay triangulation, determines the circumcircle centers of its triangles, and connects these points according to the neighbourhood relations between the triangles.
voronoi.mosaic(x,y=NULL,duplicate="error")
voronoi.mosaic(x,y=NULL,duplicate="error")
x |
vector containing x coordinates of the data. If |
y |
vector containing y coordinates of the data. |
duplicate |
flag indicating how to handle duplicate elements.
Possible values are: |
An object of class voronoi
.
A. Gebhardt
voronoi
,voronoi.mosaic
, print.voronoi
, plot.voronoi
# example from TRIPACK: data(tritest) tritest.vm<-voronoi.mosaic(tritest$x,tritest$y) tritest.vm # use a part of the quakes data set: data(quakes) quakes.part<-quakes[(quakes[,1]<=-17 & quakes[,1]>=-19.0 & quakes[,2]<=182.0 & quakes[,2]>=180.0),] quakes.vm<-voronoi.mosaic(quakes.part$lon, quakes.part$lat, duplicate="remove") quakes.vm
# example from TRIPACK: data(tritest) tritest.vm<-voronoi.mosaic(tritest$x,tritest$y) tritest.vm # use a part of the quakes data set: data(quakes) quakes.part<-quakes[(quakes[,1]<=-17 & quakes[,1]>=-19.0 & quakes[,2]<=182.0 & quakes[,2]>=180.0),] quakes.vm<-voronoi.mosaic(quakes.part$lon, quakes.part$lat, duplicate="remove") quakes.vm
This functions extracts polygons from a voronoi.mosaic
object.
voronoi.polygons(voronoi.obj)
voronoi.polygons(voronoi.obj)
voronoi.obj |
object of class |
Returns an object of class voronoi.polygons
with unamed list
elements for each polygon. These list
elements are matrices with columns x
and y
.
Denis White
plot.voronoi.polygons
,voronoi.mosaic
##---- Should be DIRECTLY executable !! ---- ##-- ==> Define data, use random, ##-- or do help(data=index) for the standard data sets. data(tritest) tritest.vm <- voronoi.mosaic(tritest$x,tritest$y) tritest.vp <- voronoi.polygons(tritest.vm) tritest.vp
##---- Should be DIRECTLY executable !! ---- ##-- ==> Define data, use random, ##-- or do help(data=index) for the standard data sets. data(tritest) tritest.vm <- voronoi.mosaic(tritest$x,tritest$y) tritest.vp <- voronoi.polygons(tritest.vm) tritest.vp