#' @title Cut raster based on a shape
#' @name cut.raster
#'
#' @description A function to cut a raster based on an informed shapefile.
#'
#' @param abio an object with the rasters to be cut. Accepts an object of type RasterStack, generated by the function \code{\link[raster]{stack}}
#' @param shape.dir here are two options. The first is to inform an object in the shapefile format. The second option is to inform a directory containing the shapefile.
#' @param extension extent of the cut rasters. Default is .asc, see\code{\link[raster]{write.raster}} for more output format possibilities.
#' @param br logical. If TRUE, use the Brazil shape provided by \code{\link[maptools]{wrld_simpl}}. If FALSE (default), use the informed shape.
#' @param plot logical. Plots one of the cut rasters.
#' @param trim logical. If TRUE, exclude pixels with NA (very slow). If FALSE (default), keep the NAs.
#'
#' @details When the shapeir argument is a folder where it contains the reference shape, this folder must contain only the shape and its auxiliary files.
#'
#' @return Raster files in a folder "Cortados".
#'
#' @author Diogo S. B. Rocha
#'
#' @seealso \code{\link[raster]{crop}}, \code{\link[raster]{mask}}, \code{\link[raster]{stack}}
#'
#' @examples
#'
#' fnames <- list.files(path=paste(system.file(package="dismo"), '/ex', sep=''), pattern='grd', full.names=TRUE )
#' predictors <- raster::stack(fnames)
#' cut.raster(abio = predictors, br=T, extension = ".tif")
#'
#' @import raster
#'
#' @export
cut.raster =
function(abio, shape.dir, extension = ".asc", br = FALSE, plot = TRUE, diretorio = "cortados", trim = FALSE){
if (dir.exists(paste0(diretorio)) == F) {
dir.create(paste0(diretorio))
}
# Definir shape para cortar
if(missing(shape.dir)){
# Possibilidade de cortar para o Brasil
#importando shape do brasil
if(br==T){
data(wrld_simpl, package = "maptools")
shape.dir=subset(wrld_simpl, wrld_simpl$NAME=="Brazil")
shape=shape.dir
}else(stop("The directory containing the shape to cut the raster was not specified"))
}
if(class(shape.dir)== "character") {
shape = rgdal::readOGR(list.files(shape.dir, pattern = ".shp", full.names = T)[1])
} else (shape = shape.dir)
crs(shape) = crs(abio)
# loop para cortar todos os rasters
# sem trim
if (trim == F) {
ini = Sys.time()
for (i in 1:length(names(abio))) {
ini1 = Sys.time()
mask(crop(abio[[i]], extent(shape)), shape, filename = paste("./Cortados", "/", names(abio)[i],
sep = "", extension), overwrite = TRUE)
print(Sys.time())
cat("\n", paste("Tá indo", i))
fim1 = Sys.time()
cat(paste("\n", round(as.numeric(fim1 - ini1), 2), units(fim1 - ini1)))
if (i == length(names(abio))) {
cat("\n", "Acabou!")
fim = Sys.time()
fim - ini
}
}
}
# com trim
if (trim == T) {
ini = Sys.time()
for (i in 1:length(names(abio))) {
ini1 = Sys.time()
writeRaster(trim(mask(crop(abio[[i]], extent(shape)), shape)),
filename = paste(".\\Cortados", "\\", names(abio)[i], ".tif", sep = ""),
format = "GTiff", overwrite = TRUE, NAflag = -9999)
cat("\r")
print(Sys.time())
cat("\n", paste("T? indo", i))
fim1 = Sys.time()
cat(paste("\n", round(as.numeric(fim1 - ini1), 2), units(fim1 - ini1)))
if (i == length(names(abio))) {
cat("\n", "Acabou!", "\n")
fim = Sys.time()
fim - ini
}
}
}
# plotando a primeira variavel cortada
if (plot == T) {
plot(raster(list.files("./cortados", pattern = "", full.names = TRUE)[1], native = TRUE))
plot(shape, add = T)
}
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.