move POINT geometries out of a list column, into a matrix

NAMESPACE

@@ -8,6 +8,7 @@ S3method("$<-",sf)
 S3method("[",sf)
 S3method("[",sfc)
 S3method("[<-",sfc)
+S3method("[[",sfc)
 S3method("[[<-",sf)
 S3method("st_agr<-",sf)
 S3method("st_crs<-",bbox)

R/arith.R

@@ -128,6 +128,11 @@ Ops.sfc <- function(e1, e2) {
 	if (length(e1) == 0) # empty set
 		return(e1)
 
+	if (inherits(e1, "sfc")) # realize
+		e1 = e1[]
+	if (inherits(e2, "sfc"))
+		e2 = e2[]
+
 	if (is.numeric(e2) && !is.matrix(e2) && length(e2) <= 2 && .Generic %in% c("+", "-")) {
 		if (.Generic == "-")
 			e2 <- -e2

R/bbox.R

@@ -9,6 +9,10 @@ bb_wrap = function(bb) {
 	structure(as.double(bb), names = c("xmin", "ymin", "xmax", "ymax"), class = "bbox")
 }
 
+bbox.pointmatrix = function(obj, ...) {
+	bb_wrap(as.vector(t(apply(obj[,1:2,drop=FALSE], 2, range))))
+}
+
 bbox.Set = function(obj, ...) {
 	sel = !sfc_is_empty(obj)
 	if (! any(sel))
@@ -134,7 +138,9 @@ print.bbox = function(x, ...) {
 }
 
 compute_bbox = function(obj) {
-	switch(class(obj)[1],
+	if (!is.null(pts <- attr(obj, "points")))
+		bbox.pointmatrix(pts)
+	else switch(class(obj)[1],
 		sfc_POINT = bb_wrap(bbox.Set(obj)),
 		sfc_MULTIPOINT = bb_wrap(bbox.MtrxSet(obj)),
 		sfc_LINESTRING = bb_wrap(bbox.MtrxSet(obj)),

R/bind.R

@@ -29,7 +29,11 @@ rbind.sf = function(..., deparse.level = 1) {
 		if (!all(equal_crs))
 			stop("arguments have different crs", call. = FALSE)
 	}
-	ret = st_sf(rbind.data.frame(...), crs = crs0, sf_column_name = sf_column)
+	for (i in seq_along(dots)) {
+		if (all(sapply(unclass(st_geometry(dots[[i]])), is.null)))
+			st_geometry(dots[[i]]) = st_geometry(dots[[i]])[] # realize
+	}
+	ret = st_sf(do.call(rbind.data.frame, dots), crs = crs0, sf_column_name = sf_column)
 	st_geometry(ret) = st_sfc(st_geometry(ret)) # might need to reclass to GEOMETRY
 	bb = do.call(rbind, lapply(dots, st_bbox))
 	bb = bb_wrap(c(min(bb[,1L], na.rm = TRUE), min(bb[,2L], na.rm = TRUE),

R/geom-transformers.R

@@ -384,7 +384,7 @@ st_minimum_rotated_rectangle.sf = function(x, dTolerance, ...) {
 #'  n = 100
 #'  pts = st_as_sf(data.frame(matrix(runif(n), , 2), id = 1:(n/2)), coords = c("X1", "X2"))
 #'  # compute Voronoi polygons:
-#'  pols = st_collection_extract(st_voronoi(do.call(c, st_geometry(pts))))
+#'  pols = st_collection_extract(st_voronoi(st_combine(pts)))
 #'  # match them to points:
 #'  pts$pols = pols[unlist(st_intersects(pts, pols))]
 #'  plot(pts["id"], pch = 16) # ID is color
@@ -651,8 +651,14 @@ st_segmentize.sf = function(x, dfMaxLength, ...) {
 #' @examples
 #' nc = st_read(system.file("shape/nc.shp", package="sf"))
 #' st_combine(nc)
-st_combine = function(x)
-	st_sfc(do.call(c, st_geometry(x)), crs = st_crs(x)) # flatten/merge
+st_combine = function(x) {
+	x = st_geometry(x)
+	if (inherits(x, "sfc_POINT") && !is.null(pts <- attr(x, "points"))) {
+		mp = structure(list(st_multipoint(pts, attr(x, "point_dim"))), bbox = st_bbox(x))
+		st_sfc(mp, crs = st_crs(x))
+	} else
+		st_sfc(do.call(c, x), crs = st_crs(x)) # flatten/merge
+}
 
 # x: object of class sf
 # y: object of class sf or sfc

R/m_range.R

@@ -9,6 +9,10 @@ mb_wrap = function(mb) {
 	structure(mb, names = c("mmin", "mmax"), class = "m_range")
 }
 
+m_range.pointmatrix = function(obj, ...) {
+	mb_wrap(range(obj[,ncol(obj)]))
+}
+
 m_range.Set = function(obj, ...) {
 	sel = vapply(obj, function(x) { length(x) && !all(is.na(x)) }, TRUE)
 	if (! any(sel))
@@ -130,7 +134,9 @@ print.m_range = function(x, ...) {
 }
 
 compute_m_range = function(obj) {
-	switch(class(obj)[1],
+	if (!is.null(pts <- attr(obj, "points")))
+		m_range.pointmatrix(pts)
+	else switch(class(obj)[1],
 		   sfc_POINT = mb_wrap(m_range.Set(obj)),
 		   sfc_MULTIPOINT = mb_wrap(m_range.MtrxSet(obj)),
 		   sfc_LINESTRING = mb_wrap(m_range.MtrxSet(obj)),

R/plot.R

@@ -312,7 +312,8 @@ plot.sfc_POINT = function(x, y, ..., pch = 1, cex = 1, col = 1, bg = 0, lwd = 1,
 	col = rep(col, length.out = npts)
 	bg = rep(bg, length.out = npts)
 	cex = rep(cex, length.out = npts)
-	mat = t(matrix(unlist(x, use.names = FALSE), ncol = length(x))) #933
+	#mat = t(matrix(unlist(x, use.names = FALSE), ncol = length(x))) #933
+	mat = st_coordinates(x)
 	if (!is.null(mat)) {
 		ne = !is.na(rowMeans(mat))  ## faster than apply; #933
 		points(mat[ne,, drop = FALSE], pch = pch[ne], col = col[ne], bg = bg[ne],

R/sf.R

@@ -48,14 +48,22 @@ st_as_sf.data.frame = function(x, ..., agr = NA_agr_, coords, wkt,
 		if (na.fail && any(is.na(cc)))
 			stop("missing values in coordinates not allowed")
 		# classdim = getClassDim(rep(0, length(coords)), length(coords), dim, "POINT")
-		x$geometry = structure( points_rcpp(as.matrix(cc), dim),
+		# x$geometry = structure( points_rcpp(attr(x, "points"), dim),
+		if (length(coords) == 2)
+			dim = "XY"
+		stopifnot(length(coords) == nchar(dim), dim %in% c("XY", "XYZ", "XYZM", "XYM"))
+		points = as.matrix(cc)
+		dimnames(points) = NULL
+		x$geometry = structure(vector("list", length = nrow(cc)),
+			points = points,
+			points_dim = dim,
 			n_empty = 0L, precision = 0, crs = NA_crs_,
 			bbox = structure(
 				c(xmin = min(cc[[1]], na.rm = TRUE),
 				ymin = min(cc[[2]], na.rm = TRUE),
 				xmax = max(cc[[1]], na.rm = TRUE),
 				ymax = max(cc[[2]], na.rm = TRUE)), class = "bbox"),
-			class =  c("sfc_POINT", "sfc" ), names = NULL)
+			class =  c("sfc_POINT", "sfc"), names = NULL)
 
 		if (is.character(coords))
 			coords = match(coords, names(x))
@@ -406,9 +414,10 @@ print.sf = function(x, ..., n = getOption("sf_max_print", default = 10)) {
 		app = paste0(app, "\n", "Active geometry column: ", attr(x, "sf_column"))
 	print(st_geometry(x), n = 0, what = "Simple feature collection with", append = app)
 	if (n > 0) {
-		if (inherits(x, "tbl_df"))
+		if (inherits(x, "tbl_df")) {
+			st_geometry(x) = x[[attr(x, "sf_column")]][] # note the extra []: this reloads points
 			NextMethod()
-		else {
+		} else {
 			y <- x
 			if (nrow(y) > n) {
 				cat(paste("First", n, "features:\n"))

R/sfc.R

@@ -48,14 +48,24 @@ st_sfc = function(..., crs = NA_crs_, precision = 0.0, check_ring_dir = FALSE, d
 		lst = lst[[1]]
 	stopifnot(is.numeric(crs) || is.character(crs) || inherits(crs, "crs"))
 
+	points_in_attr <- !is.null(attr(lst, "points"))
+
 	# check for NULLs:
 	a = attributes(lst)
-	is_null = sfc_is_null(lst)
+	is_null = if (points_in_attr)
+			rep(FALSE, length(lst))
+		else
+			sfc_is_null(lst)
 	lst = unclass(lst)
-	lst = lst[! is_null]
+	if (!points_in_attr)
+		lst = lst[! is_null]
+
 	attributes(lst) = a
 
-	dims_and_types = sfc_unique_sfg_dims_and_types(lst)
+	dims_and_types = if (points_in_attr)
+			list(class_dim = attr(lst, "points_dim"), class_type = "POINT")
+		else
+			sfc_unique_sfg_dims_and_types(lst)
 
 	cls = if (length(lst) == 0) # empty set: no geometries read
 		c("sfc_GEOMETRY", "sfc")
@@ -140,8 +150,12 @@ st_sfc = function(..., crs = NA_crs_, precision = 0.0, check_ring_dir = FALSE, d
 "[.sfc" = function(x, i, j, ..., op = st_intersects) {
 	if (!missing(i) && (inherits(i, "sf") || inherits(i, "sfc") || inherits(i, "sfg")))
 		i = lengths(op(x, i, ...)) != 0
-	st_sfc(NextMethod(), crs = st_crs(x), precision = st_precision(x),
-		dim = if(length(x)) class(x[[1]])[1] else "XY")
+	if (inherits(x, "sfc_POINT") && !is.null(attr(x, "points")))
+		st_sfc(restore_points(x, i), crs = st_crs(x), precision = st_precision(x),
+			dim = if(length(x)) class(x[[1]])[1] else "XY")
+	else
+		st_sfc(NextMethod() , crs = st_crs(x), precision = st_precision(x),
+			dim = if(length(x)) class(x[[1]])[1] else "XY")
 }
 
 
@@ -487,16 +501,22 @@ st_coordinates.sfc = function(x, ...) {
 		return(matrix(nrow = 0, ncol = 2))
 
 	ret = switch(class(x)[1],
-		sfc_POINT = matrix(unlist(x, use.names = FALSE), nrow = length(x), byrow = TRUE,
-		     dimnames = NULL),
+		sfc_POINT = if (is.null(attr(x, "points"))) {
+					matrix(unlist(x, use.names = FALSE), nrow = length(x), byrow = TRUE, dimnames = NULL)
+				} else {
+					attr(x, "points")
+				},
 		sfc_MULTIPOINT = ,
 		sfc_LINESTRING = coord_2(x),
 		sfc_MULTILINESTRING = ,
 		sfc_POLYGON = coord_3(x),
 		sfc_MULTIPOLYGON = coord_4(x),
 		stop(paste("not implemented for objects of class", class(x)[1]))
 	)
-	Dims = class(x[[1]])[1]
+	Dims = if (!is.null(attr(x, "points_dim")))
+			attr(x, "points_dim")
+		else
+			class(x[[1]])[1]
 	ncd = nchar(Dims)
 	colnames(ret)[1:ncd] = vapply(seq_len(ncd), function(i) substr(Dims, i, i), "")
 	ret
@@ -579,3 +599,24 @@ st_as_sfc.bbox = function(x, ...) {
 	box = st_polygon(list(matrix(x[c(1, 2, 3, 2, 3, 4, 1, 4, 1, 2)], ncol = 2, byrow = TRUE)))
 	st_sfc(box, crs = st_crs(x))
 }
+
+#' @export
+`[[.sfc` = function(x, i, j, ..., exaxt = TRUE) {
+	if (inherits(x, "sfc_POINT") && !is.null(attr(x, "points")))
+		restore_point(x, i)
+	else
+		NextMethod()
+}
+
+restore_point = function(x, i = TRUE) {
+	restore_points(x, i)[[1]]
+}
+
+restore_points = function(x, i = TRUE) {
+	a = attributes(x)
+	points = a$points[i, , drop=FALSE]
+	structure(points_rcpp(points, a$points_dim),
+		n_empty = 0L, precision = a$precision, crs = a$crs,
+		bbox = bbox.pointmatrix(points), class = a$class, 
+		points = NULL, points_dim = NULL)
+}

R/sfg.R

@@ -258,6 +258,8 @@ c.sfg = function(..., recursive = FALSE, flatten = TRUE) {
 	Paste0 = function(lst) lapply(lst, unclass)
 	Paste1 = function(lst) do.call(c, lapply(lst, unclass))
 	lst = list(...)
+	if (length(lst) && is.null(lst[[1]]))
+		stop("to combine POINTs into MULTIPOINT, use st_combine(), or realize them first using x[]")
 	if (flatten) {
 		cls = vapply(lst, function(x) class(x)[2], "")
 		ucls = unique(cls)

R/tidyverse.R

@@ -497,6 +497,8 @@ nest.sf = function (.data, ...) {
 	if (!requireNamespace("tidyr", quietly = TRUE))
 		stop("tidyr required: install first?")
 
+	if (inherits(g <- st_geometry(.data), "sfc_POINT") && !is.null(attr(g, "points")))
+		st_geometry(.data) = g[] # realize
 	class(.data) <- setdiff(class(.data), "sf")
 	ret = tidyr::nest(.data, ...)
 	lst = which(sapply(ret, inherits, "list"))[1]
@@ -582,7 +584,7 @@ type_sum.sfc <- function(x, ...) {
 #' Summarize simple feature item for tibble
 #' @name tibble
 obj_sum.sfc <- function(x) {
-	vapply(x, function(sfg) format(sfg, width = 15L), "")
+	vapply(x[], function(sfg) format(sfg, width = 15L), "")
 }
 
 #' @name tibble

R/z_range.R

@@ -9,6 +9,10 @@ zb_wrap = function(zb) {
 	structure(zb, names = c("zmin", "zmax"), class = "z_range")
 }
 
+z_range.pointmatrix = function(obj, ...) {
+	zb_wrap(range(obj[,3]))
+}
+
 z_range.Set = function(obj, ...) {
 	sel = vapply(obj, function(x) { length(x) && !all(is.na(x)) }, TRUE)
 	if (! any(sel))
@@ -129,7 +133,9 @@ print.z_range = function(x, ...) {
 }
 
 compute_z_range = function(obj) {
-	switch(class(obj)[1],
+	if (!is.null(pts <- attr(obj, "points")))
+		z_range.pointmatrix(pts)
+	else switch(class(obj)[1],
 		   sfc_POINT = zb_wrap(z_range.Set(obj)),
 		   sfc_MULTIPOINT = zb_wrap(z_range.MtrxSet(obj)),
 		   sfc_LINESTRING = zb_wrap(z_range.MtrxSet(obj)),