0

How can I calculate the 3D surface area (including topography) for an elevation raster map in R, preferably with package terra? Something like the r.surf.area module of GRASS GIS, or the surfaceArea function of the R sp package, but (if possible) also accounting for latitude/longitude:

library(terra)
library(sp)

r <- terra::rast(system.file("ex/elev.tif", package="terra"))
terra::expanse(r)  # computes 2D surface area, considering the relationship between longlat degrees and meters across the globe

m <- sp::SpatialPixelsDataFrame(points = crds(r), data = data.frame(r))
sp::surfaceArea(m)  # computes 3D area, considering the sloping nature of the surface
AMBarbosa
  • 250
  • 1
  • 8

1 Answers1

1

I suppose you could adjust the area with (multiplier) 1/cos(slope). In that case, for a slope of 45° the adjustment would be sqrt(2)

1/cos(pi/4)
#[1] 1.414214

With that you can do

library(terra)
r <- rast(system.file("ex/elev.tif", package="terra"))

area <- cellSize(r)  
slope <- terrain(r, "slope", unit="radians") 
adjarea <- area / cos(slope)

I see that in GRASS and sp another (perhaps better) approach is used. They create triangles between neighboring cell centers and add up the areas. I could implement that in terra as well.

Robert Hijmans
  • 40,301
  • 4
  • 55
  • 63