Skip to content

Instantly share code, notes, and snippets.

Embed
What would you like to do?
library(sf)
library(cartography)
library(osrm)
# set osrm to my own server
options(osrm.server = "http://address.of.my.server/", osrm.profile = "driving")
# destination point (useR2019 conference in Toulouse)
dst <- data.frame(id="dst", x = 1.4344, y = 43.6113)
dst <- st_as_sf(dst, coords = c('x','y'), crs = 4326)
dst <- st_transform(dst, 2154)
# origin points (regular grid within a 100km radius circle around destination)
dst_buff <- st_buffer(x = dst, dist = 100000)
grid <- st_sf(geometry = st_make_grid(x = dst_buff, n = c(30,30),
what = "centers",
square = FALSE))
# selection of points within the 100km radius circle
grid <- grid[dst_buff,]
# Query the routes
l <- vector(mode = "list", length = nrow(grid))
for (i in seq_len(nrow(grid))){
l[[i]] <- osrmRoute(
src = grid[i,],
dst = dst,
returnclass = "sf",
overview = "full"
)
}
roads <- do.call(rbind,l)
# tranform routes to polygons ==> better display
roads <- roads[order(roads$duration, decreasing = T),]
roads <- st_buffer(roads, dist = 200)
# Display the routes
png("tolose0.png", width = 500, height = 595, res = 100)
par(mar = c(0,0,0,0), bg = "black")
choroLayer(x = roads, var = "duration", border = NA,
col = hcl.colors(n = 13, palette = "Oslo", alpha = 1),
breaks = c(8,seq(10,120,10),197), legend.pos = "bottomleft",
legend.frame = TRUE, legend.horiz = TRUE,
legend.title.txt = "Road distance by car, in minutes")
mtext(text = " No Time Toulouse", side = 3, line = -2, adj = 0.5, cex = 1.5, col = "white")
mtext(text = "Routing: OSRM - http://project-osrm.org\n© OpenStreetMap contributors. www.openstreetmap.org/copyright.", side = 4, line = -1,
col = "white", adj = 0, cex = .6)
dev.off()
# Just a nice picture
png("tolose2.png", width = 1600, height = 1200, res = 190)
par(mar = c(0,0,0,0), bg="black")
plot(roads["duration"], border = NA, bg="black",
pal = hcl.colors(n = 13, palette = "Oslo", alpha = .35),
breaks = c(8,seq(10,120,10),197),key.pos = NULL)
dev.off()
# Compute and display isochrones
iso <- osrmIsochrone(loc=c(x = 1.4344, y = 43.6113), returnclass = "sf",
breaks = seq(0,60,5), res = 100)
osm <- getTiles(iso, crop=TRUE, zoom = 9)
png("tolose3.png", width = 675, height = 500, res = 100)
par(mar = c(0.1,0,0.1,0.5))
tilesLayer(osm)
choroLayer(iso, var="center", border = NA,
col = hcl.colors(n = 13, palette = "Cividis", alpha = .9),
breaks = seq(0,60,5), add=T, legend.pos = "bottomleft",
legend.frame = TRUE, legend.horiz = TRUE,
legend.title.txt = "Road distance by car, in minutes")
mtext(text = "Routing: OSRM - http://project-osrm.org\n© OpenStreetMap contributors. Tiles style under CC BY-SA, www.openstreetmap.org/copyright. ",
side = 4, line = -.5,
col = "black", adj = 0, cex = .6)
mtext(text = "Time to Toulouse", side = 3, line = -2,
adj = 0.5, cex = 1.5, col = "Black")
dev.off()
sessionInfo()
# R version 3.6.0 (2019-04-26)
# Platform: x86_64-pc-linux-gnu (64-bit)
# Running under: Debian GNU/Linux 9 (stretch)
#
# Matrix products: default
# BLAS: /usr/lib/openblas-base/libblas.so.3
# LAPACK: /usr/lib/libopenblasp-r0.2.19.so
#
# locale:
# [1] LC_CTYPE=fr_FR.UTF-8 LC_NUMERIC=C LC_TIME=fr_FR.UTF-8
# [4] LC_COLLATE=fr_FR.UTF-8 LC_MONETARY=fr_FR.UTF-8 LC_MESSAGES=fr_FR.UTF-8
# [7] LC_PAPER=fr_FR.UTF-8 LC_NAME=C LC_ADDRESS=C
# [10] LC_TELEPHONE=C LC_MEASUREMENT=fr_FR.UTF-8 LC_IDENTIFICATION=C
#
# attached base packages:
# [1] stats graphics grDevices utils datasets methods base
#
# other attached packages:
# [1] osrm_3.3.0 cartography_2.2.0 sf_0.7-4
#
# loaded via a namespace (and not attached):
# [1] Rcpp_1.0.1 raster_2.8-19 magrittr_1.5 units_0.6-2
# [5] lattice_0.20-38 plyr_1.8.4 tools_3.6.0 rgdal_1.4-3
# [9] grid_3.6.0 packrat_0.5.0 png_0.1-7 KernSmooth_2.23-15
# [13] e1071_1.7-1 DBI_1.0.0 rgeos_0.4-3 class_7.3-15
# [17] abind_1.4-5 lwgeom_0.1-6 bitops_1.0-6 codetools_0.2-16
# [21] curl_3.3 isoband_0.2.0 RCurl_1.95-4.12 rosm_0.2.4
# [25] sp_1.3-1 compiler_3.6.0 classInt_0.3-3 jsonlite_1.6
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
You can’t perform that action at this time.