Utils.R 34.1 KB
Newer Older
1
.onAttach <- function(libname, pkgname) {
2
  if (packageVersion("htmlwidgets") <= "1.5.2") {
3
    base::packageStartupMessage("\n---------------------------\n")
4
    base::packageStartupMessage("This version of 'airGRteaching' is designed to work with 'htmlwidgets' >= 1.5.2.9000 (troubles with 'dygraphs')")
5
    base::packageStartupMessage("Install the latest version of 'htmlwidgets' from GitHub with the following command lines:")
6
7
8
9
10
    base::packageStartupMessage("\tinstall.packages(\"remotes\")\n\tremotes::install_github(\"ramnathv/htmlwidgets\")")
    base::packageStartupMessage("\n---------------------------\n")
  }
}

11
12
13



14
## =================================================================================
15
## commands to avoid warnings during package checking when global variables are used
16
## =================================================================================
17
18

if (getRversion() >= "2.15.1") {
19
20
  utils::globalVariables(c(".ShinyGR.args"))
  utils::suppressForeignCheck(c(".ShinyGR.args"))
21
22
  utils::globalVariables(c(".ShinyGR.hist"))
  utils::suppressForeignCheck(c(".ShinyGR.hist"))
23
24
}

25
26


27

28
29
30
31
32
33
## =================================================================================
## function to test if a remote url exists
## =================================================================================

.CheckUrl <- function(url, timeout = 2) {
  con <- url(description = url)
34
35
  check <- suppressWarnings(try(open.connection(con = con, open = "rt", timeout = timeout),
                                silent = TRUE)[1])
36
37
38
39
40
41
42
  suppressWarnings(try(close.connection(con), silent = TRUE))
  is.null(check)
}




43
44
45
46
## =================================================================================
## function to compute the start and stop id of equal values in a vector
## =================================================================================

47
.StartStop <- function(x, FUN) {
48
49
50
  naQ_rle <- rle(FUN(x))
  naQ_ide <- cumsum(naQ_rle$lengths)[naQ_rle$values]   + 1
  naQ_ids <- naQ_ide - naQ_rle$lengths[naQ_rle$values] - 1
51
52
  idNA <- data.frame(start = naQ_ids, stop = naQ_ide)
  idNA$start <- ifelse(idNA$start < 1        , 1        , idNA$start)
53
54
55
56
57
58
59
60
61
62
63
  idNA$stop  <- ifelse(idNA$stop  > length(x), length(x), idNA$stop )
  idNA
}




## =================================================================================
## function for drawing several shadows of dygraphic regions simultaneously
## =================================================================================

64
.DyShadingMulti <- function(dygraph, ts, idStart, IdStop, ...) {
65
66
67
68
69
70
71
72
73
74
  for (i in seq_along(idStart)) {
    dygraph <- dygraphs::dyShading(dygraph = dygraph,
                                   from    = as.character(ts)[idStart[i]],
                                   to      = as.character(ts)[IdStop[i]],
                                   ...)
  }
  dygraph
}


75
## =================================================================================
76
## function to manage the model units
77
78
## =================================================================================

79
.TypeModelGR <- function(x) {
80

81
82
  if (!is.list(x)) {
    x <- list(TypeModel = x)
83
  }
84
  if (any(class(x) %in% c("PrepGR", "CalGR", "SimGR")) || names(x) %in% "TypeModel") {
85
    x <- x$TypeModel
86
  }
87

88
  StrName    <- "(.*)(GR)(\\d{1})(\\D{1})"
89
  NameModel  <- gsub(StrName, "\\2\\3\\4", x)
90
  TimeUnitFR <- gsub(StrName, "\\4", x)
91
92
93
94
  # TimeUnit   <- switch(TimeUnitFR, H = "hour", J = "day", M = "month", A = "year")
  # TimeLag    <- switch(TimeUnit, "hour" = 3600, "day" = 3600*24, "month" = 3600*24*31, "year" = 366)
  TimeUnit   <- switch(TimeUnitFR, H = "h", J = "d", M = "month", A = "y")
  TimeLag    <- switch(TimeUnit, "h" = 3600, "d" = 3600*24, "month" = 3600*24*31, "y" = 366)
95
96
  NbParam    <- gsub(StrName, "\\3", x)
  isCN       <- grepl("CemaNeige"  , x)
97
98

  res <- list(TypeModel = x, NameModel = NameModel, CemaNeige = isCN,
99
100
              NbParam = as.numeric(NbParam),
              TimeUnit = TimeUnit, TimeLag = TimeLag)
101
  return(res)
102
103
104
}


105
106


107
## =================================================================================
108
## function to plot the gr models diagrams (only GR4J and GR5J)
109
## =================================================================================
110

111
.DiagramGR <- function(OutputsModel, Param, SimPer, EventDate, HydroModel, CemaNeige, Theme = NULL) {
112
113


114
  # --------------------------------------------------------------------------------
115
  # PARAMETRES
116
  # --------------------------------------------------------------------------------
117

118
  # Parametres
119
  mgp             <- c(0, 0.75, 0)
120
  col_P           <- rgb(066, 139, 202, maxColorValue = 255) #"royalblue"
121
  col_E           <- rgb(164, 196, 000, maxColorValue = 255) #"forestgreen"
122
  col_Q           <- "orangered"
123
124
  col_SP          <- adjustcolor("cyan4"   , alpha.f = 0.60)
  col_SR          <- adjustcolor("darkblue", alpha.f = 0.60)
125
126
127
  col_R           <- rgb(066, 139, 202, maxColorValue = 255) #rgb(037, 155, 210, maxColorValue = 255)
  col_mod_bg      <- rgb(245, 245, 245, maxColorValue = 255)
  col_mod_bd      <- rgb(231, 231, 231, maxColorValue = 255)
128
  xy_E            <- c(250, 980)
129
130
  xy_PE           <- c(250, 940)
  xy_AE           <- c(250, 860)
131
  xy_P            <- c(600, 980)
132
  xy_Precip       <- c(600, 950)
133
  xy_Q            <- c(700,  30)
134
135
136
137
138
  x_Ps            <- 440
  x_PnPs          <- 700
  y_interception  <- 900
  y_rendement     <- 815
  y_percolation   <- 575
139
140
  xy_Q9           <- c(400, 310)
  xy_Q1           <- c(800, 310)
141
  y_routage       <- 100
142
  fact_res        <- 200/800/3
143
  fact_resExp     <- 1
144
145
146
147
  base_res        <- 300
  NH              <- 10
  xy_min_PROD     <- c(200, 610)
  xy_min_ROUT     <- c(250, 150)
148
  xy_min_EXPO     <- c(200, 250)
149
  y_entreeUH      <- 500
150
151
152
153
  xy_UH1          <- c(500, 420)
  xy_UH2          <- c(900, 420)
  y_Ech_Q1        <- 170 #200
  y_Ech_Q9        <- 150 #180
154
155
156
  D               <- 5/2
  xpad            <- 1.5
  ypad            <- 1.5
157
  max_triangle    <- max(unlist(OutputsModel[c("Perc", "PR", "Q9", "Q1", "QR", "QD", "Pn", "Ps",
158
                                               "AE", "Precip", "PotEvap", "AExch1", "AExch2")]))
159
  fact_var        <- 40
160
  fact_triangle   <- 100#25
161
162
  cex_max_poly    <- 2 # 0.005
  cex_tri         <- function(cex, fact = 25, max) suppressWarnings(log(abs(cex) * fact + 1) / max)
163
164
  radius1         <- 0
  radius2         <- 60
165
166
167
  tri_R           <- -0x25BA
  tri_B           <- -0x25BC
  tri_L           <- -0x25C4
168
  tri_T           <- -0x25B2
169

170
  par(col.axis = par("fg"), cex.axis = 1.3, cex.lab = 1.3, cex = 0.7, mgp = mgp)
171

172
173
174
175
176
177
178
179
180
  if (!is.null(Theme)) {
    if (Theme == "Cyborg") {
      col_mod_bg    <- rgb(255-245, 255-245, 255-245, maxColorValue = 255)
      col_mod_bd    <- rgb(255-231, 255-231, 255-231, maxColorValue = 255)
    }
    if (Theme == "Flatly") {
      col_mod_bg    <- "#ECF0F1"
      col_mod_bd    <- "#ECF0F1"
    }
181
  }
182

183
  # Pas de temps
184
185
186
  dates_deb       <- EventDate
  n_pdt           <- length(which(OutputsModel$DatesR >= EventDate & OutputsModel$DatesR <= SimPer[2L]))
  i_pdt           <- which(format(OutputsModel$DatesR, "%Y%m%d") == format(EventDate, "%Y%m%d"))
187
188


189
  # --------------------------------------------------------------------------------
190
  # UH 1 & 2
191
  # --------------------------------------------------------------------------------
192

193
  if (HydroModel %in% c("GR4J", "GR6J")) {
194
195
196
197
198
199
200
    # Calcul des ordonnees SH1 de l' "hydrogramme unitaire cumule" UH1
    SH1     <- array(NA, NH)
    for (i in 1:NH) {
      if (i <= 0)                  SH1[i] <- 0
      if (i > 0 & i < Param[4])    SH1[i] <- (i/Param[4])^(D)
      if (i >= Param[4])           SH1[i] <- 1
    }
201

202
203
204
205
206
207
208
    # Calcul des ordonnees UH1 de l' "hydrogramme unitaire discret" UH1
    UH1     <- array(NA, NH)
    for (j in 1:NH) {
      if (j == 1) {
        UH1[j] <- SH1[j]
      } else {
        UH1[j] <- SH1[j] - SH1[j-1]
209
      }
210
    }
211

212
213
    # Parametres
    max_UH1 <- log(sqrt(max(max(UH1)*OutputsModel$PR*0.9))+1)
214

215
  }
216

217
218
219
220
221
222
223
224
225
  if (HydroModel != "GR2M") {
    # Calcul des ordonnees SH2 de l' "hydrogramme unitaire cumule" UH2
    SH2     <- array(NA, 2*NH)
    for (i in 1:(2*NH)) {
      if (i <= 0)                           SH2[i] <- 0
      if (i > 0 & i < Param[4])             SH2[i] <- 0.5*(i/Param[4])^(D)
      if (i >= Param[4] & i < 2*Param[4])   SH2[i] <- 1 - (0.5*(2-i/Param[4])^(D))
      if (i >= 2*Param[4])                  SH2[i] <- 1
    }
226

227
228
229
230
231
232
233
    # Calcul des ordonnees UH2 de l' "hydrogramme unitaire discret" UH2
    UH2     <- array(NA, 2*NH)
    for (j in 1:(2*NH)) {
      if (j == 1) {
        UH2[j] <- SH2[j]
      } else {
        UH2[j] <- SH2[j] - SH2[j-1]
234
      }
235
    }
236

237
    # Parametres
238
    max_UH2 <- log(sqrt(max(max(UH2)*OutputsModel$PR*0.1))+1)
239
  }
240

241
242
243
  # --------------------------------------------------------------------------------
  # PARTITIONNEMENT FENETRE GRAPHIQUE
  # --------------------------------------------------------------------------------
244

245
  # layout(matrix(c(1:4, 4, 4), nrow = 3, ncol = 2, byrow = FALSE), widths = c(1.0, 0.6))
246

247
  # --------------------------------------------------------------------------------
248
  # PLUIE ET ETP
249
  # --------------------------------------------------------------------------------
250

251
  # P
252
253
254
255
256
  # par(mar = c(2, 4, 1, 1), mgp = mgp)
  # plot(OutputsModel$Dates, OutputsModel$Precip, type = "h", col = col_P, ylim = rev(range(OutputsModel$Precip)), xaxt = "n", ylab = "precip. [mm/d]")
  # rect(xleft = EventDate, ybottom = par("usr")[3], xright =  par("usr")[2], ytop =  par("usr")[4], col = adjustcolor(par("bg"), alpha.f = 0.75), border = NA)
  # abline(v = EventDate, col = "grey", lwd = 2, lty = 2)
  # box()
257

258
  # ETP
259
260
261
262
  # par(mar = c(2, 4, 1, 1), mgp = mgp)
  # plot(OutputsModel$Dates, OutputsModel$PotEvap, pch = 19, col = col_E, xaxt = "n", ylab = "evapo. [mm/d]")
  # rect(xleft = EventDate, ybottom = par("usr")[3], xright =  par("usr")[2], ytop =  par("usr")[4], col = adjustcolor(par("bg"), alpha.f = 0.75), border = NA)
  # abline(v = EventDate, col = "grey", lwd = 2, lty = 2)
263
264
265
  # box()


266
267
268
  # --------------------------------------------------------------------------------
  # DEBIT
  # --------------------------------------------------------------------------------
269

270
  # Q
271
272
273
274
275
  # par(mar = c(2, 4, 1, 1), mgp = mgp)
  # plot(OutputsModel$Dates, OutputsModel$Qobs, type = "l", ylab = "flow [mm/d]")
  # lines(OutputsModel$Dates[1:i_pdt], OutputsModel$Qsim[1:i_pdt], type = "l", col = "orangered")
  # rect(xleft = EventDate, ybottom = par("usr")[3], xright =  par("usr")[2], ytop =  par("usr")[4], col = adjustcolor(par("bg"), alpha.f = 0.75), border = NA)
  # abline(v = EventDate, col = "grey", lwd = 2, lty = 2)
276
277
278
  # box()


279
  # --------------------------------------------------------------------------------
280
  # SCHEMAS MODELES
281
  # --------------------------------------------------------------------------------
282

283
  # Cadre
284
  par(mar = rep(0.2, 4))
285
  par(fg = par("fg"))
286
  plot(x = 0, type = "n", xlab = "", ylab = "", axes = FALSE, ylim = c(0, 1000), xlim = c(0, 1000))
287

288
  # Le modele
289
  rect(xleft = 0, xright = 1000, ybottom = 50, ytop = 970, col = col_mod_bg, border = col_mod_bd)
290
291


292
  # --------------------------------------------------------------------------------
293
  # ENTREES / SORTIES
294
  # --------------------------------------------------------------------------------
295

296
  # Entrees P et ETP
297
298
299
  if (CemaNeige) {
    text(x = xy_P[1]*1.65, y = xy_P[2]*0.98, labels = "+ CemaNeige", adj = c(1, 1), font = 2, col = "grey40", cex = 1.6)
  }
300
301
  text(x = xy_P[1], y = xy_P[2], labels = "P", pos = 3, font = 2, col = col_P, cex = 1.8)
  text(x = xy_E[1], y = xy_E[2], labels = "E", pos = 3, font = 2, col = col_E, cex = 1.8)
302

303
  # Sorties Q
304
  text(x = xy_Q[1], y = xy_Q[2], labels = "Q", pos = 1, font = 2, col = col_Q, cex = 1.8)
305

306
307
  # Parametres
  tmp_decal   <- 20
308
309


310
  # --------------------------------------------------------------------------------
311
  # NEUTRALISATION DE P
312
  # --------------------------------------------------------------------------------
313

314
315
316
317
318
319
  if (HydroModel != "GR2M") {
    # Interception
    segments(x0 = xy_E[1]-50, x1 = xy_P[1]+50,
             y0 = y_interception+tmp_decal, y1 = y_interception+tmp_decal)
    text(x = xy_P[1]+50, y = y_interception+20, labels = "Interception", pos = 4, font = 1, cex = 1.4)
  }
320

321
322
  # E vers Es et P vers Ps ou Pn
  y_Xs <- ifelse(HydroModel == "GR2M", y_rendement+2*tmp_decal, y_interception+tmp_decal)
323

324
  # P vers Pn
325
  segments(x0 = xy_P[1], x1 = xy_P[1], y0 = xy_P[2], y1 = y_Xs)
326

327
  # Pn vers Ps
328
329
330
331
332
333
  segments(x0 = xy_P[1], x1 = xy_P[1],
           y0 = y_interception, y1 = y_rendement+2*tmp_decal)
  segments(x0 = x_Ps, x1 = xy_P[1],
           y0 = y_rendement+2*tmp_decal, y1 = y_rendement+2*tmp_decal)
  segments(x0 = x_Ps, x1 = x_Ps,
           y0 = y_rendement+2*tmp_decal, y1 = y_rendement)
334

335
  # Pn vers Pn - Ps (P vers Pn si GR2M)
336
337
338
339
  segments(x0 = xy_P[1], x1 = x_PnPs,
           y0 = y_rendement+2*tmp_decal, y1 = y_rendement+2*tmp_decal)
  segments(x0 = x_PnPs , x1 = x_PnPs,
           y0 = y_rendement+2*tmp_decal, y1 = y_rendement)
340

341
  # Pn - Ps vers Pr (Pn vers Pr si GR2M)
342
343
  segments(x0 = x_PnPs, x1 = x_PnPs,
           y0 = y_rendement, y1 = y_percolation)
344

345
  # E vers En puis Es
346
  segments(x0 = xy_E[1], x1 = xy_E[1],
347
           y0 = xy_E[2], y1 = y_Xs)
348
349
  segments(x0 = xy_E[1], x1 = xy_E[1],
           y0 = y_interception, y1 = y_rendement)
350

351
  if (HydroModel != "GR2M") {
352
  # Ecriture
353
354
355
356
  plotrix::boxed.labels(x = xy_P[1], y = y_interception, labels = "Pn",
                        bg = col_mod_bg, border = NA, xpad = xpad, ypad = ypad)
  plotrix::boxed.labels(x = xy_E[1], y = y_interception, labels = "En",
                        bg = col_mod_bg, border = NA, xpad = xpad, ypad = ypad)
357
  }
358

359
360
361
362
363
364
  # ETP
  if (OutputsModel$PotEvap[i_pdt] != 0) {
    points(x = xy_PE[1], y =  xy_PE[2],
           type = "p", pch = tri_T, col = col_E,
           cex = cex_tri(OutputsModel$PotEvap[i_pdt], fact = fact_triangle, max = cex_max_poly))
  }
365

366
367
368
369
370
371
  # Precipitation
  if (OutputsModel$Precip[i_pdt] != 0) {
    points(x = xy_Precip[1], y =  xy_Precip[2],
           type = "p", pch = tri_B, col = col_P,
           cex = cex_tri(OutputsModel$Precip[i_pdt], fact = fact_triangle, max = cex_max_poly))
  }
372

373
  # Pn et Ps
374
375
376
  points(x = x_Ps, y = y_rendement+1.2*tmp_decal,
         type = "p", pch = tri_B, col = col_P,
         cex = cex_tri(OutputsModel$Ps[i_pdt], fact = fact_triangle, max = cex_max_poly))
377
  if (HydroModel != "GR2M") {
378
379
380
381
382
    points(x = x_PnPs, y = y_rendement+1.2*tmp_decal,
           type = "p", pch = tri_B, col = col_P,
           cex = cex_tri(OutputsModel$Pn[i_pdt] - OutputsModel$Ps[i_pdt], fact = fact_triangle, max = cex_max_poly))
  } else {
    points(x = x_PnPs, y = y_rendement+1.2*tmp_decal,
383
           type = "p", pch = tri_B, col = col_P,
384
           cex = cex_tri(OutputsModel$Pn[i_pdt], fact = fact_triangle, max = cex_max_poly))
385
  }
386
387


388
  # --------------------------------------------------------------------------------
389
  # FONCTION DE RENDEMENT
390
  # --------------------------------------------------------------------------------
391

392
  # Es
393
394
  plotrix::boxed.labels(x = xy_E[1], y = y_rendement, labels = "Es",
                        bg = col_mod_bg, border = NA, xpad = xpad, ypad = ypad)
395

396
397
398
399
400
401
  # Evaporation reelle
  if (OutputsModel$AE[i_pdt] != 0) {
    points(x = xy_AE[1], y =  xy_AE[2],
           type = "p", pch = tri_T, col = col_P,
           cex = cex_tri(OutputsModel$AE[i_pdt], fact = fact_triangle, max = cex_max_poly))
  }
402

403
  # Ps et Pn - Ps (Ps et Pn si GR2M)
404
405
  plotrix::boxed.labels(x = x_Ps  , y = y_rendement, labels = "Ps"   ,
                        bg = col_mod_bg, border = NA, xpad = xpad, ypad = ypad)
406
  plotrix::boxed.labels(x = x_PnPs, y = y_rendement, labels = ifelse(HydroModel != "GR2M", "Pn - Ps", "Pn"),
407
                        bg = col_mod_bg, border = NA, xpad = xpad, ypad = ypad)
408

409
  # Reservoir de production
410
411
  rect(xleft   = xy_min_PROD[1], xright = xy_min_PROD[1]+base_res,
       ybottom = xy_min_PROD[2], ytop   = xy_min_PROD[2]+OutputsModel$Prod[i_pdt]*fact_res,
412
       col = col_SP, border = NA)
413
414
415
416
417
418
  segments(x0 = xy_min_PROD[1], x1 = xy_min_PROD[1]+base_res,
           y0 = xy_min_PROD[2], y1 = xy_min_PROD[2])
  segments(x0 = xy_min_PROD[1], x1 = xy_min_PROD[1],
           y0 = xy_min_PROD[2], y1 = xy_min_PROD[2]+Param[1]*fact_res)
  segments(x0 = xy_min_PROD[1]+base_res, x1 = xy_min_PROD[1]+base_res,
           y0 = xy_min_PROD[2], y1 = xy_min_PROD[2]+Param[1]*fact_res)
419
  text(x = 30, y = xy_min_PROD[2]+15, labels = "Prod.\nstore", cex = 1.4, pos = 4)
420
421


422
  # --------------------------------------------------------------------------------
423
  # PERCOLATION
424
  # --------------------------------------------------------------------------------
425

426
  # Reservoir de production vers Pr
427
428
429
430
  segments(x0 = xy_min_PROD[1]+base_res/2, x1 = xy_min_PROD[1]+base_res/2,
           y0 = xy_min_PROD[2], y1 = y_percolation)
  segments(x0 = xy_min_PROD[1]+base_res/2, x1 = x_PnPs,
           y0 = y_percolation, y1 = y_percolation)
431

432
  # Perc
433
434
  plotrix::boxed.labels(x = xy_min_PROD[1]+base_res/2, y = y_percolation, labels = "Perc.",
                        bg = col_mod_bg, border = NA, xpad = xpad, ypad = ypad)
435

436
  # Valeur de Perc
437
  if (OutputsModel$Perc[i_pdt] != 0) {
438
439
440
    points(x = xy_min_PROD[1]+base_res+75, y = y_percolation,
           type = "p", pch = tri_R, col = col_P,
           cex = cex_tri(OutputsModel$Perc[i_pdt], fact = fact_triangle, max = cex_max_poly))
441
  }
442

443
  # parametres
444
445
  tmp_decal   <- (y_percolation - y_entreeUH) / 2

446
447
  # Pr vers UH (Pr vers reservoir de routage si GR2M)
  k <- ifelse(HydroModel == "GR2M", 0.5, 1)
448
  segments(x0 = x_PnPs, x1 = x_PnPs,
449
           y0 = y_percolation, y1 = (y_entreeUH*k) + tmp_decal/2)
450
451


452
453
454
  if (HydroModel == "GR2M") {
    plotrix::boxed.labels(x = x_PnPs, y = y_percolation, labels = "Pr",
                          bg = col_mod_bg, border = NA, xpad = xpad, ypad = ypad)
455

456
457
458
459
460
461
    if (OutputsModel$PR[i_pdt] != 0) {
      points(x = x_PnPs[1], y = y_entreeUH+tmp_decal,
             type = "p", pch = tri_B, col = col_P,
             cex = cex_tri(OutputsModel$PR[i_pdt], fact = fact_triangle, max = cex_max_poly))
    }
  }
462
463


464
  if (HydroModel %in% c("GR4J", "GR6J")) {
465

466
467
468
    # --------------------------------------------------------------------------------
    # SEPARATION DE PR
    # --------------------------------------------------------------------------------
469

470
471
472
473
474
    # Pr vers UH1
    segments(x0 = xy_Q9[1], x1 = x_PnPs,
             y0 = y_entreeUH+tmp_decal/2, y1 = y_entreeUH+tmp_decal/2)
    segments(x0 = xy_Q9[1], x1 = xy_Q9[1],
             y0 = y_entreeUH+tmp_decal/2, y1 = xy_Q9[2])
475

476
477
478
479
480
    # Pr vers UH2
    segments(x0 = x_PnPs, x1 = xy_Q1[1],
             y0 = y_entreeUH+tmp_decal/2, y1 = y_entreeUH+tmp_decal/2)
    segments(x0 = xy_Q1[1], x1 = xy_Q1[1],
             y0 = y_entreeUH+tmp_decal/2, y1 = y_routage)
481

482
    # Pr
483
484
    plotrix::boxed.labels(x = x_PnPs, y = y_percolation, labels = "Pr",
                          bg = col_mod_bg, border = NA, xpad = xpad, ypad = ypad)
485

486
487
488
489
490
491
    # Pr
    if (OutputsModel$PR[i_pdt] != 0) {
      points(x = x_PnPs[1], y = y_entreeUH+tmp_decal,
             type = "p", pch = tri_B, col = col_P,
             cex = cex_tri(OutputsModel$PR[i_pdt], fact = fact_triangle, max = cex_max_poly))
    }
492
493


494
495
496
    # --------------------------------------------------------------------------------
    # HYDROGRAMME UNITAIRE 1
    # --------------------------------------------------------------------------------
497

498
499
500
501
502
503
    # Entree de UH1
    if (OutputsModel$PR[i_pdt] != 0) {
      points(x = xy_Q9[1], y =y_entreeUH,
             type = "p", pch = tri_B, col = col_P,
             cex =  cex_tri(OutputsModel$PR[i_pdt]*0.9, fact = fact_triangle, max = cex_max_poly))
    }
504

505
506
507
508
509
510
511
    # Remplissage de UH1
    PR_mat_UH1_lg <- ceiling(Param[4])
    PR_mat_UH1_id <- max(i_pdt-PR_mat_UH1_lg+1, 1):i_pdt
    PR_mat_UH1 <- matrix(rep(c(rep(0, times = PR_mat_UH1_lg-length(PR_mat_UH1_id)+1),
                               OutputsModel$PR[PR_mat_UH1_id]), times = PR_mat_UH1_lg),
                         ncol = PR_mat_UH1_lg+1)[, -1L]
    PR_mat_UH1[lower.tri(PR_mat_UH1)] <- 0
512
513


514
515
516
    # --------------------------------------------------------------------------------
    # HYDROGRAMME UNITAIRE 2
    # --------------------------------------------------------------------------------
517

518
519
520
521
522
523
    # Entree de UH2
    if (OutputsModel$PR[i_pdt] != 0) {
      points(x = xy_Q1[1], y = y_entreeUH,
             type = "p", pch = tri_B, col = col_P,
             cex = cex_tri(OutputsModel$PR[i_pdt]*0.1, fact = fact_triangle, max = cex_max_poly))
    }
524

525
526
527
528
529
530
531
    # Remplissage de UH2
    PR_mat_UH2_lg <- ceiling(Param[4]*2)
    PR_mat_UH2_id <- max(i_pdt-PR_mat_UH2_lg+1, 1):i_pdt
    PR_mat_UH2 <- matrix(rep(c(rep(0, times = PR_mat_UH2_lg-length(PR_mat_UH2_id)+1),
                               OutputsModel$PR[PR_mat_UH2_id]), times = PR_mat_UH2_lg),
                         ncol = PR_mat_UH2_lg+1)[, -1L]
    PR_mat_UH2[lower.tri(PR_mat_UH2)] <- 0
532
533
534

  }

535
  if (HydroModel == "GR5J") {
536

537
538
539
    # --------------------------------------------------------------------------------
    # SEPARATION DE PR
    # --------------------------------------------------------------------------------
540

541
542
543
    # sortie UH
    segments(x0 = x_PnPs, x1 = x_PnPs,
             y0 = y_entreeUH-2*tmp_decal, y1 = y_entreeUH-3*tmp_decal)
544

545
546
547
548
549
    # sortie UH vers branche 1
    segments(x0 = xy_Q9[1], x1 = x_PnPs,
             y0 = y_entreeUH-3*tmp_decal, y1 = y_entreeUH-3*tmp_decal)
    segments(x0 = xy_Q9[1], x1 = xy_Q9[1],
             y0 = y_entreeUH-3*tmp_decal, y1 = xy_Q9[2])
550

551
552
553
554
555
    # sortie UH vers branche 2
    segments(x0 = x_PnPs, x1 = xy_Q1[1],
             y0 = y_entreeUH-3*tmp_decal, y1 = y_entreeUH-3*tmp_decal)
    segments(x0 = xy_Q1[1], x1 = xy_Q1[1],
             y0 = y_entreeUH-3*tmp_decal, y1 = y_routage)
556

557
    # Pr
558
559
    plotrix::boxed.labels(x = x_PnPs, y = y_percolation, labels = "Pr",
                          bg = col_mod_bg, border = NA, xpad = xpad, ypad = ypad)
560
561


562
563
564
    # --------------------------------------------------------------------------------
    # HYDROGRAMME UNITAIRE
    # --------------------------------------------------------------------------------
565

566
567
568
569
570
571
    # Entree de UH (PR)
    if (OutputsModel$PR[i_pdt] != 0) {
      points(x = x_PnPs[1], y = y_entreeUH+tmp_decal,
             type = "p", pch = tri_B, col = col_P,
             cex = cex_tri(OutputsModel$PR[i_pdt], fact = fact_triangle, max = cex_max_poly))
    }
572
573


574
575
576
577
578
579
580
    # Remplissage de UH2
    PR_mat_UH2_lg <- ceiling(Param[4]*2)
    PR_mat_UH2_id <- max(i_pdt-PR_mat_UH2_lg+1, 1):i_pdt
    PR_mat_UH2 <- matrix(rep(c(rep(0, times = PR_mat_UH2_lg-length(PR_mat_UH2_id)+1),
                               OutputsModel$PR[PR_mat_UH2_id]), times = PR_mat_UH2_lg),
                         ncol = PR_mat_UH2_lg+1)[, -1L]
    PR_mat_UH2[lower.tri(PR_mat_UH2)] <- 0
581
582


583
584
585
586
587
588
    # Sorties de UH
    if (PR_mat_UH2[1] != 0) {
      points(x = x_PnPs[1], y =  y_entreeUH-5*tmp_decal/2,
             type = "p", pch = tri_B, col = col_P,
             cex = cex_tri(PR_mat_UH2[1], fact = fact_triangle, max = cex_max_poly))
    }
589
590
591

  }

592
593
594
595
596
597
  # sortie de UH1 vers reservoirs exponentiel et de routage
  if (HydroModel == "GR6J") {
    segments(x0 = xy_Q9[1], x1 = xy_Q9[1],
             y0 = y_entreeUH-3*tmp_decal, y1 = xy_Q9[2])
    segments(x0 = xy_Q9[1]*0.80, x1 = xy_Q9[1]*1.30,
             y0 = xy_Q9[2], y1 = xy_Q9[2])
598
    segments(x0 = xy_Q9[1]*1.30, x1 = xy_Q9[1]*1.30,
599
             y0 = xy_Q9[2], y1 = xy_Q9[2]*0.65)
600
601
    segments(x0 = xy_Q9[1]*0.80, x1 = xy_Q9[1]*0.80,
             y0 = xy_Q9[2], y1 = xy_Q9[2]*0.90)
602
603
604
605
    segments(x0 = xy_Q9[1]*0.55, x1 = xy_Q9[1]*0.55,
             y0 = xy_Q9[2]*0.70, y1 = y_routage)
    segments(x0 = xy_Q9[1]*0.55, x1 = xy_min_ROUT[1]+base_res/2,
             y0 = y_routage, y1 = y_routage)
606
  }
607

608
609
610
611
  if (HydroModel != "GR2M") {
    # Q9
    if (OutputsModel$Q9[i_pdt] != 0) {
      points(x = xy_Q9[1], y = xy_Q9[2]+tmp_decal,
612
             type = "p", pch = tri_B, col = col_P,
613
614
615
616
617
618
             cex = cex_tri(OutputsModel$Q9[i_pdt], fact = fact_triangle, max = cex_max_poly))
      if (HydroModel == "GR6J") {
        # Q9 exp
        points(x = xy_Q9[1]*0.80, y = xy_Q9[1]*0.73,
               type = "p", pch = tri_B, col = col_P,
               cex = cex_tri(OutputsModel$Q9[i_pdt]*0.4, fact = fact_triangle, max = cex_max_poly))
619
        # Q9 rout
620
621
622
623
624
625
626
627
628
        points(x = xy_Q9[1]*1.30, y = xy_Q9[1]*0.73,
               type = "p", pch = tri_B, col = col_P,
               cex = cex_tri(OutputsModel$Q9[i_pdt]*0.6, fact = fact_triangle, max = cex_max_poly))
        # QrExp
        plotrix::boxed.labels(x = xy_Q9[1]*0.55, y = y_routage, labels = "QrExp", bg = col_mod_bg, border = NA, xpad = xpad, ypad = ypad)
      }
    }
    plotrix::boxed.labels(x = xy_Q9[1], y = xy_Q9[2], labels = "Q9",
                          bg = col_mod_bg, border = NA, xpad = xpad, ypad = ypad)
629
630


631
632
633
    # Q1
    if (OutputsModel$Q1[i_pdt] != 0) {
      points(x = xy_Q1[1], y =  xy_Q1[2]+tmp_decal,
634
             type = "p", pch = tri_B, col = col_P,
635
636
             cex = cex_tri(OutputsModel$Q1[i_pdt], fact = fact_triangle, max = cex_max_poly))
      segments(x0 = xy_Q[1], x1 = xy_Q1[1], y0 = y_routage, y1 = y_routage)
637
    }
638

639
    plotrix::boxed.labels(x = xy_Q1[1], y = xy_Q1[2], labels = "Q1", bg = col_mod_bg, border = NA, xpad = xpad, ypad = ypad)
640

641
642
643
644
645
646
    # Valeur de Qd
    if (OutputsModel$QD[i_pdt] != 0) {
      points(x = xy_Q[1]+30, y =  y_routage,
             type = "p", pch = tri_L, col = col_P,
             cex = cex_tri(OutputsModel$QD[i_pdt], fact = fact_triangle, max = cex_max_poly))
    }
647

648
649
    # Qd
    plotrix::boxed.labels(x = xy_Q1[1], y = y_routage, labels = "Qd", bg = col_mod_bg, border = NA, xpad = xpad, ypad = ypad)
650

651
  }
652

653
  # --------------------------------------------------------------------------------
654
  # RESERVOIR DE ROUTAGE
655
  # --------------------------------------------------------------------------------
656

657
  # Triche pour la taille du reservoire de routage
658
  tmp_triche   <- 0#80
659

660
  # Reservoir de routage
661
  if (HydroModel == "GR2M") {
662
    xy_min_ROUT[1] <- x_PnPs - base_res/2
663
664
    Param[3] <- 600
  }
665
666
  rect(xleft = xy_min_ROUT[1], xright = xy_min_ROUT[1]+base_res,
       ybottom = xy_min_ROUT[2], ytop = xy_min_ROUT[2]+OutputsModel$Rout[i_pdt]*fact_res+tmp_triche,
667
       col = col_SR, border = NA)
668
  segments(x0 = xy_min_ROUT[1], x1 = xy_min_ROUT[1]+base_res,
669
670
671
672
673
           y0 = xy_min_ROUT[2], y1 = xy_min_ROUT[2])
  segments(x0 = xy_min_ROUT[1], x1 = xy_min_ROUT[1],
           y0 = xy_min_ROUT[2], y1 = xy_min_ROUT[2]+Param[3]*fact_res+tmp_triche)
  segments(x0 = xy_min_ROUT[1]+base_res, x1 = xy_min_ROUT[1]+base_res,
           y0 = xy_min_ROUT[2], y1 = xy_min_ROUT[2]+Param[3]*fact_res+tmp_triche)
674
  text(x = 30, y = xy_min_ROUT[2]+15, labels = "Routing\nstore", cex = 1.4, pos = 4)
675

676
  # Sorties du reservoir
677
678
679
680
  segments(x0 = xy_min_ROUT[1]+base_res/2, x1 = xy_min_ROUT[1]+base_res/2,
           y0 = xy_min_ROUT[2], y1 = y_routage)
  segments(x0 = xy_min_ROUT[1]+base_res/2, x1 = xy_Q[1],
           y0 = y_routage, y1 = y_routage)
681
682


683
684
  if (HydroModel != "GR2M") {
    # Qr
685
    if (HydroModel != "GR6J") {
686
687
688
689
690
691
692
      plotrix::boxed.labels(x = xy_min_ROUT[1]+base_res/2, y = y_routage, labels = "Qr",
                            bg = col_mod_bg, border = NA, xpad = xpad, ypad = ypad)
    }
    if (HydroModel == "GR6J") {
      plotrix::boxed.labels(x = xy_min_ROUT[1]+base_res/1.5, y = (xy_min_ROUT[2]+y_routage)/2, labels = "Qr",
                            bg = col_mod_bg, border = NA, xpad = xpad, ypad = ypad)
    }
693

694
695
696
697
698
699
700
701
702
703
704
    # Valeur de Qr
    if (OutputsModel$QR[i_pdt] != 0) {
      if (HydroModel != "GR6J") {
        points(x = xy_Q[1]-100, y = y_routage,
               type = "p", pch = tri_R, col = col_P,
               cex = cex_tri(OutputsModel$QR[i_pdt], fact = fact_triangle, max = cex_max_poly))
      } else {
        points(x = xy_min_ROUT[1]+base_res/2, y = (xy_min_ROUT[2]+y_routage)/2,
               type = "p", pch = tri_B, col = col_P,
               cex = cex_tri(OutputsModel$QR[i_pdt], fact = fact_triangle, max = cex_max_poly))
      }
705
    }
706
  }
707
708


709
710
711
712
713
  # --------------------------------------------------------------------------------
  # RESERVOIR EXPONENTIEL
  # --------------------------------------------------------------------------------

  if (HydroModel == "GR6J") {
714

715
716
    # Triche pour la taille du reservoire exponentiel
    tmp_triche   <- 0#80
717

718
719
720
721
722
    # Exp en log
    signExp <- ifelse(OutputsModel$Exp[i_pdt] > 0, +1, -1)
    logExpIpdt <- log(abs(OutputsModel$Exp[i_pdt])+1e-6) * signExp
    logExpMax  <- log(max(abs(OutputsModel$Exp  ))+1e-6)

723
724
    # Reservoir exponentiel
    rect(xleft = xy_min_EXPO[1], xright = xy_min_EXPO[1]+base_res,
725
         ybottom = xy_min_EXPO[2], ytop = xy_min_EXPO[2]+logExpIpdt*fact_resExp+tmp_triche,
726
         col = ifelse(OutputsModel$Exp[i_pdt] > 0, "#10B510", "#FF0303"), border = NA)
727
    segments(x0 = xy_min_EXPO[1], x1 = xy_min_EXPO[1]+base_res,
728
729
             y0 = xy_min_EXPO[2], y1 = xy_min_EXPO[2])
    segments(x0 = xy_min_EXPO[1], x1 = xy_min_EXPO[1],
730
             y0 = xy_min_EXPO[2], y1 = xy_min_EXPO[2]+logExpMax*fact_resExp+tmp_triche)
731
    segments(x0 = xy_min_EXPO[1]+base_res, x1 = xy_min_EXPO[1]+base_res,
732
             y0 = xy_min_EXPO[2], y1 = xy_min_EXPO[2]+logExpMax*fact_resExp+tmp_triche)
733
    segments(x0 = xy_min_EXPO[1], x1 = xy_min_EXPO[1],
734
             y0 = xy_min_EXPO[2], y1 = xy_min_EXPO[2]-logExpMax*fact_resExp-tmp_triche)
735
    segments(x0 = xy_min_EXPO[1]+base_res, x1 = xy_min_EXPO[1]+base_res,
736
             y0 = xy_min_EXPO[2], y1 = xy_min_EXPO[2]-logExpMax*fact_resExp-tmp_triche)
737
738
    text(x = 30, y = xy_min_EXPO[2]+00, labels = "Exp.\nstore", cex = 1.4, pos = 4)
    points(x = 180, y = xy_min_EXPO[2]+20, pch = 43, # +
739
         cex = 2.0, col = "#10B510")
740
    points(x = 178, y = xy_min_EXPO[2]-20, pch = 95, # -
741
742
743
744
745
746
747
748
         cex = 1.6, col = "#FF0303")

    # Valeur de QrExp
    if (OutputsModel$QR[i_pdt] != 0) {
      points(x = xy_Q[1]-350, y = y_routage,
             type = "p", pch = tri_R, col = col_P,
             cex = cex_tri(OutputsModel$QRExp[i_pdt], fact = fact_triangle, max = cex_max_poly))
    }
749

750
    # Valeur de Qr + QrExp
751
    if (OutputsModel$QR[i_pdt] != 0) {
752
753
754
755
756
      points(x = xy_Q[1]-100, y = y_routage,
             type = "p", pch = tri_R, col = col_P,
             cex = cex_tri(OutputsModel$QR[i_pdt]+OutputsModel$QRExp[i_pdt], fact = fact_triangle, max = cex_max_poly))
    }
  }
757
758


759
760
761
  # --------------------------------------------------------------------------------
  # Q FINAL
  # -------------------------------------------------------------------------------
762