füge Multiband-Raster in r hinzu

Question

Ich arbeite mit Landsat-Daten, die einige topographische Korrekturen unter Verwendung der von Hantson & Chuvieco (2010) vorgeschlagenen Methode machen, die eine Deckungstrennungsmethode basierend auf einem bestimmten NDVI-Wert beinhaltet, um zwischen vegetativ und nicht-vegetativ zu unterscheiden Bereiche.

Also gehe ich von einem Multiband-Raster-Objekt zu zwei, die nebeneinander verarbeitet und anschließend kombiniert werden, um die ursprüngliche, topographisch korrigierte Szene zu bilden.

Das Problem ist, dass ich einige NA-Werte in der Szene habe, also wenn ich sie nicht los werde, bekomme ich ein leeres Bild als Ergebnis. Ich benutze die sum Funktion mit der Option na.rm definiert als TRUE, um dieses Problem anzugehen, aber das Ergebnis ist ein Ein-Band-Raster-Objekt anstelle eines Stack-oder Rasterbrick-Objekt.

Dies ist mein 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) 

Gibt es eine Funktion, die Multiband-Rastern mit einem na.rm Tag hinzufügen können? Ich habe die do.call Funktion gefunden, aber es scheint mir, dass es mehrere Raster zu einem Einbandobjekt zusammenfasst.

Danke für Ihre Hilfe.

ps: Hantson, S. Chuvieco, E. 2010. Auswertung verschiedener topographischer Korrekturmethoden für Landsat-Bilder. Internationale Zeitschrift für Angewandte Erdbeobachtung und Geoinformation 13 (2011): 691-700p.

Answer 1

Bitte geben Sie ein reproduzierbares Beispiel

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

Jetzt können Sie

r1 <- s1 + s2 

tun Aber das nicht kümmert sie von NA Werten. Also müssen Sie etwas greifen, wie

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

oder diese

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

PS. Ich verstehe, dass sum(s1, s2) erwartet ein mehrschichtiges Objekt zurückgeben, aber ich denke, Rückkehr eine Schicht mit base::sum konsistent ist:

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