diff --git a/R/RunModel.GRiwrmInputsModel.R b/R/RunModel.GRiwrmInputsModel.R
index 624415f..db0d171 100644
--- a/R/RunModel.GRiwrmInputsModel.R
+++ b/R/RunModel.GRiwrmInputsModel.R
@@ -31,5 +31,7 @@ RunModel.GRiwrmInputsModel <- function(x, RunOptions, Param, ...) {
     )
   }
   attr(OutputsModel, "Qm3s") <- OutputsModelQsim(x, OutputsModel, RunOptions[[1]]$IndPeriod_Run)
+  attr(OutputsModel, "GRiwrm") <- attr(x, "GRiwrm")
+  attr(OutputsModel, "TimeStep") <- attr(x, "TimeStep")
   return(OutputsModel)
 }
diff --git a/R/RunModel_Reservoir.R b/R/RunModel_Reservoir.R
index 91b211a..e839402 100644
--- a/R/RunModel_Reservoir.R
+++ b/R/RunModel_Reservoir.R
@@ -2,7 +2,7 @@
 #'
 #' The reservoir model is a model combining a lag model and the calculation of
 #' the water storage time series according to the released flow time series
-#' from the `Qinf` parameter of [CreateInputsModel.GRiwrm].
+#' from the `Qrelease` parameter of [CreateInputsModel.GRiwrm].
 #'
 #' @details
 #' The simulated flow corresponds to the released flow except when the reservoir
diff --git a/R/plot.GRiwrmOutputsModel.R b/R/plot.GRiwrmOutputsModel.R
index ac96f5e..b8d3bfc 100644
--- a/R/plot.GRiwrmOutputsModel.R
+++ b/R/plot.GRiwrmOutputsModel.R
@@ -4,6 +4,7 @@
 #' @param Qobs (optional) [matrix] time series of observed flows
 #'        (for the same time steps than simulated) (mm/time step) with one column
 #'        by hydrological model output named with the node ID (See [CreateGRiwrm] for details)
+#' @param unit (optional) [character] flows unit ("m3/s" or "mm")
 #' @param ... Further arguments for [airGR::plot.OutputsModel] and [plot]
 #'
 #' @return [list] of plots.
@@ -13,8 +14,14 @@
 #'
 #' @example man-examples/RunModel.GRiwrmInputsModel.R
 #'
-plot.GRiwrmOutputsModel <- function(x, Qobs = NULL, ...) {
+plot.GRiwrmOutputsModel <- function(x, Qobs = NULL, unit = "m3/s", ...) {
 
+  # Arguments checks
+  stopifnot(is.null(Qobs) || is.matrix(Qobs) || is.data.frame(Qobs),
+            is.character(unit),
+            unit %in% c("mm", "m3/s"))
+
+  griwrm <- attr(x, "GRiwrm")
   ## define outer margins and a title inside it
   oldpar <- par(no.readonly = TRUE)
   on.exit(par(oldpar))
@@ -22,7 +29,7 @@ plot.GRiwrmOutputsModel <- function(x, Qobs = NULL, ...) {
 
   lapply(
     names(x),
-    function(id, OutputsModels) {
+    function(id) {
       Qobs_id <- NULL
       if (!is.null(Qobs)) {
         if (id %in% colnames(Qobs)) {
@@ -31,10 +38,15 @@ plot.GRiwrmOutputsModel <- function(x, Qobs = NULL, ...) {
           warning("Column \"", id, "\" not found in Qobs")
         }
       }
-      plot(OutputsModels[[id]], Qobs = Qobs_id, ...)
+
+      BasinArea <- griwrm$area[griwrm$id == id & !is.na(griwrm$model) & griwrm$model != "Diversion"]
+      if (unit == "m3/s" && length(BasinArea) == 1 && !is.na(BasinArea)) {
+        plot(x[[id]], Qobs = Qobs_id, BasinArea = BasinArea, ...)
+      } else {
+        plot(x[[id]], Qobs = Qobs_id, ...)
+      }
       title(main = id, outer = TRUE, line = 1.2, cex.main = 1.4)
-    },
-    OutputsModels = x
+    }
   )
   invisible(NULL)
 }
diff --git a/R/plot.OutputsModelReservoir.R b/R/plot.OutputsModelReservoir.R
index e27e92c..c283359 100644
--- a/R/plot.OutputsModelReservoir.R
+++ b/R/plot.OutputsModelReservoir.R
@@ -1,7 +1,7 @@
 #' Plot simulated reservoir volume, inflows and released flows time series on a reservoir node
 #'
 #' @param x Object returned by [RunModel_Reservoir]
-#' @param Qobs [numeric] time series of observed flow (for the same time steps than simulated) \[m3/time step\]
+#' @param Qobs (optional) [numeric] time series of observed released flow \[m3/time step\]
 #' @param ... Further arguments passed to [plot.Qm3s]
 #'
 #' @return Function used for side effect.
@@ -9,6 +9,14 @@
 #'
 #' @example man-examples/RunModel_Reservoir.R
 plot.OutputsModelReservoir <- function(x, Qobs = NULL, ...) {
+
+  # Argument checks
+  if (!is.null(Qobs)) {
+    stopifnot(is.numeric(Qobs),
+              length(Qobs) == length(x$Qsim_m3))
+
+  }
+
   oldpar <- par(mfrow=c(2,1),
                 mar = c(2,3.3,1.2,0.5),
                 mgp = c(2,1,0))
diff --git a/man-examples/RunModel_Reservoir.R b/man-examples/RunModel_Reservoir.R
index 85a6764..a330072 100644
--- a/man-examples/RunModel_Reservoir.R
+++ b/man-examples/RunModel_Reservoir.R
@@ -86,4 +86,4 @@ plot(OutputsModel, Qobs = Qobs)
 # N.B. "Observed releases" should be considered as "Target releases" here
 
 # The plot for the reservoir can also be plotted alone
-plot(OutputsModel$Reservoir)
+plot(OutputsModel$Reservoir, Qobs = Qobs[, "Reservoir"])
diff --git a/man/RunModel_Reservoir.Rd b/man/RunModel_Reservoir.Rd
index ee504de..48f3f77 100644
--- a/man/RunModel_Reservoir.Rd
+++ b/man/RunModel_Reservoir.Rd
@@ -21,7 +21,7 @@ completed with an item \code{Vsim} representing the water volume time series in
 \description{
 The reservoir model is a model combining a lag model and the calculation of
 the water storage time series according to the released flow time series
-from the \code{Qinf} parameter of \link{CreateInputsModel.GRiwrm}.
+from the \code{Qrelease} parameter of \link{CreateInputsModel.GRiwrm}.
 }
 \details{
 The simulated flow corresponds to the released flow except when the reservoir
@@ -143,5 +143,5 @@ plot(OutputsModel, Qobs = Qobs)
 # N.B. "Observed releases" should be considered as "Target releases" here
 
 # The plot for the reservoir can also be plotted alone
-plot(OutputsModel$Reservoir)
+plot(OutputsModel$Reservoir, Qobs = Qobs[, "Reservoir"])
 }
diff --git a/man/plot.GRiwrmOutputsModel.Rd b/man/plot.GRiwrmOutputsModel.Rd
index 2493e95..974f682 100644
--- a/man/plot.GRiwrmOutputsModel.Rd
+++ b/man/plot.GRiwrmOutputsModel.Rd
@@ -4,7 +4,7 @@
 \alias{plot.GRiwrmOutputsModel}
 \title{Function which creates screen plots giving an overview of the model outputs in the GRiwrm network}
 \usage{
-\method{plot}{GRiwrmOutputsModel}(x, Qobs = NULL, ...)
+\method{plot}{GRiwrmOutputsModel}(x, Qobs = NULL, unit = "m3/s", ...)
 }
 \arguments{
 \item{x}{[object of class \emph{GRiwrmOutputsModel}] see \link{RunModel.GRiwrmInputsModel} for details}
@@ -13,6 +13,8 @@
 (for the same time steps than simulated) (mm/time step) with one column
 by hydrological model output named with the node ID (See \link{CreateGRiwrm} for details)}
 
+\item{unit}{(optional) \link{character} flows unit ("m3/s" or "mm")}
+
 \item{...}{Further arguments for \link[airGR:plot.OutputsModel]{airGR::plot.OutputsModel} and \link{plot}}
 }
 \value{
diff --git a/man/plot.OutputsModelReservoir.Rd b/man/plot.OutputsModelReservoir.Rd
index b0796d9..7b3de6f 100644
--- a/man/plot.OutputsModelReservoir.Rd
+++ b/man/plot.OutputsModelReservoir.Rd
@@ -9,7 +9,7 @@
 \arguments{
 \item{x}{Object returned by \link{RunModel_Reservoir}}
 
-\item{Qobs}{\link{numeric} time series of observed flow (for the same time steps than simulated) [m3/time step]}
+\item{Qobs}{(optional) \link{numeric} time series of observed released flow [m3/time step]}
 
 \item{...}{Further arguments passed to \link{plot.Qm3s}}
 }
@@ -108,5 +108,5 @@ plot(OutputsModel, Qobs = Qobs)
 # N.B. "Observed releases" should be considered as "Target releases" here
 
 # The plot for the reservoir can also be plotted alone
-plot(OutputsModel$Reservoir)
+plot(OutputsModel$Reservoir, Qobs = Qobs[, "Reservoir"])
 }