---
title: "Hillshade Example"
knit: litedown:::knit
vignette: >
  %\VignetteIndexEntry{Hillshade Example}
  %\VignetteEngine{litedown::vignette}
  %\VignetteEncoding{UTF-8}
---

```{r, include = FALSE}
EVAL <- isTRUE(as.logical(Sys.getenv("R_RGEEDIM_RUN_EXAMPLES"))) &&
  requireNamespace("terra", quietly = TRUE) &&
  rgeedim::gd_is_initialized(project = Sys.getenv("GOOGLE_CLOUD_QUOTA_PROJECT", "rgeedim-demo"))

litedown::reactor(
  eval = EVAL,
  collapse = TRUE,
  fig.width = 8,
  fig.align = 'center'
)
```

```{r, eval=EVAL}
library(rgeedim)
library(terra)

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

## hillshade example
b <- gd_bbox(
  xmin = -120.296,
  xmax = -120.227,
  ymin = 37.9824,
  ymax = 38.0071
)

# download 1m DEM in AEA
# need to set resampling method in composite step
x <- "USGS/3DEP/1m" |>
  gd_collection_from_name() |>
  gd_search(region = b) |>
  gd_composite(
    resampling = "bilinear"
  ) |> 
  gd_download(
    region = b,
    bands = list("elevation"),
    crs = "EPSG:5070",
    scale = 10,
    filename = "dem.tif"
  ) |>
  rast()

# calculate slope, aspect, and hillshade with terra
slp <- terrain(x, "slope", unit = "radians")
asp <- terrain(x, "aspect", unit = "radians")
hs <- shade(slp, asp)

# compare elevation v.s. hillshade
plot(c(x, hillshade = hs))
```

```{r, include=FALSE, eval=EVAL}
unlink("dem.tif")
```