sln_add_node for multiple nodes

[richard-davey-arcadisgen]

Would it be possible to add another function or augment sln_add_node to add multiple nodes at one time? It can only add one node per function call right now.

I have a long list of nodes to add to a spatial lines network and cannot implement an lapply method due to the nature of sf objects.

[Robinlovelace]

Hi @richard-davey-arcadisgen, great to hear from you and good question. Can you provide a reproducible example to illustrate exactly what you're trying to achieve and what you have tried so far? You could also try using the sfnetworks package and could ask a question in the Discussion page there, although I suggest you need a reprex to get support here or on the https://github.com/luukvdmeer/sfnetworks repo.

[richard-davey-arcadisgen]

If you ignore the from and to arguments then I think this is a fair reprex (it was the one I had on hand to achieve something related).

library(stplanr)
from <- c(rnorm(1, -1.535181, 0.002), rnorm(1, 53.82534, 0.02))
to <- c(rnorm(1, -1.52446, 0.002), rnorm(1, 53.80949, 0.02))
from_sf <- sf::st_sfc(sf::st_point(from), crs = sf::st_crs(route_network_sf))
to_sf <- sf::st_sfc(sf::st_point(to), crs = sf::st_crs(route_network_sf))

## My workaround

#' Adds start and endpoints to a spatial lines network ready for route planning
#'
#' @param sln a spatial lines network object created by stplanr::SpatialLinesNetwork()
#' @param start a vector of sf geometries, each row representing the start of a line segment
#' @param end a vector of sf geometries, each row representing the end of a line segment
#' @export
add_endpoints <- function(sln, start, end){
  stopifnot(typeof(sln) == "S4")
  stopifnot(class(sln)[1] == "sfNetwork")
  pairs <- sf::st_union(start, end, by_feature = TRUE) %>%
    sf::st_cast("POINT")
  sln_out <- sln
  for(i in 1:length(pairs)){
    sln_out <- stplanr::sln_add_node(sln = sln_out, pairs[i])
  }
  sln_out
}

## The one-by-one method
sln <- SpatialLinesNetwork(route_network_sf)
sln_new <- sln_add_node(sln, from_sf) %>%
  sln_add_node(to_sf)

## The slow loop method
sln <- SpatialLinesNetwork(route_network_sf)
sln_new <- add_endpoints(sln, from_sf, to_sf)

[Robinlovelace]

Hi @richard-davey-arcadisgen many thanks, this is a great reprex, thanks for providing a function and a clear reproducible example. No solution for now, but here's the output of what your code on my computer (new discovery: #' text renders as text with the reprex package)!

library(stplanr)
from <- c(rnorm(1, -1.535181, 0.002), rnorm(1, 53.82534, 0.02))
to <- c(rnorm(1, -1.52446, 0.002), rnorm(1, 53.80949, 0.02))
from_sf <- sf::st_sfc(sf::st_point(from), crs = sf::st_crs(route_network_sf))
to_sf <- sf::st_sfc(sf::st_point(to), crs = sf::st_crs(route_network_sf))

## My workaround

Adds start and endpoints to a spatial lines network ready for route planning

@param sln a spatial lines network object created by stplanr::SpatialLinesNetwork()
@param start a vector of sf geometries, each row representing the start of a line segment
@param end a vector of sf geometries, each row representing the end of a line segment
@export

add_endpoints <- function(sln, start, end){
  stopifnot(typeof(sln) == "S4")
  stopifnot(class(sln)[1] == "sfNetwork")
  pairs <- sf::st_union(start, end, by_feature = TRUE) %>%
    sf::st_cast("POINT")
  sln_out <- sln
  for(i in 1:length(pairs)){
    sln_out <- stplanr::sln_add_node(sln = sln_out, pairs[i])
  }
  sln_out
}

## The one-by-one method
sln <- SpatialLinesNetwork(route_network_sf)
sln_new <- sln_add_node(sln, from_sf) %>%
  sln_add_node(to_sf)
#> Linking to GEOS 3.9.0, GDAL 3.2.1, PROJ 7.2.1
## The slow loop method
sln <- SpatialLinesNetwork(route_network_sf)
sln_new <- add_endpoints(sln, from_sf, to_sf)


plot(sln_new)

waldo::compare(sln, sln_new)
#> `names(old@sl)`: "All" "geometry" "length"  
#> `names(new@sl)`: "All" "length"   "geometry"
#> 
#> `attr(old@sl, 'row.names')` is an integer vector (1, 2, 3, 4, 5, ...)
#> `attr(new@sl, 'row.names')` is a character vector ('1', '2', '3', '4', '5', ...)
#> 
#> `old@sl$All[78:80]`: 0 0 0    
#> `new@sl$All[78:82]`: 0 0 0 0 0
#> 
#> `old@sl$geometry` is length 80
#> `new@sl$geometry` is length 82
#> 
#> `dim(old@sl$geometry[[18]])`: 24 2
#> `dim(new@sl$geometry[[18]])`: 22 2
#> 
#>      old@sl$geometry[[18]] | new@sl$geometry[[18]]                
#>  [1] -1.535906             | -1.535906             [1]            
#>  [2] -1.535878             - -1.535135             [2]            
#>  [3] -1.535838             - -1.534722             [3]            
#>  [4] -1.535788             - -1.534319             [4]            
#>  [5] -1.535733             - -1.53394              [5]            
#>  [6] -1.535443             - -1.533227             [6]            
#>  [7] -1.534933             - -1.532639             [7]            
#>  [8] -1.533518             - -1.532338             [8]            
#>  [9] -1.532966             - -1.531974             [9]            
#> [10] -1.532606             - -1.531982             [10]           
#>  ... ...                     ...                   and 38 more ...
#> 
#> `dim(old@sl$geometry[[19]])`: 22 2
#> `dim(new@sl$geometry[[19]])`:  3 2
#> 
#>      old@sl$geometry[[19]] | new@sl$geometry[[19]]                
#>  [1] -1.535906             - -1.535397             [1]            
#>  [2] -1.535135             - -1.535563             [2]            
#>  [3] -1.534722             - -1.535722             [3]            
#>  [4] -1.534319             - 53.826564             [4]            
#>  [5] -1.53394              - 53.826534             [5]            
#>  [6] -1.533227             - 53.826505             [6]            
#>  [7] -1.532639             -                                      
#>  [8] -1.532338             -                                      
#>  [9] -1.531974             -                                      
#> [10] -1.531982             -                                      
#>  ... ...                     ...                   and 34 more ...
#> 
#> `dim(old@sl$geometry[[20]])`: 3 2
#> `dim(new@sl$geometry[[20]])`: 6 2
#> 
#>     old@sl$geometry[[20]] | new@sl$geometry[[20]]               
#> [1] -1.535397             | -1.535397             [1]           
#> [2] -1.535563             - -1.534878             [2]           
#> [3] -1.535722             - -1.534453             [3]           
#> [4] 53.826564             - -1.533934             [4]           
#> [5] 53.826534             - -1.533407             [5]           
#> [6] 53.826505             - -1.532425             [6]           
#>                           - 53.826564             [7]           
#>                           - 53.82664              [8]           
#>                           - 53.826699             [9]           
#>                           - 53.826786             [10]          
#> ... ...                     ...                   and 2 more ...
#> 
#> And 189 more differences ...

^{Created on 2021-09-16 by the reprex package (v2.0.1)}

[Robinlovelace]

Just to clarify, are you proposing to add add_endpoints() as a function in stplanr?

[agila5]

Hi @richard-davey-arcadisgen! I'm not sure how to fix the problem within stplanr, but I think that, at the moment, you can use the following workaround:

Load packages

library(sf)
#> Linking to GEOS 3.9.1, GDAL 3.2.1, PROJ 7.2.1
library(stplanr)
library(sfnetworks)

Define toy data (just 3 segments)

segments <- st_as_sf(st_sfc(
  st_linestring(rbind(c(0, 0), c(1, 0))), 
  st_linestring(rbind(c(1, 0), c(1, 1))), 
  st_linestring(rbind(c(1, 1), c(0, 0)))
))
sln <- SpatialLinesNetwork(segments)
sfn <- as_sfnetwork(sln)

Sample 5 points in the bbox of the segments

set.seed(1)
pts <- st_sample(st_as_sfc(st_bbox(sln@sl)), 5L)

Split (or blend) the network at those points

sfn_split <- st_network_blend(sfn, pts)
#> Warning: st_network_blend assumes attributes are constant over geometries

Recreate the sln object

sln_split <- SpatialLinesNetwork(st_as_sf(sfn_split, "edges"))

Plot. The big grey segments are the original network, while the rainbow
segments represent the splitted network

plot(sln@sl, col = "grey", lwd = 15, main = "", reset = FALSE)
plot(pts, add = TRUE, pch = 16)
plot(st_geometry(sln_split@sl), col = rainbow(8), add = TRUE, lwd = 3)

^{Created on 2021-09-27 by the reprex package (v2.0.1)}

[richard-davey-arcadisgen]

Just to clarify, are you proposing to add add_endpoints() as a function in stplanr?

Sorry for the delayed response. No, I'm not proposing adding my example function above to stplanr, simply because it would be very slow for many nodes and I'm sure there would be a better way to do it. I would welcome a vectorised or more efficient version in stplanr though.

[richard-davey-arcadisgen]

@agila5, +1 for rainbows! I was unaware of the functionality of sfnetworks and I think that st_network_blend is exactly what I need. Thanks very much.