Retain sf point columns after st_network_blend()

252 Views Asked by At

I have an sfnetwork (net) and pts that I am combining using st_network_blend(). After calculating some new attributes for these points I would like to extract them and regain the columns of the original points. A spatial join using st_join won't work because the points are snapped during st_network_blend() (therefore no longer overlap) and cbind() won't work because not all points are within the snapping tolerance.

Is there a way of preserving an ID from pts through the blending process so it can be used for joining afterwards to out_pts? Or will left_join() give a correct answer when joining by rowid_to_column() applied to pts and afterwards on out_pts?

library(sfnetworks)
library(sf)
library(dplyr)

# Create a network and a set of points to blend.
n11 = st_point(c(0,0))
n12 = st_point(c(1,1))
e1 = st_sfc(st_linestring(c(n11, n12)), crs = 3857)

n21 = n12
n22 = st_point(c(0,2))
e2 = st_sfc(st_linestring(c(n21, n22)), crs = 3857)

n31 = n22
n32 = st_point(c(-1,1))
e3 = st_sfc(st_linestring(c(n31, n32)), crs = 3857)

net = as_sfnetwork(c(e1,e2,e3))

pts = net %>%
  st_bbox() %>%
  st_as_sfc() %>%
  st_sample(10, type = "random") %>%
  st_set_crs(3857) %>%
  st_cast('POINT') %>%
  st_as_sf() %>%
  mutate(is_point = TRUE,
         dummy_var = paste0("A", seq_along(1:nrow(.)))) %>%
  tibble::rowid_to_column("node_id")

plot(net)
plot(st_geometry(pts), col = "red", pch = 20, add = TRUE)

# Blend points into the network.
# --> Not all points get blended 
tol = units::set_units(0.2, "m")
net_blend = st_network_blend(net, pts["is_point"], tolerance = tol)

# now only interested in points that fell within tol
# but dummy_var column was not preserved from pts
# added node_id column, will this differ from pts??
out_pts <- st_as_sf(net_blend, "nodes") %>%
  tibble::rowid_to_column("node_id") %>%
  filter(is_point) %>%
  select(-is_point)

# there's now a mismatch because pts and out_pts are different lengths
# as not all pts made it through st_network_blend()
nrow(pts)
nrow(out_pts)

# try st_join anyway
# dummy_var should != NA
st_join(out_pts, pts) # nope, doesn't work    

# pts without geometry
pts_cols <-  pts %>%
  as_data_frame()%>%
  select(-is_point)
pts_cols$x <- NULL

# what about cbind()
cbind(out_pts, pts_cols) # nope

# how do I tell if dummy_var has been joined to the correct point?
left_join(out_pts, pts_cols, by = "node_id")

I've also tried

tol = units::set_units(0.2, "m")
net_blend = st_network_blend(net, pts, tolerance = tol)

sf_attr(net_blend, "agr", active = "nodes")

But this gives

  node_id  is_point dummy_var 
     <NA>      <NA>      <NA> 
Levels: constant aggregate identity
1

There are 1 best solutions below

0
On

Attributes are preserved...

I changed

net_blend = st_network_blend(net, pts[is_point], tolerance = tol)

out_pts <- st_as_sf(net_blend, "nodes") %>%
  tibble::rowid_to_column("node_id") %>%
  filter(is_point) %>%
  select(-is_point)

To

net_blend = st_network_blend(net, pts, tolerance = tol)#

out_pts <- st_as_sf(net_blend, "nodes") %>%
  filter(is_point) %>%
  select(-is_point)

To give

> out_pts
Simple feature collection with 2 features and 2 fields
Geometry type: POINT
Dimension:     XY
Bounding box:  xmin: 0.6040494 ymin: 1.27107 xmax: 0.7289301 ymax: 1.395951
Projected CRS: WGS 84 / Pseudo-Mercator
# A tibble: 2 x 3
  node_id dummy_var                    x
    <int> <chr>              <POINT [m]>
1       6 A6         (0.7289301 1.27107)
2       3 A3        (0.6040494 1.395951)