# INBO CODING CLUB

14 December 2023

Welcome!

## Share your code snippet

If you want to share your code snippet, copy paste your snippet within a section of three backticks (```):

As an **example**:

library(tidyverse)

(*you can copy paste this example and add your code further down*)

## Yellow sticky notes

No yellow sticky notes online. Put your name + " | " and add a "*" each time you solve a challenge (see below).

## Participants

Name | Challenges
--- | ---
Damiano Oldoni | ***
Pieter Huybrechts |**
Nele Mullens | **
Britt Lonneville | *
Ward Langeraert | *
Sophia Woolcock |
Dirk Maes |
Marijke Thoonen |
Cécile Herr |**
Wouter Depaepe |*
Karen Wuyts |**
Raïsa |*
Amber | *
Ward Standaert | **

## Challenge 1

### Ward

# 1.
qgis_get_argument_specs(algorithm = "native:buffer")
qgis_show_help(algorithm = "native:buffer")

result2 <- qgis_run_algorithm(
  "native:buffer",
  INPUT = input,
  DISTANCE = 1,
  DISSOLVE = FALSE
)
output_sf2 <- sf::st_as_sf(result2)

mapview(output_sf2)

#> Buffer per rij en bij DISOLVE = TRUE gaat hij eerst alles joinen tot 1
#> geometry

# 2.
qgis_algorithms() %>%
  filter(provider == "native",
         grepl("shortest", algorithm)) %>%
  select(provider, algorithm, algorithm_title)

qgis_get_argument_specs(algorithm = "native:shortestline") %>%
  select(name, description)
qgis_show_help(algorithm = "native:shortestline")

springs_to_water <- qgis_run_algorithm(
  "native:shortestline",
  SOURCE = hab_springs,
  DESTINATION = water
)
springs_to_water_sf <- sf::st_as_sf(springs_to_water)
mapview(springs_to_water_sf)

# 3.
water_to_springs <- qgis_run_algorithm(
  "native:shortestline",
  SOURCE = water,
  DESTINATION = hab_springs
)
water_to_springs_sf <- sf::st_as_sf(water_to_springs)
mapview(water_to_springs_sf)

# 4.
qgis_run_algorithm(
  "native:shortestline",
  SOURCE = water,
  DESTINATION = hab_springs
  ) |>
  st_as_sf() |>
  mapview()

### Britt

# Calculate the ‘shortest’ path (lines) between the hab_springs points and the nearest watersurface water.

algorithms_dist <- qgis_search_algorithms(algorithm = "distance")
qgis_show_help(algorithm = "qgis:distancetonearesthublinetohub")

result_12 <- qgis_run_algorithm(
  "qgis:distancetonearesthublinetohub",
  INPUT = hab_springs,
  HUBS = water,
  FIELD = "polygon_id",
  UNIT = 3,
  OUTPUT = "data/result_challenge12.geojson"
)
output_sf_12 <- sf::st_as_sf(result_12)
mapview(output_sf_12)

# Calculate the ‘shortest’ paths (lines) between the polygons water and the hab_springs points.

result_13 <- qgis_run_algorithm(
  "qgis:distancetonearesthublinetohub",
  INPUT = water,
  HUBS = hab_springs,
  FIELD = "id",
  UNIT = 3,
  OUTPUT = "data/result_challenge13.geojson"
)
output_sf_13 <- sf::st_as_sf(result_13)
mapview(output_sf_13)

# We are all pipe (%>%)* lovers, aren't we? Can you use qgis_run_algorithm() in a pipe? No. How to combine qgisprocess and pipes? Hint: check the cheat sheet.

sf::st_read("data/20231214_watersurfaces_hab.gpkg", layer = "watersurfaces_hab_polygons") %>% 
  qgis_run_algorithm_p(
    "qgis:distancetonearesthublinetohub",
    HUBS = hab_springs,
    FIELD = "id",
    UNIT = 3,
    OUTPUT = "data/result_challenge13.geojson"
    ) %>%
  st_as_sf() %>%
  mapview()

### Nele

# 1.1 Dissolve = FALSE? Show algorithm used. 
input <- sf::read_sf(system.file("shape/nc.shp", package = "sf"))
result <- qgis_run_algorithm(
  "native:buffer",
  INPUT = input,
  DISTANCE = 1,
  DISSOLVE = FALSE
)
result

output_sf <- sf::st_as_sf(result)
plot(sf::st_geometry(output_sf))

qgis_show_help(algorithm="native:buffer")

qgis_get_argument_specs(algorithm=
  "native:buffer")

#1.2
qgis_search_algorithms(algorithm = "Shortest")
#native:shortestline
qgis_show_help(algorithm="native:shortestline")


result <- qgis_run_algorithm(
  "native:shortestline",
  SOURCE = hab_springs,
  DESTINATION = water, 
  )

output_sf <- sf::st_as_sf(result)
plot(sf::st_geometry(output_sf))

mapview(output_sf)

# input needed;
# hab_springs
# water

#1.3
qgis_search_algorithms(algorithm = "Shortest")
#native:shortestline
qgis_show_help(algorithm="native:shortestline")


result <- qgis_run_algorithm(algorithm = 
  "native:shortestline",
  SOURCE = water,
  DESTINATION = hab_springs, 
)

output_sf <- sf::st_as_sf(result)
plot(sf::st_geometry(output_sf))

mapview(output_sf)

# 1.4 
hab_springs |> qgis_run_algorithm_p(algorithm = 
  "native:shortestline",
  DESTINATION = water, 
) |> st_as_sf() |> mapview()

### Raïsa

#1. In challenge 0 we specified (a little blindly) some arguments, e.g. DISSOLVE = TRUE. What does this argument do? What happens if you set DISSOLVE = FALSE? How can you get a description of the chosen algorithm? How can you retrieve information about all the arguments allowed by the algorithm? Or about the generated output? Hint: use the cheat sheet!
qgis_show_help(algorithm = "native:buffer")
qgis_get_argument_specs(algorithm = "native:buffer") %>%
  filter(name == "DISSOLVE") %>% 
  dplyr::pull(description)

#2. Calculate the ‘shortest’ path (lines) between the hab_springs points and the nearest watersurface water.
qgis_algorithms() %>%
  filter(str_detect(algorithm, "shortest")) %>%
  pull(algorithm)

qgis_show_help(algorithm = "native:shortestline")

hab_springs %>%
  qgis_run_algorithm_p(
    algorithm = "native:shortestline",
    DESTINATION = water
) %>%
  st_as_sf() %>%
  mapview() +
  mapview(hab_springs)

#3. Calculate the ‘shortest’ paths (lines) between the polygons water and the hab_springs points.
water %>%
  qgis_run_algorithm_p(
    algorithm = "native:shortestline",
    DESTINATION = hab_springs
  ) %>%
  st_as_sf() %>%
  mapview() +
  mapview(hab_springs)

### Pieter

# dissolve is false: many vectors
result <- qgis_run_algorithm(
  "native:buffer",
  INPUT = input,
  DISTANCE = 1,
  DISSOLVE = FALSE
)
output_sf <- sf::st_as_sf(result)
plot(sf::st_geometry(output_sf))
mapview(output_sf)
# get information about algorithm
qgis_show_help(algorithm = "native:buffer")
# information about the parameters
qgis_get_argument_specs(algorithm = "native:buffer")
## 1.2 Calculate the ‘shortest’ path (lines) between the hab_springs points and
#the nearest watersurface water.
qgis_search_algorithms(
  algorithm = "shortest",
  provider = "native"
)
qgis_show_help(algorithm = "native:shortestline")
qgis_get_argument_specs(algorithm = "native:shortestline")

result <- qgis_run_algorithm(
  "native:shortestline",
  SOURCE = hab_springs,
  DESTINATION = water,
  METHOD = 0
)
output_sf <- sf::st_as_sf(result)
mapview(output_sf)

## 1.3 
result <- qgis_run_algorithm(
  "native:shortestline",
  SOURCE = water,
  DESTINATION = hab_springs,
  METHOD = 0
)
output_sf <- sf::st_as_sf(result)
mapview(output_sf)

## 1.4
hab_springs %>% 
qgis_run_algorithm_p(
  algorithm = "native:shortestline",
  DESTINATION = water,
  METHOD = 0,
) |> 
  sf::st_as_sf() |>
  mapview()`

## Challenge 2

### Pieter

# Your version of QGIS might change what algorithms are available! Or even the
# results, so you might want to check the version is still the same as it was
# when you wrote the script
qgis_version_expected <- "3.32.2-Lima"
assertthat::assert_that(qgis_version() == qgis_version_expected,
                        msg = glue::glue("The QGIS version has changed, ",
                                         "it used to be {qgis_version_expected}",
                                         ", but now is {qgis_version_current}",
                                         qgis_version_expected = qgis_version_expected,
                                         qgis_version_current = qgis_version()))

## 2.1 Count the number of occurrences of Ludwigia grandiflora in each of the
#water polygons
qgis_search_algorithms("count", provider = "native")
qgis_show_help("native:countpointsinpolygon")

ludwigia_grandiflora <-
  sf::st_read("data/20231214/20231214_ludwigia_grandiflora.geojson")

result <- qgis_run_algorithm(
  algorithm = "native:countpointsinpolygon",
  POLYGONS = water,
  POINTS = ludwigia_grandiflora
)
# get a tibble out and sort by number of points in polygon, descending
sf::st_as_sf(result) %>%
  dplyr::arrange(-NUMPOINTS)

## 2.2 Count the number of individuals (column individualCount) of Ludwigia
## grandiflora in each of the water polygons
result_ind <- qgis_run_algorithm(
  algorithm = "native:countpointsinpolygon",
  POLYGONS = water,
  POINTS = ludwigia_grandiflora,
  WEIGHT = "individualCount"
)
# get a tibble out and sort by number of individuals, descending
sf::st_as_sf(result_ind) %>%
  dplyr::arrange(-NUMPOINTS)

## 2.3 Create a buffer around hab_springs points with a dynamic radius expressed
## as the square root of area_m2/pi

# See: https://r-spatial.github.io/qgisprocess/articles/qgis_expressions.html#example-applying-a-data-defined-override
result_buff <-
  qgis_run_algorithm(
    "native:buffer",
    INPUT = hab_springs,
    DISTANCE = 'expression: "area_m2" / pi'
  )
# lets have a look at these buffers then:
hab_springs$area_m2/pi

# output to mapview
sf::st_as_sf(result_buff) %>%
  mapview()
`

### Nele

## 2.1 
qgis_search_algorithms(algorithm = "count")
#native:countpointsinpolygon
qgis_show_help(algorithm="native:countpointsinpolygon")

ludwigia_grandiflora <- sf::st_read("./data/20231214/20231214_ludwigia_grandiflora.geojson")


Countpoints <- qgis_run_algorithm(
  "native:countpointsinpolygon",
  POLYGONS = water,
  POINTS = ludwigia_grandiflora
)

st_as_sf(Countpoints)

# 2.2
Countpoints_ind <- qgis_run_algorithm(
  "native:countpointsinpolygon",
  POLYGONS = water,
  POINTS = ludwigia_grandiflora, 
  WEIGHT = "individualCount"
)

st_as_sf(Countpoints_ind) %>% View()

# 2.3 
qgis_search_algorithms(algorithm = "buffer")
#native:buffer
qgis_show_help(algorithm="native:buffer")

buffer_var <- qgis_run_algorithm(
  "native:buffer",
  INPUT = hab_springs,
  DISTANCE = 'expression: "area_m2" / pi'
)

output_buf <- st_as_sf(buffer_var)
mapview(output_buf)

### Ward

# 1.
qgis_search_algorithms(algorithm = "count",
                       provider = "native") %>%
  select(provider, algorithm, algorithm_title)

qgis_get_argument_specs(algorithm = "native:countpointsinpolygon") %>%
  select(name, description)
qgis_show_help(algorithm = "native:countpointsinpolygon")

ludwigia_grandiflora <- st_read(paste0("./data/20231214/",
                                       "20231214_ludwigia_grandiflora.geojson"))

ludwigia_in_water <- qgis_run_algorithm(
  "native:countpointsinpolygon",
  POLYGONS = water,
  POINTS = ludwigia_grandiflora
  ) |>
  st_as_sf()

mapview(ludwigia_in_water)

# 2.
ludwigia_in_water_ind <- qgis_run_algorithm(
  "native:countpointsinpolygon",
  POLYGONS = water,
  POINTS = ludwigia_grandiflora,
  WEIGHT = "individualCount"
  ) |>
  st_as_sf()

mapview(ludwigia_in_water_ind)

waldo::compare(ludwigia_in_water, ludwigia_in_water_ind)

# 3.
qgis_get_argument_specs(algorithm = "native:buffer") %>%
  select(name, description)
qgis_show_help(algorithm = "native:buffer")

# following does not work:
# Error in `processx::run("cmd.exe", c("/c", "call", path, args), ...)`:
# ! System command 'cmd.exe' failed
# QGIS version too old (?)
buffer_springs <- hab_springs |>
  mutate(dynamic_buf = sqrt(area_m2 / pi)) |>
  qgis_run_algorithm_p(
    "native:buffer",
    DISTANCE = "dynamic_buf"
  ) |> 
  st_as_sf()

mapview(buffer_springs)

## Challenge 3

### Pieter

## 3.1 How can you get general statistics about the temperature difference in
## January

qgis_search_algorithms("statistics", provider = "native")
qgis_show_help("native:rasterlayerstatistics")

raster_stats <- qgis_function("native:rasterlayerstatistics")

raster_stats(INPUT = temperature$tdiff_01,
             BAND = 1)



## 3.2 How can you do the same for the other two layers (minimum/maximum
## temperature)?

qgis_run_algorithm(
  "native:rasterlayerstatistics",
  INPUT = temperature$tmin_01,
  BAND = 1
)

qgis_run_algorithm(
  "native:rasterlayerstatistics",
  INPUT = temperature$tmax_01,
  BAND = 1
)

# or in one go
raster::as.list(temperature) %>%
  purrr::map(\(x) raster_stats(INPUT = x, BAND = 1)) %>%
  purrr::set_names(names(temperature))

# or one at a time using band
 
qgis_run_algorithm(
  "native:rasterlayerstatistics",
  INPUT = temperature, 
  BAND = 2
)

qgis_run_algorithm(
  "native:rasterlayerstatistics",
  INPUT = temperature, 
  BAND = 3
)

# or using band and purrr
purrr::map(seq_along(names(temperature)),
           \(band) raster_stats(INPUT = temperature, BAND = band)) %>% 
  purrr::set_names(names(temperature))


## 3.3 Let's do some zonal statistics: calculate the most and the least
## occurring land use category in each of the Natura2000 protected areas.

# NOTE

# Notice that if you want to do the same purely in R, you need to use the R
# package extractr. See the last exercise in slide 14 (Nov 24, 2022) or the last
# exercise of slide 13 (Nov 25, 2021).

qgis_search_algorithms("zonal", provider = "native")
qgis_show_help("native:zonalstatistics")

qgis_run_algorithm(
  algorithm = "native:zonalstatistics",
  INPUT_RASTER = land_use,
  RASTER_BAND = 1,
  INPUT_VECTOR = prot_areas
)