convhulln {geometry} | R Documentation |
Returns an index matrix to the points of simplices (“triangles”) that form the smallest convex simplical complex of a set of input points in N-dimensional space. This function interfaces the Qhull library.
convhulln(p, options = "Tv", output.options = NULL, return.non.triangulated.facets = FALSE)
p |
An M-by-N matrix. The rows of |
options |
String containing extra options for the underlying Qhull command; see details below and Qhull documentation at ../doc/qhull/html/qconvex.html#synopsis. |
output.options |
String containing Qhull options to control
output. Currently |
return.non.triangulated.facets |
logical defining whether the
output facets should be triangulated; |
For silent operation, specify the option Pp
.
If return.non.triangulated.facets
is FALSE
(default), return an M-by-N index matrix of which
each row defines an N-dimensional “triangle”.
If return.non.triangulated.facets
is TRUE
then the
number of columns equals the maximum number of vertices in a
facet, and each row defines a polygon corresponding to a facet
of the convex hull with its vertices followed by NA
s
until the end of the row. The indices refer to the rows in
p
.
If the output.options
or options
argument contains
FA
or n
, return a list with class convhulln
comprising the named elements:
hull
The convex hull, represented as a matrix, as described above
area
If FA
is specified, the generalised area of
the hull. This is the surface area of a 3D hull or the length of
the perimeter of a 2D hull.
See ../doc/qhull/html/qh-optf.html#FA.
vol
If FA
is specified, the generalised volume of
the hull. This is volume of a 3D hull or the area of a 2D hull.
See ../doc/qhull/html/qh-optf.html#FA.
normals
If n
is specified, this is a matrix
hyperplane normals with offsets. See ../doc/qhull/html/qh-opto.html#n.
This is a port of the Octave's (http://www.octave.org) geometry library. The Octave source was written by Kai Habel.
See further notes in delaunayn
.
Raoul Grasman, Robert B. Gramacy, Pavlo Mozharovskyi and David Sterratt david.c.sterratt@ed.ac.uk
Barber, C.B., Dobkin, D.P., and Huhdanpaa, H.T., “The Quickhull algorithm for convex hulls,” ACM Trans. on Mathematical Software, Dec 1996.
convex.hull
, delaunayn
,
surf.tri
, distmesh2d
, intersectn
# example convhulln # ==> see also surf.tri to avoid unwanted messages printed to the console by qhull ps <- matrix(rnorm(3000), ncol=3) # generate points on a sphere ps <- sqrt(3)*ps/drop(sqrt((ps^2) %*% rep(1, 3))) ts.surf <- t(convhulln(ps)) # see the qhull documentations for the options ## Not run: rgl.triangles(ps[ts.surf,1],ps[ts.surf,2],ps[ts.surf,3],col="blue",alpha=.2) for(i in 1:(8*360)) rgl.viewpoint(i/8) ## End(Not run)