luukvdmeer / sfnetworks

Tidy Geospatial Networks in R
https://luukvdmeer.github.io/sfnetworks/
Other
347 stars 20 forks source link

Function that returns junction points #14

Closed agila5 closed 4 years ago

agila5 commented 4 years ago

Is your feature request related to a problem? Please describe. I think we should include a function that returns junction points of a street network. As far as I know, there is no clear definition of junction point but I think that a reasonable definition (for OSM data) is the following: a function point is a point that is repeated three or more times in the city network.

Describe the solution you'd like We should create a function that, given an object of class sfnetwork, returns an sf object with a POINT geometry with the coordinates of the junction points.

Describe alternatives you've considered This is a reprex of the problem with a possible solution:

# packages
library(sf)
#> Linking to GEOS 3.6.1, GDAL 2.2.3, PROJ 4.9.3
library(tmap)

# create a very simple network
network_column <- st_sfc(
  st_linestring(matrix(c(0, -2, 0, -1, 0, 0), ncol = 2, byrow = TRUE)), 
  st_linestring(matrix(c(0, 2, 0, 1, 0, 0), ncol = 2, byrow = TRUE)), 
  st_linestring(matrix(c(-2, 0, -1, 0, 0, 0), ncol = 2, byrow = TRUE)), 
  st_linestring(matrix(c(2, 0, 1, 0, 0, 0), ncol = 2, byrow = TRUE)), 
  st_linestring(matrix(c(2, -2, 2, -1, 2, 0), ncol = 2, byrow = TRUE)),
  st_linestring(matrix(c(2, 2, 2, 1, 2, 0), ncol = 2, byrow = TRUE)), 
  st_linestring(matrix(c(2, 0, 3, 0, 4, 0), ncol = 2, byrow = TRUE)), 
  crs = 3003
)

# find junction points
points_coordinates <- st_coordinates(network_column)[, c(1, 2)]
# points that occur two or more times
duplicated_points <- points_coordinates[duplicated(points_coordinates), ]
# points that occur three or more times
duplicated_duplicated_points <- duplicated_points[duplicated(duplicated_points), ]
# exclude further duplicates
duplicated_duplicated_points <- unique(duplicated_duplicated_points)

# transform back into an sf object
junction_points <- st_sfc(st_multipoint(duplicated_duplicated_points), crs = st_crs(network_column))

# plot
tm_shape(network_column) + 
  tm_lines() +
tm_shape(junction_points) + 
  tm_dots(col = "red", size = 0.5)

Created on 2020-02-17 by the reprex package (v0.3.0)

I've already discussed this problem with @Robinlovelace here: https://github.com/ropensci/stplanr/issues/357 (and maybe we can also adopt the same approach of that reprex to close that issue but I would just wait for sfnetwork) and I think this is the first step for #12.

luukvdmeer commented 4 years ago

Do you mean something like this?

library(sfnetworks)
library(sf)
library(tidygraph)

net = as_sfnetwork(roxel) %>%
  mutate(degree = centrality_degree()) %>%
  mutate(is_junction = degree > 2)
# A tbl_graph: 701 nodes and 851 edges
#
# A directed multigraph with 14 components
#
# Node Data: 701 x 3 (active)
             geometry degree is_junction
          <POINT [°]>  <dbl> <lgl>      
1 (7.533722 51.95556)      2 FALSE      
2 (7.533461 51.95576)      2 FALSE      
3 (7.532442 51.95422)      3 TRUE       
4  (7.53209 51.95328)      1 FALSE      
5 (7.532709 51.95209)      2 FALSE      
6 (7.532869 51.95257)      2 FALSE      
# … with 695 more rows
#
# Edge Data: 851 x 5
   from    to name                  type                                                         geometry
  <int> <int> <fct>                 <fct>                                                <LINESTRING [°]>
1     1     2 Havixbecker Strasse   residential                    (7.533722 51.95556, 7.533461 51.95576)
2     3     4 Pienersallee          secondary     (7.532442 51.95422, 7.53236 51.95377, 7.53209 51.95328)
3     5     6 Schulte-Bernd-Strasse residential (7.532709 51.95209, 7.532823 51.95239, 7.532869 51.95257)
# … with 848 more rows
net %>%
  filter(is_junction) %>%
  st_as_sf()
Simple feature collection with 8 features and 2 fields
geometry type:  POINT
dimension:      XY
bbox:           xmin: 7.528016 ymin: 51.94436 xmax: 7.53881 ymax: 51.95764
CRS:            EPSG:4326
# A tibble: 8 x 3
             geometry degree is_junction
          <POINT [°]>  <dbl> <lgl>      
1 (7.532442 51.95422)      3 TRUE       
2  (7.532275 51.9556)      3 TRUE       
3  (7.53881 51.95679)      3 TRUE       
4 (7.535912 51.95764)      3 TRUE       
5 (7.536379 51.95587)      3 TRUE       
6   (7.5282 51.94436)      3 TRUE       
7 (7.528016 51.95665)      3 TRUE       
8 (7.531559 51.94562)      3 TRUE  
agila5 commented 4 years ago

I think so, and it's simpler than I expected. Do you want to create an ad-hoc function that wraps centrality degree and filter for degree > 2? If not, then we can simply close the issue.

luukvdmeer commented 4 years ago

I think this would be a good example of a function that belongs in our extension to the 'Quering node types' set of functions that tidygraph has: https://tidygraph.data-imaginist.com/reference/node_types.html

That is, a function node_is_junction() returning TRUE or FALSE. The filtering is then one extra step that a user needs to do for him/herself.

However, I think we should only add it if there is some widely used definition in the literature of what a junction really is. I am not sure if there is!

luukvdmeer commented 4 years ago

Close for now, and may be reopened if we find a common definition of a junction. Agree @agila5 >