Nothing
#' @keywords plotting
#' @export
#'
#' @title Create a rasterBrick from an Esri tiled image server
#'
#' @param centerLon Map center longitude.
#' @param centerLat Map center latitude.
#' @param maptype Selects the appropriate Esri tile server. Options include:
#' \itemize{
#' \item "world_topo"
#' \item "world_imagery"
#' \item "world_terrain"
#' \item "de_Lorme"
#' \item "world_grey"
#' \item "world_streets"
#' }
#' @param zoom map Zoom level.
#' @param width Width of image, in pixels.
#' @param height Height of image, in pixels.
#' @param bbox Bounding box vector (lonLo, latLo, lonHi, latHi). If not null,
#' \code{centerLon}, \code{centerLat}, and \code{zoom} are ignored.
#' @param maxTiles Maximum number of tiles to be returned. The greater the
#' number, the slower the performance -- arbitrarily set to 20 by default.
#' @param crs Object of class CRS. The Coordinate Reference System (CRS) for the
#' returned map. If the CRS of the downloaded map does not match, it will be
#' projected to the specified CRS using \code{raster::projectRaster}.
#' @param tileCacheDir Optional location for cached tiles.
#'
#' @description Uses the input coordinates to fetch and composite a raster from
#' the tile server. Returns a \code{raster::rasterBrick} object. This
#' can then passed as the \code{rasterBrick} object to the
#' \code{staticmap_plotRasterBrick()} function for plotting.
#'
#' As of July 2019, this list is a handy reference to the freely available tile
#' servers which can be previewed at the following URL:
#'
#' \url{https://leaflet-extras.github.io/leaflet-providers/preview/}
#'
#' @note The spatial reference of the image when it is downloaded is 3857. If
#' the crs argument is different, projecting may cause the size and extent of
#' the image to differ very slightly from the input, on a scale of 1-2 pixels or
#' 10^-3 degrees.
#'
#' If bbox is specified and the bbox aspect ratio does not match the
#' width/height aspect ratio the extent is resized to prevent the map image from
#' appearing stretched, so the map extent may not match the bbox argument
#' exactly.
#'
#' @note If both \code{zoom} and \code{maxTiles} are specified, \code{maxTiles}
#' takes precedence. To get a specified zoom level, set \code{maxTiles = NULL}.
#'
#' @return A rasterBrick object which can be plotted with
#' \code{staticmap_plotRasterBrick()} or \code{raster::plotRGB()} and serve as a
#' base plot.
#'
#' @examples
#' \dontrun{
#' # Fail gracefully if any resources are not available
#' try({
#'
#' rasterBrick <- staticmap_getEsrimapBrick(-122.3318, 47.668)
#' staticmap_plotRasterBrick(rasterBrick)
#'
#' }, silent = FALSE)
#' }
#' @seealso \code{\link{staticmap_getRasterBrick}}
#' @seealso \code{\link{staticmap_plotRasterBrick}}
staticmap_getEsrimapBrick <- function(
centerLon = NULL,
centerLat = NULL,
maptype = "world_topo",
zoom = 12,
width = 640,
height = 640,
bbox = NULL,
maxTiles = 20,
crs = sp::CRS("+init=epsg:4326"),
tileCacheDir = tempdir()
) {
# ===== DEBUGGING ============================================================
if ( FALSE ) {
centerLon = -117.2554
centerLat = 47.68013
maptype <- "world_topo"
zoom <- 12
width <- 640
height <- 640
bbox <- NULL
maxTiles <- 40
crs <- sp::CRS("+init=epsg:4326")
tileCacheDir <- tempdir()
}
# ----- Validate parameters --------------------------------------------------
validMapTypes <- c("world_topo", "world_imagery", "world_terrain",
"de_Lorme", "world_gray", "world_streets")
if ( !maptype %in% validMapTypes )
stop(paste0("Required parameter maptype = '", maptype, "' is not recognized."))
if ( is.null(centerLon) || is.null(centerLat) ) {
if ( is.null(bbox) ) {
stop("Required arameters centerLat + centerLon, or bbox must be specified")
}
}
# ----- Determine bounding box -----------------------------------------------
# Projection information:
# https://spatialreference.org/
# https://en.wikipedia.org/wiki/Web_Mercator_projection
# EPSG:4326 -- lon lat pairs on the WGS84 reference ellipsoid
# EPSG:3857 -- most (all?) tiling servers have a "web Mercator" projection
# Calculate degrees per pixel from zoom to determine bbox:
# * google maps tiles are 256x256 pixels
# * zoom level 0 includes 360 degrees-EW
# * every time zoom level increases, scale halves. Thus, for degrees EW:
# pixels/degree = 256pixels*2^zoomLevel/360degrees
# degrees/pixel = 360degrees/(256pixels*2^zoomLevel)
# degrees-NS/pixel = degreesPerPixelNS*cos(latitude)
# Source: https://gis.stackexchange.com/questions/7430/what-ratio-scales-do-google-maps-zoom-levels-correspond-to
# Create lon/lat ranges
if ( !is.null(centerLon) && !is.null(centerLat) ) {
degreesPerPixelEW <- 360/(256 * 2^zoom)
degreesPerPixelNS <- degreesPerPixelEW*cos(centerLat * pi/180) # radians
xmin <- centerLon - degreesPerPixelEW*(width/2 - 0.5)
xmax <- centerLon + degreesPerPixelEW*(width/2 - 0.5)
ymin <- centerLat - degreesPerPixelNS*(height/2 - 0.5)
ymax <- centerLat + degreesPerPixelNS*(height/2 - 0.5)
bbox <-
sf::st_bbox(c(xmin = xmin,
xmax = xmax,
ymin = ymin,
ymax = ymax),
crs = crs)
} else if ( !is.null(bbox) ) {
xmin <- bbox[1]
xmax <- bbox[3]
ymin <- bbox[2]
ymax <- bbox[4]
bbox <-
sf::st_bbox(c(xmin = bbox[1],
xmax = bbox[3],
ymin = bbox[2],
ymax = bbox[4]),
crs = crs)
} else {
stop("centerLat + centerLon, or bbox must be specified")
}
# convert sf bbox to slippymath tile grid
tile_grid <- slippymath::bbox_to_tile_grid(bbox,
zoom = zoom,
max_tiles = maxTiles)
# Use slippymath zoom level (because we set max_tiles)
zoom <- tile_grid$zoom
# ----- Download tiles -------------------------------------------------------
# Examples of additional, non-Esri tile servers can be found here:
# https://leaflet-extras.github.io/leaflet-providers/preview/
esriBaseUrl <- "https://server.arcgisonline.com/ArcGIS/rest/services/"
slippyMapApi <- "/MapServer/tile/{zoom}/{y}/{x}"
world_topo <- paste0(esriBaseUrl, "World_Topo_Map", slippyMapApi)
world_imagery <- paste0(esriBaseUrl, "World_Imagery", slippyMapApi)
world_terrain <- paste0(esriBaseUrl, "World_Terrain_Base", slippyMapApi)
de_Lorme <- paste0(esriBaseUrl, "Specialty/DeLorme_World_Base_Map", slippyMapApi)
world_gray <- paste0(esriBaseUrl, "Canvas/World_Light_Gray_Base", slippyMapApi)
world_streets <- paste0(esriBaseUrl, "World_Street_Map", slippyMapApi)
tileserverQueryString <- eval(parse(text=maptype))
# curl images
images <-
purrr::pmap(tile_grid$tiles,
function(x, y, zoom) {
filename <- glue::glue("{maptype}_{zoom}_{x}_{y}.jpg")
filepath <- file.path(tileCacheDir, filename)
if ( !file.exists(filepath) ) {
curl::curl_download(url = glue::glue(tileserverQueryString),
destfile = filepath)
}
return(filepath)
},
zoom = tile_grid$zoom)
# ----- Assemble image -------------------------------------------------------
# Uncropped image
rasterBrick <- slippymath::compose_tile_grid(tile_grid, images)
# Determine the actual extent of tiles we received
e <- raster::extent(rasterBrick)
mercatorMatrix <- matrix(c(e@xmin, e@xmax, e@ymin, e@ymax), nrow = 2)
lonLatMatrix <- slippymath::merc_to_lonlat(mercatorMatrix)
# Assign a geographic coordinate system
raster::extent(rasterBrick) <- c(
lonLatMatrix[1,1],
lonLatMatrix[2,1],
lonLatMatrix[1,2],
lonLatMatrix[2,2]
)
raster::crs(rasterBrick) <- sp::CRS("+init=epsg:4326")
# ----- Crop image -----------------------------------------------------------
extent <- raster::extent(xmin, xmax, ymin, ymax)
rasterBrick <- raster::crop(rasterBrick, extent)
raster::crs(rasterBrick) <- sp::CRS("+init=epsg:4326")
# ----- Return ---------------------------------------------------------------
return(rasterBrick)
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.