Add "emissions" workflow (#28) for Beijing

DESCRIPTION

@@ -22,6 +22,8 @@ Description: Used in Project KEYS for generating inputs to
 License: MIT + file LICENSE
 URL: https://github.com/KWB-R/urbanAnnualRunoff
 BugReports: https://github.com/KWB-R/urbanAnnualRunoff/issues
+Depends: 
+    R >= 3.5.0
 Imports: 
     caret,
     dplyr,
@@ -30,12 +32,14 @@ Imports:
     fs,
     kwb.utils,
     lubridate,
+    magrittr,
     randomForest,
     raster,
     rlang,
     rgdal,
     remotes,
-    magrittr
+    stringr,
+    tidyselect
 Suggests: 
     covr,
     knitr,

LICENSE

@@ -1,2 +1,22 @@
-YEAR: 2020
-COPYRIGHT HOLDER: Kompetenzzentrum Wasser Berlin gGmbH (KWB)
+MIT License
+
+Copyright (c) 2020-2021 Kompetenzzentrum Wasser Berlin gGmbH (KWB)
+
+Permission is hereby granted, free of charge, to any person obtaining a copy
+of this software and associated documentation files (the "Software"), to deal
+in the Software without restriction, including without limitation the rights
+to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
+copies of the Software, and to permit persons to whom the Software is
+furnished to do so, subject to the following conditions:
+
+The above copyright notice and this permission notice shall be included in all
+copies or substantial portions of the Software.
+
+THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
+SOFTWARE.
+

LICENSE.md

@@ -1,6 +1,6 @@
 # MIT License
 
-Copyright (c) 2020 Kompetenzzentrum Wasser Berlin gGmbH (KWB)
+Copyright (c) 2020-2021 Kompetenzzentrum Wasser Berlin gGmbH (KWB)
 
 Permission is hereby granted, free of charge, to any person obtaining a copy
 of this software and associated documentation files (the "Software"), to deal

NAMESPACE

@@ -2,6 +2,7 @@
 
 export("%>%")
 export(buildClassMod)
+export(calculate_loads)
 export(computeABIMOclimate)
 export(makeFLGES)
 export(makeOverlay)
@@ -11,15 +12,20 @@ export(makeVG)
 export(padCODE)
 export(postProcessABIMO)
 export(predictSurfClass)
+export(read_concentrations)
 importFrom(caret,createFolds)
 importFrom(caret,train)
 importFrom(caret,trainControl)
 importFrom(doParallel,registerDoParallel)
+importFrom(dplyr,bind_rows)
+importFrom(dplyr,c_across)
 importFrom(dplyr,group_by)
 importFrom(dplyr,left_join)
+importFrom(dplyr,rowwise)
 importFrom(dplyr,summarize)
 importFrom(foreign,read.dbf)
 importFrom(kwb.utils,multiSubstitute)
+importFrom(kwb.utils,substSpecialChars)
 importFrom(lubridate,month)
 importFrom(lubridate,year)
 importFrom(magrittr,"%>%")
@@ -34,5 +40,8 @@ importFrom(raster,res)
 importFrom(raster,shapefile)
 importFrom(raster,writeRaster)
 importFrom(rlang,.data)
+importFrom(stats,setNames)
+importFrom(stringr,str_replace)
+importFrom(tidyselect,all_of)
 importFrom(utils,read.table)
 importFrom(utils,write.table)

R/emissions.R

@@ -0,0 +1,88 @@
+#' Emissions: read concentrations from OgRe database
+#' @description imports data from OgRe database and selects relevant substances
+#' for case study sites (Beijing, Jinxi) and calculates mean concentrations
+#' over all structures (column: "mean"). In addition new columns (short_name,
+#' unit_load, label_load) are created
+#' @param path path to OgRe database file "annual_mean_conc.csv"
+#' @return data frame with selected substances and column
+#' @export
+#' @importFrom kwb.utils multiSubstitute substSpecialChars
+#' @importFrom rlang .data
+#' @importFrom stringr str_replace
+#' @importFrom dplyr c_across rowwise
+#' @importFrom tidyselect all_of
+#'
+read_concentrations <- function(path) {
+
+  # grab AMC from OgRe
+  x_conc <- read.table(file = path,
+                       sep = ";",
+                       dec = ".",
+                       stringsAsFactors = FALSE,
+                       header = TRUE)
+
+  ### select substances of interest
+  shortnames_list <- list("Biologischer Sauerstoffbedarf" = "BOD" ,
+                          "Chemischer Sauerstoffbedarf" = "COD",
+                          "Abfiltrierbare Stoffe" = "TSS",
+                          "Blei$" = "Pb",
+                          "Cadmium$" = "Cd",
+                          "Chrom$" = "Cr",
+                          "Kupfer$" = "Cu",
+                          "Nickel$" = "Ni",
+                          "Titan$" = "Ti",
+                          "Vanadium$" = "Va",
+                          "Zink$" = "Zn")
+
+  x_conc$short_name <- kwb.utils::multiSubstitute(
+    strings = x_conc$VariableName,
+    replacements = shortnames_list)
+
+  short_names <- as.character(unlist(shortnames_list))
+
+  structures <- c("ALT", "NEU", "STR", "EFH", "GEW", "ANDERE")
+
+  # average across catchment types
+  concentrations  <- x_conc %>%
+    dplyr::filter(short_name %in% short_names) %>%
+    dplyr::mutate(unit_load = kwb.utils::substSpecialChars(.data$UnitsAbbreviation) %>%
+                    stringr::str_replace("L", "m2_year"),
+                  label_load = sprintf("%s.%s",
+                                       .data$short_name,
+                                       .data$unit_load)) %>%
+    dplyr::rowwise(VariableID) %>%
+    dplyr::mutate(mean = mean(dplyr::c_across(tidyselect::all_of(structures))))
+
+  concentrations
+}
+
+#' Emissions: calculate loads
+#' @description  The annual load is calculated with V x c. For for heavy metals
+#' -> l/m2-year x ug/l = ug/m2-year;  for BOD/COD/TSS -> l/m2-year x mg/l = mg/m2-year
+#' @param abimo_inpout data.frame or SpatialPolygonsDataFrame with ABIMO input and
+#' output as retrieved by \code{\link{postProcessABIMO}}
+#' @param concentrations concentrations data frame as retrieved by
+#' \code{\link{read_concentrations}}
+#' @return add calculated loads as additional colums  to abimo_inpout data.frame
+#' or SpatialPolygonsDataFrame
+#' @export
+#' @importFrom dplyr bind_rows
+#' @importFrom stats setNames
+calculate_loads <- function(abimo_inpout,
+                            concentrations) {
+  # annual load = V * c
+  # for heavy metals -> l/m2-year * ug/l = ug/m2-year
+  # for TSS -> l/m2-year * mg/l = mg/m2-year
+  loads <- dplyr::bind_rows(
+    stats::setNames(lapply(X = concentrations$mean,
+                    FUN = function(a){
+                      abimo_inpout$ROW*a/1e3
+                    }
+    ),
+    nm = concentrations$label_load
+    )
+  )
+
+  # add computed loads to ABIMO dataset
+  cbind(abimo_inpout, loads)
+}

vignettes/workflow_beijing.Rmd

@@ -41,7 +41,7 @@ path_list <- list(
   abimo_exe = "<abimo>/Abimo3_2.exe",
   gis = "<root_path>/<data>/_DataAnalysis/gis",
   climate = "<root_path>/<data>/_DataAnalysis/climate",
-  emissions_input = "<root_path>/<data>/_DataAnalysis/emissions/input",
+  emissions_input = "<root_path>/<data>/_DataAnalysis/emissions/input/annual_mean_conc.csv",
   emissions_output = "<root_path>/<data>/_DataAnalysis/emissions/output"
 )
 
@@ -297,56 +297,13 @@ abimo_inpout <- urbanAnnualRunoff::postProcessABIMO(
 
 
 # compute emissions Tongzhou with OgRe AMC (annual mean concentration) database
+conc <- urbanAnnualRunoff::read_concentrations(path = paths$emissions_input)
 
-cDataFile <- file.path(paths$emissions_input, "annual_mean_conc.csv")
-
-# grab AMC from OgRe
-x_conc <- read.table(file = cDataFile, 
-                     sep = ";", 
-                     dec = ".", 
-                     stringsAsFactors = FALSE, 
-                     header = TRUE)
-
-# select substances
-#c_heavy_metals <- x_conc[5:13, 4:ncol(x_conc)]
-#rownames(c_heavy_metals) <- x_conc$VariableName[5:13]
-
-structures <- c("ALT", "NEU", "STR", "EFH", "GEW", "ANDERE")
-heavy_metals <- c("Abfiltrierbare Stoffe", 
-                  "Blei", 
-                  "Cadmium", 
-                  "Chrom", 
-                  "Kupfer", 
-                  "Nickel", 
-                  "Titan", 
-                  "Vanadium", 
-                  "Zink") 
-
-# average across catchment types
-c_heavy_metals_means  <- x_conc %>% 
-  dplyr::filter(VariableName %in% heavy_metals) %>% 
-  dplyr::mutate(unit_load = kwb.utils::substSpecialChars(.data$UnitsAbbreviation) %>%  
-                  stringr::str_replace("L", "m2_year"),
-                label_load = sprintf("%s.%s", 
-                                     .data$VariableName, 
-                                     .data$unit_load)) %>% 
-  dplyr::rowwise(VariableID) %>% 
-  dplyr::mutate(conc_mean = mean(dplyr::c_across(tidyselect::all_of(structures))))
-
-# annual load = V * c 
-# for heavy metals -> l/m2-year * ug/l = ug/m2-year
-# for TSS -> l/m2-year * mg/l = mg/m2-year
-loads <- dplyr::bind_rows(
-  setNames(lapply(X = c_heavy_metals_means$conc_mean,
-                FUN = function(a){
-                  abimo_inpout$ROW*a/1e3
-                }),
-           nm = c_heavy_metals_means$label_load)
+abimo_inpout_emissions <- urbanAnnualRunoff::calculate_loads(
+  abimo_inpout = abimo_inpout, 
+  concentrations = conc
   )
 
-# add computed loads to ABIMO dataset
-abimo_inpout_emissions <- cbind(abimo_inpout, loads)
-
 # write out joined table output as shapefile
 outFile <- file.path(paths$emissions_output, 
                      sprintf("abimo_%s_emissions.shp", paths$site))