R/hf_hike.R

In kbvernon/hiker: Model Hiking Paths Over Costly Terrains

#' Shortest path between points
#'
#' Calculate the "shortest" or least-cost paths between multiple start and end points.
#' For each from-point, the shortest path to every to-point is calculated.
#'
#' @param x a cost `terrain` generated by `hf_terrain()`.
#' @param from an `sf` specifying start locations, must be POINT or MULTIPOINT.
#' @param to an `sf` specifying end locations, must be POINT or MULTIPOINT.
#' @param add_cost logical, whether to include a column of travel costs.
#'
#' @return An `sf` object with least-cost `LINE` features and two columns: from and to.
#'   Values in the from- and to-columns are `1:nrow(from)` and `1:nrow(to)`, respectively.
#'   If `add_cost = TRUE`, the `sf` object also includes a travel cost column.
#'
#' @export
#'
#' @examples
#' library(sf)
#' library(terra)
#'
#' fn <- system.file("extdata/red_butte_dem.tif", package = "hiker")
#' red_butte_dem <- rast(fn)
#'
#' from <- st_sf(geometry = st_sfc(st_point(c(432000, 4514000)),
#'                                 st_point(c(434000, 4518000)),
#'                                 crs = 26912))
#'
#' to <- st_sf(geometry = st_sfc(st_point(c(431000, 4515000)),
#'                               st_point(c(436500, 4518500)),
#'                               crs = 26912))
#'
#' terrain <- hf_terrain(red_butte_dem)
#'
#' short_paths <- hf_hike(terrain, from, to)
#'
#' plot(red_butte_dem)
#' plot(st_geometry(short_paths), lty = 2, add = TRUE)
#' plot(st_geometry(from), pch = 19, col = "red2", add = TRUE)
#' plot(st_geometry(to), pch = 19, col = "dodgerblue4", add = TRUE)
#'
hf_hike <- function(x, from, to, add_cost = FALSE) {

  stop_if_not_terrain(x)

  stop_if_not_sf(from, to)

  stop_if_not_point(from, to)

  stop_if_not_crs_equal(x$crs, from, to)

  from_xy <- sf::st_coordinates(from)[, 1:2, drop = FALSE]
  to_xy   <- sf::st_coordinates(to)[, 1:2, drop = FALSE]

  rr <- terra::rast(
    nrow   = x$nrow,
    ncol   = x$ncol,
    extent = terra::ext(x$bb8),
    crs    = x$crs
  )

  from_cells <- terra::cellFromXY(rr, from_xy)
  to_cells   <- terra::cellFromXY(rr, to_xy)

  graph <- igraph::graph_from_adjacency_matrix(
    x$conductance,
    mode = "directed",
    weighted = TRUE
  )

  # invert conductance to get travel cost
  igraph::E(graph)$weight <- (1 / igraph::E(graph)$weight)

  # return a list of lists of vectors of vertices on paths
  shorties <-
    lapply(
      from_cells,
      function(z) {

        res <-
          igraph::shortest_paths(
            graph,
            from = z,
            to = to_cells,
            mode = "out"
          )

        res$vpath

      })

  # collapse (list of lists of vectors) to (list of vectors)
  shorties <- unlist(shorties, recursive = FALSE)

  # now we can just make a bunch of linestrings out of them
  line_list <-
    lapply(
      shorties,
      function(z) {

        cells <- as.vector(z)

        # to handle case where from-vertex == to-vertex
        if (length(cells) == 1) cells <- rep(cells, 2)

        xy <- terra::xyFromCell(rr, cells)

        sf::st_linestring(xy)

      }
    )

  # and then an sf
  sf_col <- sf::st_sfc(line_list, crs = x$crs)

  short_paths <- sf::st_sf(geometry = sf_col)

  # add indices
  short_paths <- transform(
    short_paths,
    from = rep(1:nrow(from_xy), each = nrow(to_xy)),
    to   = rep(1:nrow(to_xy), times = nrow(from_xy))
  )

  if (add_cost) {

    cost <- igraph::distances(graph, from_cells, to_cells, mode = "out")

    short_paths <- transform(short_paths, cost = c(t(cost)))

  }

  return(short_paths)

}