R/countryCellRast.R

In pik-piam/geodata: Landuse geodata library

### get the index/cells of a raster for each country
# cache=T allows to read from cache
# writecache=T allows to write to cache in case the file didn't exist yet




#' Return cell index values per country
#' 
#' This function returns a list of countries and all the cell indices of a
#' raster to it. The output can for instance be used to access the cell values
#' of individual countries. Function is based on the function cellFromPolygon
#' from the raster package and only considers cells whose center is within the
#' country's polygon.
#' 
#' 
#' @usage countryCellRast(x, countries="DIVA GIS", rm.length0=T, cache=T,
#' writecache=T)
#' @param x Raster file where cell index values should be calculated for
#' @param countries Dataset that is used to sperate between counties. Currently
#' only polygons from 'DIVA GIS' supported
#' @param rm.length0 If TRUE, all countries with no associated cells are
#' removed from the list
#' @param cache If TRUE, will look if the output already exists in cache and
#' read from there
#' @param writecache If TRUE, will also store the newly calculated output in
#' the cache folder
#' @return List of countries named by their ISO codes with raster cell indices
#' of each country
#' @author Ulrich Kreidenweis
#' @seealso \code{\link{countryFunRast}}
#' @examples
#' 
#' \dontrun{countryCellRast()}
#' 
#' @export countryCellRast
#' @importFrom raster res yres xres
countryCellRast <- function(x, countries="DIVA GIS", rm.length0=T,  cache=T, writecache=T) {
  
  if (countries != "DIVA GIS") stop("Currently only 'DIVA GIS' allowed as countries")
  
  # read the countries information this extraction is based on
  boundaries <- readCountries(countries = countries, cache=cache)
  
  fromcache <- FALSE
  
  ## check if the cache file exists, cache is true, and y and x resolution are the same. Only if all these are fullfilled load from cache
  if ( file.exists(paste0(geodata$config$mainfolder, "cache/countryCellRast/", countries,"_", signif(res(x)[1], 6),".rda")) & cache == T & xres(x) == yres(x) ) {
    fromcache <- TRUE
    load(file=paste0(geodata$config$mainfolder, "cache/countryCellRast/", countries,"_", signif(res(x)[1], 6),".rda"))
  
  } else {
    dat <- raster::cellFromPolygon(x, boundaries, weights=FALSE)
    ## name the list with the ISO codes
    names(dat) <- sapply(slot(boundaries, "polygons"), function(i) slot(i, "ID"))
  }
    
  # write to cache if this is enabled and if if hasn't been read from there
  if (writecache==T & fromcache==F) {
    print("writing cell index to cache")
    # create the functions cache directory if it doesn't exist yet
    dir.create(file.path(geodata$config$mainfolder, "cache/countryCellRast"), showWarnings = FALSE)
    # write the result to the cache
    save(dat, file=paste0(geodata$config$mainfolder, "cache/countryCellRast/", countries,"_", signif(res(x)[1], 6),".rda"))
  }
  
  # remove elements from the list which have no values
  if (rm.length0==T) {
    dat <- dat[lapply(dat,length)>0]
  }
  
  return(dat)
}