tmsppc
tmsppc

Reputation: 117

Raster series sum

I have a loop with a series of raster images and I want to extract the values ​​equal to 150 and then add the total amount of pixels for the entire length of the loop. With the code that I have only managed to get the total values ​​for each image separately and not in total form. Thanks

   m=52419 #total pixels basin
   for(i in 1:4){
   b1<-raster(myras1[i])
   bc = b1 == 150 #Values eq 150
   nbc = cellStats(bc,stat="sum")
   print(nbc)
   [1] 34962
   [1] 38729
   [1] 52389
   [1] 52176
   pc=nbc*100/m
   }

Upvotes: 1

Views: 1235

Answers (1)

dww
dww

Reputation: 31452

In general, it is recommended not to use loops in R for things that can be easily vectorised. Rather than trying to fix the (several) problems with your loop, I show instead a better way. You can perform the whole calculation in a single vectorised line:

sum(cellStats(myras1==150, stat="sum")) * 100/m

Breaking this down: cellStats performed on a raster stack will return a vector of values, one for each layer. sum then adds these together. Then we divide by number of cells in the whole stack (all layers combined) and multiply by 100 to convert to a percnetage.

Testing this on some reproducible dummy test data:

set.seed(123)
myras1 = list(
  raster(nrows = 100, ncols = 100, vals = sample(140:150,10000,T)),
  raster(nrows = 100, ncols = 100, vals = sample(140:150,10000,T)),
  raster(nrows = 100, ncols = 100, vals = sample(140:150,10000,T)),
  raster(nrows = 100, ncols = 100, vals = sample(140:150,10000,T))
)
myras1 = stack(myras1)
m = ncol(myras1) * nrow(myras1) * nlayers(myras1)

sum(cellStats(myras1==150, stat="sum")) * 100/m
# [1] 8.815

Upvotes: 2

Related Questions