Rasters and H3 in R

# In reference to https://fosstodon.org/@TheRealJimShady/109559871547727070
# Caution! Much like h3jsr itself, just because you can do this doesn't mean you should :P

# Nonetheless, a demonstration:

library(terra)
library(sf)
library(dplyr)
library(h3jsr)

# from ?rast
f <- system.file("ex/elev.tif", package="terra")
r <- rast(f)

# get cell center coords
cds <- xyFromCell(r, seq(ncell(r)))
pts <- st_as_sf(as.data.frame(cds), coords = c('x','y'), crs = 4326)
pts$ID <- seq(nrow(pts))

# add h3 address at a few resolutions
h3tab <- h3jsr::point_to_cell(pts, res = 4:8, simple = FALSE)
pts_h3 <- dplyr::inner_join(pts, h3tab) |> 
  mutate(across(matches('h3'), as.factor))

# lets try plotting as points
plot(pts_h3['h3_resolution_6'], pch = 20)

# can't apply these to a raster as a RAT as they're linked to cell ID/loc
# rather than cell value
# However, can create a new categorical raster like so:

rhex <- rast(r) # empty grid
values(rhex) <- pts_h3$h3_resolution_6

# trim back to extent of r with overlay
rhex <- mask(rhex, r)

plot(rhex)

# oh, and: index raster data on to the H3 table
pts_h3$elevation <- values(r)

plot(pts_h3['elevation'], pch = 20)

# mean elevation at res 6 is then
pts_h3 <- pts_h3 |> 
  group_by(h3_resolution_6) |> 
  mutate(elevation_mean_6 = mean(elevation)) |> 
  ungroup()

plot(pts_h3['elevation_mean_6'], pch = 20)

# \o/