feat: hillshade

NAMESPACE

@@ -37,6 +37,7 @@ export(wbw_is_int)
 export(wbw_is_rgb)
 export(wbw_max_procs)
 export(wbw_mean_filter)
+export(wbw_multidirectional_hillshade)
 export(wbw_random_sample)
 export(wbw_read_raster)
 export(wbw_read_vector)

R/geomorphometry.R

@@ -213,7 +213,7 @@ S7::method(wbw_ruggedness_index, WhiteboxRaster) <-
 #' @return [WhiteboxRaster] object
 #'
 #' @eval rd_wbw_link("fill_missing_data")
-#' @eval rd_example("fill_missing_data")
+#' @eval rd_example("wbw_fill_missing_data")
 #'
 #' @export
 wbw_fill_missing_data <-

R/hillshade.R

@@ -0,0 +1,82 @@
+#' Multidirectional Hillshade
+#' @keywords geomorphometry
+#'
+#' @description
+#' This tool performs a hillshade operation (also called shaded relief) on 
+#' an input digital elevation model (DEM) with multiple sources of 
+#' illumination.
+#'
+#' @details
+#' The hillshade value (HS) of a DEM grid cell is calculate as:
+#' \deqn{HS = \frac{\tan(s)}{\sqrt{1 - \tan(s)^2}} \times [\frac{\sin(Alt)}{\tan(s)} -
+#'   \cos(Alt) \times \sin(Az - a)]}
+#' where \eqn{s} and \eqn{a} are the local slope gradient and aspect (orientation)
+#' respectively and \eqn{Alt} and \eqn{Az} are the illumination source altitude and azimuth
+#' respectively. Slope and aspect are calculated using Horn's (1981) 3rd-order
+#' finate difference method.
+#'
+#' @eval rd_input_raster("dem")
+#' @param altitude \code{double}, the altitude of the illumination sources.
+#' i.e. the elevation of the sun above the horizon, measured as an angle from
+#' 0 to 90 degrees
+#' @param z_factor \code{double}, Z conversion factor is only important
+#' when the vertical and horizontal units are not the same in the DEM.
+#' When this is the case, the algorithm will multiply each elevation in the
+#' DEM by the Z conversion factor
+#' @param full_360_mode \code{boolean}, default \code{FALSE}. I.e. whether or
+#' not to use full 360-degrees of illumination sources. When \code{FALSE}
+#' (default) the tool will perform a weighted summation of the hillshade
+#' images from four illumination azimuth positions at 225, 270, 315, and
+#' 360 (0) degrees, given weights of 0.1, 0.4, 0.4, and 0.1 respectively. When
+#' run in the full 360-degree mode, eight illumination source azimuths are
+#' used to calculate the output at 0, 45, 90, 135, 180, 225, 270, and 315
+#' degrees, with weights of 0.15, 0.125, 0.1, 0.05, 0.1, 0.125, 0.15,
+#' and 0.2 respectively.
+#'
+#' @return [WhiteboxRaster] object
+#'
+#' @eval rd_wbw_link("multidirectional_hillshade")
+#' @references 
+#' Horn B.K.P., 1981, Hill shading and the reflectance map, Proceedings of 
+#' the I.E.E.E. 69, 14
+#' 
+#' @eval rd_example("wbw_multidirectional_hillshade")
+#'
+#' @export
+wbw_multidirectional_hillshade <-
+  S7::new_generic(
+    name = "wbw_multidirectional_hillshade",
+    dispatch_args = "dem",
+    fun = function(dem,
+                   altitude = 30,
+                   z_factor = 1,
+                   full_360_mode = FALSE) {
+      S7::S7_dispatch()
+    }
+  )
+
+S7::method(wbw_multidirectional_hillshade, WhiteboxRaster) <-
+  function(dem,
+           altitude = 30,
+           z_factor = 1,
+           full_360_mode = FALSE) {
+    # Checks
+    check_env(wbe)
+    checkmate::assert_double(altitude, len = 1, lower = 0, upper = 90)
+    checkmate::assert_double(z_factor, len = 1)
+    checkmate::assert_logical(full_360_mode, len = 1)
+
+    # WBT
+    out <- wbe$multidirectional_hillshade(
+      dem = dem@source,
+      altitude = altitude,
+      z_factor = z_factor,
+      full_360_mode = full_360_mode
+    )
+
+    # Return
+    WhiteboxRaster(
+      name = paste0(dem@name, "(Hillshade)"),
+      source = out
+    )
+  }

_pkgdown.yml

@@ -13,15 +13,26 @@ authors:
 template:
   bootstrap: 5
   light-switch: true
+  math-rendering: mathjax
   bslib:
     base_font: {google: "Inter"}
     heading_font: {google: "Recursive"}
     code_font: {google: "Fira Code"}
 
 navbar:
+  components:
+    faq:
+      text: FAQ
+      href: articles/FAQ.html
+    tutorials:
+      text: Tutorials
+      menu:
+      - text: Tutorial 1
+      - text: Tutorial 2
+        href: articles/Tutorial-2.html
   structure:
-    left:  [reference, articles, tutorials, news]
-    right: [search, github, lightswitch]
+    left:  [reference, tutorials, news]
+    right: [faq, search, github, lightswitch]
 
 reference:
 - title: "Input/Output"
@@ -75,4 +86,4 @@ reference:
 - title: "Methods"
   desc: "S3 methods for WhiteboxRaster and WhiteboxVector objects"
   contents:
-  - has_keyword("methods")
+  - has_keyword("methods")

tests/testthat/test-geomorphometry.R

@@ -28,6 +28,13 @@ test_that(
     expect_error(wbw_fill_missing_data(x, weight = "2.5"))
     expect_error(wbw_fill_missing_data(x, exclude_edge_nodata = "YES"))
     expect_error(wbw_fill_missing_data(NULL))
+
+    # wbw_multidirectional_hillshade
+    expect_error(wbw_multidirectional_hillshade(mtcars))
+    expect_error(wbw_multidirectional_hillshade(x, altitude = 100))
+    expect_error(wbw_multidirectional_hillshade(x, z_factor = "2.5"))
+    expect_error(wbw_multidirectional_hillshade(x, full_360_mode = "YES"))
+    expect_error(wbw_multidirectional_hillshade(NULL))
   }
 )
 
@@ -63,12 +70,18 @@ test_that(
       wbw_ruggedness_index(x),
       WhiteboxRaster
     )
-    
+
     # wbw_fill_missing_data
     expect_s7_class(
       wbw_fill_missing_data(x),
       WhiteboxRaster
     )
+
+    # wbw_multidirectional_hillshade
+    expect_s7_class(
+      wbw_multidirectional_hillshade(x),
+      WhiteboxRaster
+    )
   }
 )
 

vignettes/articles/FAQ.Rmd

@@ -0,0 +1,41 @@
+---
+title: "Frequently Asked Questions"
+output: 
+  rmarkdown::html_vignette:
+    toc: false
+vignette: >
+  %\VignetteIndexEntry{FAQ}
+  %\VignetteEngine{knitr::rmarkdown}
+  %\VignetteEncoding{UTF-8}
+  %\VignetteDepends{terra}
+---
+
+```{r knitr_setup, include = FALSE}
+knitr::opts_chunk$set(
+  collapse = TRUE,
+  comment = "#>",
+  warning = FALSE,
+  message = FALSE,
+  dpi = 300,
+  retina = 2
+)
+
+requireNamespace("terra", quietly = TRUE)
+```
+
+# How can I plot `WhiteboxRaster` objects?
+
+While there exists a built-in `plot()` method for all `WhiteboxRaster` and `WhiteboxVector` objects, it is highly advisable to use `{terra}`'s plotting functionality by converting the `Whitebox...` objects to `SpatRaster` or `SpatVector` accordingly as follows:
+
+```{r terra_plot, message=FALSE, warning=FALSE}
+library(wbw)
+library(terra)
+
+raster_path <- system.file("extdata/dem.tif", package = "wbw")
+
+wbw_read_raster(raster_path) |>
+  wbw_multidirectional_hillshade() |> 
+  as_rast() |> 
+  plot()
+
+```

vignettes/articles/Tutorial-2.Rmd

@@ -14,7 +14,7 @@ vignette: >
 knitr::opts_chunk$set(
   collapse = TRUE,
   comment = "#>",
-  warn = FALSE,
+  warning = FALSE,
   message = FALSE,
   dpi = 300,
   retina = 2
@@ -23,7 +23,7 @@ knitr::opts_chunk$set(
 requireNamespace("terra", quietly = TRUE)
 ```
 
-> NB! This is R representation of the of the [Whitebox Workflows for Python (WbW) Tutorial 2: Geomorphometric Analysis](https://github.com/jblindsay/jblindsay.github.io/blob/master/WhiteboxTutorials/WbW_tutorials/WbW_tutorial2.ipynb)
+> NB! This tutorial is based on the [Whitebox Workflows for Python (WbW) Tutorial 2: Geomorphometric Analysis](https://github.com/jblindsay/jblindsay.github.io/blob/master/WhiteboxTutorials/WbW_tutorials/WbW_tutorial2.ipynb)
 
 # 1. Download sample data
 
@@ -64,4 +64,62 @@ dem_filled |>
   plot()
 ```
 
-To be continued ...
+# 3. Vizualization
+
+To properly vizualize the DEM, let's create a hillshade map of the relief:
+
+```{r hillshade}
+dem_hillshade <- wbw_multidirectional_hillshade(dem_filled)
+
+dem_hillshade |> 
+  as_rast() |> 
+  plot(col = grDevices::gray.colors(10))
+```
+
+# 4. Raster summary stats
+
+One can get `WhiteboxRaster` summary statistics by running familiar functions from base R (and `{terra}`):
+
+```{r summary_stats}
+# Raw DEM
+summary(dem)
+wbw_cols(dem) # Number of columns
+wbw_rows(dem) # Number of rows
+num_cells(dem) # Number of cells
+wbw_res(dem) # Raster resolution
+
+# NA-Filled DEM
+summary(dem_filled)
+wbw_cols(dem_filled) # Number of columns
+wbw_rows(dem_filled) # Number of rows
+num_cells(dem_filled) # Number of cells
+wbw_res(dem_filled) # Raster resolution
+```
+
+# 5. Extracting land-surface parameters (LSPs)
+Now let's extract some common land-surface parameters (LSPs), the basic building blocks of a geomorphometric analysis.
+
+Slope and aspect are two of the most common LSPs.
+
+```{r slope_and_aspect}
+# Slope
+slope <- wbw_slope(dem_filled)
+aspect <- wbw_aspect(dem_filled)
+
+# Set up a 1x2 plotting layout
+par(mfrow = c(1, 2))
+
+# Plot the first raster (Slope)
+plot(as_rast(slope),
+     main = "Slope",
+     legend = FALSE,
+     col = terrain.colors(100),
+     mar=c(1,1,1,1))
+
+# Plot the second raster (Aspect)
+plot(as_rast(aspect),
+     main = "Aspect",
+     legend = FALSE,
+     mar=c(1,1,1,1))
+
+```