CreateIniStates.R 7.54 KB
Newer Older
1
CreateIniStates <- function(FUN_MOD, InputsModel,
2
                            ProdStore = 350, RoutStore = 90, ExpStore = NULL,
3
4
                            UH1 = NULL, UH2 = NULL,
                            GCemaNeigeLayers = NULL, eTGCemaNeigeLayers = NULL,
5
                            GthrCemaNeigeLayers = NULL, GlocmaxCemaNeigeLayers = NULL,
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
                            verbose = TRUE) {
  
  
  ObjectClass <- NULL
  
  UH1n <- 20L
  UH2n <- UH1n * 2L
  
  
  ## check FUN_MOD
  BOOL <- FALSE
  if (identical(FUN_MOD, RunModel_GR4H)) {
    ObjectClass <- c(ObjectClass, "GR", "hourly")
    BOOL <- TRUE
  }
  if (identical(FUN_MOD, RunModel_GR4J) |
      identical(FUN_MOD, RunModel_GR5J) |
      identical(FUN_MOD, RunModel_GR6J)) {
    ObjectClass <- c(ObjectClass, "GR", "daily")
    BOOL <- TRUE
  }
  if (identical(FUN_MOD, RunModel_GR2M)) {
    ObjectClass <- c(ObjectClass, "GR", "monthly")
    BOOL <- TRUE
  }
  if (identical(FUN_MOD, RunModel_GR1A)) {
    stop("'RunModel_GR1A' does not require 'IniStates' object")
  }
  if (identical(FUN_MOD, RunModel_CemaNeige)) {
    ObjectClass <- c(ObjectClass, "CemaNeige", "daily")
    BOOL <- TRUE
  }
  if (identical(FUN_MOD, RunModel_CemaNeigeGR4J) |
      identical(FUN_MOD, RunModel_CemaNeigeGR5J) |
      identical(FUN_MOD, RunModel_CemaNeigeGR6J)) {
    ObjectClass <- c(ObjectClass, "GR", "CemaNeige", "daily")
    BOOL <- TRUE
  }
  if (!BOOL) {
    stop("Incorrect 'FUN_MOD' for use in 'CreateIniStates'")
  }
  
  ## check InputsModel
  if (!inherits(InputsModel, "InputsModel")) {
    stop("'InputsModel' must be of class 'InputsModel'")
  }
  if ("GR" %in% ObjectClass & !inherits(InputsModel, "GR")) {
    stop("'InputsModel' must be of class 'GR'")
  }
  if ("CemaNeige" %in% ObjectClass &
      !inherits(InputsModel, "CemaNeige")) {
57
    stop("'InputsModel' must be of class 'CemaNeige'")
58
59
60
61
  }
  
  
  ## check states
62
63
64
  if (any(eTGCemaNeigeLayers > 0)) {
    stop("Positive values are not allowed for 'eTGCemaNeigeLayers'")
  }  
65
  
66
67
68
69
  if (identical(FUN_MOD, RunModel_GR6J) | identical(FUN_MOD, RunModel_CemaNeigeGR6J)) {
    if (is.null(ExpStore)) {
      stop("'RunModel_*GR6J' need an 'ExpStore' value")
    }
70
71
  } else if (!is.null(ExpStore)) {
    if (verbose) {
72
      warning(sprintf("'%s' does not require 'ExpStore'. Value set to NA", as.character(substitute(FUN_MOD))))
73
    }
74
75
76
    ExpStore <- Inf
  }
  
77
  if (identical(FUN_MOD, RunModel_GR2M)) {
78
    if (!is.null(UH1)) {
79
      if (verbose) {
80
        warning(sprintf("'%s' does not require 'UH1'. Values set to NA", as.character(substitute(FUN_MOD))))
81
      }
82
83
84
      UH1 <- rep(Inf, UH1n)
    }
    if (!is.null(UH2)) {
85
      if (verbose) {
86
        warning(sprintf("'%s' does not require 'UH2'. Values set to NA", as.character(substitute(FUN_MOD))))
87
      }
88
89
90
91
      UH2 <- rep(Inf, UH2n)
    }
  }
  
92
93
  if ((identical(FUN_MOD, RunModel_GR5J) | identical(FUN_MOD, RunModel_CemaNeigeGR5J)) & !is.null(UH1)) {
    if (verbose) {
94
      warning(sprintf("'%s' does not require 'UH1'. Values set to NA", as.character(substitute(FUN_MOD))))
95
    }
96
97
    UH1 <- rep(Inf, UH1n)
  }
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
 
  if ("CemaNeige" %in% ObjectClass & ! "GR" %in% ObjectClass) {
    if (!is.null(ProdStore)) {
      if (verbose) {
        warning(sprintf("'%s' does not require 'ProdStore'. Values set to NA", as.character(substitute(FUN_MOD))))
      }
    }
    ProdStore <- Inf
    if (!is.null(RoutStore)) {
      if (verbose) {
        warning(sprintf("'%s' does not require 'RoutStore'. Values set to NA", as.character(substitute(FUN_MOD))))
      }
    }
    RoutStore <- Inf
    if (!is.null(ExpStore)) {
      if (verbose) {
        warning(sprintf("'%s' does not require 'ExpStore'. Values set to NA", as.character(substitute(FUN_MOD))))
      }
    }
    ExpStore <- Inf
    if (!is.null(UH1)) {
      if (verbose) {
        warning(sprintf("'%s' does not require 'UH1'. Values set to NA", as.character(substitute(FUN_MOD))))
      }
    }
    UH1 <- rep(Inf, UH1n)
    if (!is.null(UH2)) {
      if (verbose) {
        warning(sprintf("'%s' does not require 'UH2'. Values set to NA", as.character(substitute(FUN_MOD))))
      }
    }
    UH2 <- rep(Inf, UH2n)
  }
131
132
133
134
135
  if("CemaNeige" %in% ObjectClass &
     (is.null(GCemaNeigeLayers) | is.null(eTGCemaNeigeLayers))) {
    stop("'RunModel_CemaNeigeGR*' need values for 'GCemaNeigeLayers' and 'GCemaNeigeLayers'")
  }
  if(!"CemaNeige" %in% ObjectClass &
136
137
     (!is.null(GCemaNeigeLayers) | !is.null(eTGCemaNeigeLayers))) {
    if (verbose) {
138
      warning(sprintf("'%s' does not require 'GCemaNeigeLayers' and 'GCemaNeigeLayers'. Values set to NA", as.character(substitute(FUN_MOD))))
139
    }
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
    GCemaNeigeLayers   <- Inf
    eTGCemaNeigeLayers <- Inf
  }
  
  
  ## set states
  if("CemaNeige" %in% ObjectClass) {
    NLayers <- length(InputsModel$LayerPrecip)
  } else {
    NLayers <- 1
  }
  
  
  ## manage NULL values
  if (is.null(ExpStore)) {
    ExpStore <- Inf 
  }
  if (is.null(UH1)) {
    if ("hourly"  %in% ObjectClass) {
      k <- 24
    } else {
      k <- 1
    }
163
    UH1 <- rep(Inf, UH1n * k)
164
165
166
167
168
169
170
  }
  if (is.null(UH2)) {
    if ("hourly"  %in% ObjectClass) {
      k <- 24
    } else {
      k <- 1
    }
171
    UH2 <- rep(Inf, UH2n * k)
172
173
174
175
176
177
  }
  if (is.null(GCemaNeigeLayers)) {
    GCemaNeigeLayers <- rep(Inf, NLayers)
  }
  if (is.null(eTGCemaNeigeLayers)) {
    eTGCemaNeigeLayers <- rep(Inf, NLayers)
178
  }
179
180
181
182
183
184
  if (is.null(GthrCemaNeigeLayers)) {
    GthrCemaNeigeLayers <- rep(Inf, NLayers)
  }
  if (is.null(GlocmaxCemaNeigeLayers)) {
    GlocmaxCemaNeigeLayers <- rep(Inf, NLayers)
  }
185
  
186
187
188
189
190
191
192
  
  # check negative values
  if (any(ProdStore < 0) | any(RoutStore < 0) |
      any(UH1 < 0) | any(UH2 < 0) |
      any(GCemaNeigeLayers < 0)) {
    stop("Negative values are not allowed for any of 'ProdStore', 'RoutStore', 'UH1', 'UH2', 'GCemaNeigeLayers'")
  }
193
194
195
196
197
198
199
200
201
202
203
204
  
  
  ## check length
  if (!is.numeric(ProdStore) || length(ProdStore) != 1L) {
    stop("'ProdStore' must be numeric of length one")
  }
  if (!is.numeric(RoutStore) || length(RoutStore) != 1L) {
    stop("'RoutStore' must be numeric of length one")
  }
  if (!is.numeric(ExpStore) || length(ExpStore) != 1L) {
    stop("'ExpStore' must be numeric of length one")
  }
205
  if ( "hourly" %in% ObjectClass & (!is.numeric(UH1) || length(UH1) != UH1n * 24)) {
206
207
    stop(sprintf("'UH1' must be numeric of length 480 (%i * 24)", UH1n))
  }
208
  if (!"hourly" %in% ObjectClass & (!is.numeric(UH1) || length(UH1) != UH1n)) {
209
210
    stop(sprintf("'UH1' must be numeric of length %i", UH1n))
  }
211
  if ( "hourly" %in% ObjectClass & (!is.numeric(UH2) || length(UH2) != UH2n * 24)) {
212
213
    stop(sprintf("'UH2' must be numeric of length 960 (%i * 24)", UH2n))
  }
214
  if (!"hourly" %in% ObjectClass & (!is.numeric(UH2) || length(UH2) != UH2n)) {
215
216
217
218
219
220
221
222
223
    stop(sprintf("'UH2' must be numeric of length %i (2 * %i)", UH2n, UH1n))
  }
  if (!is.numeric(GCemaNeigeLayers) || length(GCemaNeigeLayers) != NLayers) {
    stop(sprintf("'GCemaNeigeLayers' must be numeric of length %i", NLayers))
  }
  if (!is.numeric(eTGCemaNeigeLayers) || length(eTGCemaNeigeLayers) != NLayers) {
    stop(sprintf("'eTGCemaNeigeLayers' must be numeric of length %i", NLayers))
  }
  
224
225

  ## format output
226
227
  IniStates <- list(Store = list(Prod = ProdStore, Rout = RoutStore, Exp = ExpStore),
                    UH = list(UH1 = UH1, UH2 = UH2),
228
229
                    CemaNeigeLayers = list(G = GCemaNeigeLayers, eTG = eTGCemaNeigeLayers,
                                           Gthr = GthrCemaNeigeLayers, Glocmax = GlocmaxCemaNeigeLayers))
230
231
232
233
234
235
236
237
238
  IniStatesNA <- unlist(IniStates)
  IniStatesNA[is.infinite(IniStatesNA)] <- NA
  IniStatesNA <- relist(IniStatesNA, skeleton = IniStates)
  
  class(IniStatesNA) <- c("IniStates", ObjectClass)
  return(IniStatesNA)
  
  
}