diff --git a/DESCRIPTION b/DESCRIPTION index 3295a7a25087c56b5854db43588eb7818005a2e5..acaf9106109f3a0506765e7a12f35903c62b933b 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,7 +1,7 @@ Package: SWATdoctR Type: Package Title: Finding the right diagnoses and treatments for SWAT+ models -Version: 0.1.21 +Version: 0.1.22 Author: c(person("Christoph", "Schürz", email = "christoph.schuerz@ufz.de", role = c("aut", "cre")), @@ -14,6 +14,8 @@ Imports: data.table, dplyr, ggplot2, + ggrepel, + grDevices, lubridate, patchwork, png, diff --git a/NAMESPACE b/NAMESPACE index 44d53bee96333ff38ae27da39d8f34a88b784a65..a91bc53b9c829e765a24ecb5a0bccfc264913658 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -69,6 +69,8 @@ importFrom(ggplot2,lims) importFrom(ggplot2,theme) importFrom(ggplot2,theme_bw) importFrom(ggplot2,theme_void) +importFrom(ggrepel,geom_text_repel) +importFrom(grDevices,palette.colors) importFrom(lubridate,ceiling_date) importFrom(lubridate,floor_date) importFrom(lubridate,int_end) @@ -117,6 +119,7 @@ importFrom(tidyr,pivot_longer) importFrom(tidyr,pivot_wider) importFrom(tidyr,unite) importFrom(tidyselect,all_of) +importFrom(tidyselect,any_of) importFrom(tidyselect,ends_with) importFrom(tidyselect,everything) importFrom(tidyselect,starts_with) diff --git a/R/globals.R b/R/globals.R index ce423d3caa3fcf0de2ee2b4beb9359b3f887b765..8785e546a943fb2af8edc6157826179ee93f22eb 100644 --- a/R/globals.R +++ b/R/globals.R @@ -8,4 +8,5 @@ utils::globalVariables(c("%&&%", "%//%", ".", "crop", "day", "ecanopy", "eplant" "wndspd", "yr", "rhum", "rm_skp", "Date", "Values", "surq_gen", "latq", "perc", "description", "flo", "..density..", "wateryld", "cn", "sw_final", "esr", "mmdd", "lu_mgt_ini", "has_suffix", "suffix_upd", "lu_mgt_upd", "schedule_upd", - "plnt_com_upd", "n_chr", "simulation")) + "plnt_com_upd", "n_chr", "simulation", "start", "end", "mid", "stress", "n_op", + "grw_group", "date_diff", "no_plntkill", "grw_group", "no_plntkill")) diff --git a/R/plot_hru_pw.R b/R/plot_hru_pw.R index 404815891964653364b4399e8edd33416ec99e7a..9aa650c0f6ecbb93c80849574b44624c5f73fee4 100644 --- a/R/plot_hru_pw.R +++ b/R/plot_hru_pw.R @@ -53,6 +53,7 @@ get_hru_id_by_attribute <- function(sim_verify, lum = NULL, mgt = NULL, soil = N #' @importFrom tibble tibble #' @importFrom tidyr pivot_longer pivot_wider #' @importFrom tidyselect all_of +#' @importFrom grDevices palette.colors #' @import ggplot2 #' @import patchwork #' diff --git a/R/plot_ps_tile.R b/R/plot_ps_tile.R index 78624c883ae71cfc28f10a2f0493b36c2b5fb1bb..0771a1c4c81ee39fab899d8324b08504911edf54 100644 --- a/R/plot_ps_tile.R +++ b/R/plot_ps_tile.R @@ -9,6 +9,8 @@ #' @importFrom dplyr %>% mutate group_by group_map select bind_rows case_when #' @importFrom ggplot2 ggplot geom_line facet_wrap labs theme_bw theme aes #' @importFrom tidyr pivot_longer +#' @importFrom ggrepel geom_text_repel +#' @importFrom tidyselect any_of #' @export #' #' @examples @@ -17,31 +19,39 @@ #' } plot_ps <- function(sim_verify, conc = FALSE){ + ##Function to calculate concentration + calc_conc <- function(df) { + pivot_longer(df, !any_of(c("yr", "name", "flo")), names_to = 'var', values_to = 'Values') %>% + mutate(Values = ifelse(var == "sed", 1000000*Values/flo, 1000*Values/flo)) %>% + select(any_of(c("yr", "name", "var", "Values"))) %>% + bind_rows(select(df, any_of(c("yr", "name", "flo"))) %>% + mutate(var = "flo") %>% + rename(Values = flo) %>% + mutate(Values = Values/(24*60*60*365.25))) %>% + mutate(var = case_when(var == 'flo' ~ "flo m3/s", + var %in% c("orgn", "no3", "nh3", "no2") ~ paste(var, "N mg/l"), + var %in% c("sedp", "solp") ~ paste(var, "P mg/l"), + TRUE ~ paste(var, "mg/l"))) + } + ##Function to calculate loads + calc_load <- function(df) { + pivot_longer(df, !any_of(c("yr", "name")), names_to = 'var', values_to = 'Values') %>% + mutate(var = case_when(var == 'flo' ~ "flo m3/d", + var == 'sed' ~ "sed t/y", + var %in% c("orgn", "no3", "nh3", "no2") ~ paste(var, "N kg/y"), + var %in% c("sedp", "solp") ~ paste(var, "P kg/y"), + TRUE ~ paste(var, "kg/y"))) + } + ##Case non constant point sources if(!is.null(sim_verify$recall_yr)){ df <- sim_verify$recall_yr[, -c(1:3,6)] %>% .[, colSums(.!= 0) > 0] %>% - mutate(name = gsub("hru00", "ps", name)) + mutate(name = gsub("hru00", "ps", name), + yr = as.factor(yr)) if(conc){ - df <- df %>% - pivot_longer(-c(yr, name, flo), names_to = 'var', values_to = 'Values') %>% - mutate(Values = ifelse(var == "sed", 1000000*Values/flo, 1000*Values/flo)) %>% - select(yr, name, var, Values) %>% - bind_rows(df[c("yr", "name", "flo")] %>% - mutate(var = "flo") %>% - rename(Values = flo) %>% - mutate(Values = Values/(24*60*60*365.25))) %>% - mutate(var = case_when(var == 'flo' ~ "flo m3/s", - var %in% c("orgn", "no3", "nh3", "no2") ~ paste(var, "N mg/l"), - var %in% c("sedp", "solp") ~ paste(var, "P mg/l"), - TRUE ~ paste(var, "mg/l"))) + df <- calc_conc(df) } else { - df <- df %>% - pivot_longer(-c(yr, name), names_to = 'var', values_to = 'Values') %>% - mutate(var = case_when(var == 'flo' ~ "flo m3/d", - var == 'sed' ~ "sed t/y", - var %in% c("orgn", "no3", "nh3", "no2") ~ paste(var, "N kg/y"), - var %in% c("sedp", "solp") ~ paste(var, "P kg/y"), - TRUE ~ paste(var, "kg/y"))) + df <- calc_load(df) } fig <- ggplot(df, aes(x=yr, y=Values, group=name, colour=name))+ geom_line(linewidth=1)+ @@ -52,10 +62,45 @@ plot_ps <- function(sim_verify, conc = FALSE){ strip.text = element_text(color = "white", face="bold"), panel.border = element_rect(colour = "grey80"), axis.text.x = element_text(angle = 25, hjust=1)) - - return(fig) + } + ##Case for constant point sources + if(!is.null(sim_verify$exco_om)){ + df <- sim_verify$exco_om %>% + .[, colSums(.!= 0) > 0] %>% + mutate(name = gsub("hru00", "ps", name)) + if(conc){ + df <- calc_conc(df) + } else { + df <- calc_load(df) + } + ##Putting into single figure + if(exists("fig")){ + fig <- fig + + geom_hline(data = df, + aes(yintercept = Values), linetype = 'dashed', linewidth=0.75) + + geom_text_repel(data = df, aes(0, Values,label = name, vjust = -1, hjust = -1), + size = 3, show.legend = FALSE, inherit.aes = FALSE, segment.color = "grey") + } else { + fig <- ggplot()+ + facet_wrap(~var, scales = "free_y")+ + geom_hline(data = df, + aes(yintercept = Values), linetype = 'dashed', linewidth=0.75) + + geom_text_repel(data = df, aes(0, Values,label = name, vjust = -1, hjust = -1), + size = 3, show.legend = FALSE, inherit.aes = FALSE, segment.color = "grey")+ + theme_bw()+ + theme(strip.background = element_rect(fill = "grey50", colour = "grey80"), + strip.text = element_text(color = "white", face="bold"), + panel.border = element_rect(colour = "grey80"), + axis.text.x = element_blank(), + axis.ticks.x = element_blank(), + axis.title.x = element_blank()) + } + } + ##If no point sources exist + if(is.null(sim_verify$recall_yr) & is.null(sim_verify$exco_om)){ + return(print("No point sources exists in this model setup!!!")) } else { - print("No point sources exists in this model setup!!!") + return(fig) } } diff --git a/R/run_swat_verify.R b/R/run_swat_verify.R index 4b6173cddc22d2c0ca8fc67904d3fd933c11138c..2133948f513f7a3be8954c64f4b88d4e4a9e585a 100644 --- a/R/run_swat_verify.R +++ b/R/run_swat_verify.R @@ -85,6 +85,12 @@ run_swat_verification <- function(project_path, outputs = c('wb', 'mgt', 'plt'), error = function(e) { model_output$recall_yr <- NULL }) + tryCatch({ + model_output$exco_om <- read_tbl('exco_om.exc', run_path, 2) %>% lwr + }, + error = function(e) { + model_output$exco_om <- NULL + }) } if ('mgt' %in% outputs) { model_output$mgt_out <- read_mgt(run_path) %>% lwr