Skip to content

Instantly share code, notes, and snippets.

@mdsumner
Last active September 23, 2021 07:00
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save mdsumner/cab31b2b3cd30768711bcad934e31917 to your computer and use it in GitHub Desktop.
Save mdsumner/cab31b2b3cd30768711bcad934e31917 to your computer and use it in GitHub Desktop.

https://twitter.com/mdsumner/status/1440933887354441731?s=20

## build mapbox WMTS data source name from a custom style
mapbox_wmts <- function(type = "", key = NULL) {
  if (nchar(type) < 1) stop("type must be a useable mapbox style name e.g. 'styles/v1/<username>/<style_id>'")
  if (is.null(key)) {
    key <- Sys.getenv("MAPBOX_API_KEY")
    if (is.null(key)) stop("provide mapbox 'key' or via MAPBOX_API_KEY env var")
  }
  base <- "WMTS:https://api.mapbox.com/%s/wmts?access_token=%s"
  sprintf(base,
          type, key)
}

## styles I have in my account, published publically
## Satellite
sat <- mapbox_wmts("styles/v1/mdsumner/cjy6m24oh1amo1cmy74uk23n6")
## Basic
bas <- mapbox_wmts("styles/v1/mdsumner/ckb4o07v00s5i1irmzw7obvbr")

## install.packages("vapour")
## remotes::install_github("hypertidy/gdalio)
library(gdalio)  
## load format-specific functions like gdalio_raster()
source(system.file("raster_format/raster_format.codeR", package = "gdalio", mustWork = TRUE))
## set desired grid: extent, dimension, projection
g <- list(extent = c(-1, 1, -.5, .5) * 3e5, dimension = 768 * c(1, 0.5), 
           projection = "+proj=laea +lon_0=138 +lat_0=-35")
gdalio_set_default_grid(g)
sat_im <- gdalio_terra(sat, bands = 1:3, resample = "cubic", band_output_type = "Int32")
bas_im <- gdalio_terra(bas, bands = 1:3, resample = "cubic", band_output_type = "Int32")
op <- par(mfrow = c(2, 1))
terra::plotRGB(sat_im)
terra::plotRGB(bas_im)
par(op)


@mdsumner
Copy link
Author

image

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment