reprex_for_som.md

# packages
library(sf)
#> Linking to GEOS 3.6.1, GDAL 2.2.3, PROJ 4.9.3
library(osmdata)
#> Data (c) OpenStreetMap contributors, ODbL 1.0. https://www.openstreetmap.org/copyright
library(tmap)
#> Warning: replacing previous import 'sf::st_make_valid' by
#> 'lwgeom::st_make_valid' when loading 'tmap'
tmap_mode("view")
#> tmap mode set to interactive viewing

# Download the polygonal bbox for the london city
london_poly <- getbb("City of London", format_out = "sf_polygon")

# The first entry correspond to the city of london
london_poly <- london_poly[1, ]

# Check the result
tm_shape(st_boundary(london_poly)) + 
  tm_lines(lwd = 2)

# Now we can extract the bbox
london_bbox <- st_bbox(london_poly)

# and expand it a little bit to include also the highways that are slightly
# outside of the bbox
london_bbox[c(1, 3)] <- london_bbox[c(1, 3)] + c(-0.05, 0.05) * diff(london_bbox[c(1, 3)])
london_bbox[c(2, 4)] <- london_bbox[c(2, 4)] + c(-0.05, 0.05) * diff(london_bbox[c(2, 4)])

# And download the OSM data for the highways
london_highway <- opq(london_bbox) %>% 
  add_osm_feature(key = "highway") %>% 
  osmdata_sf() %>% 
  osm_poly2line()

# This is the result
london_highway
#> Object of class 'osmdata' with:
#>                  $bbox : 51.506047015,-0.11587585,51.524143885,-0.07067135
#>         $overpass_call : The call submitted to the overpass API
#>                  $meta : metadata including timestamp and version numbers
#>            $osm_points : 'sf' Simple Features Collection with 21558 points
#>             $osm_lines : 'sf' Simple Features Collection with 7343 linestrings
#>          $osm_polygons : 'sf' Simple Features Collection with 319 polygons
#>        $osm_multilines : NULL
#>     $osm_multipolygons : 'sf' Simple Features Collection with 45 multipolygons

# Extract only the osm_lines object
london_higways_lines <- london_highway$osm_lines

# Plot the result
tm_shape(london_higways_lines) + 
  tm_lines() + 
tm_shape(st_boundary(london_poly)) + 
  tm_lines(lwd = 2)

# Filter only the highways strictly the polygon
london_higways_lines <- london_higways_lines[london_poly, op = st_covered_by]
#> although coordinates are longitude/latitude, st_covered_by assumes that they are planar

# repeat the plotting component
tm_shape(london_higways_lines) + 
  tm_lines() + 
  tm_shape(st_boundary(london_poly)) + 
  tm_lines(lwd = 2)

^{Created on 2020-03-14 by the reprex package (v0.3.0)}