Add logic for a smarter default z setting based on extent of request

[mdsumner]

I'd like to see the ability to get a raster for a specific extent, with some logic for a default scale. That way, you give it any Spatial/sf/raster object and it gives you a "reasonable" resolution within that exact extent.

The motivation is that extent is usually the first natural way to request, but then you might tweak the scale to pull more or less detail as an alternative. As it is the scale expands/contracts the query extent. Both modes are important to have of course, so perhaps this is a new function with slightly different logic.

Inspiration might come from dismo::gmap and/or raster::getData. I'll try to contribute if I can.

(And my apologies if it's already possible)

[jhollist]

Like both ideas!

1.) specific extent could just be an added argument. Currently is returns the full tiles that overlap the mins and maxs of the extent. That seems like a reasonable default to me. Adding the crop/mask to the exact extent would be nice.

2.) Be great to have some logic on picking the zoom level. I currently just have a default zoom. Any thoughts on how to specify?

[mdsumner]

Oh good, it did seem to be tile-aligned but I didn't quite get to checking.

I think start with a loose-ish target maximum of pixels, determine the area and work out the resolution. You can leverage raster to do all the work actually. Something like this:

library(raster)

maxpix <- 1e6  ## something

ex <- extent(0, 1000, 0, 500)  ## again, could be anything 

width <- xmax(ex) - xmin(ex)
height <- ymax(ex) - ymin(ex)

## average out the width/height for the appropriate gridded extent
raster(ex, res = sqrt(width * height / maxpix))

That res value translates into the scale somehow that I'd have to look up. Something I really want to do is extract all this nice logic from raster into a more abstract package that other packages can use. Might be a good example of how/what to do that for here. Elsewhere I also use a "buffer" on an extent, so that I can always get a tidy alignment for an exact extent, that might be useful too:

#' Whole cell buffer for an extent
#'
#' Ensure a raster extent aligns to whole parts.
#' @param e1 input \code{\link[raster]{extent}}
#' @param e2 grain size
#' @examples
#' library(raster)
#' bufext(extent(0.1, 2.2, 0, 3), 2)
#' @importFrom raster extent xmin xmax ymin ymax
#' @export
bufext <- function(e1, e2) {
  e1 <- extent(e1)
  if (e2 == 0) return(e1)
  num0 <- as.double(e1)
  extent((num0 %/% e2) * e2 + c(0, e2, 0, e2))
}

as.double.Extent <- function(x, ...) {
  c(xmin(x), xmax(x), ymin(x), ymax(x))
}

as.integer.Extent <- function(x, ...) {
  as.integer(as.double(x))
}

(Good night! I have no more for now, thanks again!)

[jhollist]

Got a new release coming soon. Did add ability to clip resultant DEM with a bounding box or the input spatial data. Doesn't really get to this issue, but at least easier way to not just get full tiles.

I do like this idea of having smarter defaults via the code here or in MilesMcBain/slippymath. Will wait till that is stable and on CRAN and then add in for a future release. Keeping this live though, so I don't forget!