Added find_snap() R wrapper function around Java's findSnapPoints()

Added `use_elevation` parameter to `setup_r5()`, defaulting to FALSE. Also added a disclaimer message that support for elevation in r5r is experimental.

r-package/NAMESPACE

@@ -3,6 +3,7 @@
 export(accessibility)
 export(detailed_itineraries)
 export(download_r5)
+export(find_snap)
 export(isochrones)
 export(setup_r5)
 export(stop_r5)

r-package/R/accessibility.R

@@ -137,13 +137,14 @@
 accessibility <- function(r5r_core,
                           origins,
                           destinations,
+                          opportunities_colname = "opportunities",
                           mode = "WALK",
                           mode_egress = "WALK",
                           departure_datetime = Sys.time(),
                           time_window = 1L,
                           percentiles = 50L,
-                          cutoffs = 30L,
                           decay_function = "step",
+                          cutoffs = 30L,
                           decay_value = 1.0,
                           max_walk_dist = Inf,
                           max_trip_duration = 120L,

r-package/R/elevation_utils.R

@@ -24,14 +24,20 @@ tobler_hiking <- function(slope) {
 #'              the terrain's slopes
 #'
 #' @param r5r_core a rJava object to connect with R5 routing engine
-#' @param raster_file string. Path to a raster file containing the study area's
-#'                    topography.
+#' @param raster_file string. Path to raster files containing the study area's
+#'                    topography. If a list is provided, all the rasters are
+#'                    automatically merged.
 #'
 #' @return No return value, called for side effects.
 #' @family elevation support functions
-apply_elevation <- function(r5r_core, raster_file) {
-  # load raster file containing elevation data
-  dem <- raster::raster(raster_file)
+apply_elevation <- function(r5r_core, raster_files) {
+  # load raster files containing elevation data
+  if (length(raster_files) == 1) {
+    dem <- raster::raster(raster_files[1])
+  } else {
+    dem_files <- lapply(raster_files, raster::raster)
+    dem <- do.call(raster::merge, dem_files)
+  }
 
   # extract street edges from r5r_core
   edges <- r5r_core$getEdges()

r-package/R/find_snap.R

@@ -0,0 +1,56 @@
+#' Find snapped locations of input points on street network
+#'
+#' @param r5r_core a rJava object to connect with R5 routing engine
+#' @param points a spatial sf POINT object, or a data.frame
+#'               containing the columns 'id', 'lon', 'lat'
+#' @param mode string. Defaults to "WALK", also allows "BICYCLE", and "CAR".
+#' @param n_threads numeric. The number of threads to use in parallel computing.
+#'                  Defaults to use all available threads (Inf).
+#'
+#' @return A data.table with the original points as well as their respective
+#'         snapped coordinates on the street network.
+#' @export
+#'
+find_snap <- function(r5r_core,
+                      points,
+                      mode = "WALK",
+                      n_threads = Inf) {
+
+  # set data.table options --------------------------------------------------
+
+  old_options <- options()
+  old_dt_threads <- data.table::getDTthreads()
+
+  on.exit({
+    options(old_options)
+    data.table::setDTthreads(old_dt_threads)
+  })
+
+  options(datatable.optimize = Inf)
+
+  # check inputs ------------------------------------------------------------
+
+  # r5r_core
+  checkmate::assert_class(r5r_core, "jobjRef")
+
+  # modes
+  if (!(mode %in% c('WALK','BICYCLE','CAR'))) {
+    stop(paste0(mode, " is not a valid 'mode'.\nPlease use one of the following: WALK, BICYCLE, CAR"))
+  }
+
+  # origins and destinations
+  points  <- assert_points_input(points, "points")
+
+  # set number of threads to be used by r5 and data.table
+  set_n_threads(r5r_core, n_threads)
+
+  # snap points to street network
+  snap_df <- r5r_core$findSnapPoints(points$id, points$lat, points$lon, mode)
+  snap_df <- jdx::convertToR(snap_df)
+  setDT(snap_df)
+
+  snap_df[found == FALSE, `:=`(snap_lat = NA, snap_lon = NA, distance = NA)]
+
+  return(snap_df)
+}
+

r-package/R/setup_r5.R

@@ -15,6 +15,9 @@
 #'                FALSE to show only eventual ERROR and WARNING messages.
 #' @param temp_dir logical, whether the R5 Jar file should be saved in temporary
 #'                 directory. Defaults to FALSE
+#' @param use_elevation boolean. If TRUE, load any tif files containing elevation
+#'                      found in the data_path folder and calculate impedances for
+#'                      walking and cycling based on street slopes.
 #'
 #' @return An rJava object to connect with R5 routing engine
 #' @family setup
@@ -33,7 +36,8 @@
 setup_r5 <- function(data_path,
                      version = "6.2.0",
                      verbose = TRUE,
-                     temp_dir = FALSE) {
+                     temp_dir = FALSE,
+                     use_elevation = FALSE) {
 
   # check inputs ------------------------------------------------------------
   checkmate::assert_logical(verbose)
@@ -118,14 +122,20 @@ setup_r5 <- function(data_path,
 
   }
 
-  # check for any elevation files in data_path (*.tif)
-  tif_files <- list.files(path = data_path, pattern = "*.tif$", full.names = TRUE)
+  # elevation
+  if (use_elevation) {
+    # check for any elevation files in data_path (*.tif)
+    tif_files <- list.files(path = data_path, pattern = "*.tif$", full.names = TRUE)
 
-  # if there are any .tif files in the data_path folder, apply elevetion to street network
-  if (length(tif_files) > 0) {
-    apply_elevation(r5r_core, tif_files[1])
+    # if there are any .tif files in the data_path folder, apply elevation to street network
+    if (length(tif_files) > 0) {
+      message(sprintf("%s TIF files found in data path. Loading elevation into street edges.\n", length(tif_files)),
+              "DISCLAIMER: this is an r5r specific feature, and it will be deprecated once native support for elevation data is added to R5.")
+      apply_elevation(r5r_core, tif_files)
+    }
   }
 
+  # finish R5's setup by pre-calculating distances between transit stops and street network
   r5r_core$buildDistanceTables()
 
   return(r5r_core)

r-package/tests_marcus/test_access.R

@@ -29,7 +29,8 @@ decays %>%
   ggplot(aes(x=secs/60, y=decay, color=fun)) +
   geom_point() +
   geom_vline(xintercept = 30) +
-  facet_wrap(~fun)
+  facet_wrap(~fun) +
+  theme(legend.position = "none")
 
 
 
@@ -216,3 +217,20 @@ ggplot(df_access) +
        subtitle = "considering different decay functions",
        caption = "cutoff = 30min, exponential decay = 0.001, linear width = 10, logistic st.dev = 20")
 
+
+######### ELEVATION
+tobler_hiking <- function(slope) {
+  C <- 1.19403
+
+  tobler_factor <- C * exp(-3.5 * abs(slope+0.05))
+
+  return(1 / tobler_factor)
+}
+
+t_factor = tobler_hiking(9L:-13L)
+plot(t_factor)
+
+slopes <- as.double(9.0:-13.0)
+alts <- as.double(rep(0.0, 23))
+b_factor <- r5r_core$bikeSpeedCoefficientOTP(slopes, alts)
+plot(b_factor)

r-package/tests_marcus/test_snap.R

@@ -1,5 +1,6 @@
 # library(tidyverse)
 options(java.parameters = '-Xmx10G')
+devtools::load_all(".")
 library(r5r)
 library(sf)
 library(tidyverse)
@@ -10,23 +11,23 @@ data_path <- system.file("extdata/poa", package = "r5r")
 r5r_core <- setup_r5(data_path, verbose = FALSE)
 
 # get regular grid at resolution 8
-grid_df <- r5r_core$getGrid(11L)
+grid_df <- r5r_core$getGrid(8L)
 grid_df <- jdx::convertToR(grid_df)
 grid_df$geometry <- st_as_sfc(grid_df$geometry)
 grid_df <- st_as_sf(grid_df, crs = 4326)
 
-grid_df %>% mutate(id = as.integer(id)) %>% mapview::mapview(zcol = "id")
-  mapview::mapview(snap_df %>% filter(point_id=="4"), xcol="lon", ycol="lat") +
+grid_df %>% mutate(id = as.integer(id)) %>% mapview::mapview(alpha.regions = 0)
+  mapview::mapview(snap_df, xcol="lon", ycol="lat")
   mapview::mapview(snap_df %>% filter(point_id=="4"), xcol="snap_lon", ycol="snap_lat")
 
 
 # snap grid to street network
-snap_df <- r5r_core$findSnapPoints(grid_df$id, grid_df$lat, grid_df$lon, "WALK")
+snap_df <- r5r_core$findSnapPoints(grid_df$id, grid_df$lat, grid_df$lon, "CAR")
 snap_df <- jdx::convertToR(snap_df)
 snap_df <- snap_df %>%
   filter(found == TRUE)
 
-mapview::mapview(snap_df, xcol="lon", ycol="lat", crs=4326) +
+mapview::mapview(snap_df, xcol="lon", ycol="lat", crs=4326)
 mapview::mapview(snap_df, xcol="snap_lon", ycol="snap_lat", zcol="found", crs=4326)
 
 leafsync::sync(mv1, mv2)
@@ -46,8 +47,8 @@ ttm_join <- left_join(ttm_orig, ttm_snap, by=c("fromId", "toId"),
 
 street_net <- street_network_to_sf(r5r_core)
 
-street_net$vertices %>% mapview() +
-street_net$edges %>% mapview()
+street_net$vertices %>% mapview()
+street_net$edges %>% ggplot() + geom_sf()
 
 hex <- read_csv(system.file("extdata/poa", "poa_hexgrid.csv", package = "r5r"))