Eddie, (I'm answering through the list, perhaps somebody else is interested)
You need one element in "indices" for each input layer, so if you have 5 layers you should use mx <- stackApply(s,indices=c(1,1,1,1,1),max) The actual value of each index indicates to which output layer that input layer is contributing. In your case, there is only one output layer, thus all values must be 1. But see the help page for an example with different index values. Agus On Tue, Oct 8, 2013 at 9:41 PM, Agustin Lobo <agustin.l...@ictja.csic.es> wrote: > ---------- Forwarded message ---------- > From: Eddie Smith <eddie...@gmail.com> > Date: Tue, Oct 8, 2013 at 7:40 PM > Subject: Re: [R-sig-Geo] Creating a maximum raster layer from RasterStack > To: agustin.l...@ictja.csic.es > > > Dear Agus, > > Thank you very much for your reply. Your script below is working > however, would you mind if I ask a basic question on what is the > purpose of [indices=c(1,1)] ? I tried this code with 5 raster layers > and the error message is >> mx <- stackApply(s,indices=c(1,1),max) > Warning message: > In stackApply(s, indices = c(1, 1), max) : > number of items to replace is not a multiple of replacement length > > I am sorry if my question is too basic. I am new to R. > > Thank you. > > library(raster) > r1 <- raster(matrix(rnorm(100,5,10),10,10)) > r2 <- raster(matrix(rnorm(100,2,30),10,10)) > s <- stack(r1,r2) > mx <- stackApply(s,indices=c(1,1),max) > miwhichmax <- function(x,na.rm=TRUE) which.max(x) > mxd <- stackApply(s,indices=c(1,1),miwhichmax) > > > On Mon, Oct 7, 2013 at 7:20 PM, Agustin Lobo <alobolis...@gmail.com> wrote: >> >> oops, copied the wrong line, >> it must be >> mx <- stackApply(s,indices=c(1,1),max) >> >> Agus >> >> >> On Mon, Oct 7, 2013 at 8:18 PM, Agustin Lobo <alobolis...@gmail.com> wrote: >> > What about this: >> > r1 <- raster(matrix(rnorm(100,5,10),10,10)) >> > r2 <- raster(matrix(rnorm(100,2,30),10,10)) >> > s <- stack(r1,r2) >> > mx <- stackApply(s,max) >> > >> > which.max() fails because raster requires the na.rm argument (I think, >> > Robert?) >> > >> > mxd <- stackApply(s,indices=c(1,1),which.max) >> > Error in FUN(newX[, i], ...) : unused argument (na.rm = TRUE) >> > >> > so make your own function icluding the na.rm argument: >> > >> > miwhichmax <- function(x,na.rm=TRUE) which.max(x) >> > mxd <- stackApply(s,indices=c(1,1),miwhichmax) >> > >> > Agus >> > >> > >> > On Thu, Oct 3, 2013 at 6:55 PM, Eddie Smith <eddie...@gmail.com> wrote: >> >> Dear list, >> >> >> >> I have a RasterLayer that contains 365 bands with a dimension of 1032, >> >> 1656, 1708992 (nrow, ncol, ncell). >> >> >> >> library(raster) >> >> img <- list.files(pattern='\\.img$') >> >> stack <- stack(img) >> >> data2002 <- writeRaster(stack, filename="2002.tif", format="GTiff", >> >> overwrite=TRUE) >> >> >> >> 1. How could I produce a raster layer that contains a maximum value from >> >> 365 bands for each pixel? >> >> I tried maxValue() but that only give me a maximum value for each band. >> >> >> >> 2. If I manage to do step 1, is it possible to know from which band is >> >> actually the maximum value originate? >> >> >> >> [[alternative HTML version deleted]] >> >> >> >> _______________________________________________ >> >> R-sig-Geo mailing list >> >> R-sig-Geo@r-project.org >> >> https://stat.ethz.ch/mailman/listinfo/r-sig-geo > > > > > -- > -- > Dr. Agustin Lobo > Institut de Ciencies de la Terra "Jaume Almera" (CSIC) > Lluis Sole Sabaris s/n > 08028 Barcelona > Spain > Tel. 34 934095410 > Fax. 34 934110012 > e-mail agustin.l...@ictja.csic.es _______________________________________________ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo
