From 12a31f4c130c551f82232d9ef7dfb608bd52c53f Mon Sep 17 00:00:00 2001 From: Johannes Ranke Date: Tue, 27 Sep 2016 23:00:48 +0200 Subject: Reorganise repository using standard package layout --- pkg/R/FOCUS_GW_scenarios_2012.R | 11 - pkg/R/GUS.R | 98 --------- pkg/R/PEC_soil.R | 215 ------------------- pkg/R/PEC_sw_drainage_UK.R | 66 ------ pkg/R/PEC_sw_drift.R | 65 ------ pkg/R/PEC_sw_sed.R | 50 ----- pkg/R/SFO_actual_twa.R | 36 ---- pkg/R/SSLRC_mobility_classification.R | 42 ---- pkg/R/TOXSWA_cwa.R | 380 ---------------------------------- pkg/R/drift_data_JKI.R | 44 ---- pkg/R/endpoint.R | 110 ---------- pkg/R/geomean.R | 38 ---- pkg/R/pfm_degradation.R | 48 ----- pkg/R/soil_scenario_data_EFSA_2015.R | 40 ---- 14 files changed, 1243 deletions(-) delete mode 100644 pkg/R/FOCUS_GW_scenarios_2012.R delete mode 100644 pkg/R/GUS.R delete mode 100644 pkg/R/PEC_soil.R delete mode 100644 pkg/R/PEC_sw_drainage_UK.R delete mode 100644 pkg/R/PEC_sw_drift.R delete mode 100644 pkg/R/PEC_sw_sed.R delete mode 100644 pkg/R/SFO_actual_twa.R delete mode 100644 pkg/R/SSLRC_mobility_classification.R delete mode 100644 pkg/R/TOXSWA_cwa.R delete mode 100644 pkg/R/drift_data_JKI.R delete mode 100644 pkg/R/endpoint.R delete mode 100644 pkg/R/geomean.R delete mode 100644 pkg/R/pfm_degradation.R delete mode 100644 pkg/R/soil_scenario_data_EFSA_2015.R (limited to 'pkg/R') diff --git a/pkg/R/FOCUS_GW_scenarios_2012.R b/pkg/R/FOCUS_GW_scenarios_2012.R deleted file mode 100644 index 49cf7dd..0000000 --- a/pkg/R/FOCUS_GW_scenarios_2012.R +++ /dev/null @@ -1,11 +0,0 @@ -#' A very small subset of the FOCUS Groundwater scenario defitions -#' -#' Currently, only a small subset of the soil definitions are provided. -#' -#' @name FOCUS_GW_scenarios_2012 -#' @references FOCUS (2012) Generic guidance for Tier 1 FOCUS ground water assessments. Version 2.1. -#' FOrum for the Co-ordination of pesticde fate models and their USe. -#' http://focus.jrc.ec.europa.eu/gw/docs/Generic_guidance_FOCV2_1.pdf -#' @examples -#' FOCUS_GW_scenarios_2012 -NULL diff --git a/pkg/R/GUS.R b/pkg/R/GUS.R deleted file mode 100644 index 4a6532d..0000000 --- a/pkg/R/GUS.R +++ /dev/null @@ -1,98 +0,0 @@ -# Copyright (C) 2015 Johannes Ranke -# Contact: jranke@uni-bremen.de -# This file is part of the R package pfm - -# This program is free software: you can redistribute it and/or modify it under -# the terms of the GNU General Public License as published by the Free Software -# Foundation, either version 3 of the License, or (at your option) any later -# version. - -# This program is distributed in the hope that it will be useful, but WITHOUT -# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS -# FOR A PARTICULAR PURPOSE. See the GNU General Public License for more -# details. - -# You should have received a copy of the GNU General Public License along with -# this program. If not, see - -#' Groundwater ubiquity score based on Gustafson (1989) -#' -#' The groundwater ubiquity score GUS is calculated according to -#' the following equation -#' \deqn{GUS = \log_{10} DT50_{soil} (4 - \log_{10} K_{oc})}{GUS = log10 DT50soil * (4 - log10 Koc)} -#' -#' @references Gustafson, David I. (1989) Groundwater ubiquity score: a simple -#' method for assessing pesticide leachability. \emph{Environmental -#' toxicology and chemistry} \bold{8}(4) 339–57. -#' @inheritParams endpoint -#' @param DT50 Half-life of the chemical in soil. Should be a field -#' half-life according to Gustafson (1989). However, leaching to the sub-soil -#' can not completely be excluded in field dissipation experiments and Gustafson -#' did not refer to any normalisation procedure, but says the field study should -#' be conducted under use conditions. -#' @param Koc The sorption constant normalised to organic carbon. Gustafson -#' does not mention the nonlinearity of the sorption constant commonly -#' found and usually described by Freundlich sorption, therefore it is -#' unclear at which reference concentration the Koc should be observed -#' (and if the reference concentration would be in soil or in porewater). -#' @param chent If a chent is given with appropriate information present in its -#' chyaml field, this information is used, with defaults specified below. -#' @param degradation_value Which of the available degradation values should -#' be used? -#' @param lab_field Should laboratory or field half-lives be used? This -#' defaults to lab in this implementation, in order to avoid -#' double-accounting for mobility. If comparability with the original GUS -#' values given by Gustafson (1989) is desired, non-normalised first-order -#' field half-lives obtained under actual use conditions should be used. -#' @param redox Aerobic or anaerobic degradation data -#' @param sorption_value Which of the available sorption values should be used? -#' Defaults to Kfoc as this is what is generally available from the European -#' pesticide peer review process. These values generally use a reference -#' concentration of 1 mg/L in porewater, that means they would be expected to -#' be Koc values at a concentration of 1 mg/L in the water phase. -#' @param degradation_aggregator Function for aggregating half-lives -#' @param sorption_aggregator Function for aggregation Koc values -#' @param ... Included in the generic to allow for further arguments later. Therefore -#' this also had to be added to the specific methods. -#' @return A list with the DT50 and Koc used as well as the resulting score -#' of class GUS_result -#' @author Johannes Ranke -#' @export -GUS <- function(...) UseMethod("GUS") - -#' @rdname GUS -#' @export -GUS.numeric <- function(DT50, Koc, ...) { - score <- log10(DT50) * (4 - log10(Koc)) - res <- list(DT50 = DT50, Koc = Koc, score = score) - class(res) <- "GUS_result" - return(res) -} - -#' @rdname GUS -#' @export -GUS.chent <- function(chent, - degradation_value = "DT50ref", - lab_field = "laboratory", - redox = "aerobic", - sorption_value = "Kfoc", - degradation_aggregator = geomean, - sorption_aggregator = geomean, - ...) -{ - DT50 = soil_DT50(chent, lab_field = lab_field, redox = redox, - value = degradation_value, - aggregator = degradation_aggregator, signif = 5) - Koc = soil_Kfoc(chent, value = sorption_value, - aggregator = sorption_aggregator, signif = 5) - GUS.numeric(DT50, Koc) -} - -#' @rdname GUS -#' @export -#' @param x An object of class GUS_result to be printed -#' @param digits The number of digits used in the print method -print.GUS_result = function(x, ..., digits = 1) { - cat("GUS: ", round(x$score, digits = 1), "\n") - cat("calculated from DT50 ", x$DT50, " and Koc ", x$Koc, "\n") -} diff --git a/pkg/R/PEC_soil.R b/pkg/R/PEC_soil.R deleted file mode 100644 index 0263e47..0000000 --- a/pkg/R/PEC_soil.R +++ /dev/null @@ -1,215 +0,0 @@ -# Copyright (C) 2015 Johannes Ranke -# Contact: jranke@uni-bremen.de -# This file is part of the R package pfm - -# This program is free software: you can redistribute it and/or modify it under -# the terms of the GNU General Public License as published by the Free Software -# Foundation, either version 3 of the License, or (at your option) any later -# version. - -# This program is distributed in the hope that it will be useful, but WITHOUT -# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS -# FOR A PARTICULAR PURPOSE. See the GNU General Public License for more -# details. - -# You should have received a copy of the GNU General Public License along with -# this program. If not, see - -# Register global variables -if(getRversion() >= '2.15.1') utils::globalVariables(c("destination", "study_type", "TP_identifier", - "soil_scenario_data_EFSA_2015")) - -#' Calculate predicted environmental concentrations in soil -#' -#' This is a basic calculation of a contaminant concentration in bulk soil -#' based on complete, instantaneous mixing. If an interval is given, an -#' attempt is made at calculating a long term maximum concentration using -#' the concepts layed out for example in the PPR panel opinion (EFSA 2012). -#' -#' This assumes that the complete load to soil during the time specified by -#' 'interval' (typically 365 days) is dosed at once. As in the PPR panel -#' opinion cited below (PPR panel 2012), only temperature correction using the -#' Arrhenius equation is performed. -#' -#' Total soil and porewater PEC values for the scenarios as defined in the EFSA -#' guidance (2015, p. 13) can easily be calculated. -#' -#' @note If temperature information is available in the selected scenarios, as -#' e.g. in the EFSA scenarios, the DT50 for groundwater modelling -#' (destination 'PECgw') is taken from the chent object, otherwise the DT50 -#' with destination 'PECsoil'. -#' @importFrom methods is -#' @param rate Application rate in units specified below -#' @param rate_units Defaults to g/ha -#' @param interception The fraction of the application rate that does not reach the soil -#' @param mixing_depth Mixing depth in cm -#' @param interval Period of the deeper mixing, defaults to 365, which is a year if -#' rate units are in days -#' @param n_periods Number of periods to be considered for long term PEC calculations -#' @param PEC_units Requested units for the calculated PEC. Only mg/kg currently supported -#' @param PEC_pw_units Only mg/L currently supported -#' @param tillage_depth Periodic (see interval) deeper mixing in cm -#' @param chent An optional chent object holding substance specific information. Can -#' also be a name for the substance as a character string -#' @param DT50 If specified, overrides soil DT50 endpoints from a chent object -#' If DT50 is not specified here and not available from the chent object, zero -#' degradation is assumed -#' @param Koc If specified, overrides Koc endpoints from a chent object -#' @param Kom Calculated from Koc by default, but can explicitly be specified -#' as Kom here -#' @param t_avg Averaging times for time weighted average concentrations -#' @param scenarios If this is 'default', the DT50 will be used without correction -#' and soil properties as specified in the REACH guidance (R.16, Table -#' R.16-9) are used for porewater PEC calculations. If this is "EFSA_2015", -#' the DT50 is taken to be a modelling half-life at 20°C and pF2 (for when -#' 'chents' is specified, the DegT50 with destination 'PECgw' will be used), -#' and corrected using an Arrhenius activation energy of 65.4 kJ/mol. Also -#' model and scenario adjustment factors from the EFSA guidance are used. -#' @param porewater Should equilibrium porewater concentrations be estimated -#' based on Kom and the organic carbon fraction of the soil instead of total -#' soil concentrations? Based on equation (7) given in the PPR panel opinion -#' (EFSA 2012, p. 24) and the scenarios specified in the EFSA guidance (2015, -#' p. 13). -#' @return The predicted concentration in soil -#' @references EFSA Panel on Plant Protection Products and their Residues (2012) -#' Scientific Opinion on the science behind the guidance for scenario -#' selection and scenario parameterisation for predicting environmental -#' concentrations of plant protection products in soil. \emph{EFSA Journal} -#' \bold{10}(2) 2562, doi:10.2903/j.efsa.2012.2562 -#' -#' EFSA (European Food Safety Authority) (2015) EFSA guidance document for -#' predicting environmental concentrations of active substances of plant -#' protection products and transformation products of these active substances -#' in soil. \emph{EFSA Journal} \bold{13}(4) 4093 -#' doi:10.2903/j.efsa.2015.4093 -#' @author Johannes Ranke -#' @export -#' @examples -#' PEC_soil(100, interception = 0.25) -#' -#' # This is example 1 starting at p. 79 of the EFSA guidance (2015) -#' PEC_soil(1000, interval = 365, DT50 = 250, t_avg = c(0, 21), -#' scenarios = "EFSA_2015") -#' PEC_soil(1000, interval = 365, DT50 = 250, t_av = c(0, 21), -#' Kom = 1000, scenarios = "EFSA_2015", porewater = TRUE) -#' -#' # The following is from example 4 starting at p. 85 of the EFSA guidance (2015) -#' # Metabolite M2 -#' # Calculate total and porewater soil concentrations for tier 1 scenarios -#' # Relative molar mass is 100/300, formation fraction is 0.7 * 1 -#' results_pfm <- PEC_soil(100/300 * 0.7 * 1 * 1000, interval = 365, DT50 = 250, t_avg = c(0, 21), -#' scenarios = "EFSA_2015") -#' results_pfm_pw <- PEC_soil(100/300 * 0.7 * 1000, interval = 365, DT50 = 250, t_av = c(0, 21), -#' Kom = 100, scenarios = "EFSA_2015", porewater = TRUE) - -PEC_soil <- function(rate, rate_units = "g/ha", interception = 0, - mixing_depth = 5, - PEC_units = "mg/kg", PEC_pw_units = "mg/L", - interval = NA, n_periods = Inf, - tillage_depth = 20, - chent = NA, - DT50 = NA, - Koc = NA, Kom = Koc / 1.724, - t_avg = 0, - scenarios = c("default", "EFSA_2015"), - porewater = FALSE) -{ - rate_to_soil = (1 - interception) * rate - rate_units = match.arg(rate_units) - PEC_units = match.arg(PEC_units) - scenarios = match.arg(scenarios) - sce <- switch(scenarios, - default = data.frame(rho = 1.5, T_arr = NA, theta_fc = 0.2, f_oc = 0.02, - f_sce = 1, f_mod = 1, row.names = "default"), - EFSA_2015 = if (porewater) soil_scenario_data_EFSA_2015[4:6, ] - else soil_scenario_data_EFSA_2015[1:3, ] - ) - n_sce = nrow(sce) - - soil_volume = 100 * 100 * (mixing_depth/100) # in m3 - soil_mass = soil_volume * sce$rho * 1000 # in kg - - # The following is C_T,ini from EFSA 2012, p. 22, but potentially with interception > 0 - PEC_soil_ini = rate_to_soil * 1000 / soil_mass # in mg/kg - - # Decide which DT50 to take, or set degradation to zero if no DT50 available - if (is.na(DT50) & is(chent, "chent")) { - if (all(is.na(sce$T_arr))) { # No temperature correction - DT50 <- subset(chent$soil_degradation_endpoints, destination == "PECsoil")$DT50 - } else { - DT50 <- subset(chent$soil_degradation_endpoints, destination == "PECgw")$DT50 - } - if (length(DT50) > 1) stop("More than one PECsoil DT50 in chent object") - if (length(DT50) == 0) DT50 <- Inf - } - k = log(2)/DT50 - - # Temperature correction of degradation (accumulation) - if (all(is.na(sce$T_arr))) { # No temperature correction - f_T = 1 - } else { - # Temperature correction as in EFSA 2012 p. 23 - f_T = ifelse(sce$T_arr == 0, - 0, - exp(- (65.4 / 0.008314) * (1/(sce$T_arr + 273.15) - 1/293.15))) - } - - # X is the fraction left after one period (EFSA guidance p. 23) - X = exp(- k * f_T * interval) - - # f_accu is the fraction left after n periods (X + X^2 + ...) - f_accu = 0 - if (!is.na(interval)) { - if (n_periods == Inf) { - f_accu = X/(1 - X) - } else { - for (i in 1:n_periods) { - f_accu = f_accu + X^i - } - } - } - - f_tillage = mixing_depth / tillage_depth - - PEC_background = f_accu * f_tillage * PEC_soil_ini - - PEC_soil = (1 + f_accu * f_tillage) * PEC_soil_ini - - # Get porewater PEC if requested - if (porewater) { - - # If Kom is not specified, try to get K(f)oc - if (is.na(Kom)) { - # If Koc not specified, try to get K(f)oc from chent - if (is.na(Koc) & is(chent, "chent")) { - Koc <- soil_Kfoc(chent) - } - Kom <- Koc / 1.724 - } - - if (is.na(Kom)) stop("No Kom information specified") - - PEC_soil = PEC_soil/((sce$theta_fc/sce$rho) + sce$f_om * Kom) - } - - # Scenario adjustment factors - PEC_soil_sce = PEC_soil * sce$f_sce - - # Model adjustment factors - PEC_soil_sce_mod = PEC_soil_sce * sce$f_mod - - result <- matrix(NA, ncol = n_sce, nrow = length(t_avg), - dimnames = list(t_avg = t_avg, scenario = rownames(sce))) - - result[1, ] <- PEC_soil_sce_mod - - for (i in seq_along(t_avg)) { - t_av_i <- t_avg[i] - if (t_av_i > 0) { - # Equation 10 from p. 24 (EFSA 2015) - result[i, ] <- PEC_soil_sce_mod/(t_av_i * f_T * k) * (1 - exp(- f_T * k * t_av_i)) - } - } - - return(result) -} diff --git a/pkg/R/PEC_sw_drainage_UK.R b/pkg/R/PEC_sw_drainage_UK.R deleted file mode 100644 index e53f179..0000000 --- a/pkg/R/PEC_sw_drainage_UK.R +++ /dev/null @@ -1,66 +0,0 @@ -# Copyright (C) 2015 Johannes Ranke -# Contact: jranke@uni-bremen.de -# This file is part of the R package pfm - -# This program is free software: you can redistribute it and/or modify it under -# the terms of the GNU General Public License as published by the Free Software -# Foundation, either version 3 of the License, or (at your option) any later -# version. - -# This program is distributed in the hope that it will be useful, but WITHOUT -# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS -# FOR A PARTICULAR PURPOSE. See the GNU General Public License for more -# details. - -# You should have received a copy of the GNU General Public License along with -# this program. If not, see - -#' Calculate initial predicted environmental concentrations in surface water due to drainage using the UK method -#' -#' This implements the method specified in the UK data requirements handbook and was checked against the spreadsheet -#' published on the CRC website -#' -#' @param rate Application rate in g/ha -#' @param interception The fraction of the application rate that does not reach the soil -#' @param Koc The sorption coefficient normalised to organic carbon in L/kg -#' @param latest_application Latest application date, formatted as e.g. "01 July" -#' @param soil_DT50 Soil degradation half-life, if SFO kinetics are to be used -#' @param model The soil degradation model to be used. Either one of "FOMC", -#' "DFOP", "HS", or "IORE", or an mkinmod object -#' @param model_parms A named numeric vector containing the model parameters -#' @return The predicted concentration in surface water in µg/L -#' @export -#' @author Johannes Ranke -#' @examples -#' PEC_sw_drainage_UK(150, Koc = 100) -PEC_sw_drainage_UK <- function(rate, interception = 0, Koc, - latest_application = NULL, soil_DT50 = NULL, - model = NULL, model_parms = NULL) -{ - percentage_lost <- SSLRC_mobility_classification(Koc)[[2]] - amount_available <- rate * (1 - interception) # g/ha - - if (!missing(latest_application)) { - lct <- Sys.getlocale("LC_TIME") - tmp <- Sys.setlocale("LC_TIME", "C") - latest <- as.Date(paste(latest_application, "1999"), "%d %b %Y") - tmp <- Sys.setlocale("LC_TIME", lct) - degradation_time <- as.numeric(difftime(as.Date("1999-10-01"), units = "days", latest)) - if (!missing(soil_DT50)) { - k = log(2)/soil_DT50 - as.Date(paste(latest_application, "1999"), "%d %B %Y") - - amount_available <- amount_available * exp(-k * degradation_time) - if (!missing(model)) stop("You already supplied a soil_DT50 value, implying SFO kinetics") - } - if (!missing(model)) { - fraction_left <- pfm_degradation(model, parms = model_parms, - times = degradation_time)[1, "parent"] - amount_available <- fraction_left * amount_available - } - } - - volume = 130000 # L/ha - PEC = 1e6 * (percentage_lost/100) * amount_available / volume - return(PEC) -} diff --git a/pkg/R/PEC_sw_drift.R b/pkg/R/PEC_sw_drift.R deleted file mode 100644 index 261706c..0000000 --- a/pkg/R/PEC_sw_drift.R +++ /dev/null @@ -1,65 +0,0 @@ -# Copyright (C) 2015 Johannes Ranke -# Contact: jranke@uni-bremen.de -# This file is part of the R package pfm - -# This program is free software: you can redistribute it and/or modify it under -# the terms of the GNU General Public License as published by the Free Software -# Foundation, either version 3 of the License, or (at your option) any later -# version. - -# This program is distributed in the hope that it will be useful, but WITHOUT -# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS -# FOR A PARTICULAR PURPOSE. See the GNU General Public License for more -# details. - -# You should have received a copy of the GNU General Public License along with -# this program. If not, see - -#' Calculate predicted environmental concentrations in surface water due to drift -#' -#' This is a basic, vectorised form of a simple calculation of a contaminant -#' concentration in surface water based on complete, instantaneous mixing -#' with input via spray drift. -#' -#' @param rate Application rate in units specified below -#' @param applications Number of applications for selection of drift percentile -#' @param drift_percentages Percentage drift values for which to calculate PECsw. -#' 'drift_data' and 'distances' if not NULL. -#' @param drift_data Source of drift percentage data -#' @param crop Crop name (use German names for JKI data), defaults to "Ackerbau" -#' @param distances The distances in m for which to get PEC values -#' @param water_depth Depth of the water body in cm -#' @param rate_units Defaults to g/ha -#' @param PEC_units Requested units for the calculated PEC. Only µg/L currently supported -#' @return The predicted concentration in surface water -#' @export -#' @author Johannes Ranke -#' @examples -#' PEC_sw_drift(100) -PEC_sw_drift <- function(rate, - applications = 1, - water_depth = 30, - drift_percentages = NULL, - drift_data = "JKI", - crop = "Ackerbau", - distances = c(1, 5, 10, 20), - rate_units = "g/ha", - PEC_units = "\u00B5g/L") -{ - rate_units <- match.arg(rate_units) - PEC_units <- match.arg(PEC_units) - drift_data <- match.arg(drift_data) - water_volume <- 100 * 100 * (water_depth/100) * 1000 # in L (for 1 ha) - PEC_sw_overspray <- rate * 1e6 / water_volume # in µg/L - dist_index <- as.character(distances) - - if (is.null(drift_percentages)) { - drift_percentages <- pfm::drift_data_JKI[[applications]][dist_index, crop] - names(drift_percentages) <- paste(dist_index, "m") - } else { - names(drift_percentages) <- paste(drift_percentages, "%") - } - - PEC_sw_drift <- PEC_sw_overspray * drift_percentages / 100 - return(PEC_sw_drift) -} diff --git a/pkg/R/PEC_sw_sed.R b/pkg/R/PEC_sw_sed.R deleted file mode 100644 index a4a3d6a..0000000 --- a/pkg/R/PEC_sw_sed.R +++ /dev/null @@ -1,50 +0,0 @@ -# Copyright (C) 2015 Johannes Ranke -# Contact: jranke@uni-bremen.de -# This file is part of the R package pfm - -# This program is free software: you can redistribute it and/or modify it under -# the terms of the GNU General Public License as published by the Free Software -# Foundation, either version 3 of the License, or (at your option) any later -# version. - -# This program is distributed in the hope that it will be useful, but WITHOUT -# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS -# FOR A PARTICULAR PURPOSE. See the GNU General Public License for more -# details. - -# You should have received a copy of the GNU General Public License along with -# this program. If not, see - -#' Calculate predicted environmental concentrations in sediment from surface -#' water concentrations -#' -#' The method 'percentage' is equivalent to what is used in the CRD spreadsheet -#' PEC calculator -#' -#' @param PEC_sw Numeric vector or matrix of surface water concentrations in µg/L for -#' which the corresponding sediment concentration is to be estimated -#' @param percentage The percentage in sediment, used for the percentage method -#' @param method The method used for the calculation -#' @param sediment_depth Depth of the sediment layer -#' @param water_depth Depth of the water body in cm -#' @param sediment_density The density of the sediment in L/kg (equivalent to -#' g/cm3) -#' @param PEC_sed_units The units of the estimated sediment PEC value -#' @return The predicted concentration in sediment -#' @export -#' @author Johannes Ranke -#' @examples -#' PEC_sw_sed(PEC_sw_drift(100, distances = 1), percentage = 50) -PEC_sw_sed <- function(PEC_sw, percentage = 100, method = "percentage", - sediment_depth = 5, water_depth = 30, - sediment_density = 1.3, - PEC_sed_units = c("\u00B5g/kg", "mg/kg")) -{ - method = match.arg(method) - PEC_sed_units = match.arg(PEC_sed_units) - if (method == "percentage") { - PEC_sed = PEC_sw * (percentage/100) * (water_depth / sediment_depth) * (1 / sediment_density) - if (PEC_sed_units == "mg/kg") PEC_sed <- PEC_sed / 1000 - } - return(PEC_sed) -} diff --git a/pkg/R/SFO_actual_twa.R b/pkg/R/SFO_actual_twa.R deleted file mode 100644 index 7facb6a..0000000 --- a/pkg/R/SFO_actual_twa.R +++ /dev/null @@ -1,36 +0,0 @@ -# Copyright (C) 2015 Johannes Ranke -# Contact: jranke@uni-bremen.de -# This file is part of the R package pfm - -# This program is free software: you can redistribute it and/or modify it under -# the terms of the GNU General Public License as published by the Free Software -# Foundation, either version 3 of the License, or (at your option) any later -# version. - -# This program is distributed in the hope that it will be useful, but WITHOUT -# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS -# FOR A PARTICULAR PURPOSE. See the GNU General Public License for more -# details. - -# You should have received a copy of the GNU General Public License along with -# this program. If not, see - -#' Actual and maximum moving window time average concentrations for SFO kinetics -#' -#' @param DT50 The half-life. -#' @param times The output times, and window sizes for time weighted average concentrations -#' @export -#' @author Johannes Ranke -#' @source FOCUS (2014) Generic Guidance for Estimating Persistence and Degradation -#' Kinetics from Environmental Fate Studies on Pesticides in EU Registratin, Version 1.1, -#' 18 December 2014, p. 251 -#' @examples -#' SFO_actual_twa(10) -SFO_actual_twa <- function(DT50 = 1000, times = c(0, 1, 2, 4, 7, 14, 21, 28, 42, 50, 100)) -{ - k = log(2)/DT50 - result <- data.frame(actual = 1 * exp(-k * times), - twa = (1 - exp(-k * times))/(k * times), - row.names = times) - return(result) -} diff --git a/pkg/R/SSLRC_mobility_classification.R b/pkg/R/SSLRC_mobility_classification.R deleted file mode 100644 index deda5cf..0000000 --- a/pkg/R/SSLRC_mobility_classification.R +++ /dev/null @@ -1,42 +0,0 @@ -# Copyright (C) 2015 Johannes Ranke -# Contact: jranke@uni-bremen.de -# This file is part of the R package pfm - -# This program is free software: you can redistribute it and/or modify it under -# the terms of the GNU General Public License as published by the Free Software -# Foundation, either version 3 of the License, or (at your option) any later -# version. - -# This program is distributed in the hope that it will be useful, but WITHOUT -# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS -# FOR A PARTICULAR PURPOSE. See the GNU General Public License for more -# details. - -# You should have received a copy of the GNU General Public License along with -# this program. If not, see - -#' Determine the SSLRC mobility classification for a chemical substance from its Koc -#' -#' This implements the method specified in the UK data requirements handbook and was -#' checked against the spreadsheet published on the CRC website -#' -#' @param Koc The sorption coefficient normalised to organic carbon in L/kg -#' @return A list containing the classification and the percentage of the -#' compound transported per 10 mm drain water -#' @export -#' @author Johannes Ranke -#' @examples -#' SSLRC_mobility_classification(100) -SSLRC_mobility_classification <- function(Koc) -{ - if (!is.numeric(Koc) | length(Koc) != 1) stop("Please give a single number") - result <- list("Non mobile", 0.01) - if (Koc < 4000) result <- list("Slightly mobile", 0.02) - if (Koc < 1000) result <- list("Slightly mobile", 0.5) - if (Koc < 500) result <- list("Moderately mobile", 0.7) - if (Koc < 75) result <- list("Mobile", 1.9) - if (Koc < 15) result <- list("Very mobile", 1.9) - names(result) <- c("Mobility classification", - "Percentage drained per mm of drain water") - return(result) -} diff --git a/pkg/R/TOXSWA_cwa.R b/pkg/R/TOXSWA_cwa.R deleted file mode 100644 index c5ddce9..0000000 --- a/pkg/R/TOXSWA_cwa.R +++ /dev/null @@ -1,380 +0,0 @@ -# Copyright (C) 2014,2015,2016 Johannes Ranke -# Contact: jranke@uni-bremen.de -# This file is part of the R package pfm - -# This program is free software: you can redistribute it and/or modify it under -# the terms of the GNU General Public License as published by the Free Software -# Foundation, either version 3 of the License, or (at your option) any later -# version. - -# This program is distributed in the hope that it will be useful, but WITHOUT -# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS -# FOR A PARTICULAR PURPOSE. See the GNU General Public License for more -# details. - -# You should have received a copy of the GNU General Public License along with -# this program. If not, see - -#' Read TOXSWA surface water concentrations -#' -#' Read TOXSWA hourly concentrations of a chemical substance in a specific -#' segment of a TOXSWA surface water body. Per default, the data for the last -#' segment are imported. As TOXSWA 4 reports the values at the end of the hour -#' (ConLiqWatLayCur) in its summary file, we use this value as well instead -#' of the hourly averages (ConLiqWatLay). -#' -#' @param filename The filename of the cwa file (TOXSWA 2.x.y or similar) or the -#' out file (FOCUS TOXSWA 4, i.e. TOXSWA 4.4.2 or similar). -#' @param basedir The path to the directory where the cwa file resides. -#' @param zipfile Optional path to a zip file containing the cwa file. -#' @param segment The segment for which the data should be read. Either "last", or -#' the segment number. -#' @param total Set this to TRUE in order to read total concentrations as well. This is -#' only necessary for .out files as generated by TOXSWA 4.4.2 or similar, not for .cwa -#' files. For .cwa files, the total concentration is always read as well. -#' @param substance For TOXSWA 4 .out files, the default value "parent" leads -#' to reading concentrations of the parent compound. Alternatively, the substance -#' of interested can be selected by its code name. -#' @param windows Numeric vector of width of moving windows in days, for calculating -#' maximum time weighted average concentrations and areas under the curve. -#' @param thresholds Numeric vector of threshold concentrations in µg/L for -#' generating event statistics. -#' @importFrom readr read_fwf fwf_empty -#' @return An instance of an R6 object of class -#' \code{\link{TOXSWA_cwa}}. -#' @export -#' @author Johannes Ranke -#' @examples -#' H_sw_D4_pond <- read.TOXSWA_cwa("00001p_pa.cwa", -#' basedir = "SwashProjects/project_H_sw/TOXSWA", -#' zipfile = system.file("testdata/SwashProjects.zip", -#' package = "pfm")) -read.TOXSWA_cwa <- function(filename, basedir = ".", zipfile = NULL, - segment = "last", substance = "parent", - total = FALSE, - windows = NULL, thresholds = NULL) -{ - if (!missing(filename)) { - cwa <- TOXSWA_cwa$new(filename, basedir, zipfile, - substance = substance, total = total) - if (!is.null(windows[1])) cwa$moving_windows(windows) - if (!is.null(thresholds[1])) cwa$get_events(thresholds) - invisible(cwa) - } else { - message("You need to specify a filename for the cwa file to be read") - } -} - -#' Plot TOXSWA surface water concentrations -#' -#' Plot TOXSWA hourly concentrations of a chemical substance in a specific -#' segment of a TOXSWA surface water body. -#' -#' @import graphics -#' @param x The TOXSWA_cwa object to be plotted. -#' @param xlab,ylab Labels for x and y axis. -#' @param time_column What should be used for the time axis. If "t_firstjan" is chosen, -#' the time is given in days relative to the first of January in the first year. -#' @param add Should we add to an existing plot? -#' @param total Should the total concentration in water be plotted, including substance sorbed -#' to suspended matter? -#' @param LC_TIME Specification of the locale used to format dates -#' @param ... Further arguments passed to \code{plot} if we are not adding to an existing plot -#' @export -#' @author Johannes Ranke -#' @examples -#' H_sw_D4_pond <- read.TOXSWA_cwa("00001p_pa.cwa", -#' basedir = "SwashProjects/project_H_sw/TOXSWA", -#' zipfile = system.file("testdata/SwashProjects.zip", -#' package = "pfm")) -#' plot(H_sw_D4_pond) -plot.TOXSWA_cwa <- function(x, time_column = c("datetime", "t", "t_firstjan", "t_rel_to_max"), - xlab = "default", ylab = "default", - add = FALSE, - total = FALSE, LC_TIME = "C", ...) -{ - time_column = match.arg(time_column) - cwa_column = ifelse(total, "cwa_tot_mug_per_L", "cwa_mug_per_L") - lct <- Sys.getlocale("LC_TIME") - tmp <- Sys.setlocale("LC_TIME", LC_TIME) - if (xlab == "default") { - xlab = switch(time_column, - datetime = "Time", - t = "Time [days]", - t_firstjan = "Time since first of January [days]", - t_rel_to_max = "Time relative to maximum concentration [days]") - } - if (ylab == "default") { - ylab = paste( ifelse(total, "Total concentration", "Concentration"), "[\u03bcg/L]") - } - if (add) { - lines(x$cwas[c(time_column, cwa_column)], xlab = xlab, ylab = ylab, ...) - } else{ - if (time_column == "datetime") { - plot(x$cwas$datetime, x$cwas[[cwa_column]], type = "l", - xlab = xlab, ylab = ylab, xaxt = "n", ...) - seq_x <- seq(min(x$cwas$datetime), max(x$cwas$datetime), by = "quarter") - axis(1, at = seq_x, labels = format(seq_x, "%b %Y")) - } else { - plot(x$cwas[c(time_column, cwa_column)], type = "l", - xlab = xlab, ylab = ylab, ...) - } - } - tmp <- Sys.setlocale("LC_TIME", lct) -} - -#' R6 class for holding TOXSWA cwa concentration data and associated statistics -#' -#' An R6 class for holding TOXSWA cwa concentration data and some associated statistics. -#' Usually, an instance of this class will be generated by \code{\link{read.TOXSWA_cwa}}. -#' -#' @docType class -#' @importFrom R6 R6Class -#' @export -#' @format An \code{\link{R6Class}} generator object. -#' @field filename Length one character vector. -#' @field basedir Length one character vector. -#' @field segment Length one integer, specifying for which segment the cwa data were read. -#' @field cwas Dataframe holding the concentrations. -#' @field events List of dataframes holding the event statistics for each threshold. -#' @field windows Matrix of maximum time weighted average concentrations (TWAC_max) -#' and areas under the curve in µg/day * h (AUC_max_h) or µg/day * d (AUC_max_d) -#' for the requested moving window sizes in days. -#' @section Methods: -#' \describe{ -#' \item{\code{get_events(threshold, total = FALSE)}}{ -#' Populate a datataframe with event information for the specified threshold value -#' in µg/L. If \code{total = TRUE}, the total concentration including the amount -#' adsorbed to suspended matter will be used. The resulting dataframe is stored in the -#' \code{events} field of the object. -#' } -#' \item{\code{moving_windows(windows, total = FALSE)}}{ -#' Add to the \code{windows} field described above. -#' Again, if \code{total = TRUE}, the total concentration including the amount -#' adsorbed to suspended matter will be used. -#' } -#' } -#' @examples -#' H_sw_R1_stream <- read.TOXSWA_cwa("00003s_pa.cwa", -#' basedir = "SwashProjects/project_H_sw/TOXSWA", -#' zipfile = system.file("testdata/SwashProjects.zip", -#' package = "pfm")) -#' H_sw_R1_stream$get_events(c(2, 10)) -#' H_sw_R1_stream$moving_windows(c(7, 21)) -#' print(H_sw_R1_stream) -#' @keywords data - -TOXSWA_cwa <- R6Class("TOXSWA_cwa", - public = list( - filename = NULL, - basedir = NULL, - zipfile = NULL, - segment = NULL, - substance = NULL, - cwas = NULL, - windows = NULL, - events = list(), - initialize = function(filename, basedir, zipfile = NULL, - segment = "last", substance = "parent", total = FALSE) { - self$filename <- filename - self$basedir <- basedir - self$zipfile <- zipfile - if (!is.null(zipfile)) { - try(file_connection <- unz(zipfile, paste0(basedir, "/", filename))) - } else { - try(file_connection <- file(file.path(basedir, filename), "rt")) - } - - - if (grepl(".cwa$", filename)) { - # cwa file from FOCUS TOXSWA 3 (TOXSWA 2.x.y) - cwa_all_segments <- try( - read.table(file_connection, - sep = "", skip = 40, - encoding = "UTF-8", - colClasses = c("character", "numeric", - "integer", rep("numeric", 5)), - col.names = c("datetime", "t", "segment", - "xcd", "cwa_tot", "cwa", "Xss", "Xmp"))) - - if (is.null(zipfile)) close(file_connection) # only needed for files - - if (!inherits(cwa_all_segments, "try-error")) { - available_segments = 1:max(cwa_all_segments$segment) - if (segment == "last") segment = max(available_segments) - if (!segment %in% available_segments) stop("Invalid segment specified") - self$segment <- segment - cwa <- subset(cwa_all_segments, segment == self$segment, - c("datetime", "t", "segment", "cwa", "cwa_tot")) - lct <- Sys.getlocale("LC_TIME"); Sys.setlocale("LC_TIME", "C") - cwa$datetime <- strptime(cwa$datetime, "%d-%b-%Y-%H:%M") - Sys.setlocale("LC_TIME", lct) - startyear = format(cwa$datetime[1], "%Y") - firstjan <- strptime(paste0(startyear, "-01-01"), "%Y-%m-%d") - cwa$t_firstjan <- as.numeric(difftime(cwa$datetime, - firstjan, units = "days")) - - t_max = cwa[which.max(cwa$cwa), "t"] - cwa$t_rel_to_max = cwa$t - t_max - cwa$cwa_mug_per_L <- cwa$cwa * 1000 - cwa$cwa_tot_mug_per_L <- cwa$cwa_tot * 1000 - self$cwas <- cwa[c("datetime", "t", "t_firstjan", - "t_rel_to_max", - "cwa_mug_per_L", - "cwa_tot_mug_per_L")] - } else { - stop("Could not read ", filename) - } - } else { - # out file from FOCUS TOXSWA 4 (TOXSWA 4.4.2 or similar) - outfile <- try(readLines(file_connection)) - - close(file_connection) # only needed for files - - if (inherits(outfile, "try-error")) { - stop("Could not read ", filename) - } else { - # Get the substance name(s) - sub_lines <- grep(".*0.000.*ConLiqWatLayCur_", outfile[1:50], value = TRUE) - substances <- gsub(".*ConLiqWatLayCur_(.*?) *[0-9].*", "\\1", sub_lines) - if (!substance %in% c("parent", substances)) { - stop("No data for substance ", substance, " present in the .out file.") - } - - # Generate field name for the concentrations at end of hour for the - # substance of interest - if (substance == "parent") { - cwa_string = paste0("ConLiqWatLayCur_", substances[1]) - } else { - cwa_string = paste0("ConLiqWatLayCur_", substance) - } - - cwa_lines <- grep(cwa_string, outfile, value = TRUE) - - cwa_all_segments <- read_fwf(paste(cwa_lines, collapse = "\n"), - fwf_empty(paste(tail(cwa_lines), collapse = "\n"))) - - # Append time "-00h00" to datetime in first row, as this is not (always?) present - # in the line ConLiqWatLayCur - if (nchar(cwa_all_segments[1, "X2"]) == 11) { - cwa_all_segments[1, "X2"] = paste0(cwa_all_segments[1, "X2"], "-00h00") - } - - available_segments = 1:(ncol(cwa_all_segments) - 3) - if (segment == "last") segment = max(available_segments) - if (!segment %in% available_segments) stop("Invalid segment specified") - self$segment <- segment - cwa <- data.frame( - datetime = as.character(cwa_all_segments$X2), - t = cwa_all_segments$X1, - cwa = cwa_all_segments[[3 + segment]] - ) - if (total) { - cwa_tot_lines <- outfile[grep("ConSysWatLay_", outfile)] # hourly total conc. - cwa_tot_all_segments <- read.table(text = cwa_lines) - cwa$cwa_tot = cwa_tot_all_segments[[3 + segment]] - } - lct <- Sys.getlocale("LC_TIME"); Sys.setlocale("LC_TIME", "C") - cwa$datetime <- strptime(cwa$datetime, "%d-%b-%Y-%Hh%M") - Sys.setlocale("LC_TIME", lct) - - startyear = format(cwa$datetime[1], "%Y") - firstjan <- strptime(paste0(startyear, "-01-01"), "%Y-%m-%d") - cwa$t_firstjan <- as.numeric(difftime(cwa$datetime, - firstjan, units = "days")) - t_max = cwa[which.max(cwa$cwa), "t"] - cwa$t_rel_to_max = cwa$t - t_max - cwa$cwa_mug_per_L <- cwa$cwa * 1000 - self$cwas <- cwa[c("datetime", "t", "t_firstjan", - "t_rel_to_max", - "cwa_mug_per_L")] - if (total) { - self$cwas$cwa_tot_mug_per_L <- cwa$cwa_tot * 1000 - } - } - } - }, - moving_windows = function(windows, total = FALSE) { - window_names = paste(windows, "days") - n = length(window_names) - self$windows <- data.frame(window = window_names, - max_TWAC = numeric(n), - max_AUC_h = numeric(n), - max_AUC_d = numeric(n)) - if (missing(windows)) { - stop("You need to specify at least one moving window size in days") - } - cwa_column = ifelse(total, "cwa_tot_mug_per_L", "cwa_mug_per_L") - for (i in seq_along(windows)) { - window_size = windows[i] - filter_size = window_size * 24 - max_TWAC = max(filter(self$cwas[cwa_column], rep(1/filter_size, filter_size), - "convolution"), na.rm = TRUE) - max_AUC_h = max_TWAC * filter_size - max_AUC_d = max_TWAC * window_size - - self$windows[i, -1] = c(max_TWAC, max_AUC_h, max_AUC_d) - } - invisible(self) - }, - get_events = function(thresholds, total = FALSE) { - if (missing(thresholds)) { - stop("You need to specify at least one threshold concentration in \u03bcg/L") - } - for (threshold in thresholds) { - events = data.frame(t_start = numeric(), cwa_max = numeric(), - duration = numeric(), pre_interval = numeric(), - AUC_h = numeric(), AUC_d = numeric()) - cwa_column = ifelse(total, "cwa_tot_mug_per_L", "cwa_mug_per_L") - - event_end = 0 - event = FALSE - event_max = 0 - event_nr = 0 - n_rows = nrow(self$cwas) - for (i in 1:n_rows) { - cwa_cur = self$cwas[i, cwa_column] - if (event == FALSE) { - if (cwa_cur > threshold) { - event_start = self$cwas[i, "t"] - pre_interval = event_start - event_end - i_start = i - event = TRUE - event_max = cwa_cur - } - } else { - if (cwa_cur > event_max) event_max = cwa_cur - if (cwa_cur < threshold || i == n_rows) { - event_nr = event_nr + 1 - i_end = i - if (i == n_rows) i_end = i - event_end = self$cwas[i_end, "t"] - event_length = event_end - event_start - event_cwas <- self$cwas[i_start:(i_end - 1), cwa_column] - event_AUC_h = sum(event_cwas) - event_AUC_d = event_AUC_h / 24 - events[event_nr, ] = c(event_start, event_max, event_length, - pre_interval, event_AUC_h, event_AUC_d) - event = FALSE - } - } - } - - self$events[[as.character(threshold)]] <- events - } - invisible(self) - }, - print = function() { - cat(" data from file", self$filename, "segment", self$segment, "\n") - print(head(self$cwas)) - cat("Moving window analysis\n") - print(self$windows) - for (threshold in names(self$events)) { - cat("Event statistics for threshold", threshold, "\n") - if (nrow(self$events[[threshold]]) == 0) cat("No events found\n") - else print(self$events[[threshold]]) - } - } - ) -) -# vim: set ts=2 sw=2 expandtab: diff --git a/pkg/R/drift_data_JKI.R b/pkg/R/drift_data_JKI.R deleted file mode 100644 index 888d3b5..0000000 --- a/pkg/R/drift_data_JKI.R +++ /dev/null @@ -1,44 +0,0 @@ -#' Deposition from spray drift expressed as percent of the applied dose as -#' published by the JKI -#' -#' Deposition from spray drift expressed as percent of the applied dose as -#' published by the German Julius-Kühn Institute (JKI). -#' -#' The data were extracted from the spreadsheet cited below using the R code -#' given in the example section. The spreadsheet is not included in the package -#' as its licence is not clear. -#' -#' -#' @name drift_data_JKI -#' @docType data -#' @format A list currently containing matrices with spray drift percentage -#' data for field crops (Ackerbau), and Pome/stone fruit, early and late -#' (Obstbau frueh, spaet). -#' @source JKI (2010) Spreadsheet 'Tabelle der Abdrifteckwerte.xls', retrieved -#' from -#' http://www.jki.bund.de/no_cache/de/startseite/institute/anwendungstechnik/abdrift-eckwerte.html -#' on 2015-06-11 -#' @keywords datasets -#' @examples -#' -#' \dontrun{ -#' # This is the code that was used to extract the data -#' library(readxl) -#' abdrift_path <- "inst/extdata/Tabelle der Abdrifteckwerte.xls" -#' JKI_crops <- c("Ackerbau", "Obstbau frueh", "Obstbau spaet") -#' names(JKI_crops) <- c("Field crops", "Pome/stone fruit, early", "Pome/stone fruit, late") -#' drift_data_JKI <- list() -#' -#' for (n in 1:8) { -#' drift_data_raw <- read_excel(abdrift_path, sheet = n + 1, skip = 2) -#' drift_data <- as.matrix(drift_data_raw[1:9, 2:4]) -#' dimnames(drift_data) <- list(distance = as.integer(drift_data_raw[1:9, 1]), -#' crop = JKI_crops) -#' drift_data_JKI[[n]] <- drift_data -#' } -#' save(drift_data_JKI, file = "data/drift_data_JKI.RData") -#' } -#' -#' # And this is the resulting data -#' drift_data_JKI -NULL diff --git a/pkg/R/endpoint.R b/pkg/R/endpoint.R deleted file mode 100644 index 314acd2..0000000 --- a/pkg/R/endpoint.R +++ /dev/null @@ -1,110 +0,0 @@ -#' Retrieve endpoint information from the chyaml field of a chent object -#' -#' R6 class objects of class \code{\link{chent}} represent chemical entities -#' and can hold a list of information loaded from a chemical yaml file in their -#' chyaml field. Such information is extracted and optionally aggregated by -#' this function. -#' -#' The functions \code{soil_*} are functions to extract soil specific endpoints. -#' For the Freundlich exponent, the capital letter \code{N} is used in order to -#' facilitate dealing with such data in R. In pesticide fate modelling, this -#' exponent is often called 1/n. -#' -#' @import chents -#' @export -#' @param chent The \code{\link{chent}} object to get the information from -#' @param medium The medium for which information is sought -#' @param type The information type -#' @param lab_field If not NA, do we want laboratory or field endpoints -#' @param redox If not NA, are we looking for aerobic or anaerobic data -#' @param value The name of the value we want. The list given in the -#' usage section is not exclusive -#' @param aggregator The aggregator function. Can be mean, -#' \code{\link{geomean}}, or identity, for example. -#' @param raw Should the number(s) be returned as stored in the chent -#' object (could be a character value) to retain original information -#' about precision? -#' @param signif How many significant digits do we want -#' @return The result from applying the aggregator function to -#' the values converted to a numeric vector, rounded to the -#' given number of significant digits, or, if raw = TRUE, -#' the values as a character value, retaining any implicit -#' information on precision that may be present. -#' -endpoint <- function(chent, - medium = "soil", - type = c("degradation", "sorption"), - lab_field = c(NA, "laboratory", "field"), - redox = c(NA, "aerobic", "anaerobic"), - value = c("DT50ref", "Kfoc", "N"), - aggregator = geomean, - raw = FALSE, - signif = 3) -{ - if (!is(chent, "chent")) { - stop("Please supply a chent object as created using the package 'chents' available from jrwb.de") - } - ep_list <- chent$chyaml[[medium]][[type]] - if (!is.na(lab_field[1])) { - ep_list <- ep_list[[lab_field]] - } - if (!is.na(redox[1])) { - ep_list <- ep_list[[redox]] - } - values <- ep_list$data[[value]] - if (raw) return(values) - else return(signif(aggregator(as.numeric(values)), signif)) -} - -#' @inheritParams endpoint -#' @rdname endpoint -#' @export -soil_DT50 <- function(chent, aggregator = geomean, signif = 3, - lab_field = "laboratory", value = "DT50ref", - redox = "aerobic", raw = FALSE) { - ep <- endpoint(chent, medium = "soil", type = "degradation", - lab_field = "laboratory", redox = redox, - signif = signif, - value = value, aggregator = aggregator, raw = raw) - return(ep) -} - -#' @inheritParams endpoint -#' @rdname endpoint -#' @export -soil_Kfoc <- function(chent, aggregator = geomean, signif = 3, - value = "Kfoc", raw = FALSE) { - ep <- endpoint(chent, medium = "soil", type = "sorption", - signif = signif, - value = value, aggregator = aggregator, raw = raw) - return(ep) -} - -#' @inheritParams endpoint -#' @rdname endpoint -#' @export -soil_N <- function(chent, aggregator = mean, signif = 3, raw = FALSE) { - ep <- endpoint(chent, medium = "soil", type = "sorption", - signif = signif, - value = "N", aggregator = aggregator, raw = raw) - return(ep) -} - -#' @inheritParams endpoint -#' @rdname endpoint -#' @param values The values to be returned -#' @param aggregators A named vector of aggregator functions to be used -#' @export -soil_sorption <- function(chent, values = c("Kfoc", "N"), - aggregators = c(Kfoc = geomean, Koc = geomean, N = mean), - signif = c(Kfoc = 3, N = 3), - raw = FALSE) { - res <- sapply(values, - function(x) { - endpoint(chent, medium = "soil", type = "sorption", - signif = signif[[x]], - value = x, aggregator = aggregators[[x]], raw = raw) - } - ) - return(res) -} diff --git a/pkg/R/geomean.R b/pkg/R/geomean.R deleted file mode 100644 index 626829b..0000000 --- a/pkg/R/geomean.R +++ /dev/null @@ -1,38 +0,0 @@ -# Copyright (C) 2015 Johannes Ranke -# Contact: jranke@uni-bremen.de -# This file is part of the R package pfm - -# This program is free software: you can redistribute it and/or modify it under -# the terms of the GNU General Public License as published by the Free Software -# Foundation, either version 3 of the License, or (at your option) any later -# version. - -# This program is distributed in the hope that it will be useful, but WITHOUT -# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS -# FOR A PARTICULAR PURPOSE. See the GNU General Public License for more -# details. - -# You should have received a copy of the GNU General Public License along with -# this program. If not, see - -#' Calculate the geometric mean -#' -#' Based on some posts in a thread on Stackoverflow -#' \url{http://stackoverflow.com/questions/2602583/geometric-mean-is-there-a-built-in} -#' This function checks for negative values, removes NA values per default and -#' returns 0 if at least one element of the vector is 0. -#' -#' @param x Vector of numbers -#' @param na.rm Should NA values be omitted? -#' @return The geometric mean -#' @export -#' @author Johannes Ranke -#' @examples -#' geomean(c(1, 3, 9)) -#' geomean(c(1, 3, NA, 9)) -#' \dontrun{geomean(c(1, -3, 9)) # returns an error} -geomean = function(x, na.rm = TRUE){ - if (any(is.na(x)) & na.rm == FALSE) stop("Removal of NA values was not requested") - if (any(x < 0, na.rm = na.rm)) stop("Only defined for positive numbers") - exp(mean(log(x), na.rm = na.rm)) -} diff --git a/pkg/R/pfm_degradation.R b/pkg/R/pfm_degradation.R deleted file mode 100644 index 6c8610e..0000000 --- a/pkg/R/pfm_degradation.R +++ /dev/null @@ -1,48 +0,0 @@ -# Copyright (C) 2015 Johannes Ranke -# Contact: jranke@uni-bremen.de -# This file is part of the R package pfm - -# This program is free software: you can redistribute it and/or modify it under -# the terms of the GNU General Public License as published by the Free Software -# Foundation, either version 3 of the License, or (at your option) any later -# version. - -# This program is distributed in the hope that it will be useful, but WITHOUT -# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS -# FOR A PARTICULAR PURPOSE. See the GNU General Public License for more -# details. - -# You should have received a copy of the GNU General Public License along with -# this program. If not, see - -#' Calculate a time course of relative concentrations based on an mkinmod model -#' -#' @import mkin -#' @param model The degradation model to be used. Either a parent only model like -#' 'SFO' or 'FOMC', or an mkinmod object -#' @param DT50 The half-life. This is only used when simple exponential decline -#' is calculated (SFO model). -#' @param parms The parameters used for the degradation model -#' @param years For how many years should the degradation be predicted? -#' @param step_days What step size in days should the output have? -#' @param times The output times -#' @export -#' @author Johannes Ranke -#' @examples -#' head(pfm_degradation("SFO", DT50 = 10)) -pfm_degradation <- function(model = "SFO", DT50 = 1000, parms = c(k_parent_sink = log(2)/DT50), - years = 1, step_days = 1, - times = seq(0, years * 365, by = step_days)) -{ - if (model %in% c("SFO", "FOMC", "DFOP", "HS", "IORE")) { - model <- mkinmod(parent = list(type = model)) - } - initial_state = c(1, rep(0, length(model$diffs) - 1)) - names(initial_state) <- names(model$diffs) - time_course <- mkinpredict(model, odeparms = parms, - odeini = initial_state, - outtimes = times, - solution_type = ifelse(length(model$spec) == 1, - "analytical", "deSolve")) - invisible(time_course) -} diff --git a/pkg/R/soil_scenario_data_EFSA_2015.R b/pkg/R/soil_scenario_data_EFSA_2015.R deleted file mode 100644 index fb096bc..0000000 --- a/pkg/R/soil_scenario_data_EFSA_2015.R +++ /dev/null @@ -1,40 +0,0 @@ -#' Properties of the predefined scenarios from the EFSA guidance from 2015 -#' -#' Properties of the predefined scenarios used at Tier 1, Tier 2A and Tier 3A for the -#' concentration in soil as given in the EFSA guidance (2015, p. 13/14). Also, the -#' scenario and model adjustment factors from p. 15 and p. 17 are included. -#' -#' @name soil_scenario_data_EFSA_2015 -#' @docType data -#' @format A data frame with one row for each scenario. Row names are the scenario codes, -#' e.g. CTN for the Northern scenario for the total concentration in soil. Columns are -#' mostly self-explanatory. \code{rho} is the dry bulk density of the top soil. -#' @source EFSA (European Food Safety Authority) (2015) -#' EFSA guidance document for predicting environmental concentrations -#' of active substances of plant protection products and transformation products of these -#' active substances in soil. \emph{EFSA Journal} \bold{13}(4) 4093 -#' doi:10.2903/j.efsa.2015.4093 -#' @keywords datasets -#' @examples -#' \dontrun{ -#' # This is the code that was used to define the data -#' soil_scenario_data_EFSA_2015 <- data.frame( -#' Zone = rep(c("North", "Central", "South"), 2), -#' Country = c("Estonia", "Germany", "France", "Denmark", "Czech Republik", "Spain"), -#' T_arit = c(4.7, 8.0, 11.0, 8.2, 9.1, 12.8), -#' T_arr = c(7.0, 10.1, 12.3, 9.8, 11.2, 14.7), -#' Texture = c("Coarse", "Coarse", "Medium fine", "Medium", "Medium", "Medium"), -#' f_om = c(0.118, 0.086, 0.048, 0.023, 0.018, 0.011), -#' theta_fc = c(0.244, 0.244, 0.385, 0.347, 0.347, 0.347), -#' rho = c(0.95, 1.05, 1.22, 1.39, 1.43, 1.51), -#' f_sce = c(3, 2, 2, 2, 1.5, 1.5), -#' f_mod = c(2, 2, 2, 4, 4, 4), -#' stringsAsFactors = FALSE, -#' row.names = c("CTN", "CTC", "CTS", "CLN", "CLC", "CLS") -#' ) -#' save(soil_scenario_data_EFSA_2015, file = '../data/soil_scenario_data_EFSA_2015.RData') -#' } -#' -#' # And this is the resulting dataframe -#' soil_scenario_data_EFSA_2015 -NULL -- cgit v1.2.1