add multiband rasters in r

Question

i'm working with landsat data doing some topographic corrections using the method proposed by Hantson & Chuvieco (2010) which involves a cover separation method based on a certain NDVI value to tell between vegetated and non-vegetated areas.

So, I go from having one multiband raster object to two, which are processed sepparately and subsequently combined to form the original, topographically corrected, scene.

The problem is that I have some NA values in the scene, so if I don't get rid of them I get a blank image as a result. I use the sum function with the na.rm option defined as TRUE to tackle this problem, but the result is a one-band raster object instead of a stack or rasterbrick object.

This is my code:

#topographic correction
VEGB1TC<- topocorr(VEGB1,slopeG,aspectG,sunelev,sunazimuth,method="ccorrection",na.value=NA)
VEGB2TC<- topocorr(VEGB2,slopeG,aspectG,sunelev,sunazimuth,method="ccorrection",na.value=NA)
....        
SDB1TC<- topocorr(SDB1,slopeG,aspectG,sunelev,sunazimuth,method="ccorrection",na.value=NA)
SDB2TC<- topocorr(SDB2,slopeG,aspectG,sunelev,sunazimuth,method="ccorrection",na.value=NA)
...

#rasterize corrected band and stack them 
VB1<-raster(VEGB1TC)
VB2<-raster(VEGB2TC)
....

ADEVEG<-stack(VB1, VB2,...)

SB1<-raster(SDB1TC)
SB2<-raster(SDB2TC)
....

ADESD<-stack(SB1, SB2,...)

#-----------combine images-------

ADE<-sum(stack(ADEVEG,ADESD),na.rm=T)

Is there a function that can add multiband rasters with a na.rm tag? I've found the do.call function but it seems to me that it sums several rasters into a single-band object.

Thanks for your help.

ps: Hantson, S. Chuvieco, E. 2010. Evaluation of different topographic correction methods for Landsat imagery. International Journal of Applied Earth Observations and Geoinformation 13(2011):691-700p.

If you sum a number of raster layers you will get one band only. Or do you want to sum stacks ADEVEG and ADESD separately and then stack them? May be you want to do is ADE<-stack(sum(ADEVEG,na.rm=T),sum(ADESD,na.rm=T)) — rar
I want to fuse the two stacks into one as they represent the barren and vegetated portions of a scene, combining the multiple bands at the same time. In the end I just replaced NA vaues with 0 manualy and then proceed to add both stacks by simple adition. Your method is not what I intend but thanks a lot for the quick answer. — Benja Sotomayor

Robert Hijmans Robert Hijmans · Accepted Answer · 2016-07-08T01:19:23

Please provide a reproducible example

library(raster)
s1 <- stack(system.file("external/rlogo.grd", package="raster")) 
s2 <- flip(s1, 'x')
s2[1000:4000] <- NA

Now you can do

r1 <- s1 + s2

But that does not take care of NA values. So do you would need to resort to something like

r2 <- overlay(s1, s2, fun=function(x) sum(x, na.rm=TRUE))

or this

 s2r <- reclassify(s2, cbind(NA, 0))
 r3 <- s1 + s2r

PS. I understand that expected sum(s1, s2) to return a multi-layered object, but I think returning one layer is consistent with base::sum:

sum(1:4, 1:4)
# [1] 20

add multiband rasters in r

1 Answers