[GIS] How to make Maximum values of a stacked raster to form a new raster in R using raster package

rraster

I am calculating the Maximum GPP of a year based on a 8-day GPP raster layers. I want to know how to calculate the maximum of the stacked raster at the position of each pixel and to output a new raster. I know Cell statistics in Arcmap can solve this problem,but I have to calculate a lot of data. max() and maxValue() in "raster" package are not worked.

My code is :

library(raster)
GPP2000049=raster("GPP.2000049.h23v03.tif") GPP2000057=raster("GPP.2000057.h23v03.tif")
GPP2000065=raster("GPP.2000065.h23v03.tif") GPP2000073=raster("GPP.2000073.h23v03.tif")
GPP2000081=raster("GPP.2000081.h23v03.tif") GPP2000089=raster("GPP.2000089.h23v03.tif")
GPP2000097=raster("GPP.2000097.h23v03.tif") GPP2000105=raster("GPP.2000105.h23v03.tif")
GPP2000113=raster("GPP.2000113.h23v03.tif") GPP2000121=raster("GPP.2000121.h23v03.tif")
GPP2000129=raster("GPP.2000129.h23v03.tif") GPP2000137=raster("GPP.2000137.h23v03.tif")
GPP2000145=raster("GPP.2000145.h23v03.tif") GPP2000153=raster("GPP.2000153.h23v03.tif")
GPP2000161=raster("GPP.2000161.h23v03.tif") GPP2000169=raster("GPP.2000169.h23v03.tif")
GPP2000177=raster("GPP.2000177.h23v03.tif") GPP2000185=raster("GPP.2000185.h23v03.tif")
GPP2000193=raster("GPP.2000193.h23v03.tif") GPP2000201=raster("GPP.2000201.h23v03.tif")
GPP2000209=raster("GPP.2000209.h23v03.tif") GPP2000217=raster("GPP.2000217.h23v03.tif")
GPP2000225=raster("GPP.2000225.h23v03.tif") GPP2000233=raster("GPP.2000233.h23v03.tif")
GPP2000241=raster("GPP.2000241.h23v03.tif") GPP2000249=raster("GPP.2000249.h23v03.tif")
GPP2000257=raster("GPP.2000257.h23v03.tif") GPP2000265=raster("GPP.2000265.h23v03.tif")
GPP2000273=raster("GPP.2000273.h23v03.tif") GPP2000281=raster("GPP.2000281.h23v03.tif")
GPP2000289=raster("GPP.2000289.h23v03.tif") GPP2000297=raster("GPP.2000297.h23v03.tif")
GPP2000305=raster("GPP.2000305.h23v03.tif") GPP2000313=raster("GPP.2000313.h23v03.tif")
GPP2000321=raster("GPP.2000321.h23v03.tif") GPP2000329=raster("GPP.2000329.h23v03.tif")
GPP2000337=raster("GPP.2000337.h23v03.tif") GPP2000345=raster("GPP.2000345.h23v03.tif")
GPP2000353=raster("GPP.2000353.h23v03.tif") GPP2000361=raster("GPP.2000361.h23v03.tif")

MaxGPP2000h23v03=stack(GPP2000049,GPP2000057,GPP2000065,GPP2000073,
                       GPP2000081,GPP2000089,GPP2000097,GPP2000105,
                       GPP2000113,GPP2000121,GPP2000129,GPP2000137,
                       GPP2000145,GPP2000153,GPP2000161,GPP2000169,
                       GPP2000177,GPP2000185,GPP2000193,GPP2000201,
                       GPP2000209,GPP2000217,GPP2000225,GPP2000233,
                       GPP2000241,GPP2000249,GPP2000257,GPP2000265,
                       GPP2000273,GPP2000281,GPP2000289,GPP2000297,
                       GPP2000305,GPP2000313, GPP2000321,GPP2000329,
                       GPP2000337,GPP2000345,GPP2000353,GPP2000361)

MaxGPP2000h23v03=max(MaxGPP2000h23v03)
writeRaster(MaxGPP2000h23v03,filename = "MaxGPP2000h23v03.tif",overwrite=TRUE)

Best Answer

You can either use stackApply or calc to reduce a rasterStack or RasterBrick to the maximum value of each pixel.

library(raster)

# Create a test RasterBrick
r <- raster(ncol=10, nrow=10)
r[]=1:ncell(r)
b <- brick(r,r,r,r,r,r)

stackApply(b, indices = rep(1, nlayers(b)), fun = max)

calc(b, function(x){max(x)})
Related Question