diff --git a/DESCRIPTION b/DESCRIPTION
index 28ae94d3a0e00b109594a6ce5eaed7accbae6152..b4f18fc58f7383911324d8099b49c44be93e4230 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,6 +1,6 @@
Package: griwrm
Title: GR based Integrated Water Resource Management
-Version: 0.2.0
+Version: 0.2.1
Authors@R:
person(given = "David",
family = "Dorchies",
diff --git a/R/CreateInputsCrit.GriwrmInputsModel.R b/R/CreateInputsCrit.GriwrmInputsModel.R
index db1769fd458dfcddfdbbac1f8cae75023cc8e282..d1cb9b57655b41edf243eb8c4ab440560e65ff7f 100644
--- a/R/CreateInputsCrit.GriwrmInputsModel.R
+++ b/R/CreateInputsCrit.GriwrmInputsModel.R
@@ -2,7 +2,7 @@
#' @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 Qobs matrix or data frame containing observed flows. Column names correspond to nodes ID
#' @param ... further arguments passed to \code{\link[airGR]{CreateInputsCrit}}.
#'
#' @return Object of class \emph{GriwrmInputsCrit}
@@ -10,7 +10,7 @@
CreateInputsCrit.GriwrmInputsModel <- function(InputsModel,
FUN_CRIT = airGR::ErrorCrit_NSE,
RunOptions,
- gits,
+ Qobs,
...) {
InputsCrit <- list()
class(InputsCrit) <- append(class(InputsCrit), "GriwrmInputsCrit")
@@ -20,7 +20,7 @@ CreateInputsCrit.GriwrmInputsModel <- function(InputsModel,
InputsModel = IM,
FUN_CRIT = FUN_CRIT,
RunOptions = RunOptions[[IM$id]],
- Obs = gits[[IM$id]]$Qobs[RunOptions[[IM$id]]$IndPeriod_Run],
+ Obs = Qobs[RunOptions[[IM$id]]$IndPeriod_Run, IM$id],
...
)
}
diff --git a/R/CreateInputsModel.Griwrm.R b/R/CreateInputsModel.Griwrm.R
index 13c80194bfca9fa2e9b67a0d29db29b39e535333..d0efe8b707ea47844bfc7f5003faea2a9f0fa451 100644
--- a/R/CreateInputsModel.Griwrm.R
+++ b/R/CreateInputsModel.Griwrm.R
@@ -2,19 +2,25 @@
#'
#' @param x Ginet object describing the diagram of the semi-distributed model, see \code{[Ginet]}.
#' @param girop Girop object giving the run-off model parameters, see \code{[Girop]}.
-#' @param gits Gits object giving the observation time series, see \code{[Gits]}.
+#' @param DateR Vector of POSIXlt observation time steps.
+#' @param Precip Matrix or data frame of numeric containing precipitation in mm. Column names correspond to node IDs.
+#' @param PotEvap Matrix or data frame of numeric containing potential evaporation in mm. Column names correspond to node IDs.
+#' @param Qobs Matrix or data frame of numeric containing potential observed flow in mm. Column names correspond to node IDs.
#' @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]{CreateInputsModel}}.
#'
#' @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, DatesR, Precip, PotEvap, Qobs, verbose = TRUE, ...) {
InputsModel <- CreateEmptyGriwrmInputsModel()
+ Qobs[is.na(Qobs)] <- -99 # airGRCreateInputsModel doesn't accept NA values
for(id in getNodeRanking(x)) {
if(verbose) cat("CreateInputsModel.griwrm: Treating sub-basin", id, "...\n")
- InputsModel[[id]] <- CreateOneGriwrmInputsModel(id, x, girop, gits, ...)
+ InputsModel[[id]] <- CreateOneGriwrmInputsModel(
+ id, x, girop, DatesR,Precip[,id], PotEvap[,id], Qobs, ...
+ )
}
return(InputsModel)
}
@@ -35,10 +41,13 @@ CreateEmptyGriwrmInputsModel <- function() {
#' @param id string of the node identifier
#' @param ginet See \code{[Ginet]}.
#' @param girop See \code{[Girop]}.
-#' @param gits See \code{[Gits]}.
-#'
+#' @param DatesR vector of dates required to create the GR model and CemaNeige module inputs.
+#' @param Precip time series of potential evapotranspiration (catchment average) (mm/time step).
+#' @param PotEvap time series of potential evapotranspiration (catchment average) (mm/time step).
+#' @param Qobs Matrix or data frame of numeric containing observed flow (mm/time step). Column names correspond to node IDs.
+##'
#' @return \emph{InputsModel} object for one.
-CreateOneGriwrmInputsModel <- function(id, ginet, girop, gits) {
+CreateOneGriwrmInputsModel <- function(id, ginet, girop, DatesR, Precip, PotEvap, Qobs) {
node <- ginet[ginet$id == id,]
FUN_MOD <- girop$model[girop$id == id]
@@ -51,7 +60,7 @@ CreateOneGriwrmInputsModel <- function(id, ginet, girop, gits) {
if(length(UpstreamNodes) > 0) {
# Sub-basin with hydraulic routing
for(idUpstrNode in UpstreamNodes) {
- Qupstream1 <- matrix(gits[[idUpstrNode]]$Qobs, ncol = 1)
+ Qupstream1 <- matrix(Qobs[,idUpstrNode], ncol = 1)
if(is.null(Qupstream)) {
Qupstream <- Qupstream1
} else {
@@ -68,9 +77,9 @@ CreateOneGriwrmInputsModel <- function(id, ginet, girop, gits) {
# Set model inputs with the airGR function
InputsModel <- CreateInputsModel(
FUN_MOD,
- DatesR = gits$date,
- Precip = gits[[id]]$Precip,
- PotEvap = gits[[id]]$PotEvap,
+ DatesR = DatesR,
+ Precip = Precip,
+ PotEvap = PotEvap,
Qupstream = Qupstream,
LengthHydro = LengthHydro,
BasinAreas = BasinAreas
diff --git a/R/gits.R b/R/gits.R
deleted file mode 100644
index 929fbc38b2e52eb1b3ba1eed7f891ca09e8f7e53..0000000000000000000000000000000000000000
--- a/R/gits.R
+++ /dev/null
@@ -1,46 +0,0 @@
-#' Title
-#'
-#' @param id string of the id of the node
-#' @param ts numeric matrix or data frame containing 3 columns for precipitation, evaporation, and observed flow
-#' @param cols named list or vector used for matching the columns of ts with the required columns names which are "Precip", "PotEvap", and "Qobs".
-#'
-#' @return \emph{Gits} class object which is a list containing a `date` element (Vector of PosiXlt timestamps) and an element named the id of the node containing a dataframe with observed data.
-#' @export
-Gits <- function(id, ts,
- cols = list(date = "date", Precip = "Precip", PotEvap = "PotEvap", Qobs = "Qobs")) {
-
- cols <- as.list(cols)
- ts <- dplyr::rename(ts, unlist(cols))
-
- if(any(is.na(ts$Qobs))) {
- stop("Qobs should not contain any NA")
- }
- if(any(ts$Qobs < 0)) {
- stop("Qobs should be strictly positive")
- }
-
- gitsOut <- list(date = ts$date)
- cols$date <- NULL
- gitsOut[[id]] <- dplyr::select(ts, names(cols))
- class(gitsOut) <- append(class(gitsOut), "Gits")
- gitsOut
-}
-
-#' Merge two gits objects with identical date time series.
-#'
-#' @param x Gits object to merge (See [Gits]).
-#' @param y Gits object to merge (See [Gits]).
-#' @param ... For merge generic function compatibility.
-#'
-#' @return Gits object merged with one item `Date` and Items corresponding to each node.
-#' @export
-merge.Gits <- function(x, y, ...) {
- if(!is(y, "Gits")) {
- stop("A Gits class object can only be merged with a Gits class object")
- }
- if(! identical(x$date, y$date)) {
- stop("Time series dates are not identical")
- }
- y$date <- NULL
- x <- utils::modifyList(x, y)
-}
diff --git a/inst/seine_data/network_gauging_stations.txt b/inst/seine_data/network_gauging_stations.txt
new file mode 100644
index 0000000000000000000000000000000000000000..83e76b49a61b27f4f2bdef23067613c659dacdb8
--- /dev/null
+++ b/inst/seine_data/network_gauging_stations.txt
@@ -0,0 +1,26 @@
+id_sgl id_hydro lambert2.x lambert2.y area description id_aval distance_aval
+TRANN_01 766626.06 2369152.25 1557.06 l'Aube à Trannes ARCIS_24 68100
+GURGY_02 H2221010 689713.00 2320549.00 3819.77 L'Yonne à Gurgy COURL_23 83612
+BRIEN_03 H2482010 695313.00 2332549.00 2979.77 L'Armançon à Brienon-sur-Armançon COURL_23 84653
+STDIZ_04 H5071010 791113.00 2407349.00 2347.53 La Marne à Saint-Dizier CHALO_21 85570
+PARIS_05 H5920010 602213.00 2427449.00 43824.66 La Seine à Paris
+BAR-S_06 H0400010 751913.00 2348349.00 2340.37 La Seine à Bar-sur-Seine MERY-_22 79766
+CHAUM_07 716518.38 2241747.00 216.50 L'Yonne à Chaumard GURGY_02 153074
+CUSSY_08 H2172310 726013.00 2275149.00 247.99 Le Cousin à Cussy-les-Forges GURGY_02 91378
+STGER_09 718512.88 2266649.25 402.74 La Cure à St-Germain GURGY_02 94152
+GUILL_10 H2322020 731013.00 2282349.00 488.71 Le Serein à Guillon CHABL_12 66026
+AISY-_11 H2452020 742413.00 2298149.00 1349.51 L'Armançon à Aisy-sur-Armançon BRIEN_03 102428
+CHABL_12 H2342010 709613.00 2314149.00 1116.27 Le Serein à Chablis COURL_23 111781
+NOGEN_13 685912.88 2389349.25 9182.39 La Seine à Nogent-sur-Seine MONTE_15 63215
+EPISY_14 H3621010 633413.00 2371049.00 3916.71 Le Loing à Épisy ALFOR_16 89196
+MONTE_15 645112.88 2376849.25 21199.39 La Seine à Montereau ALFOR_16 94475
+ALFOR_16 H4340020 606013.00 2420349.00 30784.71 La Seine à Alfortville PARIS_05 9263
+NOISI_17 H5841010 620913.00 2428949.00 12547.72 La Marne à Noisiel PARIS_05 39384
+MONTR_18 H5752020 638013.00 2431849.00 1184.81 Le Grand Morin à Montry NOISI_17 37915
+LOUVE_19 H5083050 791613.00 2393949.00 461.74 La Blaise à Louvemont CHALO_21 86165
+LASSI_20 H1362010 759513.00 2385549.00 876.53 La Voire à Lassicourt ARCIS_24 43618
+CHALO_21 H5201010 747713.00 2441349.00 6291.55 La Marne à Châlons-sur-Marne NOISI_17 237937
+MERY-_22 H0810010 714913.00 2390949.00 3899.62 La Seine à Méry-sur-Seine NOGEN_13 49933
+COURL_23 H2721010 660813.00 2370449.00 10687.35 L'Yonne à Courlon-sur-Yonne MONTE_15 26159
+ARCIS_24 H1501010 733313.00 2394749.00 3594.60 L'Aube à Arcis-sur-Aube NOGEN_13 70926
+VITRY_25 H5172010 768513.00 2418849.00 2109.14 La Saulx à Vitry-en-Perthois CHALO_21 38047
diff --git a/vignettes/V01_First_network.Rmd b/vignettes/V01_First_network.Rmd
index f834388dfa8a949ef07cbc104ed7faedb4e1cde4..5b9bb0ad5d2f126f4c77bc4710a629df2b7fb89f 100644
--- a/vignettes/V01_First_network.Rmd
+++ b/vignettes/V01_First_network.Rmd
@@ -22,8 +22,8 @@ List of nodes
```{r}
seine_nodes <- readr::read_delim(
- file = "https://stratus.irstea.fr/d/0b18e688851a45478f7a/files/?p=/climaware_hydro/griwrm_network.txt&dl=1",
- delim = ";"
+ file = system.file("seine_data", "network_gauging_stations.txt", package = "griwrm"),
+ delim = "\t"
)
seine_nodes
```
@@ -77,40 +77,45 @@ urls <-
)
names(urls) <- ginet$id
-load_ts <- function(x) {
- ts <- readr::read_delim(file = x,
- delim = ";", skip = 16, trim_ws = TRUE)
- ts$Date <- as.POSIXlt(lubridate::ymd(ts$Date))
- # Interpolation of data gap in the discharge time serie
- ts[ts$Qnat < 0, "Qnat"] <- NA
- if(is.na(ts$Qnat[nrow(ts)])) {
- ts$Qnat[nrow(ts)] <- 0 # No value at the end: time serie converge to zero
+Precip <- NULL
+PotEvap <- NULL
+Qobs <- NULL
+
+MergeTS <- function(dfOld, id, dfNew) {
+ names(dfNew) <- c("DatesR", id) # Renaming columns of the new input into date and sub-basin ID
+ if(is.null(dfOld)) {
+ dfOut <- dfNew # Creation of the first column
+ } else {
+ dfOut <- merge(dfOld, dfNew, by = "DatesR", all = TRUE) # Merge the new columns respecting to date column
}
- ts$Qnat <- zoo::na.approx(ts$Qnat)
- ts
+ return(dfOut)
}
-l <- lapply(urls, load_ts)
-
-```
-
-`Gits` object is a list containing a item named `date` with a timestamp vector of the time series and items named by the identifier of each node. These items contain a dataframe with the observations.
-
-The Gits function creates a `Gits` object
-
-```{r}
-gits <- Gits(ginet$id[1], l[[ginet$id[1]]], cols = list(date = "Date", Precip = "Ptot", PotEvap = "ETP", Qobs = "Qnat"))
-```
-
-Copy the observations for each node the ginet network:
-
-```{r}
for(id in ginet$id) {
- l[[id]]$Qnat
- l[[id]]$Qnat <- l[[id]]$Qnat * 86.4 / girop$area[girop$id == id]
- l[[id]]$Qnat[l[[id]]$Qnat < 0] <- NA
- gits <- merge(gits, Gits(id, l[[id]], cols = list(date = "Date", Precip = "Ptot", PotEvap = "ETP", Qobs = "Qnat")))
+ url <-
+ file.path(
+ "https://stratus.irstea.fr/d/0b18e688851a45478f7a/files/?p=/climaware_hydro/Q_OBS_NAT",
+ paste0(id, "_NAT.txt&dl=1")
+ )
+ ts <- readr::read_delim(file = url,
+ delim = ";", skip = 16, trim_ws = TRUE)
+ # Date conversion to POSIX
+ ts$Date <- as.POSIXlt(lubridate::ymd(ts$Date))
+ # Ptot column is merged into Precip dataframe
+ Precip <- MergeTS(Precip, id, ts[,c("Date", "Ptot")])
+ # ETP column is merged into PotEvap dataframe
+ PotEvap <- MergeTS(PotEvap, id, ts[,c("Date", "ETP")])
+ # Convert Qobs from m3/s to mm
+ ts$Qnat <- ts$Qnat * 86.4 / girop$area[girop$id == id]
+ # Setting data gaps to NA
+ ts$Qnat[ts$Qnat <= 0] <- NA
+ # Qnat column is merged into Qobs dataframe
+ Qobs <- MergeTS(Qobs, id, ts[,c("Date", "Qnat")])
}
+DatesR <- Precip$DatesR
+Precip$DateR <- NULL
+PotEvap$DateR <- NULL
+Qobs$DateR <- NULL
```
## Generate the GRIWRM InputsModel object
@@ -121,7 +126,7 @@ The airGR CreateInputsModel function is extended in order to handle the ginet ob
```{r}
-InputsModel <- CreateInputsModel(ginet, girop, gits)
+InputsModel <- CreateInputsModel(ginet, girop, DatesR, Precip, PotEvap, Qobs)
```
@@ -129,7 +134,8 @@ InputsModel <- CreateInputsModel(ginet, girop, gits)
## Save data for next vignettes
```{r}
+
dir.create("_cache", showWarnings = FALSE)
-save(ginet, girop, gits, InputsModel, file = "_cache/V01.RData")
+save(ginet, girop, Qobs, InputsModel, file = "_cache/V01.RData")
```
diff --git a/vignettes/V03_First_Calibration.Rmd b/vignettes/V03_First_Calibration.Rmd
index d31b28e54fe1aee2c03c09e03f8ba34a80782f20..919dd3116cfda8903d5e6a66b8a3155de374c0c0 100644
--- a/vignettes/V03_First_Calibration.Rmd
+++ b/vignettes/V03_First_Calibration.Rmd
@@ -27,7 +27,7 @@ load("_cache/V02.RData")
InputsCrit <- CreateInputsCrit(
InputsModel = InputsModel,
FUN_CRIT = airGR::ErrorCrit_KGE2,
- RunOptions = RunOptions, gits = gits
+ RunOptions = RunOptions, Qobs = Qobs
)
str(InputsCrit)
```
@@ -65,13 +65,20 @@ OutputsModels <- RunModel(
)
```
+## Save calibration data for next vignettes
+
+```{r}
+save(giropMichel, file = "_cache/V03.RData")
+```
+
+
## Plot the result for each basin
```{r, fig.height = 5, fig.width = 8}
htmltools::tagList(lapply(
names(OutputsModels),
function(x) {
- plot(OutputsModels[[x]], Qobs = gits[[x]]$Qobs[RunOptions[[x]]$IndPeriod_Run] , main = x)
+ plot(OutputsModels[[x]], Qobs = Qobs[RunOptions[[x]]$IndPeriod_Run,x] , main = x)
}
))
```
diff --git a/vignettes/v02_First_run.Rmd b/vignettes/v02_First_run.Rmd
index 67ca6f03cfbfd7795051b7cce3ca43853744da5a..6a6a2d877289264c278cf2c2f8a68c68248ceb4d 100644
--- a/vignettes/v02_First_run.Rmd
+++ b/vignettes/v02_First_run.Rmd
@@ -123,7 +123,7 @@ save(RunOptions, file = "_cache/V02.RData")
htmltools::tagList(lapply(
names(OutputsModelsClimAware),
function(x) {
- plot(OutputsModelsClimAware[[x]], Qobs = gits[[x]]$Qobs[IndPeriod_Run] , main = x)
+ plot(OutputsModelsClimAware[[x]], Qobs = Qobs[IndPeriod_Run,x] , main = x)
}
))
```