From b8a2f0a796297845103574f20d1842ed036b6928 Mon Sep 17 00:00:00 2001
From: hansvancalster <hans.vancalster@inbo.be>
Date: Fri, 24 Jan 2025 16:29:08 +0100
Subject: [PATCH] combine water and settlements

---
 source/pipelines/R/flea_functions.R         | 47 +++++++++++++++++++--
 source/pipelines/run_pipeline.R             | 21 +++++----
 source/pipelines/script_validation_sample.R | 28 +++++++++++-
 3 files changed, 84 insertions(+), 12 deletions(-)

diff --git a/source/pipelines/R/flea_functions.R b/source/pipelines/R/flea_functions.R
index f87fe76..9026ee8 100644
--- a/source/pipelines/R/flea_functions.R
+++ b/source/pipelines/R/flea_functions.R
@@ -473,21 +473,62 @@ download_watersurfaces <- function(path_flea_data, meta) {
   return(path)
 }
 
-get_watersurfaces <- function(path_version, polygons) {
+get_watersurfaces <- function(path_version, polygons, meta) {
   #https://inbo.github.io/n2khab/reference/read_watersurfaces.html
+
+  file_version <- switch(
+    meta$version,
+    "v1.0" = file.path(path_version, "watersurfaces.shp"),
+    "v1.1" = file.path(path_version, "watersurfaces.gpkg"),
+    "v1.2" = file.path(path_version, "watersurfaces.gpkg"),
+    "v2024" = file.path(path_version, "watersurfaces.gpkg")
+  )
+
   ws <- n2khab::read_watersurfaces(
-    path_version,
+    file = file_version,
     version = basename(path_version),
     fix_geom = TRUE
     )
   ws <- vect(ws)
 
   ws <- spatvector_crop(x = ws, y = polygons)
-
+  ws$year_flea <- meta$year_flea
   ws$layer <- basename(path_version)
   ws$value <- NA
+  ws$area_name <- NULL
+  ws$wfd_type_certain <- NULL
 
   return(ws)
 }
 
 
+combine_grb_inbo_water <- function(grb_water, inbo_water, meta) {
+  inbo_water <- vect(inbo_water)
+  inbo_water <- inbo_water[
+    inbo_water$layer == paste0("watersurfaces_",meta$version), ]
+
+  # cover: values of x that overlap with y are replaced by y
+  grb_water <- aggregate(x = grb_water, by = names(grb_water))
+  water <- cover(x = grb_water, y = inbo_water)
+  return(water)
+}
+
+combine_water_settlements <- function(water, settlements, polygons) {
+  vp <- vect(polygons)
+
+  vplist <- vector("list", nrow(vp))
+  for (i in seq_along(vp)) {
+    vp_ <- vp[i]
+    w_ <- water[vp_]
+    s_ <- settlements[vp_]
+    out <- cover(vp_, cover(w_, s_))
+    out$grts_rank <- vp_$grts_rank
+    out$cell <- vp_$cell
+    out$stratum_name <- vp_$stratum_name
+    out$changecat <- vp_$changecat
+    vplist[[i]] <- out
+  }
+  vp_wa_se <- vect(vplist)
+
+  return(vp_wa_se)
+}
diff --git a/source/pipelines/run_pipeline.R b/source/pipelines/run_pipeline.R
index 5507a82..01b1bd4 100644
--- a/source/pipelines/run_pipeline.R
+++ b/source/pipelines/run_pipeline.R
@@ -167,14 +167,19 @@ test <- get_grb_by_row(
   polygons = tar_read(validation_polygons_6e7d3123e950eb4d)[1:2,]
 )
 targets::tar_load_globals()
-targets::tar_workspace("grb_settlements_b23f12e3c6cbbb94")
-debugonce(get_grb_by_row)
-test <- get_grb_by_row(
-  layer = lyrs_settlements,
-  polygons = validation_polygons
+targets::tar_workspace("vp_water_settlements_bc312eaba6d3a034")
+debugonce(combine_water_settlements)
+test <- combine_water_settlements(
+  water = tar_read(vp_water, branches = 1)[[1]],
+  settlements = tar_read(grb_settlements_processed),
+  polygons = tar_read(validation_polygons)
 )
 
 targets::tar_load_globals()
-targets::tar_workspace("grb_parcels_processed")
-debugonce(process_parcels)
-test <- process_parcels(grb_parcels)
+targets::tar_workspace("vp_water_settlements_bc312eaba6d3a034")
+debugonce(get_watersurfaces)
+test <- get_watersurfaces(
+  path_version = zenodo_watersurface,
+  polygons = validation_polygons,
+  meta = watersurfaces_meta
+)
diff --git a/source/pipelines/script_validation_sample.R b/source/pipelines/script_validation_sample.R
index 6357e29..91cd24d 100644
--- a/source/pipelines/script_validation_sample.R
+++ b/source/pipelines/script_validation_sample.R
@@ -55,6 +55,7 @@ flea_data <- gsub(
   pattern = "flea-extent", replacement = "flea-data", x = git_root
 )
 input_names <- c("reclass_bwk2016", "reclass_bwk2020", "reclass_bwk2023")
+input_years <- c(2016, 2019, 2022)
 path_to_gdb <- "Z:/Projects/PRJ_FLEA/flea_data.gdb"
 path_to_lbg <- "Z:/Projects/PRJ_FLEA/landbouwdata.gdb"
 path_to_lyr <- "Z:/Projects/PRJ_FLEA/reclass_bwk2016.lyr"
@@ -296,6 +297,7 @@ list(
         "ortho_2019_2021",
         "ortho_2021_2023"
       ),
+      year_flea = c(2016, 2019, NA, 2022),
       version = c(
         "v1.0", "v1.1", "v1.2", "v2024"
       )
@@ -309,13 +311,37 @@ list(
     ),
     pattern = map(watersurfaces_meta)
   ),
+  # read INBO watersurfaces maps and crop with validation polygons
   geotargets::tar_terra_vect(
     name = watersurfaces_processed,
     command = get_watersurfaces(
       path_version = zenodo_watersurface,
+      polygons = validation_polygons,
+      meta = watersurfaces_meta
+    ),
+    pattern = cross(
+      map(zenodo_watersurface, watersurfaces_meta),
+      validation_polygons
+    )
+  ),
+  # combine the GRB water layer with the INBO watersurfaces
+  geotargets::tar_terra_vect(
+    name = vp_water,
+    command = combine_grb_inbo_water(
+      grb_water = grb_water_wtz_processed,
+      inbo_water = watersurfaces_processed,
+      meta = watersurfaces_meta
+    ),
+    pattern = map(watersurfaces_meta)
+  ),
+  geotargets::tar_terra_vect(
+    name = vp_water_settlements,
+    command = combine_water_settlements(
+      water = vp_water,
+      settlements = grb_settlements_processed,
       polygons = validation_polygons
     ),
-    pattern = cross(zenodo_watersurface, validation_polygons)
+    pattern = map(vp_water)
   )
   #,
   #geotargets::tar_terra_vect(