Converting sf_linestrings to sfnetwork causes some intersections to shift

[dpriss]

Hello again, First of all, I'm sorry for not using the reprex package, but I couldn't install it for some reason. I have a data set of spatial lines of which I want to create an sfnetwork. First I thought that everything works well until I discovered that some nodes were isolates although they shouldn't be. After some digging, I found that some of the intersections in the sf_linestring object, i.e. points that are shared by two or more lines, were slightly of after using as_sfnetwork(). Hence, digging even deeper, I extracted two such lines that intersect and tried to reproduce the erronous creation of nodes. However, just with the two lines, it workes perfectly again. I don't know if it is an error on my side or a bug but I would really appreciate any ideas how to solve this.

This is the code I'm using to produce the network using the whole data set hw_cast (attached as .zip):

remotes::install_github("luukvdmeer/sfnetworks")
library(sf) 
library(dplyr) 
library(sfheaders)
library(plotly)
library(sfnetworks)

## import data 
hw_cast <- st_read("hw_cast.shp")

##round geometries to avoid having too many nodes/edges because they don't coincide exactly
st_geometry(hw_cast) = st_geometry(hw_cast) %>%
  lapply(function(x) round(x, 0)) %>%
  st_sfc(crs = st_crs(hw_cast))

##create sfnetwork
net <- as_sfnetwork(hw_cast, directed = FALSE)
autoplot(net)

The result is shown in the first image with the offset intersection points in the red circle. Instead of one point at the intersection, there are now two points slightly off.

However, reproducing the same thing with just those two lines like this:

remotes::install_github("luukvdmeer/sfnetworks")
library(sf) 
library(dplyr) 
library(sfheaders)
library(plotly)
library(sfnetworks)

hw_cast <- st_read("hw_cast.shp")

## create line 1
ls1 <- data.frame(st_coordinates(hw_cast[5,]))
ls1.x <- ls1[,1]
ls1.y <- ls1[,2]

df1 <- data.frame(x = ls1.x, y = ls1.y) 
l1 <- sfheaders::sf_linestring(df1, x = "x", y = "y")

## create line 2
ls2 <- data.frame(st_coordinates(hw_cast[88,]))
ls2.x <- ls2[,1]
ls2.y <- ls2[,2]

df2 <- data.frame(x = ls2.x, y = ls2.y) 
l2 <- sfheaders::sf_linestring(df2, x = "x", y = "y")

## combine lines into one data set
lines <- st_sf(rbind(l1, l2))

## round geometries to avoid having too many nodes/edges because they don't coincide exactly
st_geometry(lines) = st_geometry(lines) %>%
  lapply(function(x) round(x, 0)) %>%
  st_sfc(crs = st_crs(lines))

## create sfnetwork
net <- as_sfnetwork(lines, directed = FALSE)
autoplot(net)

## add subdivisions
subdivision = convert(net, to_spatial_subdivision)
autoplot(subdivision)

gives the expected result with no node at the intersection but the possibility to add it using to_spatial_subdivision (second and third image). Here, the node is at the intersection as it should be.

Any ideas about why this happens and how to fix it would be highly appreciated! Thank you!

R Session Info R version 4.1.1 (2021-08-10) Platform: x86_64-w64-mingw32/x64 (64-bit) Running under: Windows 10 x64 (build 19044)

Matrix products: default

locale: [1] LC_COLLATE=English_United Kingdom.1252 LC_CTYPE=English_United Kingdom.1252 LC_MONETARY=English_United Kingdom.1252 [4] LC_NUMERIC=C LC_TIME=English_United Kingdom.1252

attached base packages: [1] stats graphics grDevices utils datasets methods base

other attached packages: [1] sfnetworks_0.6.2 plotly_4.10.1 ggplot2_3.4.0 sfheaders_0.4.0 dplyr_1.0.10 sf_1.0-6

loaded via a namespace (and not attached): [1] Rcpp_1.0.7 pillar_1.8.1 compiler_4.1.1 remotes_2.4.2 class_7.3-19 tools_4.1.1 digest_0.6.28
[8] evaluate_0.14 viridisLite_0.4.0 jsonlite_1.7.2 lifecycle_1.0.3 tibble_3.1.8 gtable_0.3.0 pkgconfig_2.0.3
[15] rlang_1.0.6 tidygraph_1.2.3 igraph_1.4.0 DBI_1.1.3 cli_3.6.0 rstudioapi_0.13 yaml_2.2.1
[22] xfun_0.26 fastmap_1.1.0 e1071_1.7-12 knitr_1.34 httr_1.4.2 withr_2.5.0 htmlwidgets_1.5.4 [29] generics_0.1.0 vctrs_0.5.1 classInt_0.4-8 grid_4.1.1 tidyselect_1.1.1 glue_1.4.2 data.table_1.14.4 [36] R6_2.5.1 fansi_0.5.0 rmarkdown_2.11 farver_2.1.0 tidyr_1.2.1 purrr_0.3.4 magrittr_2.0.1
[43] htmltools_0.5.2 scales_1.2.1 units_0.8-0 assertthat_0.2.1 colorspace_2.0-2 utf8_1.2.2 KernSmooth_2.23-20 [50] proxy_0.4-27 lazyeval_0.2.2 munsell_0.5.0 lwgeom_0.2-11 crayon_1.5.2

hw_cast.zip

[dpriss]

After some digging, it turned out that there were overlapping lines which I thought would have been removed by st_line_merge. However, they weren't, so after using hw_cast$covered = hw_cast %>% st_covered_by(.) %>% lengths > 1, I got finally rid of them and now all looks good.