# Calculate raster cell area as function of latitude?

Raster used available here: https://www.dropbox.com/s/xn7hdll2op5zcc9/MAP_global.tif?dl=0.

INTRODUCTION: The resolution of a raster is commonly related in degrees, e.g.:

``````   library(raster)
r <- raster("/CALL/IN/RASTER/FILE")
print(r)

class      : RasterLayer
dimensions : 720, 1440, 1036800  (nrow, ncol, ncell)
*resolution : 0.25, 0.25  (x, y)*
extent     : -180, 180, -90, 90  (xmin, xmax, ymin, ymax)
crs        : +proj=longlat +datum=WGS84 +no_defs
source     : MAP_global.tif
names      : MAP_global
values     : 2.561833, 10341.26  (min, max)
``````

At the equator, x = 0.25 degrees is ~27.75km, so, a 0.25 x 0.25 cell would have a horizontal area of ~770.06km^2. While the resolution (in units of degree) will remain constant as a function of latitude, the grid cell's area will shrink with distance from the equator.

QUESTION:

Can the area of a 0.25 x 0.25 degree grid cell be calculated as a function of latitude, and then be stored as a variable in that raster? If so, could someone please provide that script? Or provide me the information needed to proceed with coding it myself?

It's the raster packages "area" function!

Make a sample lat-long raster (1 degree here but any lat-long raster should work):

``````> r = raster()
> r
class      : RasterLayer
dimensions : 180, 360, 64800  (nrow, ncol, ncell)
resolution : 1, 1  (x, y)
extent     : -180, 180, -90, 90  (xmin, xmax, ymin, ymax)
crs        : +proj=longlat +datum=WGS84 +no_defs
``````

and `area` gives us:

``````> area(r)
class      : RasterLayer
dimensions : 180, 360, 64800  (nrow, ncol, ncell)
resolution : 1, 1  (x, y)
extent     : -180, 180, -90, 90  (xmin, xmax, ymin, ymax)
crs        : +proj=longlat +datum=WGS84 +no_defs
source     : memory
names      : layer
values     : 107.7756, 12308.62  (min, max)
``````

a raster where each cell is the area in km^2, approximately. See the help for details of the approximation.

Just a little remark: `{terra}` seems to handle things a little bit differently in comparison to `{raster}` since `area()` was removed from the package. `expanse()` is now what you are looking for:

``````library(terra)
#> terra 1.6.33

# init raster
r <- rast(resolution = 0.25)

# insert some unique values to remove NAs
values(r) <- 1:ncell(r)

r
#> class       : SpatRaster
#> dimensions  : 720, 1440, 1  (nrow, ncol, nlyr)
#> resolution  : 0.25, 0.25  (x, y)
#> extent      : -180, 180, -90, 90  (xmin, xmax, ymin, ymax)
#> coord. ref. : lon/lat WGS 84
#> source(s)   : memory
#> name        :   lyr.1
#> min value   :       1
#> max value   : 1036800

# area per grid cell
res <- expanse(r, unit = "km", byValue = TRUE)

summary(res[, 3])
#>    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.
#>   1.701 298.578 547.647 491.961 711.836 769.316
``````
• Thank you for the solution employing terra. It is an interesting alternative to the one suggested by SpacedMan. Unfortunately, I can only choose one solution as the 'answer', else I would flag both your's and SpacedMan's solutions. Because SpacedMan provided a solution first, I'm flagging that solution the answer. That said, I do appreciate your time and have certainly learned something - thanks! Commented Oct 30, 2022 at 14:49
• No worries, just wanted to provide some additional info to the existing answer. Commented Oct 30, 2022 at 16:00