Conditional command from Raster Calculator transfe

So far I used the following statement in the Raster Calculator of ArcGIS:

Con(("Land_use.rst" ==  -20), "Export.rst")

This calculates a new Raster which only contains the Data from Export where Land_use equals -20. That is exactly what I want. But I want to automatise this in R, as I have to do it a lot of times.

So far I got something like this:

for (catch_dir in Dir_List) {

r1 <- raster(paste0(catch_dir, '/Export.rst'))

r2 <- raster(paste0(catch_dir,'/LAND_use.rst'))

### statement that output export_streams.rst contains the data of export.rst where LAND_use.rst equals -20.
  x <- overlay(r2, r1, fun=function(x,y){ y[x!=-20] <- 0; y})

 writeRaster(x, filename = paste0(catch_dir, "/export_streams.rst"), format="IDRISI", NAflag = 0, datatype = "FLT4S",
              overwrite=TRUE)
}

That code does the following: It produces a raster which contains the values of r1 where r2 = -20. So that is good, but there are also leftovers at the border of the raster which do not equal -20 in r2. The extent of the two rasters are the same though. So the overlay probably doesnt work for my task. Any other ideas than overlay? Maybe some if(r2 == -20){ } command?

回答1:

I don't habe your data, so I cant reproduce your error. However this is how I would do it with dummy data:

#create dummy matrix
a <- matrix(c(10,10,10,10,10,10,10,10,10),nrow=3)
b <- matrix(c(-20,-20,-20,3,3,3,4,4,4),nrow=3)

#create raster from matrix
r1 <- raster(a)
r2 <- raster(b)

#Create raster with same extent and nrows/ncols as original rasters
x <- raster(r1)

#produce raster x which contains the values of r1 where r2 = -20
x[r2[]==-20] <- r1[r2[]==-20]

回答2:

You can do this:

m <- subs(r2, cbind(-20, 1))
x <- mask(r1, m)