Skip to content
GitLab
Select Git revision
  • dev default
  • master
  • CRAN_v0.3.5
  • CRAN_v0.3.4
  • CRAN_v0.3.3
  • CRAN_v0.3.2
  • CRAN_v0.3.1
  • CRAN_v0.3.0
  • CRAN_v0.2.13
  • CRAN_v0.2.12
  • CRAN_v0.2.11
  • CRAN_v0.2.10.112
  • CRAN_v0.2.9.25
  • CRAN_v0.2.9.23
  • CRAN_v0.2.8.69
  • CRAN_v0.2.6.29
  • CRAN_v0.2.6.27
  • CRAN_V0.2.2.0
  • CRAN_v0.2.0.9
  • CRAN_v0.2.0.8
  • CRAN_v0.1.11.25
Find file BlameHistoryPermalink
plot.CalGR.R 4.55 KiB
1 
plot.CalGR <- function(x, xlab = NULL, ylab =  NULL, main = NULL, which = c("perf", "iter", "ts"), ...) {
2
3 
  if (! any(class(x) %in% "CalGR")) {
4 
    stop("Non convenient data for x argument. Must be of class \"CalGR\"")
5 
  }
6 
  if (! any(which %in% c("perf", "iter", "ts"))) {
7 
    stop("Non convenient data for which argument. Must be of class \"perf\", \"iter\" or \"ts\"")
8 
  }
9
10 
  nbParamX <- .TypeModelGR(x)$NbParam #as.numeric(gsub("\\D", "", x$TypeModel))
11 
  nbParamC <- ifelse(.TypeModelGR(x)$CemaNeige, 2, 0)
12 
  nbParam <- nbParamX + nbParamC
13 
  nmParam <- c(sprintf("X%i", 1:nbParamX), sprintf("C%i", seq_len(nbParamC)))
14
15 
  opar <- par(no.readonly = TRUE)
16 
  on.exit(par(opar))
17
18 
  if (any(which[1L] %in% c("perf"))) {
19 
    plot.OutputsModel(x$OutputsModel, Qobs = x$Qobs, ...)
20 
  }
21
22 
  if (any(which[1L] %in% c("iter"))) {
23 
    layout.list <- list(matrix(c(1:2), ncol = 2),
24 
                        matrix(c(1:3, 3), ncol = 2),
25 
                        matrix(NA),
26 
                        matrix(c(1:5, 5), ncol = 3),
27 
                        matrix(c(1:5, 0, 6, 6), ncol = 4),
28 
                        matrix(c(1:7, 7), ncol = 4),
29 
                        matrix(c(1:5, 0, 6:8, 8), ncol = 5),
30 
                        matrix(c(1:9, 9), ncol = 5))
31 
    layout.mat <- layout.list[[nbParam]]
32 
    layout.w <- c(rep(1, ncol(layout.mat)-1), 2)
33 
    layout.h <- rep(1, nrow(layout.mat))
34 
    layout(mat = layout.mat, widths = layout.w, heights = layout.h)
35 
  }
36 
  # if (which[1L] %in%  c("both")) {
37 
  #   layout.mat <- rbind(layout.mat, max(layout.mat)+1, max(layout.mat)+2)
38 
  #   layout(mat = layout.mat, widths = layout.w, heights = layout.h)
39 
  # }
40 
  if (which[1L] %in% c("ts")) {
41 
    layout(mat = matrix(1:2), widths = c(1, 2), heights = c(1, 2))
42 
  }
43
44 
  if (any(which[1L] %in% c("iter"))) {
45 
    ParamLab <- data.frame(Name  = c(sprintf("X%i", 1:6), sprintf("C%i", 1:2)),
46 
                           Label = c("prod. store capacity [mm]",
47 
                                     "intercatchment exch. coef. [mm/TimeUnit]",
48 
                                     "routing store capacity [mm]",
49 
                                     "UH time constant [TimeUnit]",
50 
                                     "intercatchment exch. threshold [-]",
51 
                                     "coef. for emptying exp. store [mm]",
52 
                                     "weight for snowpack thermal state [-]",
53 
                                     "degree-day melt coef. [mm/degC/TimeUnit]"))
54 
    ParamLab$Label <- gsub("TimeUnit", substr(.TypeModelGR(x)$TimeUnit, 1, 1), ParamLab$Label)
55
56
57 
    par(mar = c(2.5, 3.5, 2.0, 1.0), mgp = c(2.0, 0.6, 0.0), oma = c(0, 0, 4, 0))
58 
    for (i in seq_len(ncol(x$OutputsCalib$HistParamR))) {
59 
      plot(x$OutputsCalib$HistParamR[, i],
60 
           col = ifelse(grepl("X", nmParam[i]), "darkorchid3", "darkolivegreen3"),
61 
           type = "b", pch = 19, cex = 0.8,
62 
           xlab = "", ylab = ParamLab[ParamLab$Name == nmParam[i], "Label"],
63 
           main = nmParam[i])
64 
    }
65 
    plot(x$OutputsCalib$HistCrit,
66 
         col = "tomato",
67 
         type = "b", pch = 19, cex = 1.0,
68 
         xlab = "", ylab = "",
69 
         main = gsub(".*_", "", x$CalCrit))
70
7172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110
mainIter <- "Evolution of parameters and efficiency criterion\nduring the iterations of the steepest-descent step" title(main = ifelse(is.null(main), mainIter, main), outer = TRUE, line = 0.6, cex.main = 1.6) } if (any(which[1L] %in% c("ts"))) { if (is.null(xlab)) { xlab <- "Time" } if (is.null(ylab)) { ylab <- paste0(c("flow [mm/", "precip. [mm/"), .TypeModelGR(x)$TimeUnit, "]") } else { if (length(ylab) < 2) { ylab <- c(ylab, "") } } data <- data.frame(DatesR = x$OutputsModel$DatesR, Precip = x$OutputsModel$Precip, Qobs = x$Qobs, Qsim = x$OutputsModel$Qsim) par(mar = c(0.1, 4, 4, 1), xaxt = "n") plot(Precip ~ DatesR, data = data, type = "h", col = "royalblue", xlab = "", ylab = ylab[2L], main = main, ylim = rev(range(data$Precip))) par(mar = c(5, 4, 0.1, 1), xaxt = "s") matplot(data$DatesR, data[, c("Qobs", "Qsim")], type = "l", col = c("black", "orangered"), lty = 1, xlab = xlab, ylab = ylab[1L], main = "", xaxt = "n") axis.POSIXct(side = 1, x = data$DatesR) # if (plot.na) { # axis(side = 1, at = as.POSIXct(data$DatesR[is.na(data$Qobs)]), # labels = FALSE, lwd.ticks = 3, col.ticks = col.na, tck = 0.025, lend = "butt") # # legend("topright", legend = "NA", pch = 15, col = col.na, bty = "n", cex = 0.8) # } } # box() }