rterra
On a categorical SpatRast object, how do I efficiently return a count of cells per factor level? Equivalent to e.g. GRASS r.stats with flag 'c'
summary(x, size = ncell(x))
returns a table object with the first few categories, but then lumps the rest into 'other'.
zonal(x, x, length)
works on small example rasters but has memory limits that mean it won't work on large ones. Is there another approach I'm missing?
It is best to avoid using x[] with SpatRasters as that forces all data into memory. That works fine on these examples, but may fail with your real data. In this case, you can use zonal in stead.
You example data:
library(terra)
f = system.file("ex/elev.tif", package="terra")
r = rast(f)
x = classify(r, cbind(-Inf, 400, NA))
y = patches(x)
To count the number of cells by patch
f = freq(y)
head(f)
# layer value count
#[1,] 1 1 866
#[2,] 1 2 3
#[3,] 1 3 1
#[4,] 1 4 7
#[5,] 1 5 1
#[6,] 1 6 1
If you want the area covered by these cells, first use cellSize
z = cellSize(y,unit="km") |> zonal(y, sum)
head(z) |> round(2)
# patches area
#1 1 479.33
#2 2 1.66
#3 3 0.55
#4 4 3.88
#5 5 0.55
#6 6 0.55
The filtered variable is a table of patch numbers and areas. The names of the table are the patch numbers.
To plot a map of TRUE/FALSE for those patches, do a boolean (logical) operator to see which values of y are in that set:
plot(y %in% as.numeric(names(filtered)))
or if you want a raster with those patches and their patch values (first I've made a copy so I don't stomp on y), set all that arent in the set to NA:
yf = y
yf[!yf %in% as.numeric(names(filtered))] = NA
plot(yf)
Best Answer
If you call
help(package="terra")you can bring up the help documents and read the function summaries. The function you are looking for isterra::freq