I am working with daily observation of climate data organized in .nc files. I read them using the stack command of the raster package. Each file (corresponding to a year) is a RasterStack element with the following characteristics:
class : RasterStack
dimensions : 360, 720, 259200, 365 (nrow, ncol, ncell, nlayers)
resolution : 0.5, 0.5 (x, y)
extent : -180, 180, -90, 90 (xmin, xmax, ymin, ymax)
each layer is the raster of values of a day.
I would like to sum the layers in order to calculate the monthly values. I believe the solution should be using calc or stackApply {raster}, but I couldn't find a way to sum from layer x to layer y or a way to subset the RasterStack before summing.
I prepared an example file with only 12 layers (to reduce the size).
I don't exactly know how to propose a code, sorry, but it should be something like:
library(raster)
setwd("myfolder")
data<-stack(mydata.nc)
datasum<- stackApply(data, ??? ,fun=sum)
Thank you
Answer
You can use stackApply to do this. Using your example data, it looks like the name of each raster layer is the date. You can use that to build the indices that you need to pass to stackApply.
The list of indices needs to have 31 1s for the January days etc.
You could do:
#get the date from the names of the layers and extract the month
indices <- format(as.Date(names(data), format = "X%Y.%m.%d"), format = "%m")
indices <- as.numeric(indices)
#sum the layers
datasum<- stackApply(data, indices, fun = sum)
The result will be a raster stack of 12 layers.
To subset raster layers from the stack, you can do data[[c(1,2]]