server.R 23.5 KB
Newer Older
1
2
# server.R

3
shinyServer(function(input, output, session) {
4
5
  
  
6
7
  ## --------------- Data preparation
  
8
  getPrep <- reactive({
9
    
10
    TMGR  <- .TypeModelGR(input$HydroModel)
11
    PARAM <- c(input$X1, input$X2, input$X3, input$X4, input$X5, input$X6)[seq_len(TMGR$NbParam)]
12
    
13
    if (input$SnowModel == "CemaNeige") {
14
15
      PARAM <- c(PARAM, input$C1, input$C2)
    }
16
    
17
18
    OBS <- ObsGR(ObsBV = get(input$Dataset), HydroModel = input$HydroModel,
                 CemaNeige = input$SnowModel == "CemaNeige",
19
20
                 Precip = .ShinyGR.args$Precip, PotEvap = .ShinyGR.args$PotEvap,
                 Qobs = get(input$Dataset), TempMean = .ShinyGR.args$TempMean, 
21
22
                 ZInputs = .ShinyGR.args$ZInputs, HypsoData = .ShinyGR.args$HypsoData,
                 NLayers = .ShinyGR.args$NLayers)
23
    
24
25
26
27
28
29
    WUPPER <- c(OBS$InputsModel$DatesR[1L], input$Period[1]-.TypeModelGR(OBS)$TimeLag)
    if (WUPPER[2] < WUPPER[1]) {
      WUPPER[2] <- WUPPER[1]
    }
    
    return(list(TMGR = TMGR, OBS = OBS, WUPPER = WUPPER))
30
    
31
  })
32
  
33
  
34
35
36
  
  ## --------------- Calibration
  
37
  ## If the user calibrate the model
38
  CAL_click <- reactiveValues(valueButton = 0)
39
  
unknown's avatar
unknown committed
40
  
41
  ## Automatic calibration
42
  observeEvent(input$CalButton, {
43

44
45
46
    CAL_opt <- list(Crit    = gsub(" .*", "", input$TypeCrit),
                    Transfo = gsub("(\\D{3} \\[)(\\w{0,4})(\\W*Q\\W*\\])", "\\2", input$TypeCrit))
    
47
    CAL     <- CalGR(ObsGR = getPrep()$OBS, CalCrit = CAL_opt$Crit, transfo = CAL_opt$Transfo,
48
49
                     WupPer = substr(getPrep()$WUPPER, 1, 10),
                     CalPer = substr(c(input$Period[1], input$Period[2]), 1, 10), verbose = FALSE)
50
    PARAM   <- CAL$OutputsCalib$ParamFinalR
51

52
53
54
55
    updateSliderInput(session, inputId = "X1", value = PARAM[1L])
    updateSliderInput(session, inputId = "X2", value = PARAM[2L])
    updateSliderInput(session, inputId = "X3", value = PARAM[3L])
    updateSliderInput(session, inputId = "X4", value = PARAM[4L])
56
    if (getPrep()$TMGR$NbParam >= 5) {
57
58
      updateSliderInput(session, inputId = "X5", value = PARAM[5L])
    }
59
    if (getPrep()$TMGR$NbParam >= 6) {
60
61
      updateSliderInput(session, inputId = "X6", value = PARAM[6L])
    }
62
    if (input$SnowModel == "CemaNeige") {
63
64
65
      updateSliderInput(session, inputId = "C1", value = PARAM[length(PARAM)-1])
      updateSliderInput(session, inputId = "C2", value = PARAM[length(PARAM)])
    }
66
67
    updateActionButton(session, inputId = "CalButton", label = "Model calibrated", icon = icon("check"))
    CAL_click$valueButton <- 1
68
  })
69
  
unknown's avatar
unknown committed
70
  
71
  ## Manual calibration
72
  observeEvent({input$Dataset ; input$HydroModel ; input$SnowModel ;
73
74
75
76
77
78
79
80
81
82
83
84
    input$X1 ; input$X2 ; input$X3 ; input$X4 ; input$X5 ; input$X6 ;
    input$TypeCrit ; input$Period}, {
      
      CAL_click$valueButton <-  CAL_click$valueButton - 1
      CAL_click$valueButton <- ifelse(CAL_click$valueButton < -1, -1, CAL_click$valueButton)
      if (CAL_click$valueButton < 0) {
        updateActionButton(session, inputId = "CalButton", label = "Run", icon = icon("refresh"))
      }
      
    })
  
  
85
  
86
  ## --------------- Simulation
87
88
  
  getRES <- reactive({
89

90
    PARAM <- c(input$X1, input$X2, input$X3, input$X4, input$X5, input$X6)[seq_len(getPrep()$TMGR$NbParam)]
91
    if (input$SnowModel == "CemaNeige") {
92
      PARAM <- c(PARAM, input$C1, input$C2)
93
    }
94
    
95
96
97
98
99
    ## Simulated flows computation
    SIM_opt <- list(Crit    = c("ErrorCrit_NSE", "ErrorCrit_KGE"),
                    Transfo = c("NO", "sqrt", "log"))
    SIM <- lapply(SIM_opt$Crit, function(iCRIT) {
      SIM_transfo <- lapply(SIM_opt$Transfo, function(iTRSF) {
100
        iTRSF <- gsub("NO", "", iTRSF)
101
102
        iSIM  <- SimGR(ObsGR = getPrep()$OBS, Param = PARAM,
                       WupPer = substr(getPrep()$WUPPER, 1, 10),
103
104
105
                       SimPer = substr(c(input$Period[1], input$Period[2]), 1, 10),
                       transfo = iTRSF, verbose = FALSE)
        iCRIT <- ErrorCrit(InputsCrit = iSIM$OptionsCrit, OutputsModel = iSIM$OutputsModel, FUN_CRIT = get(iCRIT), verbose = FALSE)
106
107
108
        iCRIT <- iCRIT[c("CritName", "CritValue")]
        return(list(SIM = iSIM, CRIT = iCRIT))
      })
109
      names(SIM_transfo) <- SIM_opt$Transfo
110
      return(SIM_transfo)
111
    })
112
    names(SIM) <- SIM_opt$Crit
113

114
    ## Criteria computation
115
    CRIT <- lapply(SIM, function(iCRIT) {
116
      lapply(SIM_opt$Transfo, function(iTRSF) {
117
118
        iCRIT[[iTRSF]][["CRIT"]]
      })
119
    })
120
    CRIT <- as.data.frame(matrix(na.omit(unlist(CRIT)), ncol = 2, byrow = TRUE), stringsAsFactors = FALSE)
121
    colnames(CRIT) <- c("Criterion", "Value")
122
    rownames(CRIT) <- NULL    
123
124
    CRIT$Value     <- as.numeric(CRIT$Value)
    CRIT$Criterion <- gsub("\\[", " [", CRIT$Criterion)
125
126
127
128
129
    
    .GlobalEnv$.ShinyGR.hist[[length(.GlobalEnv$.ShinyGR.hist)+1]] <- list(Qsim = SIM$ErrorCrit_KGE$NO$SIM$OutputsModel$Qsim,
                                                                           Param = PARAM,
                                                                           TypeModel = SIM$ErrorCrit_KGE$NO$SIM$TypeModel)

130
131
    .GlobalEnv$.ShinyGR.hist <- .GlobalEnv$.ShinyGR.hist[!(duplicated(sapply(.GlobalEnv$.ShinyGR.hist, function(x) sum(x$Param)), fromLast = TRUE) & 
                                                           duplicated(sapply(.GlobalEnv$.ShinyGR.hist, function(x) x$TypeModel ), fromLast = TRUE))]
132
    .GlobalEnv$.ShinyGR.hist <- tail(.GlobalEnv$.ShinyGR.hist, n = 2)
133
    if (length(.GlobalEnv$.ShinyGR.hist) == 2 &  is.null(names(.GlobalEnv$.ShinyGR.hist[[1]]))) {
134
      .GlobalEnv$.ShinyGR.hist[[1]] <- NULL
135
    }
136
137
    if (length(.GlobalEnv$.ShinyGR.hist) == 2 & !is.null(names(.GlobalEnv$.ShinyGR.hist[[1]]))) {
      if (length(.GlobalEnv$.ShinyGR.hist[[1]]$Qsim) != length(.GlobalEnv$.ShinyGR.hist[[2]]$Qsim)) {
138
139
        OBSold <- getPrep()$OBS
        OBSold$TypeModel <- .GlobalEnv$.ShinyGR.hist[[1]]$TypeModel
140
141
142
143
144
145
146
147
148
        if (.TypeModelGR(OBSold)$CemaNeige & !.TypeModelGR(getPrep()$OBS)$CemaNeige) { # present: No CemaNeige ; old: CemaNeige
          OBSold <- ObsGR(ObsBV = get(input$Dataset), HydroModel = .TypeModelGR(OBSold)$NameModel,
                       CemaNeige = .TypeModelGR(OBSold)$CemaNeige,
                       Precip = .ShinyGR.args$Precip, PotEvap = .ShinyGR.args$PotEvap,
                       Qobs = get(input$Dataset), TempMean = .ShinyGR.args$TempMean, 
                       ZInputs = .ShinyGR.args$ZInputs, HypsoData = .ShinyGR.args$HypsoData,
                       NLayers = .ShinyGR.args$NLayers)
        }

149
        SIMold <- SimGR(ObsGR = OBSold,
150
151
152
153
154
                        Param = .GlobalEnv$.ShinyGR.hist[[1]]$Param,
                        WupPer = substr(getPrep()$WUPPER, 1, 10),
                        SimPer = substr(c(input$Period[1], input$Period[2]), 1, 10),
                        verbose = FALSE)
        .GlobalEnv$.ShinyGR.hist[[1]]$Qsim <- SIMold$OutputsModel$Qsim
155
      }
156
157
    }
    
158
    return(list(PARAM = PARAM, SIM = SIM$ErrorCrit_KGE$NO$SIM, SIMold = .GlobalEnv$.ShinyGR.hist, Crit = CRIT))
159
160
    
  })
161
  
162
  
163
164
165
166
  
  ## --------------- Plot
  
  ## Choice
167
168
169
170
171
172
173
  getPlotType <- reactive({
    switch(input$PlotType,
           "Model performance" = 1,
           "Flow time series"  = 2,
           "State variables"   = 3,
           "Model diagram"     = 4)
  })
174
  
175
176
177
178
  
  ## Models available considering the plot type
  observe({
  if (getPlotType() == 4) {
179
      updateSelectInput(session, inputId = "HydroModel", choice = c("GR4J", "GR5J"), selected = input$HydroModel)
180
      updateSelectInput(session, inputId = "SnowModel" , choice = c("None"))
181
  } else {
182
      updateSelectInput(session, inputId = "HydroModel", choice = c("GR4J", "GR5J", "GR6J"), selected = input$HydroModel)
183
      updateSelectInput(session, inputId = "SnowModel" , choice = c("None", "CemaNeige")   , selected = input$SnowModel)
184
    }
185
186
  })
  
unknown's avatar
unknown committed
187
  
188
189
  ## Plots available considering the model type
  observe({
190
    if (input$HydroModel == "GR6J") {
191
192
193
194
195
196
197
      updateSelectInput(session, inputId = "PlotType",
                        choice = c("Flow time series", "Model performance", "State variables"),
                        selected = input$PlotType)
    } else {
      updateSelectInput(session, inputId = "PlotType",
                        choice = c("Flow time series", "Model performance", "State variables", "Model diagram"),
                        selected = input$PlotType)
198
    }
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
  })
  
  
  ## Period slider responds to changes in the selected/zoomed dateWindow 
  observeEvent({input$dyPlotTS_date_window ; input$dyPlotSVs_date_window ; input$dyPlotMDp_date_window}, {
    if (!is.null(input$dyPlotTS_date_window) && getPlotType() == 2) {
      dateWindow <- as.POSIXct(strftime(input$dyPlotTS_date_window, "%Y-%m-%d %H:%M:%S"), tz = "UTC")
    }
    if (!is.null(input$dyPlotSVs_date_window) && getPlotType() == 3) {
      dateWindow <- as.POSIXct(strftime(input$dyPlotSVs_date_window, "%Y-%m-%d %H:%M:%S"), tz = "UTC")
    }
    if (!is.null(input$dyPlotMDp_date_window) && getPlotType() == 4) {
      dateWindow <- as.POSIXct(strftime(input$dyPlotMDp_date_window, "%Y-%m-%d %H:%M:%S"), tz = "UTC")
    }
    if (exists("dateWindow")) {
      updateSliderInput(session, inputId = "Period",
215
                        value = dateWindow + .TypeModelGR(input$HydroModel)$TimeLag)
216
    }
unknown's avatar
unknown committed
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
  
  ## Reset period slider responds to dygraphs to mouse clicks
  observeEvent({input$dyPlotTS_click}, {
    updateSliderInput(session, inputId = "Period",
                      value = as.POSIXct(.ShinyGR.args$SimPer, tz = "UTC"))
  }, priority = +10)
  observeEvent({input$dyPlotSVs_click}, {
    updateSliderInput(session, inputId = "Period",
                      value = as.POSIXct(.ShinyGR.args$SimPer, tz = "UTC"))
  }, priority = +10)
  observeEvent({input$dyPlotSVq_click}, {
    updateSliderInput(session, inputId = "Period",
                      value = as.POSIXct(.ShinyGR.args$SimPer, tz = "UTC"))
  }, priority = +10)
  observeEvent({input$dyPlotMDp_click}, {
    updateSliderInput(session, inputId = "Period",
                      value = as.POSIXct(.ShinyGR.args$SimPer, tz = "UTC"))
  }, priority = +10)
  observeEvent({input$dyPlotMDe_click}, {
    updateSliderInput(session, inputId = "Period",
                      value = as.POSIXct(.ShinyGR.args$SimPer, tz = "UTC"))
  }, priority = +10)
  observeEvent({input$dyPlotMDq_click}, {
    updateSliderInput(session, inputId = "Period",
                      value = as.POSIXct(.ShinyGR.args$SimPer, tz = "UTC"))
  }, priority = +10)
  
  
unknown's avatar
unknown committed
247
248
249
  ## Target date slider
  observe({
    updateSliderInput(session, inputId = "Event",
250
                      min = input$Period[1L] + .TypeModelGR(input$HydroModel)$TimeLag,
unknown's avatar
unknown committed
251
252
253
254
                      max = input$Period[2L])
  })  
  
  
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
  ## Graphical parameters
  getPlotPar <- reactive({
    if (.GlobalEnv$.ShinyGR.args$theme == "Cyborg") {
      col_bg <- "black"
      col_fg <- "white"
      par(bg = col_bg, fg = col_fg, col.axis = col_fg, col.lab = col_fg)
    } else if (.GlobalEnv$.ShinyGR.args$theme == "Flatly") {
      col_bg <- "#2C3E50"
      col_fg <- "white"
      par(bg = col_bg, fg = col_fg, col.axis = col_bg, col.lab = col_bg)
    } else {
      col_bg <- "white"
      col_fg <- "black"
      par(bg = col_bg , fg = col_fg)
    }
    return(list(col_bg = col_bg, col_fg = col_fg, par = par(no.readonly = TRUE)))
  })
unknown's avatar
unknown committed
272

273
274
275
276
  
  ## Plot model performance
  output$stPlotMP <- renderPlot({
    OutputsModel <- getRES()$SIM$OutputsModel
277
    IndPlot <- which(OutputsModel$DatesR >= input$Period[1L] & OutputsModel$DatesR <= input$Period[2L])
278
279
    par(getPlotPar()$par)
    par(cex.axis = 1.2)
280
281
282
    if (input$SnowModel != "CemaNeige") {
      par(oma = c(20, 0, 0, 0))
    }
283
    plot(OutputsModel, Qobs = getRES()$SIM$Qobs, IndPeriod_Plot = IndPlot, cex.lab = 1.2, cex.axis = 1.4, cex.leg = 1.4)
284
  }, bg = "transparent")
285
  
unknown's avatar
unknown committed
286
  
287
288
289
290
291
292
293
  ## Plot flow time series
  output$dyPlotTS <- renderDygraph({
    op <- getPlotPar()$par
    dg1 <- dyplot(getRES()$SIM, RangeSelector = FALSE, LegendShow = "auto",
                  col.Q = c(op$fg, "orangered"), col.Precip = "#428BCA")
    dg1 <- dyOptions(dg1, axisLineColor = op$fg, axisLabelColor = op$fg, retainDateWindow = FALSE)
    dg1 <- dyLegend(dg1, show = "follow", width = 325)
294
295
  })
  
unknown's avatar
unknown committed
296
  
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
  ## Plot state variables stores
  output$dyPlotSVs <- renderDygraph({
    OutputsModel <- getRES()$SIM$OutputsModel
    data <- data.frame(DatesR = OutputsModel$DatesR,
                       prod.  = OutputsModel$Prod,
                       rout.  = OutputsModel$Rout)
    data.xts <- xts(data[, -1L], order.by = data$DatesR)
    
    op <- getPlotPar()$par
    dg2 <- dygraph(data.xts, group = "state_var", ylab = "store [mm]")
    dg2 <- dyOptions(dg2, colors = c("#00008B", "#008B8B"),
                     fillGraph = TRUE, fillAlpha = 0.3,
                     drawXAxis = FALSE, axisLineColor = op$fg, axisLabelColor = op$fg, retainDateWindow = FALSE)
    dg2 <- dyLegend(dg2, show = "always", width = 325)
    dg2 <- dyCrosshair(dg2, direction = "vertical")
312
313
  })
  
unknown's avatar
unknown committed
314
  
315
316
317
318
319
320
  ## Plot state variables Q
  output$dyPlotSVq <- renderDygraph({
    OutputsModel <- getRES()$SIM$OutputsModel
    IndPlot <- which(OutputsModel$DatesR >= input$Period[1L] & OutputsModel$DatesR <= input$Period[2L])
    OutputsModel2 <- sapply(OutputsModel[seq_len(which(names(OutputsModel) == "Qsim"))], function(x) x[IndPlot])
    OutputsModel2 <- c(OutputsModel2, Qobs = list(getRES()$SIM$Qobs[IndPlot]))
321
322
323

    data <- data.frame(DatesR = OutputsModel2$DatesR,
                       Qr     = OutputsModel2$QR,
324
                       Qd     = OutputsModel2$QD,
325
326
                       Qsim   = OutputsModel2$Qsim,
                       Qobs   = OutputsModel2$Qobs)
327
    if (input$HydroModel == "GR6J") {
328
      data$QrExp <- OutputsModel2$QRExp
329
    } else {
330
      data$QrExp <- NA
331
    }
332
333
334
335
    data.xts <- xts(data[, -1L], order.by = data$DatesR)
    
    op <- getPlotPar()$par
    dg3 <- dygraph(data.xts, group = "state_var", ylab = "flow [mm/d]", main = " ")
336
    dg3 <- dyOptions(dg3, fillAlpha = 1.0,
337
338
                     axisLineColor = op$fg, axisLabelColor = op$fg,
                     titleHeight = 10, retainDateWindow = FALSE)
339
    dg3 <- dyStackedRibbonGroup(dg3, name = c("Qd", "Qr", "QrExp"),
340
                                color = c("#FFD700", "#EE6300", "brown"), strokeBorderColor = "black")
341
342
    dg3 <- dySeries(dg3, name = "Qobs", fillGraph = FALSE, drawPoints = TRUE, color = op$fg)
    dg3 <- dySeries(dg3, name = "Qsim", fillGraph = FALSE, color = "orangered")
343
    dg3 <- dyCrosshair(dg3, direction = "vertical")
344
    dg3 <- dyLegend(dg3, show = "always", width = 325)
345
  })
346
  
unknown's avatar
unknown committed
347
  
348
349
350
351
352
353
354
355
356
357
  ## Plot model diagram precipitation
  output$dyPlotMDp <- renderDygraph({
    barChartPrecip <- scan(file = system.file("plugins/barChartPrecip.js", package = "airGRteaching"),
                           what = "character", quiet = TRUE)
    data <- data.frame(DatesR  = getRES()$SIM$OutputsModel$DatesR,
                       precip. = getRES()$SIM$OutputsModel$Precip)
    data.xts <- xts(data[, -1L, drop = FALSE], order.by = data$DatesR)
    
    dg4 <- dygraph(data.xts, group = "mod_diag", ylab = "precip. [mm/d]")
    dg4 <- dyOptions(dg4, colors = "#428BCA", drawXAxis = FALSE, plotter = barChartPrecip, retainDateWindow = FALSE)
358
    dg4 <- dyAxis(dg4, name = "y", valueRange = c(max(data.xts[, "precip."], na.rm = TRUE), -1e-3))
unknown's avatar
unknown committed
359
    dg4 <- dyEvent(dg4, input$Event, color = "orangered")
360
361
    dg4 <- dyLegend(dg4, show = "onmouseover", width = 225)
    dg4 <- dyCrosshair(dg4, direction = "vertical")
362
  })
363
  
unknown's avatar
unknown committed
364
  
365
366
367
368
  ## Plot model diagram ETP
  output$dyPlotMDe <- renderDygraph({
    op <- getPlotPar()$par
    data <- data.frame(DatesR = getRES()$SIM$OutputsModel$DatesR,
369
                       PET    = getRES()$SIM$OutputsModel$PotEvap)
370
371
    data.xts <- xts(data[, -1L, drop = FALSE], order.by = data$DatesR)
    
372
    dg5 <- dygraph(data.xts, group = "mod_diag", ylab = "PET [mm/d]", main = " ")
373
374
375
376
    dg5 <- dyOptions(dg5, colors = "#A4C400", drawPoints = TRUE,
                     strokeWidth = 0, pointSize = 2, drawXAxis = FALSE,
                     axisLineColor = op$fg, axisLabelColor = op$fg,
                     titleHeight = 10, retainDateWindow = FALSE)
unknown's avatar
unknown committed
377
    dg5 <- dyEvent(dg5, input$Event, color = "orangered")
378
379
    dg5 <- dyLegend(dg5, show = "onmouseover", width = 225)
    dg5 <- dyCrosshair(dg5, direction = "vertical")
380
  })
381
  
unknown's avatar
unknown committed
382
  
383
384
  ## Plot model diagram flow
  output$dyPlotMDq <- renderDygraph({
385
386
    if (length(.GlobalEnv$.ShinyGR.hist) == 2 & input$ShowOldQsim == "Yes") {
      QsimOld <- getRES()$SIMold[[1]]$Qsim
387
388
389
    } else {
      QsimOld <- NA
    }
390
391
392
393
    OutputsModel <- getRES()$SIM$OutputsModel
    IndPlot <- which(OutputsModel$DatesR >= input$Period[1L] & OutputsModel$DatesR <= input$Period[2L])
    OutputsModel2 <- sapply(OutputsModel[seq_len(which(names(OutputsModel) == "Qsim"))], function(x) x[IndPlot])
    OutputsModel2 <- c(OutputsModel2, Qobs = list(getRES()$SIM$Qobs[IndPlot]))
unknown's avatar
unknown committed
394
    OutputsModel2$Qsim <- ifelse(format(OutputsModel2$DatesR, "%Y%m%d") > format(input$Event, "%Y%m%d"), NA, OutputsModel2$Qsim)
395
    OutputsModel2$Qold <- ifelse(format(OutputsModel2$DatesR, "%Y%m%d") > format(input$Event, "%Y%m%d"), NA, QsimOld[IndPlot])
396

397
398
399
400
    data <- data.frame(DatesR  = OutputsModel2$DatesR,
                       Qobs    = OutputsModel2$Qobs,
                       Qsim    = OutputsModel2$Qsim,
                       QsimOld = OutputsModel2$Qold)
401
402
403
404
    data.xts <- xts(data[, -1L, drop = FALSE], order.by = data$DatesR)
    
    op <- getPlotPar()$par
    dg6 <- dygraph(data.xts, group = "mod_diag", ylab = "flow [mm/d]", main = " ")
405
    dg6 <- dyOptions(dg6, colors = c(op$fg, "grey", "orangered"), drawPoints = TRUE,
406
407
                     axisLineColor = op$fg, axisLabelColor = op$fg,
                     titleHeight = 10, retainDateWindow = FALSE)
408
409
    dg6 <- dySeries(dg6, name = "QsimOld", drawPoints = FALSE, strokePattern = "dashed")
    dg6 <- dySeries(dg6, name = "Qsim"   , drawPoints = FALSE)
unknown's avatar
unknown committed
410
    dg6 <- dyEvent(dg6, input$Event, color = "orangered")
411
412
    dg6 <- dyLegend(dg6, show = "onmouseover", width = 225)
    dg6 <- dyCrosshair(dg6, direction = "vertical")
413
  })
414
  
unknown's avatar
unknown committed
415
  
416
417
418
419
420
421
  ## Plot model diagram chart
  output$stPlotMD <- renderPlot({
    OutputsModel <- getRES()$SIM$OutputsModel
    IndPlot <- which(OutputsModel$DatesR >= input$Period[1L] & OutputsModel$DatesR <= input$Period[2L])
    OutputsModel2 <- sapply(OutputsModel[seq_len(which(names(OutputsModel) == "Qsim"))], function(x) x[IndPlot])
    OutputsModel2 <- c(OutputsModel2, Qobs = list(getRES()$SIM$Qobs[IndPlot]))
422

423
424
425
    par(getPlotPar()$par)
    airGRteaching:::DiagramGR(OutputsModel = OutputsModel2, Param = getRES()$PARAM,
                              SimPer = input$Period, EventDate = input$Event,
426
                              HydroModel = input$HydroModel)
427
  }, bg = "transparent")
428
429
  
  
430
  
431
432
  ## --------------- Criteria table
  
433
  output$Criteria <- renderTable({
434
435
    
    ## Table created in order to choose order the criteria in the table output
436
    tabCrit_gauge <- data.frame(Criterion = c("NSE [Q]", "NSE [sqrt(Q)]", "NSE [log(Q)]",
437
                                              "KGE [Q]", "KGE [sqrt(Q)]", "KGE [log(Q)]"),
438
                                ID        = 1:6, stringsAsFactors = FALSE)
439
    
440
    tabCrit_out <- merge(tabCrit_gauge, getRES()$Crit, by = "Criterion", all.x = TRUE)
441
442
    tabCrit_out <- tabCrit_out[order(tabCrit_out$ID), ]
    tabCrit_out <- tabCrit_out[, !colnames(tabCrit_out) %in% "ID"]
443
444
445
    
    ## Color the cell of the crietaia uses during the calibration
    if (CAL_click$valueButton >= 0) {
446
      CellCol <- '<div style="color: #FFFFFF; background-color: #A4C400; border: 5px solid #A4C400; position:relative; top: 0px; left: 5px; padding: 0px; margin: -5px -0px -8px -10px;">
447
<span>9999</span></div>'
448
      CellCol_id  <- which(tabCrit_out$Criterion == input$TypeCrit)
449
450
451
452
453
      tabCrit_out[CellCol_id, 1] <- gsub("9999", tabCrit_out[CellCol_id, 1], CellCol)
    }
    
    return(tabCrit_out)
  }, sanitize.text.function = function(x) x)
454
455
456
457
458
459
460
461
462
  
  
  
  ## --------------- Download buttons
  
  ## simulation table
  output$DownloadTab <- downloadHandler(
    filename = function() {
      filename <- "TabSim"
463
      filename <- sprintf("airGR_%s_%s.csv", filename, gsub("(.*)( )(\\d{2})(:)(\\d{2})(:)(\\d{2})", "\\1_\\3h\\5m\\7s", Sys.time()))
464
465
466
467
468
469
470
    },
    content = function(file) {
      OBS <- getPrep()$OBS
      SIM <- getRES()$SIM
      if (input$SnowModel != "CemaNeige") {
        PrecipSim <- NA
        FracSolid <- NA
471
        TempMean  <- NA
472
473
474
475
476
477
478
      } else {
        PrecipSol <- rowMeans(as.data.frame(OBS$InputsModel$LayerPrecip) * as.data.frame(OBS$InputsModel$LayerFracSolidPrecip), na.rm = TRUE)
        PrecipSim <- rowMeans(as.data.frame(OBS$InputsModel$LayerPrecip), na.rm = TRUE)
        FracSolid <- PrecipSol / PrecipSim
        FracSolid <- ifelse(is.na(FracSolid)  & PrecipSol == 0 & PrecipSim == 0, 0, FracSolid)
        PrecipSim <- PrecipSim[SIM$OptionsSimul$IndPeriod_Run]
        FracSolid <- FracSolid[SIM$OptionsSimul$IndPeriod_Run]
479
480
481
        FracSolid <- round(FracSolid, digits = 3)
        TempMean  <- rowMeans(as.data.frame(OBS$InputsModel$LayerTempMean), na.rm = TRUE)
        TempMean  <- TempMean[SIM$OptionsSimul$IndPeriod_Run]
482
      }
483
484
485
486
487
488
489
490
491
492
493
494
      TabSim <- data.frame(Dates                     = SIM$OutputsModel$DatesR,
                           PotEvap                   = SIM$OutputsModel$PotEvap,
                           PrecipObs                 = SIM$OutputsModel$Precip,
                           PrecipSim_CemaNeige       = PrecipSim,
                           PrecipFracSolid_CemaNeige = FracSolid,
                           TempMeanSim_CemaNeige     = TempMean,
                           Qobs                      = SIM$OptionsCrit$Qobs,
                           Qsim                      = SIM$OutputsModel$Qsim)
      colnames(TabSim) <- sprintf("%s [%s]", colnames(TabSim), c("-", rep("mm",3), "-", "°C", rep("mm", 2)))
      colnames(TabSim) <- ifelse(grepl("mm", colnames(TabSim)),
                                 gsub("mm", paste0("mm/", .TypeModelGR(OBS)$TimeUnit),colnames(TabSim)),
                                 colnames(TabSim))
495
496
497
498
499
500
501
502
503
504
505
506
      write.table(TabSim, file = file, row.names = FALSE, sep = ";")
    }
  )
  
  ## plots
  output$DownloadPlot <- downloadHandler(
    filename = function() {
      filename <- switch(input$PlotType,
                            "Model performance" = "PlotModelPerf",
                            "Flow time series"  = "PlotFlowTimeSeries",
                            "State variables"   = "PlotStateVar",
                            "Model diagram"     = "PlotModelDiag")
507
      filename <- sprintf("airGR_%s_%s.png", filename, gsub("(.*)( )(\\d{2})(:)(\\d{2})(:)(\\d{2})", "\\1_\\3h\\5m\\7s", Sys.time()))
508
509
510
    },
    content = function(file) {
      k <- 1.75
511
512
513
514
515
516
517
518
519
      ParamTitle <- c("X1", "X2", "X3", "X4", "X5", "X6")[seq_len(getPrep()$TMGR$NbParam)]
      if (input$SnowModel == "CemaNeige") {
        ParamTitle <- c(ParamTitle, "C1", "C2")
      }
      ParamTitle <- paste(ParamTitle, getRES()$PARAM, sep = " = ", collapse = ", ")
      PngTitle <- sprintf("%s - %s/%s\n%s\n%s", input$Dataset,
                          input$HydroModel, ifelse(input$SnowModel == "CemaNeige", "CemaNeige", "No snow model"),
                          paste0(input$Period, collapse = " - "),
                          ParamTitle)
520
521
      if (getPlotType() == 1) {
        png(filename = file, width = 1000*k,  height = ifelse(input$SnowModel != "CemaNeige", 700*k, 1100*k), pointsize = 14, res = 150)
522
        par(oma = c(0, 0, 4, 0))
523
        plot(getRES()$SIM)
524
        mtext(text = PngTitle, side = 3, outer = TRUE, cex = 0.8, line = 1.2)
525
526
527
528
        dev.off()
      }
      if (getPlotType() == 2) {
        png(filename = file, width = 1000*k, height = 600*k, pointsize = 14, res = 150)
529
        par(oma = c(0, 0, 4, 0))
530
        plot(getRES()$SIM, which = c( "Precip", "Flows"))
531
        mtext(text = PngTitle, side = 3, outer = TRUE, cex = 0.8, line = 1.2)
532
533
534
535
536
537
        dev.off()
      }
    }
  )
  
  
538
539
})