From 054ec7dfc15d8502c8c382e84f090bdcb7bec5cb Mon Sep 17 00:00:00 2001
From: Dorchies David <david.dorchies@irstea.fr>
Date: Fri, 29 May 2020 20:20:50 +0200
Subject: [PATCH] feat: Calibration of GR-IWRM natural model

- Add all functions for calibration
- vignette for calibration in context of natural catchment

Refs #3
---
 NAMESPACE                                | 10 ++++-
 R/Calibration.GriwrmInputsModel.R        | 56 ++++++++++++++++++++++++
 R/Calibration.InputsModel.R              | 10 +++++
 R/Calibration.R                          | 10 +++++
 R/CreateCalibOptions.GriwrmInputsModel.R | 19 ++++++++
 R/CreateCalibOptions.InputsModel.R       | 15 +++++++
 R/CreateCalibOptions.R                   | 10 +++++
 R/CreateInputsCrit.GriwrmInputsModel.R   | 29 ++++++++++++
 R/CreateInputsCrit.InputsModel.R         | 16 +++++++
 R/CreateInputsCrit.R                     | 10 +++++
 R/CreateInputsModel.Griwrm.R             |  2 +-
 R/UpdateQsimUpstream.R                   |  6 +--
 vignettes/V01_First_network.Rmd          |  2 +-
 vignettes/V03_First_Calibration.Rmd      | 46 +++++++++++++++++++
 vignettes/v02_First_run.Rmd              |  8 +++-
 15 files changed, 241 insertions(+), 8 deletions(-)
 create mode 100644 R/Calibration.GriwrmInputsModel.R
 create mode 100644 R/Calibration.InputsModel.R
 create mode 100644 R/Calibration.R
 create mode 100644 R/CreateCalibOptions.GriwrmInputsModel.R
 create mode 100644 R/CreateCalibOptions.InputsModel.R
 create mode 100644 R/CreateCalibOptions.R
 create mode 100644 R/CreateInputsCrit.GriwrmInputsModel.R
 create mode 100644 R/CreateInputsCrit.InputsModel.R
 create mode 100644 R/CreateInputsCrit.R
 create mode 100644 vignettes/V03_First_Calibration.Rmd

diff --git a/NAMESPACE b/NAMESPACE
index 4a6bcd8..32655fd 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -1,7 +1,11 @@
 # Generated by roxygen2: do not edit by hand
 
-S3method(Calibration,default)
-S3method(Calibration,griwrm)
+S3method(Calibration,GriwrmInputsModel)
+S3method(Calibration,InputsModel)
+S3method(CreateCalibOptions,GriwrmInputsModel)
+S3method(CreateCalibOptions,InputsModel)
+S3method(CreateInputsCrit,GriwrmInputsModel)
+S3method(CreateInputsCrit,InputsModel)
 S3method(CreateInputsModel,Griwrm)
 S3method(CreateInputsModel,default)
 S3method(CreateRunOptions,GriwrmInputsModel)
@@ -10,6 +14,8 @@ S3method(RunModel,GriwrmInputsModel)
 S3method(RunModel,InputsModel)
 S3method(merge,Gits)
 export(Calibration)
+export(CreateCalibOptions)
+export(CreateInputsCrit)
 export(CreateInputsModel)
 export(CreateRunOptions)
 export(Ginet)
diff --git a/R/Calibration.GriwrmInputsModel.R b/R/Calibration.GriwrmInputsModel.R
new file mode 100644
index 0000000..e580fe8
--- /dev/null
+++ b/R/Calibration.GriwrmInputsModel.R
@@ -0,0 +1,56 @@
+#' Calibration of a semi-distributed run-off model
+#'
+#' @param InputsModel object of class \emph{GriwrmInputsModel}, see \code{\link{CreateInputsModel.Griwrm}} for details.
+#' @param RunOptions object of class \emph{GriwrmRunOptions}, see \code{\link{CreateRunOptiosn.Griwrm}} for details.
+#' @param InputsCrit object of class \emph{GriwrmInputsCrit}, see \code{\link{CreateInputsCrit.Griwrm}} for details.
+#' @param CalibOptions object of class \emph{GriwrmCalibOptions}, see \code{\link{CreateCalibOptions.Griwrm}} for details.
+#' @param useUpstreamQsim boolean describing if simulated (\code{TRUE}) or observed (\code{FALSE}) flows are used for calibration. Default is \code{TRUE}.
+#' @param verbose (optional) boolean indicating if the function is run in verbose mode or not, default = \code{TRUE}
+#' @param ... further arguments passed to \code{\link[airGR]{Calibration}}.
+#'
+#' @return GriwrmOutputsCalib object which is a list of OutputsCalib (See \code{\link[airGR]{Calibration}}) for each node of the semi-distributed model.
+#' @export
+Calibration.GriwrmInputsModel <- function(InputsModel,
+                                          RunOptions,
+                                          InputsCrit,
+                                          CalibOptions,
+                                          useUpstreamQsim = TRUE,
+                                          verbose = TRUE,
+                                          ...) {
+
+  OutputsCalib <- list()
+  class(OutputsCalib) <- append(class(OutputsCalib), "GriwrmOutputsCalib")
+
+  OutputsModel <- list()
+  class(OutputsModel) <- append(class(OutputsModel), "GriwrmOutputsModel")
+
+  for(IM in InputsModel) {
+    if(verbose) cat("Calibration.GriwrmInputsModel: Treating sub-basin", IM$id, "...\n")
+
+    if(useUpstreamQsim) {
+      # Update InputsModel$QobsUpstr with simulated upstream flows
+      IM <- UpdateQsimUpstream(IM, OutputsModel)
+    }
+
+    OutputsCalib[[IM$id]] <- Calibration.InputsModel(
+      InputsModel = IM,
+      RunOptions = RunOptions[[IM$id]],
+      InputsCrit = InputsCrit[[IM$id]],
+      CalibOptions = CalibOptions[[IM$id]],
+      ...
+    )
+
+    if(useUpstreamQsim) {
+      # Run the model for the sub-basin
+      OutputsModel[[IM$id]] <- RunModel(
+        InputsModel = IM,
+        RunOptions = RunOptions[[IM$id]],
+        Param = OutputsCalib[[IM$id]]$ParamFinalR
+      )
+    }
+
+  }
+
+  return(OutputsCalib)
+
+}
diff --git a/R/Calibration.InputsModel.R b/R/Calibration.InputsModel.R
new file mode 100644
index 0000000..6520ab8
--- /dev/null
+++ b/R/Calibration.InputsModel.R
@@ -0,0 +1,10 @@
+#' Wrapper to \code{\link[airGR]{Calibration}}.
+#'
+#' @param InputsModel object of class \emph{InputsModel}, see \code{\link[airGR]{CreateInputsModel}} for details.
+#' @param ... further arguments passed to \code{\link[airGR]{Calibration}}.
+#'
+#' @return \emph{CalibOutput} object.
+#' @export
+Calibration.InputsModel <- function(InputsModel, ...) {
+  airGR::Calibration(InputsModel, FUN_MOD = InputsModel$FUN_MOD, ...)
+}
diff --git a/R/Calibration.R b/R/Calibration.R
new file mode 100644
index 0000000..04b76e1
--- /dev/null
+++ b/R/Calibration.R
@@ -0,0 +1,10 @@
+#' Calibration of either airGR model and GRIWRM semi-distributive model
+#'
+#' @param InputsModel the class of the first parameter determine which calibration is used
+#' @param ... further arguments passed to or from other methods.
+#'
+#' @return \emph{OutputsCalib} or \emph{GriwrmOutputsCalib} object
+#' @export
+Calibration <- function(InputsModel, ...) {
+  UseMethod("Calibration", InputsModel)
+}
diff --git a/R/CreateCalibOptions.GriwrmInputsModel.R b/R/CreateCalibOptions.GriwrmInputsModel.R
new file mode 100644
index 0000000..e92e48e
--- /dev/null
+++ b/R/CreateCalibOptions.GriwrmInputsModel.R
@@ -0,0 +1,19 @@
+#' Title
+#'
+#' @param InputsModel object of class \emph{GriwrmInputsModel}, see \code{\link{CreateInputsModel.Griwrm}} for details.
+#' @param ... further arguments passed to \code{\link[airGR]{CreateCalibOptions}}.
+#'
+#' @return \emph{GriwrmCalibOptions} object.
+#' @export
+CreateCalibOptions.GriwrmInputsModel <- function(InputsModel, ...) {
+
+  CalibOptions <- list()
+
+  for(IM in InputsModel) {
+    CalibOptions[[IM$id]] <- CreateCalibOptions.InputsModel(
+      InputsModel = IM,
+      ...
+    )
+  }
+  return(CalibOptions)
+}
diff --git a/R/CreateCalibOptions.InputsModel.R b/R/CreateCalibOptions.InputsModel.R
new file mode 100644
index 0000000..0099cf7
--- /dev/null
+++ b/R/CreateCalibOptions.InputsModel.R
@@ -0,0 +1,15 @@
+#' Wrapper to \code{\link[airGR]{CreateCalibOptions}}
+#'
+#' @param InputsModel object of class \emph{InputsModel}, see \code{\link[airGR]{CreateInputsModel}} for details.
+#' @param ... further arguments passed to \code{\link[airGR]{CreateCalibOptions}}.
+#'
+#' @return \emph{CalibOptions} object.
+#' @export
+CreateCalibOptions.InputsModel <- function(InputsModel,
+                               ...) {
+  airGR::CreateCalibOptions(
+    FUN_MOD = InputsModel$FUN_MOD,
+    IsSD = !is.null(InputsModel$QobsUpstr),
+    ...
+  )
+}
diff --git a/R/CreateCalibOptions.R b/R/CreateCalibOptions.R
new file mode 100644
index 0000000..4ef9311
--- /dev/null
+++ b/R/CreateCalibOptions.R
@@ -0,0 +1,10 @@
+#' CreateCalibOptions both available for \emph{InputsModel} and \emph{GrwirmInputsModel} objects
+#'
+#' @param InputsModel object of class \emph{InputsModel}, see \code{\link[airGR]{CreateInputsModel}} for details.
+#' @param ... further arguments passed to or from other methods.
+#'
+#' @return
+#' @export
+CreateCalibOptions <- function(InputsModel, ...) {
+  UseMethod("CreateCalibOptions", InputsModel)
+}
diff --git a/R/CreateInputsCrit.GriwrmInputsModel.R b/R/CreateInputsCrit.GriwrmInputsModel.R
new file mode 100644
index 0000000..db1769f
--- /dev/null
+++ b/R/CreateInputsCrit.GriwrmInputsModel.R
@@ -0,0 +1,29 @@
+#' Create \emph{GriwrmInputsCrit} object for GR-IWRM.
+#' @param InputsModel  object of class \emph{GriwrmInputsModel}, see \code{\link{CreateInputsModel.Griwrm}} for details.
+#' @param FUN_CRIT \[function (atomic or list)\] error criterion function (e.g. \code{\link[airGR]{ErrorCrit_RMSE}}, \code{\link[airGR]{ErrorCrit_NSE}})
+#' @param RunOptions object of class \emph{GriwrmRunOptions}, see \code{[CreateRunOptions.Griwrm]} for details.
+#' @param gits object of class \emph{Gits}, see [Gits].
+#' @param ... further arguments passed to \code{\link[airGR]{CreateInputsCrit}}.
+#'
+#' @return Object of class \emph{GriwrmInputsCrit}
+#' @export
+CreateInputsCrit.GriwrmInputsModel <- function(InputsModel,
+                                               FUN_CRIT = airGR::ErrorCrit_NSE,
+                                               RunOptions,
+                                               gits,
+                                               ...) {
+  InputsCrit <- list()
+  class(InputsCrit) <- append(class(InputsCrit), "GriwrmInputsCrit")
+
+  for(IM in InputsModel) {
+    InputsCrit[[IM$id]] <- CreateInputsCrit.InputsModel(
+      InputsModel = IM,
+      FUN_CRIT = FUN_CRIT,
+      RunOptions = RunOptions[[IM$id]],
+      Obs = gits[[IM$id]]$Qobs[RunOptions[[IM$id]]$IndPeriod_Run],
+      ...
+    )
+  }
+
+  return(InputsCrit)
+}
diff --git a/R/CreateInputsCrit.InputsModel.R b/R/CreateInputsCrit.InputsModel.R
new file mode 100644
index 0000000..f1a337b
--- /dev/null
+++ b/R/CreateInputsCrit.InputsModel.R
@@ -0,0 +1,16 @@
+#' Wrapper to \code{\link[airGR]{CreateInputsCrit}}
+#'
+#' @param InputsModel object of class \emph{InputsModel}, see \code{\link[airGR]{CreateInputsModel}} for details.
+#' @param FUN_CRIT \[function (atomic or list)\] error criterion function (e.g. \code{\link[airGR]{ErrorCrit_RMSE}}, \code{\link[airGR]{ErrorCrit_NSE}})
+#' @param ... further arguments passed to \code{\link[airGR]{CreateInputsCrit}}
+#'
+#' @return object of class \emph{InputsCrit} containing the data required to evaluate the model outputs. See \code{\link[airGR]{CreateInputsCrit}}
+#' @export
+CreateInputsCrit.InputsModel <- function(InputsModel,
+                                         FUN_CRIT,
+                                         ...) {
+
+  airGR::CreateInputsCrit(FUN_CRIT = FUN_CRIT,
+                          InputsModel = InputsModel,
+                          ...)
+}
diff --git a/R/CreateInputsCrit.R b/R/CreateInputsCrit.R
new file mode 100644
index 0000000..bce0108
--- /dev/null
+++ b/R/CreateInputsCrit.R
@@ -0,0 +1,10 @@
+#' Title
+#'
+#' @param InputsModel InputsModel for GR-IWRM (See \code{[CreateInputsModel.Griwrm]}) or AirGR (See \code{\link[airGR]{CreateInputsModel}})
+#' @param ... further arguments passed to or from other methods.
+#'
+#' @return
+#' @export
+CreateInputsCrit <- function(InputsModel, ...) {
+  UseMethod("CreateInputsCrit", InputsModel)
+}
diff --git a/R/CreateInputsModel.Griwrm.R b/R/CreateInputsModel.Griwrm.R
index 96e0c85..52ddce9 100644
--- a/R/CreateInputsModel.Griwrm.R
+++ b/R/CreateInputsModel.Griwrm.R
@@ -8,7 +8,7 @@
 #'
 #' @return GriwrmInputsModel object equivalent to airGR InputsModel object for a semi-distributed model (See \code{\link[airGR]{CreateInputsModel}})
 #' @export
-CreateInputsModel.Griwrm <- function(x, girop, gits, verbose = TRUE,...) {
+CreateInputsModel.Griwrm <- function(x, girop, gits, verbose = TRUE, ...) {
 
   InputsModel <- CreateEmptyGriwrmInputsModel()
 
diff --git a/R/UpdateQsimUpstream.R b/R/UpdateQsimUpstream.R
index 560984a..2bccce8 100644
--- a/R/UpdateQsimUpstream.R
+++ b/R/UpdateQsimUpstream.R
@@ -1,19 +1,19 @@
 #' Update InputsModel$QobsUpstr with simulated upstream flows provided by GriwrmOutputsModels object.
 #'
 #' @param InputsModel \emph{GriwrmInputsModel} object. See \code{[CreateInputsModel.Griwrm]}.
-#' @param OutputsModels \emph{GriwrmOutputsModel} object provided by \code{[RunModel.GriwrmInputsModel]}.
+#' @param OutputsModel \emph{GriwrmOutputsModel} object provided by \code{[RunModel.GriwrmInputsModel]}.
 #'
 #' @description This function is used by \code{\link{RunModel.GriwrmInputsModel}} and \code{\link{Calibration.GriwrmInputsModel}} in order to provide upstream simulated flows to a node.
 #'
 #' @return InputsModel object with updated QobsUpsr
 #'
-UpdateQsimUpstream <- function(InputsModel, OutputsModels) {
+UpdateQsimUpstream <- function(InputsModel, OutputsModel) {
   if(length(InputsModel$UpstreamNodes) > 0) {
     for(i in 1:length(InputsModel$UpstreamNodes)) {
       QobsUpstr1 <- matrix(
         c(
           rep(0, length(RunOptions[[InputsModel$id]]$IndPeriod_WarmUp)),
-          OutputsModels[[InputsModel$UpstreamNodes[i]]]$Qsim
+          OutputsModel[[InputsModel$UpstreamNodes[i]]]$Qsim
         ), ncol = 1
       )
       if(i == 1) {
diff --git a/vignettes/V01_First_network.Rmd b/vignettes/V01_First_network.Rmd
index d2fc99f..f834388 100644
--- a/vignettes/V01_First_network.Rmd
+++ b/vignettes/V01_First_network.Rmd
@@ -130,6 +130,6 @@ InputsModel <- CreateInputsModel(ginet, girop, gits)
 
 ```{r}
 dir.create("_cache", showWarnings = FALSE)
-save(ginet, girop, gits, InputsModel, file = "_cache/seine.RData")
+save(ginet, girop, gits, InputsModel, file = "_cache/V01.RData")
 ```
 
diff --git a/vignettes/V03_First_Calibration.Rmd b/vignettes/V03_First_Calibration.Rmd
new file mode 100644
index 0000000..7095a11
--- /dev/null
+++ b/vignettes/V03_First_Calibration.Rmd
@@ -0,0 +1,46 @@
+---
+title: "03_First_Calibration"
+author: "David Dorchies"
+vignette: >
+  %\VignetteIndexEntry{Calibration of naturalised semi-distributive model}
+  %\VignetteEngine{knitr::rmarkdown}
+  %\VignetteEncoding{UTF-8}
+---
+
+```{r setup, include=FALSE}
+knitr::opts_chunk$set(echo = TRUE)
+```
+
+## Loading network and time series data
+
+Run `vignette("01_First_network", package = "griwrm")` and `vignette("02_First_run", package = "griwrm")` before this one in order to create the Rdata files loaded below:
+
+```{r}
+load("_cache/V01.RData")
+load("_cache/V02.RData")
+```
+
+## InputsCrit object
+
+```{r cars}
+InputsCrit <- CreateInputsCrit(
+  InputsModel = InputsModel, 
+  RunOptions = RunOptions, gits = gits
+)
+str(InputsCrit)
+```
+
+## GriwrmCalibOptions object
+
+```{r, eval=FALSE}
+CalibOptions <- CreateCalibOptions(InputsModel)
+str(CalibOptions)
+```
+
+## Calibration
+
+```{r}
+OutputsCalib <- Calibration(InputsModel, RunOptions, InputsCrit, CalibOptions)
+```
+
+
diff --git a/vignettes/v02_First_run.Rmd b/vignettes/v02_First_run.Rmd
index d65efb7..ae3050f 100644
--- a/vignettes/v02_First_run.Rmd
+++ b/vignettes/v02_First_run.Rmd
@@ -25,7 +25,7 @@ library(griwrm)
 Run `vignette("01_First_network", package = "griwrm")` before this one in order to create the Rdata file loaded below:
 
 ```{r}
-load("_cache/seine.RData")
+load("_cache/V01.RData")
 ```
 
 ### Loading 
@@ -111,6 +111,12 @@ OutputsModels <- RunModel(
 )
 ```
 
+## Save data for next vignettes
+
+```{r}
+save(RunOptions, file = "_cache/V02.RData")
+```
+
 ## Plot the result for each basin
 
 ```{r, fig.height = 5, fig.width = 8}
-- 
GitLab