title: "Severn_02: Calibration of a GR4J semi-distributed model network"
author: "David Dorchies"
vignette: >
%\VignetteIndexEntry{Severn_02: Calibration of a GR4J semi-distributive model network}
%\VignetteIndexEntry{Severn_02: Calibration of a GR4J semi-distributed model network}
%\VignetteEngine{knitr::rmarkdown}
%\VignetteEncoding{UTF-8}
---
...
...
@@ -32,7 +32,7 @@ To facilitate the use of the package, there are several functions dedicated to t
### GRiwrmInputsModel object
The production method of the `GRiwrmInputsModel` object is detailed in the vignette "V01_Structure_SD_model" of the package. The following code chunk resumes all the steps of this vignette:
The method used for producing the `GRiwrmInputsModel` object is detailed in the vignette "V01_Structure_SD_model" of the package. The following code chunk resumes all the steps of this vignette:
The `CreateRunOptions()` function allows to prepare the options required for the `RunModel()` function.
...
...
@@ -60,19 +60,14 @@ The user must at least define the following arguments:
* `InputsModel`: the associated input data
* `IndPeriod_Run`: the period on which the model is run
Below, we start the run period one year after the beginning of the time series.
Below, we define a one-year warm up period and we start the run period just after the warmup period.
```{r}
IndPeriod_Run <- seq(
which(InputsModel[[1]]$DatesR == (InputsModel[[1]]$DatesR[1] + 365*24*60*60)), # Set aside warm-up period
length(InputsModel[[1]]$DatesR) # Until the end of the time series
)
```
As a consequence, the warmup period is defined as:
```{r}
IndPeriod_WarmUp = seq(1,IndPeriod_Run[1]-1)
IndPeriod_WarmUp <- seq(1, IndPeriod_Run[1] - 1)
```
Arguments of the `CreateRunOptions` function for **airGRiwrm** are the same as for the function in **airGR** and are copied for each node running a rainfall-runoff model.
The `CreateInputsCrit()` function allows to prepare the input in order to calculate a criterion. We use composed criterion with a parameter regularisation based on @delavenneRegularizationApproachImprove2019.
The `CreateInputsCrit()` function allows to prepare the input in order to calculate a criterion.
We use composed criterion with a parameter regularisation based on @delavenneRegularizationApproachImprove2019.
It needs the following arguments:
* `InputsModel`: the inputs of the `GRiwrm` network previously prepared by the `CreateInputsModel()` function
* `FUN_CRIT`: the name of the error criterion function (see the available functions description in the **airGR** package)
* `RunOptions`: the options of the `GRiwrm` network previously prepared by the `CreateRunOptions()` function
* `Qobs`: the observed variable time serie (e.g. the discharge expressed in *mm/time step*)
* `Qobs`: the observed variable time series (e.g. the discharge expressed in *mm/time step*)
* `AprioriIds`: the list of the sub-catchments IDs where to apply a parameter regularisation based on the parameters of an upstream sub-catchment (e.g. here below the parameters of the sub-catchment "54057" is regulated by the parameters of the sub-catchment "54032")
* `transfo`: a transformation function applied on the flow before criterion calculation (square-root transformation is recommended for the De Lavenne regularisation)
* `transfo`: a transformation function applied on the flow before calculation of the criterion (square-root transformation is recommended for the De Lavenne regularisation)
* `k`: coefficient used for the weighted average between the performance criterion and the gap between the optimised parameter set and an a priori parameter set (a value equal to 0.15 is recommended for the De Lavenne regularisation)
```{r InputsCrit}
InputsCrit <- CreateInputsCrit(
InputsModel = InputsModel,
FUN_CRIT = airGR::ErrorCrit_KGE2,
RunOptions = RunOptions,
Obs = Qobs[IndPeriod_Run,],
RunOptions = RunOptions,
Obs = Qobs[IndPeriod_Run,],
AprioriIds = c(
"54057" = "54032",
"54032" = "54001",
"54001" = "54095"
),
transfo = "sqrt"
transfo = "sqrt",
k = 0.15
)
str(InputsCrit)
```
...
...
@@ -153,7 +151,7 @@ OutputsModels <- RunModel(
plot(OutputsModels, Qobs = Qobs[IndPeriod_Run,])
```
The resulting flows of each node in m<sup>3</sup>/s is directly available and can be plotted with these commands:
The resulting flows of each node in m<sup>3</sup>/s are directly available and can be plotted with these commands:
Calibration options and criteria procedures are detailed in vignette "V02_Calibration_SD_model". We also apply a parameter regularisation here but only where an upstream simulated catchment is available.
Calibration options is detailed in vignette "V02_Calibration_SD_model".
We also apply a parameter regularisation here but only where an upstream simulated catchment is available.
As can be seen below, compared to results of vignette "V02_Calibration_SD_model", the use of measured flows on upstream influenced basins improves largely the model performance at downstream stations (better lowflow simulations).
As can be seen below, compared to results of vignette "V02_Calibration_SD_model", the use of measured flows on upstream influenced basins largely improves the model performance at downstream stations (better low-flow simulations).
```{r, fig.height = 5, fig.width = 8}
plot(OM_OL, Qobs = Qobs[IndPeriod_Run,], which = "Regime")
plot(OM_OL, Qobs = Qobs[IndPeriod_Run,], which = "Regime")
```
The resulting flows of each node in m<sup>3</sup>/s is directly available and can be plotted with these commands:
The resulting flows of each node in m<sup>3</sup>/s are directly available and can be plotted with these commands:
The example below is inspired by the semi-distributed model developed in the ClimAware project [@theobaldClimAwareImpactsClimate2014a].
The example below is inspired by the semi-distributed model developed in the ClimAware project [@theobaldClimAwareImpactsClimate2014a], in which we used the daily rainfall-runoff model GR4J.
## Semi-distributive network description
The aim of this vignette is to create the semi-distributed network necessary for the hydrological model.
## Semi-distributed network description
The model is distributed according to the gauging stations described in @dorchiesClimateChangeImpacts2014.
List of nodes
First, we must read the list of nodes and the associated metadata:
```{r seine_nodes}
seine_nodes <- read.table(
...
...
@@ -32,22 +34,25 @@ seine_nodes <- read.table(
seine_nodes
```
Create the GRiwrm object which lists the nodes and describes the network diagram. It's a dataframe of class `GRiwrm` and `GRiwrm` with specific column names:
Using that information, we must create the GRiwrm object that lists the nodes and describes the network diagram. It is a dataframe of class `GRiwrm` with specific column names:
- `id`: the identifier of the node in the network.
- `down`: the identifier of the next hydrological node downstream.
- `length`: hydraulic distance to the next hydrological downstream node.
- `model`: Name of the hydrological model used (E.g. "RunModel_GR4J"). `NA` for other type of node.
- `area`: Area of the sub-catchment (km<sup>2</sup>). Used for hydrological model such as GR models. `NA` if not used.
- `model`: Name of the hydrological model used (E.g. "RunModel_GR4J"). `NA` for other types of nodes.
- `area`: Area of the sub-catchment (km<sup>2</sup>). Used for hydrological models such as GR models. `NA` if not used.
`GRiwrm` function helps to rename the columns of the dataframe and assign the variable classes.
The `CreateGRiwrm` function helps to rename the columns of the dataframe and assign the variable classes.
```{r griwrm}
seine_nodes$id_aval[seine_nodes$id_aval == ""] <- NA
@@ -65,18 +70,18 @@ The daily mean precipitation and potential evaporation at the scale of the inter
The daily naturalised flow is provided by @hydratecActualisationBaseDonnees2011a.
These data are embedded in a the R package 'seinebasin' which is not publicly available.
These data are embedded in the R package 'seinebasin', which is not publicly available.
```{r QNAT, warning=FALSE, message=FALSE}
library(seinebasin)
data(QNAT)
```
## Generate the GRIWRM InputsModel object
## Generate the GRiwrm InputsModel object
The GRIWRM InputsModel object is a list of **airGR** InputsModel. The identifier of the sub-basin is used as key in the list which is ordered from upstream to downstream.
The GRiwrm InputsModel object is a list of **airGR** InputsModel objects. The identifier of the sub-basin is used as a key in the list, which is ordered from upstream to downstream.
The **airGR** CreateInputsModel function is extended in order to handle the griwrm object which describe the basin diagram:
The **airGR** CreateInputsModel function is extended in order to handle the GRiwrm object that describes the basin diagram:
The lag $\tau_0$ and route $K_0$ parameters of TGR are expressed as time delay in hours corresponding to the delay time between the farest upstream inlet and the outlet of the sub-basin.
Almost all sub basin has only a lag parameter. The only exception is for La Marne à Noisiel (NOISI_17) that has a routing parameter which can be approximated to a single lag parameter equals to $\tau_0 + K_0$.
The lag $\tau_0$ and routing $K_0$ parameters of TGR are expressed as time delay in hours corresponding to the delay time between the farthest upstream inlet and the outlet of the sub-basin.
Almost all sub basins have a routing parameter equal to 0. The only exception is for La Marne à Noisiel (NOISI_17) that has a routing parameter that can be approximated to a single lag parameter equal to $\tau_0 + K_0$.
This lag parameter has to be converted in a speed in m/s used in the **airGR** lag model:
title: "Seine_03: Calibration of a semi-distributive GR model network"
title: "Seine_03: Calibration of a semi-distributed GR model network"
author: "David Dorchies"
vignette: >
%\VignetteIndexEntry{Seine_03: Calibration of a semi-distributive GR model network}
%\VignetteIndexEntry{Seine_03: Calibration of a semi-distributed GR model network}
%\VignetteEngine{knitr::rmarkdown}
%\VignetteEncoding{UTF-8}
---
...
...
@@ -11,6 +11,8 @@ vignette: >
knitr::opts_chunk$set(echo = TRUE)
```
The aim of this vignette is to perform a first calibration of the semi-distributed hydrological model, using naturalised discharge.
## Loading libraries, network and time series data
Run `vignette("01_First_network", package = "airGRiwrm")` and `vignette("02_First_run", package = "airGRiwrm")` before this one in order to create the Rdata files loaded below:
...
...
@@ -25,11 +27,14 @@ data(QNAT)
## InputsCrit object
We need then to prepare the InputsCrit object that is necessary to define the calibration objective function. We chose here the KGE' criterion:
## Integration of the reservoir connections into the model
The goal of this vignette is to add reservoirs connections to the Seine River model previously set up. A complete description of the reservoir system of the Seine River can be found in @dehayEtudeImpactChangement2012
The aim of this vignette is to add reservoir connections to the Seine River model previously set up. A complete description of the reservoir system of the Seine River can be found in @dehayEtudeImpactChangement2012
### Add connections to the reservoirs in the gauging station network
Because of the **airGR** SD model structure, the topology of the network will differ between the physical one and the modelled one.
Because of the **airGR** SD model structure, the topology of the network will differ between the physical one and the modelled one. In the following, we provide details for each of the 4 lakes, and then we present the complete SD network.
For example, the physical topology of the Aube lake is represented below:
First, the physical topology of the Aube lake is represented below:
In the SD model, we do not model intermediate flows as well as the reservoir in the catchment. As a result, an equivalent topology compatible with **airGR** will be the one below:
In the SD model, we do not simulate intermediate flows as well as the reservoir in the catchment. As a result, an equivalent topology compatible with **airGR** will be the one below:
```{r mmd_aube_model, echo = FALSE}
mmd("
...
...
@@ -65,9 +65,9 @@ style AUBE_R3 fill:#fcc
", width = "100%")
```
`AUBE_P2` will propagate negative flow to downstream (reservoir inflows) while `AUBE_R3` will propagate positive flow to downstream (reservoir releases).
`AUBE_P2` will propagate negative flows downstream (reservoir inflows) while `AUBE_R3` will propagate positive flows downstream (reservoir releases).
With $Q_{PANNEC\_P} = - Q_{CHAUM\_07}$ as all the upstream flow fills the reservoir.
With $Q_{PANNEC\_P} = - Q_{CHAUM\_07}$ as all the upstream flows fill the reservoir.
Marne lake can be mapped as:
The Marne lake can be mapped as:
```{r mmd_marne_true, echo = FALSE}
mmd("
...
...
@@ -142,7 +142,8 @@ style MARNE_R25 fill:#fcc
style M fill:#ccf
", width = "100%")
```
And can be modeled as:
And can be modelled as:
```{r mmd_marne_model, echo = FALSE}
mmd("
...
...
@@ -197,17 +198,20 @@ data(QNAT)
data(Qreservoirs)
```
## How to handle online reservoir? The Pannecière lake case.
## How to handle an inline reservoir? The Pannecière lake case.
There are two possibilities:
- truncate the system and only simulate the release of the reservoir if we only care about downstream flows
- add a fictive inlet to the reservoir that will store all upstream discharge into the reservoir
If we know in advance the flow released by the reservoir, upstream flow informations is not usefull for the simulation. But reservoir management simulation would need upstream flow informations for simulating the reservoir state during simulation. The second alternative will be useful for the next phases.
If we know in advance the flows released by the reservoir, upstream flow informations is not useful for the simulation. However, reservoir management simulation would require upstream flow informations for simulating the reservoir state during simulation. The second alternative will therefore be useful for the next phases.
