Search code examples
rasterterrain

NA handling in raster terrain function

In calculating the slope from an elevation raster using the terrain function of the raster package, there is a border effect where NAs are returned for cells having one or more NA neighbors.

library(raster)

elevation <- getData('alt', country='ITA')
x <- terrain(elevation, 'slope', neighbors = 8)
e <- elevation
e[!is.na(e)] <- 1
e[is.na(e)] <- 2
x[!is.na(x)] <- 1
x[is.na(x)] <- 2
y <- e-x
plot(y)

I'm looking for possible ways (or alternative functions/packages) for overriding this border effect and calculating the slope for all non-NA cells based on the number of available neighbors?

Whether I see this effect appropriate where the border is artificially created due to the extent of the raster (e.g. northern Italy disconnected from Austria, Switzerland...), in other cases the border is legitimate (e.g. coastal cells).

Passing na.rm = TRUE to terrain does not change the result.

Many thanks!

Solution

  • This would be an easy workaround:

    First, download the elevation data unmask (you can mask it later if needed):

    elevation <- getData('alt', country='ITA',mask=F)

    Now you can assume, that all your NA elevation is sea/ocean surface, and therefore has a value of 0.

    So you can set your NAs to 0:

    elevation[is.na(elevation)] <- 0

    This should remove all border issues due to NA values.