Earth Engine Export Example

library(rgeedim)
library(terra)

project_id <- Sys.getenv("GOOGLE_CLOUD_QUOTA_PROJECT", "rgeedim-demo")
gd_initialize(project = project_id)

ee <- earthengine()

# use the earth engine API and reticulate to access
# World Database of Protected Areas <https://developers.google.com/earth-engine/datasets/catalog/WCMC_WDPA_current_polygons>
# filter to a specific area (yellowstone)
b <- ee$FeatureCollection("WCMC/WDPA/current/polygons")$filter(
  ee$Filter$And(
    ee$Filter$eq("NAME", "Yellowstone"),
    ee$Filter$eq("DESIG_ENG", "National Park")
  )
)
g <- gd_bbox(b)

# create a terra SpatVector from the GeoJSON-like list from gd_region()
yellowstone <- gd_region(b) |>
  gd_region_to_vect() |>
  project("EPSG:5070")

# inspect
plot(yellowstone)

# create an Earth Engine asset you can use again in the future
# 250m resolution for Yellowstone extent
"USGS/SRTMGL1_003" |>
  gd_image_from_id() |>
  gd_export(
    region = g,
    scale = 250,
    bands = list("elevation"),
    resampling = "bilinear",
    type = "asset",
    filename = "Yellowstone",
    folder = project_id,
    overwrite = TRUE,
    crs = "EPSG:5070"
  )

# create a local GeoTIFF from the above created asset
x_path <- sprintf("projects/%s/assets/Yellowstone", project_id) |>
  gd_image_from_id() |>
  gd_download(
    bands = list("elevation"),
    resampling = "bilinear",
    region = g,
    scale = 250,
    crs = "EPSG:5070"
  )

if (!inherits(x_path, "try-error")) {
  x <- rast(x_path)
  plot(x)

  # calculate hillshade locally
  slp <- terrain(x, unit = "radians")
  asp <- terrain(x, "aspect", unit = "radians")
  shd <- shade(slp, asp)

  # plot with boundary
  plot(shd)
  plot(yellowstone, add = TRUE)
  
  # clean up
  unlink(sources(x))
}