1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
|
#' Calculate FOCUS Step 1 PEC surface water
#'
#' This is reimplementation of Step 1 of the FOCUS Step 1 and 2 calculator
#' version 3.2, authored by Michael Klein. Note that results for multiple
#' applications should be compared to the corresponding results for a
#' single application. At current, this is not done automatically in
#' this implementation.
#'
#' @importFrom utils read.table
#' @references FOCUS (2014) Generic guidance for Surface Water Scenarios (version 1.4).
#' FOrum for the Co-ordination of pesticde fate models and their USe.
#' http://esdac.jrc.ec.europa.eu/public_path/projects_data/focus/sw/docs/Generic%20FOCUS_SWS_vc1.4.pdf
#' @references Website of the Steps 1 and 2 calculator at the Joint Research
#' Center of the European Union:
#' http://esdac.jrc.ec.europa.eu/projects/stepsonetwo
#' @note The formulas for input to the waterbody via runoff/drainage of the
#' parent and subsequent formation of the metabolite in water is not
#' documented in the model description coming with the calculator. As one would
#' expect, this appears to be (as we get the same results) calculated by
#' multiplying the application rate with the molar weight
#' correction and the formation fraction in water/sediment systems.
#' @note Step 2 is not implemented.
#' @export
#' @param parent A list containing substance specific parameters, e.g.
#' conveniently generated by \code{\link{chent_focus_sw}}.
#' @param rate The application rate in g/ha. Overriden when
#' applications are given explicitly
#' @param n The number of applications
#' @param i The application interval
#' @param comment A comment for the input file
#' @param met A list containing metabolite specific parameters. e.g.
#' conveniently generated by \code{\link{chent_focus_sw}}. If not NULL,
#' the PEC is calculated for this compound, not the parent.
#' @param f_drift The fraction of the application rate reaching the waterbody
#' via drift. If NA, this is derived from the scenario name and the number
#' of applications via the drift data defined by the
#' \code{\link{FOCUS_Step_12_scenarios}}
#' @param f_rd The fraction of the amount applied reaching the waterbody via
#' runoff/drainage. At Step 1, it is assumed to be 10%, be it the
#' parent or a metabolite
#' @param scenario The name of the scenario. Must be one of the scenario
#' names given in \code{\link{FOCUS_Step_12_scenarios}}
#' @param txt_file the name, and potentially the full path to the
#' Steps.12 input text file to which the specification of the run(s)
#' should be written
#' @param overwrite Should an existing file a the location specified in
#' \code{txt_file} be overwritten?
#' @param append Should the input text file be appended?
#' @examples
#' # Parent only
#' dummy_1 <- chent_focus_sw(cwsat = 6000, DT50_ws = 6, Koc = 344.8)
#' PEC_sw_focus(dummy_1, 3000, f_drift = 0)
#'
#' # Metabolite
#' new_dummy <- chent_focus_sw(mw = 250, Koc = 100)
#' M1 <- chent_focus_sw(mw = 100, cwsat = 100, DT50_ws = 100, Koc = 50, max_ws = 0, max_soil = 0.5)
#' PEC_sw_focus(new_dummy, 1000, scenario = "cereals, winter", met = M1)
PEC_sw_focus <- function(parent, rate, n = 1, i = NA,
comment = "",
met = NULL,
f_drift = NA, f_rd = 0.1,
scenario = FOCUS_Step_12_scenarios$names,
txt_file = "pesticide.txt", overwrite = FALSE, append = TRUE)
{
if (n > 1 & is.na(i)) stop("Please specify the interval i if n > 1")
if (is.na(f_drift)) {
scenario = match.arg(scenario)
f_drift = FOCUS_Step_12_scenarios$drift[scenario, "1"] / 100
# For Step 2 we would/will select the reduced percentiles for multiple apps:
# if (n <= 8) {
# f_drift = FOCUS_Step_12_scenarios$drift[scenario, as.character(n)] / 100
# } else {
# f_drift = FOCUS_Step_12_scenarios$drift[scenario, ">8"] / 100
# }
}
# Write to txt file if requested
add_line <- function(x) cat(paste0(x, "\r\n"), file = txt, append = TRUE)
add <- function(x) cat(paste(x, "\t"), file = txt, append = TRUE)
if (file.exists(txt_file)) {
if (append) {
txt <- file(txt_file, "a")
} else {
if (overwrite) {
txt <- file(txt_file, "w")
} else {
stop("The file", txt_file, "already exists, and you did not request",
"appending or overwriting it")
}
}
}
on.exit(close(txt))
# Write header to txt file
header <- c("Active Substance", "Compound", "Comment", "Mol mass a.i.",
"Mol mass met.", "Water solubility", "KOC assessed compound",
"KOC parent compound", "DT50", "Max. in Water",
"Max. in Soil asessed compound", # we reproduce the typo...
"App. Rate", "Number of App.", "Time between app.",
"App. Type", "DT50 soil parent compound", "DT50 soil",
"DT50 water", "DT50 sediment", "Region / Season",
"Interception class")
add_line(paste(header, collapse = "\t"))
if (is.null(met)) {
compound = parent$name
cwsat = parent$cwsat
mw_ratio = 1
max_soil = 1
max_ws = 1
Koc = parent$Koc
DT50_ws = parent$DT50_ws
} else {
compound = met$name
cwsat = met$cwsat
mw_ratio = met$mw / parent$mw
max_soil = met$max_soil
max_ws = met$max_ws
Koc = met$Koc
DT50_ws = met$DT50_ws
}
add(parent$name)
add(compound)
add(comment)
if (is.na(parent$mw)) add("-99.00")
else add(sprintf("%.2f", parent$mw))
if (is.na(met$mw)) add("-99.00")
else add(sprintf("%.2f", met$mw))
add(sprintf("%.2f", cwsat))
add(sprintf("%.2f", Koc))
if (is.null(met)) add("0.00E+00")
else add(sprintf("%.2f", parent$Koc))
# Rates for a single application
eq_rate_drift_s = mw_ratio * max_ws * rate
# Parent only, or metabolite formed in soil:
eq_rate_rd_s = mw_ratio * max_soil * rate
# Metabolite formed in water (this part is not documented in the Help files
# of the Steps 1/2 calculator):
if (!is.null(met)) {
eq_rate_rd_parent_s = mw_ratio * max_ws * rate
eq_rate_rd_s_tot = eq_rate_rd_s + eq_rate_rd_parent_s
} else {
eq_rate_rd_parent_s = NA
eq_rate_rd_s_tot = eq_rate_rd_s
}
# Drift input
input_drift_s = f_drift * eq_rate_drift_s / 10 # mg/m2
# Runoff/drainage input
ratio_field_wb = 10 # 10 m2 of field for each m2 of the waterbody
input_rd_s = f_rd * eq_rate_rd_s_tot * ratio_field_wb / 10
input_rd = n * input_rd_s
# No accumulation between multiple applications if 3 * DT50 < i
if (n > 1 && (3 * DT50_ws < i)) {
input_drift = input_drift_s
input_rd = input_rd_s
} else {
input_drift = n * input_drift_s
input_rd = n * input_rd_s
}
# Fraction of compound entering the water phase via runoff/drainage
depth_sw = 0.3 # m
depth_sed = 0.05 # m
depth_sed_eff = 0.01 # m, only relevant for sorption
rho_sed = 0.8 # kg/L
f_OC = 0.05 # 5% organic carbon in sediment
f_rd_sw = depth_sw / (depth_sw + (depth_sed_eff * rho_sed * f_OC * Koc))
f_rd_sed = 1 - f_rd_sw
# Initial PECs
PEC_sw_0 = (input_drift + input_rd * f_rd_sw) / depth_sw # µg/L
PEC_sed_0 = (input_rd * f_rd_sed) / (depth_sed * rho_sed) # µg/kg
# Initial PECs when assuming partitioning also of drift input
PEC_sw_0_part = (input_drift + input_rd) * f_rd_sw / depth_sw # µg/L
PEC_sed_0_part = (input_drift + input_rd) * f_rd_sed / (depth_sed * rho_sed) # µg/kg
t_out = c(0, 1, 2, 4, 7, 14, 21, 28, 42, 50, 100)
PEC = matrix(NA, nrow = length(t_out), ncol = 4,
dimnames = list(Time = t_out, type = c("PECsw", "TWAECsw", "PECsed", "TWAECsed")))
PEC["0", "PECsw"] = PEC_sw_0
PEC["0", "PECsed"] = PEC_sed_0
# Degradation after partitioning according to Koc
k_ws = log(2)/DT50_ws
PEC[as.character(t_out[-1]), "PECsw"] = PEC_sw_0_part * exp( - k_ws * t_out[-1])
PEC[as.character(t_out[-1]), "PECsed"] = PEC_sed_0_part * exp( - k_ws * t_out[-1])
# TWA concentrations
PEC_sw_1 = PEC["1", "PECsw"]
PEC_sed_1 = PEC["1", "PECsed"]
TWAEC_sw_1 = (PEC_sw_0 + PEC_sw_1) / 2
TWAEC_sed_1 = (PEC_sed_0 + PEC_sed_1) / 2
TWAEC_after = function(TWAEC_1, PEC_1, t) {
(TWAEC_1 + PEC_1 * (1 - exp(- k_ws * (t - 1))) / k_ws) / t
}
PEC["1", "TWAECsw"] = TWAEC_sw_1
PEC["1", "TWAECsed"] = TWAEC_sed_1
PEC[as.character(t_out[-1]), "TWAECsw"] = TWAEC_after(TWAEC_sw_1, PEC_sw_1, t_out[-1])
PEC[as.character(t_out[-1]), "TWAECsed"] = TWAEC_after(TWAEC_sed_1, PEC_sed_1, t_out[-1])
# Check if PEC_sw_max is above water solubility
PEC_sw_max = max(PEC[, "PECsw"])
if (PEC_sw_max > 1000 * cwsat) {
warning("The maximum PEC surface water exceeds the water solubility")
}
PEC_sed_max = max(PEC[, "PECsed"])
list(f_drift = f_drift,
eq_rate_drift_s = eq_rate_drift_s,
eq_rate_rd_s = eq_rate_rd_s,
eq_rate_rd_parent_s = eq_rate_rd_parent_s,
input_drift_s = input_drift_s,
input_rd_s = input_rd_s,
f_rd_sw = f_rd_sw, f_rd_sed = f_rd_sed,
PEC = PEC,
PEC_sw_max = PEC_sw_max,
PEC_sed_max = PEC_sed_max
)
}
#' Create a chemical compound object for FOCUS Step 1 calculations
#'
#' @export
#' @param name Length one character vector containing the name
#' @param cwsat Water solubility in mg/L
#' @param DT50_ws Half-life in water/sediment systems in days
#' @param Koc Partition coefficient between organic carbon and water
#' in L/kg.
#' @param mw Molar weight in g/mol.
#' @param max_soil Maximum observed fraction (dimensionless) in soil
#' @param max_ws Maximum observed fraction (dimensionless) in water/sediment
#' systems
#' @return A list with the substance specific properties
chent_focus_sw <- function(name, Koc, DT50_ws = NA, cwsat = 1000, mw = NA, max_soil = 1, max_ws = 1)
{
list(Koc = Koc, DT50_ws = DT50_ws, cwsat = cwsat,
mw = mw, max_soil = max_soil, max_ws = max_ws)
}
|