In RobinKohrs/soilmoistr:

knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)

Load the necessary Packages

library(soilmoistr)
library(here)
library(iffitoR)
library(tidyverse)
library(raster)
library(stars)
library(sf)

What does the function `get_sm_data` do?

The idea of this function is to extract the soilmoisture data from a large stack of reprojected, cropped, resampled and compressed GeoTIFFs.
Only Landslides with a known date are considered
The date of the landslide, plus some days ( as specified in days_before_window) and minus (days_after_window) are then compared to the available dates of the soilmoisture-"acquisitions"
The spatial object for which to extract the values can either be
points
points with a bufffer
polygons

When using points and no buffer

When we are using only points (e.g. landslide initiation points) one only needs to call for example
The extraction here is done using the st_extract-method for stars objects which extracts values at point locations for stars-objects

res = get_sm_data(<path_to_soilmoisture_data>,
                  landsld = <landslide-object>)

This will use the default values for some of the parameters:

  +  days_before_window = 5,
  +  days_after_window = 0,
  +  point_buffer = NULL,
  +  aggre_fun = NULL

This mean that we look 5 days before the landslide (or any other point like object) if we find soilmoisture information

When using points and a buffer arround them

When you want to use a buffer around your point you can specify this in the point_buffer-parameter
You directly pass the distance in units of your landsld-object. As we are mainly working in EPSG: 32632 this should be meters
As we are now working with a polygon, we can pass a vector with aggregation functions to the aggre_fun-parameter. If nothing is passed you will get back the values for all intersection cells. Look into the documentation of exactextracr to know how this works in detail.

res = get_sm_data(<path_to_soilmoisture_data>,
                  landsld = <landslide-object>,
                  point_buffer = 300,
                  aggre_fun = c("mean"))

Example

lets load the landslide data from the iffitoR-package and subset it to only some movements of interest

data = landsld[grepl("translational|rotational", landsld$second_level), ]

check its dimensions

dim(data)

lets create a copy and remove the column called date in order to see what happens

copy = data.frame(data) %>% st_as_sf()
copy[["date"]] = NULL # remove the date column

Now lets call the get_sm_data-function on that dataframe

res = get_sm_data(
  landsld = copy,
  path_sm = path_sm
)

It does not work as the dataframe really needs to have a column called date
Now lets use the original data to get the soilmoisture data

path_sm = "/mnt/CEPH_PROJECTS/Proslide/soilmoisture/32632/"

# data = data[1:50, ]

res = get_sm_data(
  landsld = data,
  path_sm = path_sm,
  quiet = F
)

Work with the extracted data

So now we extracted the values of the soilmoisture rasters at pixel location of each slide.
We considered the date of the slide as well as 5 anterior days
The result is stored as a list in the dataframe that the function get_sm_datareturns. The list-column is called sm_values.
And its all still a spatial dataframe of class sf

res %>% dplyr::select(c(date, PIFF_ID, second_level, sm_values)) %>% glimpse()

Lets check for how much landslide data we also have soilmoisture data
Remember: That even if we have intersecting days, it still is very probable that the lanslide-point overlays with a NA-value in the soilmoisture raster.

# add fonts
library(showtext)
library(ggtext)

# where does the package look for fonts?
# font_paths()
# show local fonts
# font_files()
# add the font by providing the family and the name of the truetypefont (.ttf)
font_add(family = "DejaVu Serif", regular = "DejaVuSerif-Bold.ttf")

# create color scale for true and false
pal = c("TRUE" = "#b9cb99", "FALSE" = "#c23b22")

res %>% 
  mutate(
    intersection = case_when(
      is.na(sm_values) ~ FALSE,
      TRUE ~ TRUE
    )
  ) %>% 
  count(intersection) %>% 
  ggplot() +
  geom_col(
    aes(
      x = n,
      y = intersection,
      fill = intersection
    )
  ) +
  scale_fill_manual(values = pal) +
  theme_light() + 
  labs(
    y = "",
    x = "#",
    title = "Number of Soilmoisture rasters intersecting in time with the landslide data",
    caption = "Considered were only the translational and rotational slides.
               A temoral buffer of 5 days was applied to the slides."
  ) +
  theme(
    legend.position = c(0.8,0.8),
    title = element_text(family = "mono"),
    legend.text = element_text(family = "mono"),
    legend.background = element_blank(),
    axis.text.y = element_blank()
  )