st_network_join join two networks by their nearest nodes within tolerance

[whipson]

When joining two networks, st_network_join will join only identical nodes. I'm wondering if it's possible to extend it to join by nearest node within some tolerance. The tolerance would be the maximum distance between the nodes for them to join.

Perhaps st_network_join could take a predicate function like sf::st_nearest_point?

Maybe this is already possible but there isn't much documentation on st_network_join.

Thanks

[agila5]

Hi @whipson!

Did you already check the rounding coordinates Section in the "Pre processing tasks" vignette? I'm not 100% sure, but I think that including a tolerance argument in st_network_join is not a trivial task (unless we implicitly run the same pre-processing task as detailed in the vignette). Anyway I think that, for the moment, a manual rounding to the desired precision might be the easiest and quickiest solution. A simple reprex:

Load packages

library(sf)
#> Linking to GEOS 3.9.1, GDAL 3.2.1, PROJ 7.2.1
library(sfnetworks)
Sys.setenv(`_R_S3_METHOD_REGISTRATION_NOTE_OVERWRITES_` = "false")

Simulate data

edge1 <- st_sfc(
  st_linestring(rbind(c(0, 0), c(2, 2)))
)
(net1 <- as_sfnetwork(edge1))  
#> # A sfnetwork with 2 nodes and 1 edges
#> #
#> # CRS:  NA 
#> #
#> # A rooted tree with spatially explicit edges
#> #
#> # Node Data:     2 x 1 (active)
#> # Geometry type: POINT
#> # Dimension:     XY
#> # Bounding box:  xmin: 0 ymin: 0 xmax: 2 ymax: 2
#>         x
#>   <POINT>
#> 1   (0 0)
#> 2   (2 2)
#> #
#> # Edge Data:     1 x 3
#> # Geometry type: LINESTRING
#> # Dimension:     XY
#> # Bounding box:  xmin: 0 ymin: 0 xmax: 2 ymax: 2
#>    from    to            x
#>   <int> <int> <LINESTRING>
#> 1     1     2   (0 0, 2 2)

edge2 <- st_sfc(
  st_linestring(rbind(c(sqrt(2)^2, sqrt(2)^2), c(3, 3)))
)
(net2 <- as_sfnetwork(edge2))
#> # A sfnetwork with 2 nodes and 1 edges
#> #
#> # CRS:  NA 
#> #
#> # A rooted tree with spatially explicit edges
#> #
#> # Node Data:     2 x 1 (active)
#> # Geometry type: POINT
#> # Dimension:     XY
#> # Bounding box:  xmin: 2 ymin: 2 xmax: 3 ymax: 3
#>         x
#>   <POINT>
#> 1   (2 2)
#> 2   (3 3)
#> #
#> # Edge Data:     1 x 3
#> # Geometry type: LINESTRING
#> # Dimension:     XY
#> # Bounding box:  xmin: 2 ymin: 2 xmax: 3 ymax: 3
#>    from    to            x
#>   <int> <int> <LINESTRING>
#> 1     1     2   (2 2, 3 3)

No spatial join since we notice identical and duplicated nodes:

st_network_join(net1, net2)
#> # A sfnetwork with 4 nodes and 2 edges
#> #
#> # CRS:  NA 
#> #
#> # A rooted forest with 2 trees with spatially explicit edges
#> #
#> # Node Data:     4 x 1 (active)
#> # Geometry type: POINT
#> # Dimension:     XY
#> # Bounding box:  xmin: 0 ymin: 0 xmax: 2 ymax: 2
#>         x
#>   <POINT>
#> 1   (0 0)
#> 2   (2 2)
#> 3   (2 2)
#> 4   (3 3)
#> #
#> # Edge Data:     2 x 3
#> # Geometry type: LINESTRING
#> # Dimension:     XY
#> # Bounding box:  xmin: 0 ymin: 0 xmax: 3 ymax: 3
#>    from    to            x
#>   <int> <int> <LINESTRING>
#> 1     1     2   (0 0, 2 2)
#> 2     3     4   (2 2, 3 3)

Round coords for edge2 with 7 decimal places

edge2 = st_geometry(edge2) %>%
  lapply(function(x) round(x, 7)) %>%
  st_sfc(crs = st_crs(edge2))
(net2 <- as_sfnetwork(edge2))
#> # A sfnetwork with 2 nodes and 1 edges
#> #
#> # CRS:  NA 
#> #
#> # A rooted tree with spatially explicit edges
#> #
#> # Node Data:     2 x 1 (active)
#> # Geometry type: POINT
#> # Dimension:     XY
#> # Bounding box:  xmin: 2 ymin: 2 xmax: 3 ymax: 3
#>         x
#>   <POINT>
#> 1   (2 2)
#> 2   (3 3)
#> #
#> # Edge Data:     1 x 3
#> # Geometry type: LINESTRING
#> # Dimension:     XY
#> # Bounding box:  xmin: 2 ymin: 2 xmax: 3 ymax: 3
#>    from    to            x
#>   <int> <int> <LINESTRING>
#> 1     1     2   (2 2, 3 3)

and now

st_network_join(net1, net2)
#> # A sfnetwork with 3 nodes and 2 edges
#> #
#> # CRS:  NA 
#> #
#> # A rooted tree with spatially explicit edges
#> #
#> # Node Data:     3 x 1 (active)
#> # Geometry type: POINT
#> # Dimension:     XY
#> # Bounding box:  xmin: 0 ymin: 0 xmax: 2 ymax: 2
#>         x
#>   <POINT>
#> 1   (0 0)
#> 2   (2 2)
#> 3   (3 3)
#> #
#> # Edge Data:     2 x 3
#> # Geometry type: LINESTRING
#> # Dimension:     XY
#> # Bounding box:  xmin: 0 ymin: 0 xmax: 3 ymax: 3
#>    from    to            x
#>   <int> <int> <LINESTRING>
#> 1     1     2   (0 0, 2 2)
#> 2     2     3   (2 2, 3 3)

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

[whipson]

Thanks for the suggestion and reprex @agila5. I see where you're coming from on not wanting to extend st_network_join into overly complex territory.

I did follow section 2 on network preprocessing and attempting some rounding to address my problem. The issue I ran into was that I only want to join the networks at a thin border between the two (imagine the networks sitting side by side in space and the goal is to fuse them together at a border). I tried to identify those nodes at the border and round only those nodes, but I couldn't get it working for me. Maybe the problem just requires some more thought on my end.

[luukvdmeer]

Thanks for the question @whipson and for the reply @agila5. Indeed this is currently not possible. I don't want to make the joining function too complicated with too many arguments, but am thinking of a separate function (morpher?) to update node coordinates of a network. Such a function would accept the network and an sfc object of the same length as the current node table, containing the updated node coordinates. The function would then automatically also update the affected edges to preserve the valid spatial network structure. This would allow a workflow in which you first update the coordinates of the nodes in network 2 to be equal to their nearest node in network 1 (optionally with a tolerance), and only then call st_network_join().