From ae70dd4bde96e0ec9360b1786cbecc028ea0e374 Mon Sep 17 00:00:00 2001
From: Dorchies David <david.dorchies@inrae.fr>
Date: Sun, 4 Apr 2021 12:20:52 +0200
Subject: [PATCH] fix: Allow to calibrate RunModel_Lag alone

Refs #108
---
 R/Calibration_Michel.R             |  9 ++++----
 R/CreateCalibOptions.R             | 16 +++++++++++++-
 tests/testthat/test-RunModel_Lag.R | 34 +++++++++++++++++++++++-------
 3 files changed, 46 insertions(+), 13 deletions(-)

diff --git a/R/Calibration_Michel.R b/R/Calibration_Michel.R
index 69a3da22..0dd1eec5 100644
--- a/R/Calibration_Michel.R
+++ b/R/Calibration_Michel.R
@@ -5,7 +5,8 @@ Calibration_Michel <- function(InputsModel,
                                FUN_MOD,
                                FUN_CRIT,           # deprecated
                                FUN_TRANSFO = NULL,
-                               verbose = TRUE) {
+                               verbose = TRUE,
+                               ...) {
 
 
   FUN_MOD  <- match.fun(FUN_MOD)
@@ -145,7 +146,7 @@ Calibration_Michel <- function(InputsModel,
     }
     ##Model_run
     Param <- CandidatesParamR[iNew, ]
-    OutputsModel <- RunModel(InputsModel, RunOptions, Param, FUN_MOD = FUN_MOD)
+    OutputsModel <- RunModel(InputsModel, RunOptions, Param, FUN_MOD = FUN_MOD, ...)
 
     ##Calibration_criterion_computation
     OutputsCrit <- ErrorCrit(InputsCrit, OutputsModel, verbose = FALSE)
@@ -294,7 +295,7 @@ Calibration_Michel <- function(InputsModel,
     for (iNew in 1:nrow(CandidatesParamR)) {
       ##Model_run
       Param <- CandidatesParamR[iNew, ]
-      OutputsModel <- RunModel(InputsModel, RunOptions, Param, FUN_MOD = FUN_MOD)
+      OutputsModel <- RunModel(InputsModel, RunOptions, Param, FUN_MOD = FUN_MOD, ...)
       ##Calibration_criterion_computation
       OutputsCrit <- ErrorCrit(InputsCrit, OutputsModel, verbose = FALSE)
       if (!is.na(OutputsCrit$CritValue)) {
@@ -355,7 +356,7 @@ Calibration_Michel <- function(InputsModel,
       CandidatesParamR <- FUN_TRANSFO(CandidatesParamT, "TR")
       ##Model_run
       Param <- CandidatesParamR[iNew, ]
-      OutputsModel <- RunModel(InputsModel, RunOptions, Param, FUN_MOD = FUN_MOD)
+      OutputsModel <- RunModel(InputsModel, RunOptions, Param, FUN_MOD = FUN_MOD, ...)
       ##Calibration_criterion_computation
       OutputsCrit <- ErrorCrit(InputsCrit, OutputsModel, verbose = FALSE)
       if (OutputsCrit$CritValue * OutputsCrit$Multiplier < CritOptim) {
diff --git a/R/CreateCalibOptions.R b/R/CreateCalibOptions.R
index 78302745..1d13f2ae 100644
--- a/R/CreateCalibOptions.R
+++ b/R/CreateCalibOptions.R
@@ -76,6 +76,13 @@ CreateCalibOptions <- function(FUN_MOD,
     ObjectClass <- c(ObjectClass, "CemaNeigeGR6J")
     BOOL <- TRUE
   }
+  if (identical(FUN_MOD, RunModel_Lag)) {
+    ObjectClass <- c(ObjectClass, "Lag")
+    if (IsSD) {
+      stop("RunModel_Lag should not be used with 'isSD=TRUE'")
+    }
+    BOOL <- TRUE
+  }
   if (IsHyst) {
     ObjectClass <- c(ObjectClass, "hysteresis")
   }
@@ -136,6 +143,9 @@ CreateCalibOptions <- function(FUN_MOD,
         FUN_GR <- TransfoParam_CemaNeige
       }
     }
+    if (identical(FUN_MOD, RunModel_Lag)) {
+      FUN_GR <- TransfoParam_Lag
+    }
     if (is.null(FUN_GR)) {
       stop("'FUN_GR' was not found")
       return(NULL)
@@ -151,7 +161,7 @@ CreateCalibOptions <- function(FUN_MOD,
       FUN_LAG <- TransfoParam_Lag
     }
     ## set FUN_TRANSFO
-    if (sum(ObjectClass %in% c("GR4H", "GR5H", "GR4J", "GR5J", "GR6J", "GR2M", "GR1A", "CemaNeige")) > 0) {
+    if (sum(ObjectClass %in% c("GR4H", "GR5H", "GR4J", "GR5J", "GR6J", "GR2M", "GR1A", "CemaNeige", "Lag")) > 0) {
       if (!IsSD) {
         FUN_TRANSFO <- FUN_GR
       } else {
@@ -292,6 +302,10 @@ CreateCalibOptions <- function(FUN_MOD,
   if ("CemaNeigeGR6J" %in% ObjectClass) {
     NParam <- 8
   }
+  if ("Lag" %in% ObjectClass) {
+    NParam <- 1
+  }
+
   if (IsHyst) {
     NParam <- NParam + 2
   }
diff --git a/tests/testthat/test-RunModel_Lag.R b/tests/testthat/test-RunModel_Lag.R
index cd291693..8b83924c 100644
--- a/tests/testthat/test-RunModel_Lag.R
+++ b/tests/testthat/test-RunModel_Lag.R
@@ -164,15 +164,16 @@ test_that("1 input with lag of 0.5 time step delay out gives an output delayed o
   expect_equal(QlsSdSim - QlsGR4Only, QlsUpstLagObs)
 })
 
+InputsCrit <- CreateInputsCrit(
+  FUN_CRIT = ErrorCrit_NSE,
+  InputsModel = InputsModel,
+  RunOptions = RunOptions,
+  VarObs = "Q",
+  Obs = (c(0, Qupstream[Ind_Run[1:(length(Ind_Run) - 1)]]) * BasinAreas[1L] +
+           BasinObs$Qmm[Ind_Run] * BasinAreas[2L]) / sum(BasinAreas)
+)
+
 test_that("Params from calibration with simulated data should be similar to initial params", {
-  InputsCrit <- CreateInputsCrit(
-    FUN_CRIT = ErrorCrit_NSE,
-    InputsModel = InputsModel,
-    RunOptions = RunOptions,
-    VarObs = "Q",
-    Obs = (c(0, Qupstream[Ind_Run[1:(length(Ind_Run) - 1)]]) * BasinAreas[1L] +
-             BasinObs$Qmm[Ind_Run] * BasinAreas[2L]) / sum(BasinAreas)
-  )
   CalibOptions <- CreateCalibOptions(
     FUN_MOD = RunModel_GR4J,
     FUN_CALIB = Calibration_Michel,
@@ -189,6 +190,23 @@ test_that("Params from calibration with simulated data should be similar to init
   expect_equal(OutputsCalib$ParamFinalR[1L], ParamSD[1L], tolerance = 2e-3)
 })
 
+test_that("Params from calibration with simulated data should be similar to initial params", {
+  CalibOptions <- CreateCalibOptions(
+    FUN_MOD = RunModel_Lag,
+    FUN_CALIB = Calibration_Michel,
+    IsSD = FALSE
+  )
+  OutputsCalib <- Calibration_Michel(
+    InputsModel = InputsModel,
+    RunOptions = RunOptions,
+    InputsCrit = InputsCrit,
+    CalibOptions = CalibOptions,
+    FUN_MOD = RunModel_Lag,
+    QcontribDown = OutputsGR4JOnly
+  )
+  expect_equal(OutputsCalib$ParamFinalR[1L], ParamSD[1L], tolerance = 2e-3)
+})
+
 test_that("1 no area input with lag of 1 time step delay out gives an output delayed of one time step converted to mm", {
   Qm3GR4Only <- OutputsGR4JOnly$Qsim * BasinAreas[2L] * 1e3
   # Specify that upstream flow is not related to an area
-- 
GitLab