update pipeline

- add grb
- add lbg

source/pipelines/R/flea_functions.R

@@ -288,6 +288,177 @@ get_grb_by_row <- function(layer, polygons) {
     out[[as.character(namesvec[i])]] <- grb
   }
   out <- terra::vect(out)
+  out$layer <- layer
+  # convert date(time) fields to ISO-8601 format
+  time <- terra::datatype(out) == "time"
+  datetimecols <- terra::names(out)[time]
+  if (length(datetimecols) > 0) {
+    for (i in datetimecols) {
+      out[[i]] <-
+        format(out[[i]], format = "%Y-%m-%dT%H:%M:%S.000Z", tz = "UTC")
+    }
+  }
+
   return(out)
 }
 
+process_water_wtz <- function(grb_wtz) {
+  grb <- vect(grb_wtz)
+  grb <- sf::st_as_sf(grb)
+  grb <- grb |>
+    dplyr::mutate(
+      jaar = pmin(lubridate::year(BEGINDATUM),
+                  lubridate::year(OPNDATUM), na.rm = TRUE),
+      layer = "GRB:WTZ",
+      value = NA
+    )
+
+
+  # select columns
+  grb <- grb |>
+    select(
+      gml_id,
+      grts_rank,
+      jaar,
+      layer,
+      value
+    )
+
+  grb <- vect(grb)
+
+  return(grb)
+}
+
+process_parcels <- function(grb) {
+  grb <- vect(grb)
+  grb <- grb[, c("gml_id", "grts_rank", "BEGINDATUM", "FISCDATUM")]
+  grb <- sf::st_as_sf(grb) |>
+    dplyr::mutate(
+      jaar = pmin(lubridate::year(BEGINDATUM),
+                  lubridate::year(FISCDATUM), na.rm = TRUE),
+      layer = "GRB:ADP",
+      value = NA
+    ) |>
+    dplyr::select(gml_id, grts_rank, jaar, layer, value)
+  # cast to lines
+  grb <- terra::vect(grb) |> terra::as.lines()
+  return(grb)
+}
+
+process_settlement <- function(grb) {
+  grb <- vect(grb) # this should combine grb settlement layers
+  grb <- sf::st_as_sf(grb)
+
+  # enkel TRN met bodembedekking verhard
+  grb <- grb |>
+    dplyr::filter(is.na(LBLBDMBD) | LBLBDMBD == "verhard") |>
+    tidyr::unite(
+      col = lbl,
+      c(LBLTYPE, LBLFNCT),
+      na.rm = TRUE,
+      remove = FALSE) |>
+    dplyr::mutate(
+      jaar = pmin(lubridate::year(BEGINDATUM),
+                  lubridate::year(OPNDATUM), na.rm = TRUE)
+    )
+
+  # assign values
+  # 101 1.1 Settlements - buildings
+  # 102 1.2 Settlements - sealed soil
+  # 105 1.5 Settlements - water
+  # 106 1.6 Settlements - unknown land cover
+  grb <- grb |>
+    dplyr::mutate(
+      value = case_when(
+        layer %in% c("GRB:GBG", "GRB:GBA") ~ 101,
+        layer %in% c(
+          "GRB:WBN",
+          "GRB:SBN",
+          "GRB:KNW") |
+          (layer == "GRB:TRN" & LBLBDMBD == "verhard") ~ 102,
+        TRUE ~ NA
+      )
+    )
+
+  # select columns
+  grb <- grb |>
+    select(
+      gml_id,
+      grts_rank,
+      layer,
+      jaar,
+      lbl,
+      value
+    )
+
+  grb <- vect(grb)
+  return(grb)
+}
+
+
+get_lbg_layernames <- function(path_to_lbg) {
+
+  lyrs <- terra::vector_layers(path_to_lbg)
+  lyrs <- lyrs[grepl(".+20\\d\\d.+", lyrs)]
+  return(lyrs)
+}
+
+get_lbg <- function(
+    path_to_lbg, layer, from_fields, where_field, where_values, flea_value) {
+
+  assertthat::assert_that(assertthat::is.string(layer))
+  assertthat::assert_that(is.character(from_fields))
+  assertthat::assert_that(assertthat::is.string(where_field))
+  assertthat::assert_that(is.numeric(where_values))
+
+  query <- paste0(
+    "SELECT ",
+    paste0(from_fields, collapse = ","),
+    " FROM ",
+    layer,
+    " WHERE ",
+    where_field,
+    " IN (",
+    paste0("'", where_values, "'", collapse = ","),
+    ")"
+  )
+
+  lbg <- vect(x = path_to_lbg,
+              query = query,
+              crs = "EPSG:31370")
+  lbg$layer <- layer
+  lbg$value <- flea_value
+  # convert date(time) fields to ISO-8601 format
+  time <- terra::datatype(lbg) == "time"
+  datetimecols <- terra::names(lbg)[time]
+  if (length(datetimecols) > 0) {
+    for (i in datetimecols) {
+      lbg[[i]] <-
+        format(lbg[[i]], format = "%Y-%m-%dT%H:%M:%S.000Z", tz = "UTC")
+    }
+  }
+
+  return(lbg)
+}
+
+spatvector_crop <- function(x, y) {
+  assertthat::assert_that(inherits(x, "SpatVector"))
+  assertthat::assert_that(inherits(y, "SpatVector"))
+
+  # spatial subset
+  x <- x[y]
+  # topology fix if needed
+  x <- terra::makeValid(x)
+  # if y contains overlapping polygons dissolve them
+  y <- terra::aggregate(y)
+  # tryCatch?
+  out <- terra::crop(x, y)
+  return(out)
+}
+
+
+get_watervlakken <- function() {
+  #https://inbo.github.io/n2khab/reference/read_watersurfaces.html
+}
+
+

source/pipelines/R/process_rasters.R

@@ -232,6 +232,9 @@ get_changecat_columns <- function(tempstrat) {
 }
 
 
-
+calc_mask <- function(maps, values) {
+  my_mask <- maps %in% values
+  return(my_mask)
+}
 
 

source/pipelines/run_pipeline.R

@@ -14,10 +14,12 @@ tar_make()
 
 #targets::tar_prune()
 
-targets::tar_meta(
+mt <- targets::tar_meta(
   fields = error,
   complete_only = TRUE
 )
+mt
+#View(mt)
 targets::tar_meta(fields = warnings, complete_only = TRUE)
 targets::tar_visnetwork(label = c("description", "time", "size"))
 
@@ -35,8 +37,18 @@ tar_read(mapnames)
 ct <- tar_read(catstable)
 ml <- tar_read(maps)
 
-ct
+
 ml[[1]]
+terra::cats(ml[[1]])
+settlement_mask <- app(
+  ml[[1]],
+  fun = function(x) {
+    x[!x %in% c(101, 102, 105, 106)] <- NA
+    return(x)
+  }
+)
+plot(settlement_mask, colNA = "snow4")
+
 terra::plot(ml[[1]], colNA = "orange")
 terra::values(ml[[1]], row = 5000, nrows = 1)
 terra::coltab(ml[[1]])
@@ -81,6 +93,10 @@ vs <- targets::tar_read(validation_sample)
 terra::vect(vs) |> sf::st_as_sf(crs = 31370) |>
   sf::st_drop_geometry() |>
   dplyr::count(grts_rank) |> dplyr::count(n)
+terra::vect(vs) |> sf::st_as_sf(crs = 31370) |>
+  sf::st_drop_geometry() |>
+  dplyr::count(stratum_name, changecat) |>
+  tidyr::pivot_wider(names_from = changecat, values_from = n)
 
 vp <- targets::tar_read(validation_polygons)
 lapply(vp, nrow) |> unlist() |> sum()
@@ -99,23 +115,59 @@ terra::vect(vs) |>
   st_drop_geometry() |>
   dplyr::count(stratum_name, changecat)
 
+grb_parc <- tar_read(grb_parcels)
+grb_parc <- terra::vect(grb_parc) |>
+  st_as_sf() |>
+  dplyr::mutate(source = "parcels")
+grb_set <- tar_read(grb_settlements)
+grb_set <- terra::vect(grb_set) |>
+  st_as_sf() |>
+  dplyr::mutate(source = "settlements")
+grb_water <- tar_read(grb_waterways)
+grb_water <- terra::vect(grb_water) |>
+  st_as_sf() |>
+  dplyr::mutate(source = "water")
+
+dplyr::bind_rows(
+  grb_water,
+  grb_set) |>
+  mapview::mapview(zcol = "source", alpha.regions = 0.2) +
+  mapview::mapview(grb_parc, alpha.region = 0)
+
+lbg_101 <- tar_read(lbg_101_cropped)
+lbg_104 <- tar_read(lbg_104_cropped)
+
+mapview::mapview(terra::vect(lbg_101), alpha.regions = 0.2
+                 , col.regions = "orange") +
+  mapview::mapview(terra::vect(lbg_104), alpha.regions = 0.2,
+                   col.regions = "yellow") +
+  mapview::mapview(terra::vect(vp), alpha.regions = 0)
+
 ##################
 # debug pipeline #
 ##################
 
 targets::tar_load_globals()
-tar_load(names = c(separate_grts))
-debugonce(get_changecats)
-get_changecats(separate_grts)
-
+debugonce(add_changecats_tempstrat)
+add_changecats_tempstrat(
+  tempstrat = temporal_map,
+  cats = catstable,
+  mapnames = mapnames
+)
+debugonce(get_grb_by_row)
+test <- get_grb_by_row(
+  layer = "GRB:ADP",
+  polygons = tar_read(validation_polygons_6e7d3123e950eb4d)[1:2,]
+)
 targets::tar_load_globals()
-targets::tar_workspace("grb_waterways_14e96c73491bb69c")
+targets::tar_workspace("grb_settlements_b23f12e3c6cbbb94")
 debugonce(get_grb_by_row)
 test <- get_grb_by_row(
-  layer = lyrs_waterways,
+  layer = lyrs_settlements,
   polygons = validation_polygons
 )
 
-
-
-
+targets::tar_load_globals()
+targets::tar_workspace("grb_parcels_processed")
+debugonce(process_parcels)
+test <- process_parcels(grb_parcels)