Code for my map for the 30 Day Map Challenge day 28 (November 28 2021) "The world is not flat!"

20211128.R

library(terrainr)
library(raster)
library(rayshader)
library(sf)

baker <- tmaptools::geocode_OSM("Mount Baker")
baker <- st_as_sf(as.data.frame(t(as.data.frame(baker$coords))), coords = c("x", "y"), crs = 4326)
zscale <- 5
baker <- set_bbox_side_length(
  baker,
  zscale * 3000
) |> 
  get_tiles(resolution = zscale, services = c("elevation")) |> 
  lapply(merge_rasters)

surface_matrix <- raster_to_matrix(raster(baker[[1]]))

color_pal <- function(n = 255, bias = 1) {
  pal <- colorRampPalette(
    rev(c("#e63946","#f1faee","#a8dadc","#457b9d","#1d3557")),
    bias = bias
  )
  return(pal(n))
}

surface_hillshade <- surface_matrix %>% 
  height_shade(texture = color_pal()) %>% 
  add_shadow(texture_shade(surface_matrix, detail = 0.9, brightness = 15), 0.7) %>% 
  add_shadow(ray_shade(surface_matrix, sunaltitude = 60, zscale = zscale, multicore = TRUE), 0) |> 
  add_shadow(ambient_shade(surface_matrix, zscale = zscale))

plot_3d(surface_hillshade,
        surface_matrix,
        windowsize = c(2160, 1440),
        zscale = zscale,
        zoom = 0.5,
        solid = FALSE,
        theta = -35,
        phi = 22.5)

render_highquality(
  "20211128.png",
  lightintensity = 1000,
  cache_filename = "tk.obj"
)