isoband stuff

[mdsumner]

library(raster)
#' Resample a raster, fast
#' 
#' Equivalent to `raster::aggregate` where `fun` returns a central value, but 
#' without taking any care for missing values or summarizing any values. 
#' 
#' @param x raster
#' @param fact sub-sampling factor, or up-sampling factor if < 1
#' @param method simple or bilinear, as per `raster::extract`
#' @examples 
#' decimate(raster(volcano))
#' 
#' ## also samples up
#' decimate(raster(volcano), .2)
decimate <- function(x, fact = 10, method = "simple") {
     if (fact == 1) return(x)
     out <- raster(x)
     dm <- dim(x)
     dim(out) <- c(dm[1L]/fact, ncols = dm[2L]/fact)
     setValues(out, raster::extract(x, coordinates(out), method = method))
}
#' Isoband wrapper for raster and sf
#' 
#' Input a raster and return sf polygons. Set low and high values of the band/s
#' with arguments `lo` and `hi`. Use `auto` to generate these automatically. 
#' @param x raster
#' @param lo vector of low band values
#' @param hi vector of high band values
#' @param auto if `TRUE` ignore hi and lo arguments and find values automatically from data range
#' @examples
#' ib <- isoband_raster(raster(volcano), c(94, 100, 120, 150), c(100, 120, 150, 195))
isoband_raster <- function(x, lo, hi, auto = FALSE) {
    if (auto) {
        breaks <- pretty(na.omit(values(x)))
        lo <- head(breaks, -1)
        hi <- tail(breaks, -1)
    }
    if (raster::nlayers(x) > 1) warning("multi-layer raster given, only first layer used")
    b <- isoband::isobands(xFromCol(x), yFromRow(x), as.matrix(x[[1]]), levels_low = lo, levels_hi = hi)
    sf::st_sf(lo = lo, hi = hi, geometry  = sf::st_sfc(isoband::iso_to_sfg(b), crs = projection(x)))
}
## examples

## sea surface temperature (C)
sst <- readRDS("sst.rds")
## sea surface ice concentation (%, southern)
ice <- readRDS("ice.rds")
## sea surface height (cm)
ssh <- readRDS("ssh.rds")
topo <- readRDS("topo.rds")

isoband_raster(decimate(ssh, 2), -1.3, 0.75)
isoband_raster(raster::focal(ssh, matrix(1, 3, 3)), -1.3, 0.75)

## ok
isoband_raster(sst, auto = TRUE)
## ok
iso <- isoband_raster(sst, 0, 5)