**airGR** is a package that brings into the **R software** the
hydrological modelling tools used and developed at the Catchment Hydrology Research
Group at INRAE (France),
including the **GR
rainfall-runoff models** that can be applied either on a
**lumped** or **semi-distributed** way. A snow
accumulation and melt model (**CemaNeige**)
and the associated functions for the calibration and evaluation of
models are also included. Each model core is coded in
**Fortran** to ensure low computational time. The other
package functions (i.e. mainly the calibration algorithm and the
efficiency criteria calculation) are coded in **R**.

The **airGR** package has been designed to fulfill two
major requirements: to facilitate the use by non-expert users and to
allow flexibility regarding the addition of external criteria, models or
calibration algorithms. The names of the functions and their arguments
were chosen to this end. **airGR** also contains basics
plotting facilities.

Seven hydrological models and one snow melt and accumulation model
are implemented in **airGR**. The hydrological models can
be applied either on a lumped way or on a semi-distributed way (on
sub-catchments). The snow model can either be used alone or with the
daily or hourly hydrological models. Naturally each hydrological model
can also be used alone.

The models can be called within **airGR** using the
following functions:

`RunModel_GR4H()`

: four-parameter hourly lumped hydrological model (Mathevet 2005)`RunModel_GR5H()`

: five-parameter hourly lumped hydrological model (Ficchi 2017; Ficchì, Perrin, and Andréassian 2019)`RunModel_GR4J()`

: four-parameter daily lumped hydrological model (Perrin, Michel, and Andréassian 2003)`RunModel_GR5J()`

: five-parameter daily lumped hydrological model (Le Moine 2008)`RunModel_GR6J()`

: six-parameter daily lumped hydrological model (Pushpalatha et al. 2011)`RunModel_GR2M()`

: two-parameter monthly lumped hydrological model (Mouelhi 2003; Mouelhi et al. 2006a)`RunModel_GR1A()`

: one-parameter yearly lumped hydrological model (Mouelhi 2003; Mouelhi et al. 2006b)`RunModel_CemaNeige()`

: two-parameter degree-day snowmelt and accumulation model (Valéry, Andréassian, and Perrin 2014; Riboust et al. 2019)`RunModel_CemaNeigeGR4H()`

: combined use of**GR4H**and**CemaNeige**`RunModel_CemaNeigeGR5H()`

: combined use of**GR5H**and**CemaNeige**`RunModel_CemaNeigeGR4J()`

: combined use of**GR4J**and**CemaNeige**`RunModel_CemaNeigeGR5J()`

: combined use of**GR5J**and**CemaNeige**`RunModel_CemaNeigeGR6J()`

: combined use of**GR6J**and**CemaNeige**

The **GRP**
forecasting model and the **Otamin**
predictive uncertainty tool are not available in
**airGR**.

In this vignette, we show how to prepare and run a calibration and a simulation with airGR hydrological models.

In the following example, we use a data sample contained in the
package. For real applications, the user has to import its data into
**R** and to prepare it with an adequate data.frame format
as described below.

First, it is necessary to load the **airGR**
package:

Below is presented an example of a `data.frame`

of daily
hydrometeorological observations time series for a fictional catchment
included in the **airGR** package that contains:

*DatesR*: dates in the POSIXt format*P*: average precipitation [mm/day]*T*: catchment average air temperature [℃]*E*: catchment average potential evapotranspiration [mm/day]*Qls*: outlet discharge [l/s]*Qmm*: outlet discharge [mm/day]

```
## DatesR P T E
## Min. :1984-01-01 Min. : 0.0 Min. :-18.7 Min. :0.0
## 1st Qu.:1991-04-02 1st Qu.: 0.0 1st Qu.: 4.1 1st Qu.:0.6
## Median :1998-07-02 Median : 0.3 Median : 9.1 Median :1.4
## Mean :1998-07-02 Mean : 2.9 Mean : 9.1 Mean :1.8
## 3rd Qu.:2005-10-01 3rd Qu.: 3.6 3rd Qu.: 14.5 3rd Qu.:2.9
## Max. :2012-12-31 Max. :66.8 Max. : 28.4 Max. :5.5
##
## Qls Qmm
## Min. : 70 Min. : 0.017
## 1st Qu.: 1640 1st Qu.: 0.394
## Median : 4050 Median : 0.977
## Mean : 6130 Mean : 1.473
## 3rd Qu.: 7850 3rd Qu.: 1.884
## Max. :99500 Max. :23.880
## NA's :772 NA's :802
```

The usual functions (e.g. `read.table()`

) can be used to
load real-case data sets.

To run a model, the functions of the **airGR** package
(e.g. the models, calibration and criteria calculation functions)
require data and options with specific formats.

To facilitate the use of the package, there are several functions dedicated to the creation of these objects:

`CreateInputsModel()`

: prepares the inputs for the different hydrological models (times series of dates, precipitation, observed discharge, etc.)`CreateRunOptions()`

: prepares the options for the hydrological model run (warm up period, calibration period, etc.)`CreateInputsCrit()`

: prepares the options in order to compute the efficiency criterion (choice of the criterion, choice of the transformation on discharge: “log”, “sqrt”, etc.)`CreateCalibOptions()`

: prepares the options for the hydrological model calibration algorithm (choice of parameters to optimize, predefined values for uncalibrated parameters, etc.)

To run a GR hydrological model or CemaNeige, the user has to prepare
the input data with the `CreateInputsModel()`

function. As
arguments, this function needs the function name corresponding to the
model the user wants to run, a vector of dates, a vector of
precipitation and a vector of potential evapotranspiration.

In the example below, we already have the potential
evapotranspiration. If the user does not have these data, it is possible
to compute it with the Oudin’s
formula with the `PE_Oudin()`

function (this function
only needs Julian days, daily average air temperature and latitude).

Missing values (`NA`

) of precipitation (or potential
evapotranspiration) are **not allowed**.

```
InputsModel <- CreateInputsModel(FUN_MOD = RunModel_GR4J, DatesR = BasinObs$DatesR,
Precip = BasinObs$P, PotEvap = BasinObs$E)
str(InputsModel)
```

```
## List of 3
## $ DatesR : POSIXlt[1:10593], format: "1984-01-01" "1984-01-02" ...
## $ Precip : num [1:10593] 4.1 15.9 0.8 0 0 0 0 0 2.9 0 ...
## $ PotEvap: num [1:10593] 0.2 0.2 0.3 0.3 0.1 0.3 0.4 0.4 0.5 0.5 ...
## - attr(*, "class")= chr [1:3] "InputsModel" "daily" "GR"
```

The `CreateRunOptions()`

function allows to prepare the
options required to the `RunModel*()`

functions, which are
the actual models functions.

The user must at least define the following arguments:

`FUN_MOD`

: the name of the model function to run`InputsModel`

: the associated input data`IndPeriod_Run`

: the period on which the model is run

To select a period for which the user wants to run the model, select the corresponding indexes for different time periods (not the POSIXt dates), as follows:

```
Ind_Run <- seq(which(format(BasinObs$DatesR, format = "%Y-%m-%d") == "1990-01-01"),
which(format(BasinObs$DatesR, format = "%Y-%m-%d") == "1999-12-31"))
str(Ind_Run)
```

`## int [1:3652] 2193 2194 2195 2196 2197 2198 2199 2200 2201 2202 ...`

The initialization of hydrological models is of the utmost
importance. Indeed, an inaccurate initialization causes poor quality
discharge simulations during the earliest stages of the running period.
For example, in the GR models, by default, the production and the
routing store levels store level are respectively set to 30 % and 50 %
of their capacity, which may be far from their ideal value. Two
solutions are offered to accurately initialize the GR models in
**airGR**: manually predefining the initial states
(e.g. from a previous run) or running the models during a warm up period
before the actual running period. It is generally advised to set up this
warm up period to be equal or superior to one year.

As a consequence, it is possible to define in
`CreateRunOptions()`

the following arguments:

`IniStates`

: the initial states of the 2 unit hydrographs (20 + 40 = 60 units)`IniResLevels`

: the initial levels of the production and routing stores`IndPeriod_WarmUp`

: the warm up period used to run the model, to be defined in the same format as`IndPeriod_Run`

```
RunOptions <- CreateRunOptions(FUN_MOD = RunModel_GR4J,
InputsModel = InputsModel, IndPeriod_Run = Ind_Run,
IniStates = NULL, IniResLevels = NULL, IndPeriod_WarmUp = NULL)
```

```
## Warning in CreateRunOptions(FUN_MOD = RunModel_GR4J, InputsModel = InputsModel, : model warm up period not defined: default configuration used
## the year preceding the run period is used
```

```
## List of 8
## $ IndPeriod_WarmUp: int [1:365] 1828 1829 1830 1831 1832 1833 1834 1835 1836 1837 ...
## $ IndPeriod_Run : int [1:3652] 2193 2194 2195 2196 2197 2198 2199 2200 2201 2202 ...
## $ IniStates : num [1:67] 0 0 0 0 0 0 0 0 0 0 ...
## $ IniResLevels : num [1:4] 0.3 0.5 NA NA
## $ Outputs_Cal : chr [1:2] "Qsim" "Param"
## $ Outputs_Sim : Named chr [1:22] "DatesR" "PotEvap" "Precip" "Prod" ...
## ..- attr(*, "names")= chr [1:22] "" "GR1" "GR2" "GR3" ...
## $ FortranOutputs :List of 2
## ..$ GR: chr [1:18] "PotEvap" "Precip" "Prod" "Pn" ...
## ..$ CN: NULL
## $ FeatFUN_MOD :List of 12
## ..$ CodeMod : chr "GR4J"
## ..$ NameMod : chr "GR4J"
## ..$ NbParam : int 4
## ..$ TimeUnit : chr "daily"
## ..$ Id : logi NA
## ..$ Class : chr [1:2] "daily" "GR"
## ..$ Pkg : chr "airGR"
## ..$ NameFunMod : chr "RunModel_GR4J"
## ..$ TimeStep : num 86400
## ..$ TimeStepMean: int 86400
## ..$ CodeModHydro: chr "GR4J"
## ..$ IsSD : logi FALSE
## - attr(*, "class")= chr [1:3] "RunOptions" "daily" "GR"
```

The `CreateRunOptions()`

function returns warnings if the
default initialization options are used:

`IniStates`

and`IniResLevels`

are automatically set to initialize all the model states at 0, except for the production and routing stores, which are initialized at respectively 30 % and 50 % of their capacity`IndPeriod_WarmUp`

default setting ensures a one-year warm up using the time steps preceding the`IndPeriod_Run`

, if available

The `CreateInputsCrit()`

function allows to prepare the
input in order to calculate a criterion. It is possible to define the
following arguments:

`FUN_CRIT`

: the name of the error criterion function (the available functions are introduced later on)`InputsModel`

: the inputs of the hydrological model previously prepared by the`CreateInputsModel()`

function`RunOptions`

: the options of the hydrological model previously prepared by the`CreateRunOptions()`

function`VarObs`

: the name of the considered variable (by default`"Q"`

for the discharge)`Obs`

: the observed variable time serie (e.g. the discharge expressed in*mm/time step*)

Missing values (`NA`

) are **allowed** for
observed discharge.

It is possible to compute a composite criterion (e.g. the average
between NSE computed on discharge and NSE computed on log of discharge).
In this case, users have to provide lists to the following arguments
(some of the are optional): `FUN_CRIT`

, `Obs`

,
`VarObs`

, `BoolCrit`

, `transfo`

,
`Weights.`

```
InputsCrit <- CreateInputsCrit(FUN_CRIT = ErrorCrit_NSE, InputsModel = InputsModel,
RunOptions = RunOptions, VarObs = "Q", Obs = BasinObs$Qmm[Ind_Run])
str(InputsCrit)
```

```
## List of 8
## $ FUN_CRIT:function (InputsCrit, OutputsModel, warnings = TRUE, verbose = TRUE)
## ..- attr(*, "class")= chr [1:2] "FUN_CRIT" "function"
## $ Obs : num [1:3652] 1.99 1.8 2.86 2.4 3.31 ...
## $ VarObs : chr "Q"
## $ BoolCrit: logi [1:3652] TRUE TRUE TRUE TRUE TRUE TRUE ...
## $ idLayer : logi NA
## $ transfo : chr ""
## $ epsilon : NULL
## $ Weights : NULL
## - attr(*, "class")= chr [1:2] "Single" "InputsCrit"
```

Before using the automatic calibration tool, the user needs to
prepare the calibration options with the
`CreateCalibOptions()`

function. For that, it is necessary to
define the following arguments:

`FUN_MOD`

: the name of the model function`FUN_CALIB`

: the name of the calibration algorithm

```
CalibOptions <- CreateCalibOptions(FUN_MOD = RunModel_GR4J, FUN_CALIB = Calibration_Michel)
str(CalibOptions)
```

```
## List of 4
## $ FixedParam : logi [1:4] NA NA NA NA
## $ SearchRanges : num [1:2, 1:4] 4.59e-05 2.18e+04 -1.09e+04 1.09e+04 4.59e-05 ...
## $ FUN_TRANSFO :function (ParamIn, Direction)
## $ StartParamDistrib: num [1:3, 1:4] 169.017 247.151 432.681 -2.376 -0.649 ...
## - attr(*, "class")= chr [1:4] "CalibOptions" "daily" "GR" "HBAN"
```

The evaluation of the quality of a simulation is estimated through the calculation of criteria. These criteria can be used both as objective-functions during the calibration of the model, or as a measure for evaluating its performance on a control period.

The package offers the possibility to use different criteria:

`ErrorCrit_RMSE()`

: Root mean square error (RMSE)`ErrorCrit_NSE()`

: Nash-Sutcliffe model efficiency coefficient (NSE)`ErrorCrit_KGE()`

: Kling-Gupta efficiency criterion (KGE)`ErrorCrit_KGE2()`

: modified Kling-Gupta efficiency criterion (KGE’)

It is also possible to create user-defined criteria. For doing that,
it is only necessary to define the function in **R**
following the same syntax as the criteria functions included in
**airGR**.

The objective of the calibration algorithm is to identify the model parameters: by comparing the model outputs with observed data, this algorithm determines the combination of parameters that represents the best the behavior of the watershed.

In the **airGR** package, a function called
`Calibration_Michel()`

is implemented. This functions allows
running a calibration with the package models. The calibration algorithm
optimizes the error criterion selected as objective-function. This
algorithm works in two steps:

- a screening of the parameters space is performed using either a rough predefined grid or a user-defined list of parameter sets
- a simple steepest descent local search algorithm is performed from the best set of parameters found at the first step

```
OutputsCalib <- Calibration_Michel(InputsModel = InputsModel, RunOptions = RunOptions,
InputsCrit = InputsCrit, CalibOptions = CalibOptions,
FUN_MOD = RunModel_GR4J)
```

```
## Grid-Screening in progress (0% 20% 40% 60% 80% 100%)
## Screening completed (81 runs)
## Param = 247.151, -0.020, 83.096, 2.384
## Crit. NSE[Q] = 0.7688
## Steepest-descent local search in progress
## Calibration completed (21 iterations, 234 runs)
## Param = 257.238, 1.012, 88.235, 2.208
## Crit. NSE[Q] = 0.7988
```

`## [1] 257.237556 1.012237 88.234673 2.207958`

The `Calibration_Michel()`

function is the only one
implemented in the **airGR** package to calibrate the
model, but the user can implement its own calibration function. Two
vignettes explain how it can be done (2.1 Plugging in new calibration and 2.2 MCMC parameter estimation).

The `Calibration_Michel()`

function returns a vector with
the parameters of the chosen model, which means that the number of
values can differ depending on the model that is used. It is possible to
use the `Calibration_Michel()`

function with user-implemented
hydrological models.

This step assesses the predictive capacity of the model. Control is defined as the estimation of the accuracy of the model on data sets that are not used in its construction, and in particular its calibration. The classical way to perform a control is to keep data from a period separated from the calibration period. If possible, this control period should correspond to climatic situations that differ from those of the calibration period in order to better point out the qualities and weaknesses of the model. This exercise is necessary for assessing the robustness of the model, that is to say its ability to keep stable performances outside of the calibration conditions.

Performing a model control with **airGR** is similar to
running a simulation (see below), followed by the computation of one or
several performance criteria.

To run a model, the user has to use the `RunModel*()`

functions (`InputsModel`

, `RunOptions`

and
parameters). All the data needed have already been prepared in the
previous steps defined in this guide.

```
OutputsModel <- RunModel_GR4J(InputsModel = InputsModel, RunOptions = RunOptions, Param = Param)
str(OutputsModel)
```

```
## List of 21
## $ DatesR : POSIXlt[1:3652], format: "1990-01-01" "1990-01-02" ...
## $ PotEvap : num [1:3652] 0.3 0.4 0.4 0.3 0.1 0.1 0.1 0.2 0.2 0.3 ...
## $ Precip : num [1:3652] 0 9.3 3.2 7.3 0 0 0 0 0.1 0.2 ...
## $ Prod : num [1:3652] 196 199 199 201 200 ...
## $ Pn : num [1:3652] 0 8.9 2.8 7 0 0 0 0 0 0 ...
## $ Ps : num [1:3652] 0 3.65 1.12 2.75 0 ...
## $ AE : num [1:3652] 0.2833 0.4 0.4 0.3 0.0952 ...
## $ Perc : num [1:3652] 0.645 0.696 0.703 0.74 0.725 ...
## $ PR : num [1:3652] 0.645 5.946 2.383 4.992 0.725 ...
## $ Q9 : num [1:3652] 1.78 1.52 3.86 3.17 3.45 ...
## $ Q1 : num [1:3652] 0.2 0.195 0.271 0.387 0.365 ...
## $ Rout : num [1:3652] 53.9 53.6 55.3 56.1 56.9 ...
## $ Exch : num [1:3652] 0.181 0.18 0.177 0.197 0.207 ...
## $ AExch1 : num [1:3652] 0.181 0.18 0.177 0.197 0.207 ...
## $ AExch2 : num [1:3652] 0.181 0.18 0.177 0.197 0.207 ...
## $ AExch : num [1:3652] 0.362 0.36 0.353 0.393 0.414 ...
## $ QR : num [1:3652] 2.05 1.99 2.36 2.55 2.78 ...
## $ QD : num [1:3652] 0.381 0.375 0.447 0.584 0.572 ...
## $ Qsim : num [1:3652] 2.43 2.37 2.8 3.14 3.35 ...
## $ RunOptions:List of 2
## ..$ WarmUpQsim: num [1:365] 0.779 0.815 0.807 0.731 0.675 ...
## ..$ Param : num [1:4] 257.24 1.01 88.23 2.21
## $ StateEnd :List of 3
## ..$ Store :List of 4
## .. ..$ Prod: num 189
## .. ..$ Rout: num 48.9
## .. ..$ Exp : num NA
## .. ..$ Int : num NA
## ..$ UH :List of 2
## .. ..$ UH1: num [1:20] 0.514 0.54 0.148 0 0 ...
## .. ..$ UH2: num [1:40] 0.056306 0.057176 0.042253 0.012187 0.000578 ...
## ..$ CemaNeigeLayers:List of 4
## .. ..$ G : num NA
## .. ..$ eTG : num NA
## .. ..$ Gthr : num NA
## .. ..$ Glocmax: num NA
## ..- attr(*, "class")= chr [1:3] "IniStates" "daily" "GR"
## - attr(*, "class")= chr [1:3] "OutputsModel" "daily" "GR"
```

Although it is possible for the user to design its own graphics from
the outputs of the `RunModel*()`

functions, the
**airGR** package offers the possibility to make use of the
`plot()`

function. This function returns a dashboard of
results including various graphs (depending on the model used):

- time series of total precipitation and simulated discharge (and observed discharge if provided)
- interannual average daily simulated discharge (and daily observed discharge if provided) and interannual average monthly precipitation
- cumulative frequency plot for simulated discharge (and for observed discharge if provided)
- correlation plot between simulated and observed discharge (if observed discharge provided)

Moreover, if the CemaNeige model is used, the air temperature and the simulated snowpack water equivalent time series are plotted.

To evaluate the efficiency of the model, it is possible to use the same criterion as defined at the calibration step or to use another criterion.

`## Crit. NSE[Q] = 0.7988`

```
## List of 5
## $ CritValue : num 0.799
## $ CritName : chr "NSE[Q]"
## $ CritBestValue : num 1
## $ Multiplier : num -1
## $ Ind_notcomputed: int [1:57] 2405 2406 2407 2408 2409 2410 2411 2412 2413 2414 ...
## - attr(*, "class")= chr [1:2] "NSE" "ErrorCrit"
```

`## Crit. KGE[Q] = 0.7854`

```
## SubCrit. KGE[Q] cor(sim, obs, "pearson") = 0.8985
## SubCrit. KGE[Q] sd(sim)/sd(obs) = 0.8161
## SubCrit. KGE[Q] mean(sim)/mean(obs) = 1.0437
```

```
## List of 7
## $ CritValue : num 0.785
## $ CritName : chr "KGE[Q]"
## $ SubCritValues : num [1:3] 0.898 0.816 1.044
## $ SubCritNames : chr [1:3] "r" "alpha" "beta"
## $ CritBestValue : num 1
## $ Multiplier : num -1
## $ Ind_notcomputed: int [1:57] 2405 2406 2407 2408 2409 2410 2411 2412 2413 2414 ...
## - attr(*, "class")= chr [1:2] "KGE" "ErrorCrit"
```

Ficchi, Andrea. 2017. “An Adaptive Hydrological Model for Multiple
Time-Steps: Diagnostics and Improvements Based on Fluxes
Consistency.” PhD thesis, Université Pierre et Marie Curie, Paris
6. https://www.theses.fr/2017PA066097.

Ficchì, Andrea, Charles Perrin, and Vazken Andréassian. 2019.
“Hydrological Modelling at Multiple Sub-Daily Time Steps: Model
Improvement via Flux-Matching.” *Journal of Hydrology*,
June. https://doi.org/10.1016/j.jhydrol.2019.05.084.

Le Moine, Nicolas. 2008. “Le Bassin Versant de Surface Vu Par Le
Souterrain : Une Voie d’amélioration Des Performances Et Du Réalisme Des
Modèles Pluie-Débit ?” PhD thesis, Université Pierre et Marie
Curie, Paris 6. http://webgr.irstea.fr/wp-content/uploads/2012/07/2008-LE_MOINE-THESE.pdf.

Mathevet, Thibault. 2005. “Quels Modèles Pluie-Débit Globaux Au
Pas de Temps Horaire ? Développements Empiriques Et
Comparaison de Modèles Sur Un Large Échantillon de Bassins
Versants.” PhD thesis, Paris: ENGREF. http://webgr.irstea.fr/wp-content/uploads/2012/07/2005-MATHEVET-THESE.pdf.

Mouelhi, Safouane. 2003. “Vers Une Chaîne Cohérente de Modèles
Pluie-Débit Conceptuels Globaux Aux Pas de Temps Pluriannuel, Annuel,
Mensuel Et Journalier.” PhD thesis, Paris, ENGREF. http://webgr.irstea.fr/wp-content/uploads/2012/07/2003-MOUELHI-THESE.pdf.

Mouelhi, Safouane, Claude Michel, Charles Perrin, and Vazken
Andréassian. 2006a. “Stepwise Development of a Two-Parameter
Monthly Water Balance Model.” *Journal of Hydrology* 318
(1-4): 200–214. https://doi.org/10.1016/j.jhydrol.2005.06.014.

———. 2006b. “Linking Stream Flow to Rainfall at the Annual Time
Step: The Manabe Bucket Model
Revisited.” *Journal of Hydrology* 328 (1-2): 283–96. https://doi.org/10.1016/j.jhydrol.2005.12.022.

Perrin, Charles, Claude Michel, and Vazken Andréassian. 2003.
“Improvement of a Parsimonious Model for Streamflow
Simulation.” *Journal of Hydrology* 279 (1-4): 275–89. https://doi.org/10.1016/S0022-1694(03)00225-7.

Pushpalatha, Raji, Charles Perrin, Nicolas Le Moine, Thibault Mathevet,
and Vazken Andréassian. 2011. “A Downward Structural Sensitivity
Analysis of Hydrological Models to Improve Low-Flow Simulation.”
*Journal of Hydrology* 411 (1–2): 66–76. https://doi.org/10.1016/j.jhydrol.2011.09.034.

Riboust, Philippe, Guillaume Thirel, Nicolas Le Moine, and Pierre
Ribstein. 2019. “Revisiting a Simple
Degree-Day Model for
Integrating Satellite Data:
Implementation of Swe-Sca
Hystereses.” *Journal of Hydrology and
Hydromechanics* 67 (1): 70–81. https://doi.org/10.2478/johh-2018-0004.

Valéry, Audrey, Vazken Andréassian, and Charles Perrin. 2014.
“’As Simple as Possible but Not Simpler’: What Is
Useful in a Temperature-Based Snow-Accounting Routine? Part
2 - Sensitivity Analysis of the Cemaneige Snow
Accounting Routine on 380 Catchments.” *Journal of
Hydrology*, no. 517(0): 1176–87. https://doi.org/10.1016/j.jhydrol.2014.04.058.