plot.OutputsModel.R

plot.OutputsModel <- function(x, Qobs = NULL, IndPeriod_Plot = NULL, BasinArea = NULL, which = "all", log_scale = FALSE,
                              cex.axis = 1, cex.lab = 0.9, cex.leg = 0.9, lwd = 1, verbose = TRUE, ...) {
  
  OutputsModel <- x

  if (!inherits(OutputsModel, "GR") & !inherits(OutputsModel, "CemaNeige")) {
    stop(paste("OutputsModel not in the correct format for default plotting", sep = ""))
  }
  
  
  BOOL_Dates <- FALSE; 
      if ("DatesR" %in% names(OutputsModel)) { BOOL_Dates <- TRUE; }
  BOOL_Pobs <- FALSE; 
      if ("Precip" %in% names(OutputsModel)) { BOOL_Pobs <- TRUE; }
  BOOL_Qsim <- FALSE; 
      if ("Qsim"   %in% names(OutputsModel)) { BOOL_Qsim <- TRUE; }
  BOOL_Qobs <- FALSE;
      if (BOOL_Qsim & length(Qobs) == length(OutputsModel$Qsim)) {
        if (sum(is.na(Qobs)) != length(Qobs)) {
          BOOL_Qobs <- TRUE
        }
      } else if (inherits(OutputsModel, "GR")) {
        warning("Incorrect length of Qobs. Time series of observed flow not drawn.")
      }
  BOOL_Snow <- FALSE;
      if ("CemaNeigeLayers" %in% names(OutputsModel)) { if ("SnowPack" %in% names(OutputsModel$CemaNeigeLayers[[1]])) { BOOL_Snow <- TRUE; } }
  BOOL_Psol <- FALSE;
      if ("CemaNeigeLayers" %in% names(OutputsModel)) { if ("Psol"     %in% names(OutputsModel$CemaNeigeLayers[[1]])) { BOOL_Psol <- TRUE; } }


  if ( is.null(     which)) { stop("which must be a vector of character") } 
  if (!is.vector(   which)) { stop("which must be a vector of character") } 
  if (!is.character(which)) { stop("which must be a vector of character") } 
  if (any(!which %in% c("all", "Precip", 'Temp', "SnowPack", "Flows", "Regime", "CumFreq", "CorQQ"))) {
    stop("Incorrect element found in argument which:\nit can only contain 'all', 'Precip', 'Temp', 'SnowPack', 'Flows', 'Regime', 'CumFreq' or 'CorQQ'")
  }
  if (all(which %in% c("Temp", "SnowPack")) & !inherits(OutputsModel, "CemaNeige")) {
    stop("Incorrect element found in argument which:\nwithout CemaNeige it can only contain 'all', 'Precip', 'Flows', 'Regime', 'CumFreq' or 'CorQQ'")
  }
  if (length(unique(which %in% c("Temp", "SnowPack"))) == 2 & !inherits(OutputsModel, "CemaNeige")) {
    warning("Incorrect element found in argument which:\nit can only contain 'all', 'Precip', 'Flows', 'Regime', 'CumFreq' or 'CorQQ'\nwithout CemaNeige 'Temp' and 'SnowPack' are not available")
  }  
  
  if ("all" %in% which) {
    which <- c("Precip", "Temp", "SnowPack", "Flows", "Regime", "CumFreq", "CorQQ")
  }


  if (!BOOL_Dates) {
    stop(paste("OutputsModel must contain at least DatesR to allow plotting", sep = "")) }
  if (inherits(OutputsModel, "GR") & !BOOL_Qsim) {
    stop(paste("OutputsModel must contain at least Qsim to allow plotting", sep = "")) }

  if (BOOL_Dates) {
    MyRollMean1 <- function(x, n) {
      return(filter(x, rep(1/n, n), sides = 2)); }
    MyRollMean2 <- function(x, n) {
      return(filter(c(tail(x, n%/%2), x, x[1:(n%/%2)]), rep(1/n, n), sides = 2)[(n%/%2+1):(length(x)+n%/%2)]); }
    MyRollMean3 <- function(x, n) {
      return(filter(x, filter = rep(1/n, n), sides = 2, circular = TRUE))
    }    
    BOOL_TS  <- FALSE;
    TimeStep <- difftime(tail(OutputsModel$DatesR, 1), tail(OutputsModel$DatesR, 2), units = "secs")[[1]];
    if (inherits(OutputsModel, "hourly" ) & TimeStep %in% (                  60*60)) { BOOL_TS <- TRUE; NameTS <- "hour" ; plotunit <- "[mm/h]";     formatAxis <- "%m/%Y"; }
    if (inherits(OutputsModel, "daily"  ) & TimeStep %in% (               24*60*60)) { BOOL_TS <- TRUE; NameTS <- "day"  ; plotunit <- "[mm/d]";     formatAxis <- "%m/%Y"; }
    if (inherits(OutputsModel, "monthly") & TimeStep %in% (c(28, 29, 30, 31)*24*60*60)) { BOOL_TS <- TRUE; NameTS <- "month"; plotunit <- "[mm/month]"; formatAxis <- "%m/%Y"; }
    if (inherits(OutputsModel, "yearly" ) & TimeStep %in% (    c(365, 366)*24*60*60)) { BOOL_TS <- TRUE; NameTS <- "year" ; plotunit <- "[mm/y]";     formatAxis <- "%Y"   ; }
    if (!BOOL_TS) { stop(paste("the time step of the model inputs could not be found", sep = "")) } 
  }
  if (length(IndPeriod_Plot) == 0) { IndPeriod_Plot <- 1:length(OutputsModel$DatesR); }