[R-sig-Geo] calculate selected area from stacked raster object

Francesc Montserrat f.montserrat at gmail.com
Wed Feb 22 17:45:33 CET 2017


Dear List,

As I have only recently started using R for spatial analysis, and I am 
not a geographer or spatial data specialist by any means, I have a -what 
I presume to be- relatively simple question. I am trying to calculate 
the area of part of a stacked raster object that meets certain 
conditions. More specifically, from a dataset from the deep sea in the 
south Atlantic, I have stacked two raster objects (depth and slope) that 
are further identical in coordinate system (WGS84), x-y (Lat-Long) 
position and projection. From the stacked raster object, I would like to 
extract the part that sits between (say) 1000 and 4000 m depth, with a 
slope of more than 10 degrees. I would like to know what the areal 
extent is in square km and I would like to add it to a previously 
plotted map. Below is a reproducible example:

# Raster object containing depth values
dpt <-  raster(ncol=623, nrow=815, xmx=-31.72083, xmn=-38.50417, 
ymn=-33.8875, ymx=-28.70417)
values(dpt) <- sample(-200:-5000, size=(nrow(dpt)*ncol(dpt)), replace=T)

# Raster object containing slope values
slp <- raster(ncol=623, nrow=815, xmx=-31.72083, xmn=-38.50417, 
ymn=-33.8875, ymx=-28.70417)
values(slp) <- sample(0:30, size=(nrow(slp)*ncol(slp)), replace=T)

# Stack raster objects
stk <- stack(dpt,slp)


# Colour palette
colrs <- 
colorRampPalette(c("navyblue","dodgerblue3","cyan2","green2","darkgoldenrod1"))
# Plot raster map; does not look like ocean floor because of "sample"
plot(dpt, xlab="Longitude", ylab="Latitude", col=colrs(100), font.lab=2, 
cex.lab=1.5, las=1)


# Create a blank copy of previous raster plot
selectAtt <- raster(dpt)
# Fill in cells where Attribute(s) meet(s) conditions
selectAtt[stk$layer.1 <= -1000 & stk$layer.1 >= -4000 & stk$layer.2 >= 
10] <- 90
# Set object projection
projection(selectAtt) <- CRS("+proj=longlat +ellps=WGS84")
# Plot selection in previous raster
plot(selectAtt, col="red", add=TRUE, legend=F, proj4string=crswgs84)

#####>>> What is the area (within total area) with both depth (layer.1) 
and slope (layer.2) meeting given conditions ?

a <- stk[stk$layer.1 <= -1000 & stk$layer.1 >= -4000 & stk$layer.2 >= 10]
area(a, na.rm=T)

# This gives me the error message:
Error in (function (classes, fdef, mtable)  :
             unable to find an inherited method for function ‘area’ for 
signature ‘"matrix"’

I have tried to find what this actually means, but I cannot find a way 
to solve this. What am I missing here ? Any help is greatly appreciated !



More information about the R-sig-Geo mailing list