# dim: what does the third dimension, if present, refer to? (XYZ or XYM)
getClassDim = function(x, d, dim = "XYZ", type) {
type = toupper(type)
if (d == 2)
c("XY", type, "sfg")
else if (d == 3) {
stopifnot(dim %in% c("XYZ", "XYM"))
c(dim, type, "sfg")
} else if (d == 4)
c("XYZM", type, "sfg")
else stop(paste(d, "is an illegal number of columns for a", type))
}
valid_numeric_matrix = function(x) {
stopifnot(is.numeric(x), is.matrix(x), !anyNA(x))
}
Mtrx = function(x, dim = "XYZ", type) {
valid_numeric_matrix(x)
structure(x, class = getClassDim(x, ncol(x), dim, type))
}
# creates object of class c(dim, type, "sfg") from list x, possibly checking rings are closed
MtrxSet = function(x, dim = "XYZ", type, needClosed = FALSE) {
stopifnot(is.list(x))
if (length(x) > 0) { # list()
lapply(x, valid_numeric_matrix)
nc = unique(vapply(x, ncol, 0L))
if (length(nc) != 1)
stop("matrices have unequal numbers of columns")
NotClosed = function(y) any(y[1, ] != y[nrow(y), ])
if (needClosed && any(vapply(x, NotClosed, TRUE)))
stop("polygons not (all) closed")
structure(x, class = getClassDim(x, nc, dim, type))
} else
structure(x, class = getClassDim(x, nchar(dim), dim, type))
}
# creates object of class c(dim, type, "sfg") from list x, d, possibly checking rings are closed
MtrxSetSet = function(x, dim = "XYZ", type, needClosed = FALSE) {
stopifnot(is.list(x), vapply(x, is.list, TRUE))
if (length(x)) {
nc = unique(unlist(lapply(x, function(y) vapply(y, ncol, 0L))))
if (length(nc) != 1)
stop("matrices have unequal numbers of columns")
lapply(x, function(y) lapply(y, valid_numeric_matrix))
NotClosed = function(y) any(y[1, ] != y[nrow(y), ])
if (needClosed && any(unlist(lapply(x, function(y) vapply(y, NotClosed, TRUE)))))
stop("polygons not (all) closed")
structure(x, class = getClassDim(x, nc, dim, type))
} else
structure(x, class = getClassDim(x, nchar(dim), dim, type))
}
#return "XY", "XYZ", "XYM", or "XYZM"
Dimension = function(x) {
stopifnot(inherits(x, "sfg"))
class(x)[1]
}
## internal function to get a list of sfg POINT for st_as_sf(, coords = ...)
## src/sfg.cpp
## https://github.com/r-spatial/sf/issues/700
points_rcpp <- function(pts, gdim = "XY", ...) {
stopifnot(gdim %in% c("XY", "XYZ", "XYZM", "XYM"))
if (dim(pts)[2L] == 2L && nchar(gdim) > 2L) gdim = "XY"
stopifnot(dim(pts)[2] == nchar(gdim))
points_cpp(pts, gdim)
}
#' Create simple feature from a numeric vector, matrix or list
#'
#' Create simple feature from a numeric vector, matrix or list
#' @param x for \code{st_point}, numeric vector (or one-row-matrix) of length 2, 3 or 4; for \code{st_linestring} and \code{st_multipoint}, numeric matrix with points in rows; for \code{st_polygon} and \code{st_multilinestring}, list with numeric matrices with points in rows; for \code{st_multipolygon}, list of lists with numeric matrices; for \code{st_geometrycollection} list with (non-geometrycollection) simple feature geometry (sfg) objects; see examples below
#' @param dim character, indicating dimensions: "XY", "XYZ", "XYM", or "XYZM"; only really needed for three-dimensional points (which can be either XYZ or XYM) or empty geometries; see details
#' @name st
#' @details "XYZ" refers to coordinates where the third dimension represents altitude, "XYM" refers to three-dimensional coordinates where the third dimension refers to something else ("M" for measure); checking of the sanity of \code{x} may be only partial.
#' @return object of the same nature as \code{x}, but with appropriate class attribute set
#' @examples
#' (p1 = st_point(c(1,2)))
#' class(p1)
#' st_bbox(p1)
#' (p2 = st_point(c(1,2,3)))
#' class(p2)
#' (p3 = st_point(c(1,2,3), "XYM"))
#' pts = matrix(1:10, , 2)
#' (mp1 = st_multipoint(pts))
#' pts = matrix(1:15, , 3)
#' (mp2 = st_multipoint(pts))
#' (mp3 = st_multipoint(pts, "XYM"))
#' pts = matrix(1:20, , 4)
#' (mp4 = st_multipoint(pts))
#' pts = matrix(1:10, , 2)
#' (ls1 = st_linestring(pts))
#' pts = matrix(1:15, , 3)
#' (ls2 = st_linestring(pts))
#' (ls3 = st_linestring(pts, "XYM"))
#' pts = matrix(1:20, , 4)
#' (ls4 = st_linestring(pts))
#' outer = matrix(c(0,0,10,0,10,10,0,10,0,0),ncol=2, byrow=TRUE)
#' hole1 = matrix(c(1,1,1,2,2,2,2,1,1,1),ncol=2, byrow=TRUE)
#' hole2 = matrix(c(5,5,5,6,6,6,6,5,5,5),ncol=2, byrow=TRUE)
#' pts = list(outer, hole1, hole2)
#' (ml1 = st_multilinestring(pts))
#' pts3 = lapply(pts, function(x) cbind(x, 0))
#' (ml2 = st_multilinestring(pts3))
#' (ml3 = st_multilinestring(pts3, "XYM"))
#' pts4 = lapply(pts3, function(x) cbind(x, 0))
#' (ml4 = st_multilinestring(pts4))
#' outer = matrix(c(0,0,10,0,10,10,0,10,0,0),ncol=2, byrow=TRUE)
#' hole1 = matrix(c(1,1,1,2,2,2,2,1,1,1),ncol=2, byrow=TRUE)
#' hole2 = matrix(c(5,5,5,6,6,6,6,5,5,5),ncol=2, byrow=TRUE)
#' pts = list(outer, hole1, hole2)
#' (pl1 = st_polygon(pts))
#' pts3 = lapply(pts, function(x) cbind(x, 0))
#' (pl2 = st_polygon(pts3))
#' (pl3 = st_polygon(pts3, "XYM"))
#' pts4 = lapply(pts3, function(x) cbind(x, 0))
#' (pl4 = st_polygon(pts4))
#' pol1 = list(outer, hole1, hole2)
#' pol2 = list(outer + 12, hole1 + 12)
#' pol3 = list(outer + 24)
#' mp = list(pol1,pol2,pol3)
#' (mp1 = st_multipolygon(mp))
#' pts3 = lapply(mp, function(x) lapply(x, function(y) cbind(y, 0)))
#' (mp2 = st_multipolygon(pts3))
#' (mp3 = st_multipolygon(pts3, "XYM"))
#' pts4 = lapply(mp2, function(x) lapply(x, function(y) cbind(y, 0)))
#' (mp4 = st_multipolygon(pts4))
#' (gc = st_geometrycollection(list(p1, ls1, pl1, mp1)))
#' st_geometrycollection() # empty geometry
#' @export
st_point = function(x = c(NA_real_, NA_real_), dim = "XYZ") {
stopifnot(is.numeric(x))
if (is.matrix(x))
stopifnot(nrow(x) == 1) # because we want to be able to call rbind on points
structure(x, class = getClassDim(x, length(x), dim, "POINT"))
}
#' @name st
#' @export
st_multipoint = function(x = matrix(numeric(0), 0, 2), dim = "XYZ") Mtrx(x, dim, type = "MULTIPOINT")
#' @name st
#' @export
st_linestring = function(x = matrix(numeric(0), 0, 2), dim = "XYZ") Mtrx(x, dim, type = "LINESTRING")
#' @name st
#' @export
st_polygon = function(x = list(), dim = if(length(x)) "XYZ" else "XY") {
MtrxSet(x, dim, type = "POLYGON", needClosed = TRUE)
}
#' @name st
#' @export
st_multilinestring = function(x = list(), dim = if (length(x)) "XYZ" else "XY")
MtrxSet(x, dim, type = "MULTILINESTRING", needClosed = FALSE)
#' @name st
#' @export
st_multipolygon = function(x = list(), dim = if (length(x)) "XYZ" else "XY")
MtrxSetSet(x, dim, type = "MULTIPOLYGON", needClosed = TRUE)
#' @name st
#' @param dims character; specify dimensionality in case of an empty (NULL) geometrycollection, in which case \code{x} is the empty \code{list()}.
#' @export
st_geometrycollection = function(x = list(), dims = "XY") {
cls = vapply(x, class, rep("", 3))
if (length(cls)) {
if (!is.matrix(cls) || !is.character(cls) || nrow(cls) != 3)
stop("st_geometrycollection parameter x error: list elements should be simple features")
stopifnot(all(cls[3,] == "sfg"))
stopifnot(all(cls[2,] != "GEOMETRYCOLLECTION")) # can't recurse!
# check all dimensions are equal:
dims = unique(cls[1,])
if (length(dims) > 1)
stop(paste("multiple dimensions found:", paste(dims, collapse = ", ")))
}
structure(x, class = c(dims, "GEOMETRYCOLLECTION", "sfg")) # TODO: no Z/M/ZM modifier here??
}
POINT2MULTIPOINT = function(x, dim = "XYZ") {
if (length(x) == 3) # disambiguate Z/M:
dim = class(x)[1]
st_multipoint(matrix(unclass(x), 1), dim = dim)
}
LINESTRING2MULTILINESTRING = function(x, dim = "XYZ") {
if (ncol(x) == 3) # disambiguate Z/M:
dim = class(x)[1]
st_multilinestring(list(unclass(x)), dim = dim)
}
POLYGON2MULTIPOLYGON = function(x, dim = "XYZ") {
if (st_is_empty(x)) {
return(st_multipolygon(dim = class(x)[1]))
}
if (ncol(x[[1]]) == 3) # disambiguate Z/M:
dim = class(x)[1]
st_multipolygon(list(unclass(x)), dim = dim)
}
#' @name st
#' @param width integer; number of characters to be printed (max 30; 0 means print everything)
#' @export
print.sfg = function(x, ..., width = 0) { # avoids having to write print methods for 68 classes:
f = format(x, ..., width = width)
message(f)
invisible(f)
}
#' @name st
#' @param n integer; number of elements to be selected
#' @export
head.sfg = function(x, n = 10L, ...) {
structure(head(unclass(x), n = n, ...), class = class(x))
}
#
get_start = function(x, n = 30) {
if (is.list(x)) # recurse into first element:
structure(lapply(x, get_start, n = n), class = class(x))
else # matrix:
head(x, round(n/3))
}
#' @name st
#' @export
format.sfg = function(x, ..., width = 30) {
if (is.null(width))
width = 30
if (object.size(x) > 1000 && width > 0)
x = get_start(x, n = width)
pr = st_as_text(x, ...)
if (width > 0 && nchar(pr) > width)
paste0(substr(pr, 1, width - 3), "...")
else
pr
}
#' @export
#' @name st
#' @param ... objects to be pasted together into a single simple feature
#' @param recursive logical; ignored
#' @param flatten logical; if `TRUE`, try to simplify results; if `FALSE`, return geometrycollection containing all objects
#' @examples
#' c(st_point(1:2), st_point(5:6))
#' c(st_point(1:2), st_multipoint(matrix(5:8,2)))
#' c(st_multipoint(matrix(1:4,2)), st_multipoint(matrix(5:8,2)))
#' c(st_linestring(matrix(1:6,3)), st_linestring(matrix(11:16,3)))
#' c(st_multilinestring(list(matrix(1:6,3))), st_multilinestring(list(matrix(11:16,3))))
#' pl = list(rbind(c(0,0), c(1,0), c(1,1), c(0,1), c(0,0)))
#' c(st_polygon(pl), st_polygon(pl))
#' c(st_polygon(pl), st_multipolygon(list(pl)))
#' c(st_linestring(matrix(1:6,3)), st_point(1:2))
#' c(st_geometrycollection(list(st_point(1:2), st_linestring(matrix(1:6,3)))),
#' st_geometrycollection(list(st_multilinestring(list(matrix(11:16,3))))))
#' c(st_geometrycollection(list(st_point(1:2), st_linestring(matrix(1:6,3)))),
#' st_multilinestring(list(matrix(11:16,3))), st_point(5:6),
#' st_geometrycollection(list(st_point(10:11))))
#' @details When \code{flatten=TRUE}, this method may merge points into a multipoint structure, and may not preserve order, and hence cannot be reverted. When given fish, it returns fish soup.
c.sfg = function(..., recursive = FALSE, flatten = TRUE) {
stopifnot(! recursive)
Paste0 = function(lst) lapply(lst, unclass)
Paste1 = function(lst) do.call(c, lapply(lst, unclass))
lst = list(...)
if (flatten) {
cls = vapply(lst, function(x) class(x)[2], "")
ucls = unique(cls)
if (length(ucls) == 1) {
switch(ucls,
POINT = st_multipoint(do.call(rbind, lst)),
# CURVE = st_multicurve(Paste0(lst))
# CIRCULARSTRING = st_geometrycollection(lst), # FIXME??
LINESTRING = st_multilinestring(Paste0(lst)),
# SURFACE = st_multisurface(Paste0(lst)),
POLYGON = st_multipolygon(Paste0(lst)),
# TRIANGLE = st_geometrycollection(lst),
MULTIPOINT = st_multipoint(do.call(rbind, lst)),
MULTILINESTRING = st_multilinestring(Paste1(lst)),
# MULTICURVE = st_multicurve(Paste1(lst)),
MULTIPOLYGON = st_multipolygon(Paste1(lst)),
# MULTISURFACE = st_multisurface(Paste1(lst)),
# POLYHEDRALSURFACE = st_polyhedralsurface(Paste1(lst)),
# TIN = st_tin(Paste1(lst)),
GEOMETRYCOLLECTION = st_geometrycollection(Paste1(lst)),
stop(paste("type", cls, "not supported"))
)
} else if (all(ucls %in% c("POINT", "MULTIPOINT")))
st_multipoint(do.call(rbind, lst))
else if (all(cls %in% c("LINESTRING", "MULTILINESTRING"))) {
ls = which(cls == "LINESTRING")
mls = st_multilinestring(lst[ls])
st_multilinestring(c(unlist(lst[-ls], FALSE), unclass(mls)))
} else if (all(cls %in% c("POLYGON", "MULTIPOLYGON"))) {
po = which(cls == "POLYGON")
mpo = st_multipolygon(lst[po])
st_multipolygon(c(unlist(lst[-po], FALSE), unclass(mpo)))
} else {
# unfold GC objects first, then
gc = (cls == "GEOMETRYCOLLECTION")
ret = lst[!gc]
if (any(gc)) { # append the _contents_ of GC's to the non-GC elements:
wgc = which(gc)
for (i in seq_len(length(wgc)))
ret = append(ret, lst[[wgc[i]]])
}
st_geometrycollection(ret)
}
} else # !flatten:
st_geometrycollection(lst) # breaks if one of them is a GC
}
#' @name st
#' @method as.matrix sfg
#' @export
#' @return as.matrix returns the set of points that form a geometry as a single matrix, where each point is a row; use \code{unlist(x, recursive = FALSE)} to get sets of matrices.
as.matrix.sfg = function(x, ...) {
switch(class(x)[2],
POINT = matrix(x, 1),
MULTIPOINT = as.matrix(unclass(x)),
LINESTRING = as.matrix(unclass(x)),
POLYGON = do.call(rbind, x),
MULTILINESTRING = do.call(rbind, x),
MULTIPOLYGON = do.call(rbind, lapply(x, function(y) do.call(rbind, y))),
GEOMETRYCOLLECTION = do.call(rbind, lapply(x, as.matrix)),
NextMethod()
)
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.