From e42dfeede5180b1c1e4bad98fa1b8c5340eebd8b Mon Sep 17 00:00:00 2001 From: unknown <olivier.delaigue@ANPI1430.antony.irstea.priv> Date: Wed, 5 Apr 2017 11:00:45 +0200 Subject: [PATCH] v1.0.6.1 plot.OutputsModel gains graphical parameters arguments (cex.axis, cex.lab, cex.leg, lwd) and some parts of the code were cleaned --- DESCRIPTION | 4 +- R/plot.OutputsModel.R | 433 +++++++++++++++++++++------------------ man/plot.OutputsModel.Rd | 15 +- 3 files changed, 244 insertions(+), 208 deletions(-) diff --git a/DESCRIPTION b/DESCRIPTION index 1a14feec..d593fd3c 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,8 +1,8 @@ Package: airGR Type: Package Title: Suite of GR Hydrological Models for Precipitation-Runoff Modelling -Version: 1.0.6.0 -Date: 2017-03-31 +Version: 1.0.6.1 +Date: 2017-04-05 Authors@R: c( person("Laurent", "Coron", role = c("aut", "trl")), person("Charles", "Perrin", role = c("aut", "ths")), diff --git a/R/plot.OutputsModel.R b/R/plot.OutputsModel.R index 462903b6..50677f3d 100644 --- a/R/plot.OutputsModel.R +++ b/R/plot.OutputsModel.R @@ -1,29 +1,31 @@ -plot.OutputsModel <- function(x, Qobs = NULL, IndPeriod_Plot = NULL, BasinArea = NULL, which = "all", log_scale = FALSE, verbose = TRUE, ...){ +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")){ + if (!inherits(OutputsModel, "GR") & !inherits(OutputsModel, "CemaNeige")) { stop(paste("OutputsModel not in the correct format for default plotting \n", sep = "")) return(NULL) } - + + BOOL_Dates <- FALSE; - if("DatesR" %in% names(OutputsModel)){ BOOL_Dates <- TRUE; } + if ("DatesR" %in% names(OutputsModel)) { BOOL_Dates <- TRUE; } BOOL_Pobs <- FALSE; - if("Precip" %in% names(OutputsModel)){ BOOL_Pobs <- TRUE; } + if ("Precip" %in% names(OutputsModel)) { BOOL_Pobs <- TRUE; } BOOL_Qsim <- FALSE; - if("Qsim" %in% names(OutputsModel)){ BOOL_Qsim <- TRUE; } + 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; } } + if (BOOL_Qsim & length(Qobs) == length(OutputsModel$Qsim)) { if (sum(is.na(Qobs)) != length(Qobs)) { BOOL_Qobs <- TRUE; } } BOOL_Snow <- FALSE; - if("CemaNeigeLayers" %in% names(OutputsModel)){ if("SnowPack" %in% names(OutputsModel$CemaNeigeLayers[[1]])){ BOOL_Snow <- TRUE; } } + 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 ("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 \n"); return(NULL); } - if(!is.vector( which)){ stop("which must be a vector of character \n"); return(NULL); } - if(!is.character(which)){ stop("which must be a vector of character \n"); return(NULL); } + if ( is.null( which)) { stop("which must be a vector of character \n"); return(NULL); } + if (!is.vector( which)) { stop("which must be a vector of character \n"); return(NULL); } + if (!is.character(which)) { stop("which must be a vector of character \n"); return(NULL); } 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'") return(NULL) @@ -41,37 +43,37 @@ plot.OutputsModel <- function(x, Qobs = NULL, IndPeriod_Plot = NULL, BasinArea = } - if(!BOOL_Dates){ + if (!BOOL_Dates) { stop(paste("OutputsModel must contain at least DatesR to allow plotting \n", sep = "")); return(NULL); } - if(inherits(OutputsModel, "GR") & !BOOL_Qsim){ + if (inherits(OutputsModel, "GR") & !BOOL_Qsim) { stop(paste("OutputsModel must contain at least Qsim to allow plotting \n", sep = "")); return(NULL); } - if(BOOL_Dates){ - MyRollMean1 <- function(x, n){ + if (BOOL_Dates) { + MyRollMean1 <- function(x, n) { return(filter(x, rep(1/n, n), sides = 2)); } - MyRollMean2 <- function(x, n){ + 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){ + 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 \n", sep = "")); return(NULL); } + 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 \n", sep = "")); return(NULL); } } - if(length(IndPeriod_Plot) == 0){ IndPeriod_Plot <- 1:length(OutputsModel$DatesR); } - if(inherits(OutputsModel, "CemaNeige")){ NLayers <- length(OutputsModel$CemaNeigeLayers); } - PsolLayerMean <- NULL; if(BOOL_Psol){ - for(iLayer in 1:NLayers){ - if(iLayer == 1){ PsolLayerMean <- OutputsModel$CemaNeigeLayers[[iLayer]]$Psol/NLayers; + if (length(IndPeriod_Plot) == 0) { IndPeriod_Plot <- 1:length(OutputsModel$DatesR); } + if (inherits(OutputsModel, "CemaNeige")) { NLayers <- length(OutputsModel$CemaNeigeLayers); } + PsolLayerMean <- NULL; if (BOOL_Psol) { + for(iLayer in 1:NLayers) { + if (iLayer == 1) { PsolLayerMean <- OutputsModel$CemaNeigeLayers[[iLayer]]$Psol/NLayers; } else { PsolLayerMean <- PsolLayerMean + OutputsModel$CemaNeigeLayers[[iLayer]]$Psol/NLayers; } } } - BOOL_QobsZero <- FALSE; if(BOOL_Qobs){ SelectQobsNotZero <- (round(Qobs[IndPeriod_Plot] , 4) != 0); BOOL_QobsZero <- sum(!SelectQobsNotZero, na.rm = TRUE)>0; } - BOOL_QsimZero <- FALSE; if(BOOL_Qsim){ SelectQsimNotZero <- (round(OutputsModel$Qsim[IndPeriod_Plot], 4) != 0); BOOL_QsimZero <- sum(!SelectQsimNotZero, na.rm = TRUE)>0; } - if(BOOL_QobsZero & verbose){ warning("\t zeroes detected in Qobs -> some plots in the log space will not be created using all time-steps \n"); } - if(BOOL_QsimZero & verbose){ warning("\t zeroes detected in Qsim -> some plots in the log space will not be created using all time-steps \n"); } + BOOL_QobsZero <- FALSE; if (BOOL_Qobs) { SelectQobsNotZero <- (round(Qobs[IndPeriod_Plot] , 4) != 0); BOOL_QobsZero <- sum(!SelectQobsNotZero, na.rm = TRUE)>0; } + BOOL_QsimZero <- FALSE; if (BOOL_Qsim) { SelectQsimNotZero <- (round(OutputsModel$Qsim[IndPeriod_Plot], 4) != 0); BOOL_QsimZero <- sum(!SelectQsimNotZero, na.rm = TRUE)>0; } + if (BOOL_QobsZero & verbose) { warning("\t zeroes detected in Qobs -> some plots in the log space will not be created using all time-steps \n"); } + if (BOOL_QsimZero & verbose) { warning("\t zeroes detected in Qsim -> some plots in the log space will not be created using all time-steps \n"); } BOOL_FilterZero <- TRUE; ##Plots_choices @@ -85,66 +87,66 @@ plot.OutputsModel <- function(x, Qobs = NULL, IndPeriod_Plot = NULL, BasinArea = ##Options - BLOC <- TRUE; if(BLOC){ - cexaxis <- 1.0; cexlab <- 0.9; cexleg = 1.0; lwdLine = 1.8; lineX = 2.6; lineY = 2.6; bgleg <- NA + BLOC <- TRUE + if (BLOC) { + lwdk <- 1.8 + line <- 2.6 + bg <- NA matlayout <- NULL; iPlot <- 0; Sum1 <- sum(c(BOOLPLOT_Precip, BOOLPLOT_SnowPack, BOOLPLOT_Flows)) Sum2 <- sum(c(BOOLPLOT_Regime, BOOLPLOT_CumFreq, BOOLPLOT_CorQQ)) - if(BOOLPLOT_Precip){ + if (BOOLPLOT_Precip) { matlayout <- rbind(matlayout, c(iPlot+1, iPlot+1, iPlot+1)); iPlot <- iPlot+1; } - if(BOOLPLOT_Temp){ + if (BOOLPLOT_Temp) { matlayout <- rbind(matlayout, c(iPlot+1, iPlot+1, iPlot+1), c(iPlot+1, iPlot+1, iPlot+1)); iPlot <- iPlot+1; } - if(BOOLPLOT_SnowPack){ + if (BOOLPLOT_SnowPack) { matlayout <- rbind(matlayout, c(iPlot+1, iPlot+1, iPlot+1), c(iPlot+1, iPlot+1, iPlot+1)); iPlot <- iPlot+1; } - if(BOOLPLOT_Flows){ + if (BOOLPLOT_Flows) { matlayout <- rbind(matlayout, c(iPlot+1, iPlot+1, iPlot+1), c(iPlot+1, iPlot+1, iPlot+1)); iPlot <- iPlot+1; } - if((Sum1 >= 1 & Sum2 != 0) | (Sum1 == 0 & Sum2 == 3)){ + if ((Sum1 >= 1 & Sum2 != 0) | (Sum1 == 0 & Sum2 == 3)) { matlayout <- rbind(matlayout, c(iPlot+1, iPlot+2, iPlot+3), c(iPlot+1, iPlot+2, iPlot+3)); iPlot <- iPlot+3; } - if(Sum1 == 0 & Sum2 == 2){ + if (Sum1 == 0 & Sum2 == 2) { matlayout <- rbind(matlayout, c(iPlot+1, iPlot+2)); iPlot <- iPlot+2; } - if(Sum1 == 0 & Sum2 == 1){ + if (Sum1 == 0 & Sum2 == 1) { matlayout <- rbind(matlayout, iPlot+1); iPlot <- iPlot+1; } iPlotMax <- iPlot; # isRStudio <- Sys.getenv("RSTUDIO") == "1"; - # if(!isRStudio){ - # if(Sum1 == 1 & Sum2 == 0){ width = 10; height = 05; } - # if(Sum1 == 1 & Sum2 != 0){ width = 10; height = 07; } - # if(Sum1 == 2 & Sum2 == 0){ width = 10; height = 05; } - # if(Sum1 == 2 & Sum2 != 0){ width = 10; height = 07; } - # if(Sum1 == 3 & Sum2 == 0){ width = 10; height = 07; } - # if(Sum1 == 3 & Sum2 != 0){ width = 10; height = 10; } - # if(Sum1 == 0 & Sum2 == 1){ width = 05; height = 05; } - # if(Sum1 == 0 & Sum2 == 2){ width = 10; height = 04; } - # if(Sum1 == 0 & Sum2 == 3){ width = 10; height = 03; } + # if (!isRStudio) { + # if (Sum1 == 1 & Sum2 == 0) { width = 10; height = 05; } + # if (Sum1 == 1 & Sum2 != 0) { width = 10; height = 07; } + # if (Sum1 == 2 & Sum2 == 0) { width = 10; height = 05; } + # if (Sum1 == 2 & Sum2 != 0) { width = 10; height = 07; } + # if (Sum1 == 3 & Sum2 == 0) { width = 10; height = 07; } + # if (Sum1 == 3 & Sum2 != 0) { width = 10; height = 10; } + # if (Sum1 == 0 & Sum2 == 1) { width = 05; height = 05; } + # if (Sum1 == 0 & Sum2 == 2) { width = 10; height = 04; } + # if (Sum1 == 0 & Sum2 == 3) { width = 10; height = 03; } # dev.new(width = width, height = height) # } layout(matlayout); Xaxis <- 1:length(IndPeriod_Plot); - if(BOOL_Dates){ - if(NameTS %in% c("hour", "day", "month")){ - Seq1 <- which(OutputsModel$DatesR[IndPeriod_Plot]$mday == 1 & OutputsModel$DatesR[IndPeriod_Plot]$mon %in% c(0, 3, 6, 9)); - Seq2 <- which(OutputsModel$DatesR[IndPeriod_Plot]$mday == 1 & OutputsModel$DatesR[IndPeriod_Plot]$mon == 0); - Labels2 <- format(OutputsModel$DatesR[IndPeriod_Plot], format = formatAxis)[Seq2]; + if (BOOL_Dates) { + if (NameTS %in% c("hour", "day", "month")) { + Seq1 <- which(OutputsModel$DatesR[IndPeriod_Plot]$mday == 1 & OutputsModel$DatesR[IndPeriod_Plot]$mon %in% c(0, 3, 6, 9)) + Seq2 <- which(OutputsModel$DatesR[IndPeriod_Plot]$mday == 1 & OutputsModel$DatesR[IndPeriod_Plot]$mon == 0) + Labels2 <- format(OutputsModel$DatesR[IndPeriod_Plot], format = formatAxis)[Seq2] } - if(NameTS %in% c("year")){ - Seq1 <- 1:length(OutputsModel$DatesR[IndPeriod_Plot]); - Seq2 <- 1:length(OutputsModel$DatesR[IndPeriod_Plot]); - Labels2 <- format(OutputsModel$DatesR[IndPeriod_Plot], format = formatAxis)[Seq2]; + if (NameTS %in% c("year")) { + Seq1 <- 1:length(OutputsModel$DatesR[IndPeriod_Plot]) + Seq2 <- 1:length(OutputsModel$DatesR[IndPeriod_Plot]) + Labels2 <- format(OutputsModel$DatesR[IndPeriod_Plot], format = formatAxis)[Seq2] } } - if(!is.null(BasinArea)){ - Factor_MMH_M3S <- BasinArea/( 60*60/1000); - Factor_MMD_M3S <- BasinArea/( 24*60*60/1000); - Factor_MMM_M3S <- BasinArea/(365.25/12*24*60*60/1000); - Factor_MMY_M3S <- BasinArea/( 365.25*24*60*60/1000); - if(NameTS == "hour" ){ Factor_UNIT_M3S <- Factor_MMH_M3S; } - if(NameTS == "day" ){ Factor_UNIT_M3S <- Factor_MMD_M3S; } - if(NameTS == "month"){ Factor_UNIT_M3S <- Factor_MMM_M3S; } - if(NameTS == "year" ){ Factor_UNIT_M3S <- Factor_MMY_M3S; } + if (!is.null(BasinArea)) { + # function to use switch avoiding a NOTE when the package is check + UNIT_M3S <- function(x) { + return(switch(x, hour = 60*60, day = 24*60*60, month = 365.25/12*24*60*60, year = 365.25*24*60*60)) + } + Factor_UNIT_M3S <- BasinArea/(UNIT_M3S(NameTS)/1e3) } } @@ -155,78 +157,97 @@ plot.OutputsModel <- function(x, Qobs = NULL, IndPeriod_Plot = NULL, BasinArea = seqDATA2 <- exp(seqDATA1) ##Precip - if(BOOLPLOT_Precip){ - kPlot <- kPlot+1; mar <- c(3, 5, 1, 5); + if (BOOLPLOT_Precip) { + kPlot <- kPlot + 1 + mar <- c(3, 5, 1, 5) + par(new = FALSE, mar = mar, las = 0) - ylim1 <- range(OutputsModel$Precip[IndPeriod_Plot], na.rm = TRUE); ylim2 <- ylim1 * c(1.0, 1.1); ylim2 <- rev(ylim2); - lwdP <- 0.7; if(NameTS %in% c("month", "year")){ lwdP <- 2; } - plot(Xaxis, OutputsModel$Precip[IndPeriod_Plot], type = "h", ylim = ylim2, col = "royalblue", lwd = lwdP, xaxt = "n", yaxt = "n", xlab = "", ylab = "", yaxs = "i", ...); - axis(side = 2, at = pretty(ylim1), labels = pretty(ylim1), cex.axis = cexaxis) - par(las = 0); mtext(side = 2, paste("precip.", plotunit, sep = " "), line = lineY, cex = cexlab, adj = 1); par(las = 0); - if(BOOL_Psol){ - legend("bottomright", c("solid","liquid"), col = c("lightblue", "royalblue"), lty = c(1, 1), lwd = c(lwdLine, lwdLine), bty = "o", bg = bgleg, box.col = bgleg, cex = cexleg) - par(new = TRUE); - plot(Xaxis, PsolLayerMean[IndPeriod_Plot], type = "h", ylim = ylim2, col = "lightblue", lwd = lwdP, xaxt = "n", yaxt = "n", xlab = "", ylab = "", yaxs = "i", ...); + ylim1 <- range(OutputsModel$Precip[IndPeriod_Plot], na.rm = TRUE) + ylim2 <- ylim1 * c(1.0, 1.1) + ylim2 <- rev(ylim2) + + lwdP <- lwd * 0.7 + if (NameTS %in% c("month", "year")) { + lwdP <- lwd * 2 + } + plot(Xaxis, OutputsModel$Precip[IndPeriod_Plot], + type = "h", xaxt = "n", yaxt = "n", yaxs = "i", ylim = ylim2, + col = "royalblue", lwd = lwdP * lwdk, + xlab = "", ylab = "", ...) + axis(side = 2, at = pretty(ylim1), labels = pretty(ylim1), cex.axis = cex.axis, ...) + par(las = 0) + mtext(side = 2, paste0("precip. ", plotunit), cex = cex.lab, adj = 1, line = line) + par(las = 0) + + if (BOOL_Psol) { + legend("bottomright", legend = c("solid","liquid"), + col = c("lightblue", "royalblue"), lty = c(1, 1), lwd = c(lwd, lwd), + bty = "o", bg = bg, box.col = bg, cex = cex.leg) + par(new = TRUE) + plot(Xaxis, PsolLayerMean[IndPeriod_Plot], type = "h", xaxt = "n", yaxt = "n", yaxs = "i", + ylim = ylim2, col = "lightblue", lwd = lwdP * lwdk, xlab = "", ylab = "", ...) + } + if (BOOL_Dates) { + axis(side = 1, at = Seq1, labels = FALSE, cex.axis = cex.axis, ...) + axis(side = 1, at = Seq2, labels = Labels2, lwd.ticks = 1.5, cex.axis = cex.axis, ...) + } else { + axis(side = 1, at = pretty(Xaxis), labels = pretty(Xaxis), cex.axis = cex.axis, ...) } - if(BOOL_Dates){ - axis(side = 1, at = Seq1, labels = FALSE, cex.axis = cexaxis); - axis(side = 1, at = Seq2, labels = Labels2, lwd.ticks = 1.5, cex.axis = cexaxis); - } else { axis(side = 1, at = pretty(Xaxis), labels = pretty(Xaxis), cex.axis = cexaxis); } } ##Temp - if(BOOLPLOT_Temp){ + if (BOOLPLOT_Temp) { kPlot <- kPlot+1; mar <- c(3, 5, 1, 5); par(new = FALSE, mar = mar, las = 0) ylim1 <- c(+99999, -99999) - for(iLayer in 1:NLayers){ + for(iLayer in 1:NLayers) { ylim1[1] <- min(ylim1[1], OutputsModel$CemaNeigeLayers[[iLayer]]$Temp); ylim1[2] <- max(ylim1[2], OutputsModel$CemaNeigeLayers[[iLayer]]$Temp); - if(iLayer == 1){ SnowPackLayerMean <- OutputsModel$CemaNeigeLayers[[iLayer]]$Temp/NLayers; + if (iLayer == 1) { SnowPackLayerMean <- OutputsModel$CemaNeigeLayers[[iLayer]]$Temp/NLayers; } else { SnowPackLayerMean <- SnowPackLayerMean + OutputsModel$CemaNeigeLayers[[iLayer]]$Temp/NLayers; } } plot(SnowPackLayerMean[IndPeriod_Plot], type = "n", ylim = ylim1, xlab = "", ylab = "", xaxt = "n", yaxt = "n", ...) - for(iLayer in 1:NLayers){ lines(OutputsModel$CemaNeigeLayers[[iLayer]]$Temp[IndPeriod_Plot], lty = 3, col = "orchid", lwd = lwdLine*0.8); } + for(iLayer in 1:NLayers) { lines(OutputsModel$CemaNeigeLayers[[iLayer]]$Temp[IndPeriod_Plot], lty = 3, col = "orchid", lwd = lwd * lwdk * 0.8); } abline(h = 0, col = "grey", lty = 2) - lines(SnowPackLayerMean[IndPeriod_Plot], type = "l", lwd = lwdLine*1.0, col = "darkorchid4") - axis(side = 2, at = pretty(ylim1), labels = pretty(ylim1), cex.axis = cexaxis) - par(las = 0); mtext(side = 2, expression(paste("temp. [", degree, "C]", sep = "")), padj = 0.2, line = lineY, cex = cexlab); par(las = 0); - legend("topright", c("mean", "layers"), col = c("darkorchid4", "orchid"), lty = c(1, 3), lwd = c(lwdLine*1.0, lwdLine*0.8), bty = "o", bg = bgleg, box.col = bgleg, cex = cexleg) + lines(SnowPackLayerMean[IndPeriod_Plot], type = "l", lwd = lwd * lwdk *1.0, col = "darkorchid4") + axis(side = 2, at = pretty(ylim1), labels = pretty(ylim1), cex.axis = cex.axis, ...) + par(las = 0); mtext(side = 2, expression(paste0("temp. [", degree, "C]")), padj = 0.2, line = line, cex = cex.lab); par(las = 0); + legend("topright", legend = c("mean", "layers"), col = c("darkorchid4", "orchid"), lty = c(1, 3), lwd = c(lwd*1.0, lwd*0.8), bty = "o", bg = bg, box.col = bg, cex = cex.leg) box() - if(BOOL_Dates){ - axis(side = 1, at = Seq1, labels = FALSE, cex.axis = cexaxis); - axis(side = 1, at = Seq2, labels = Labels2, lwd.ticks = 1.5, cex.axis = cexaxis); - } else { axis(side = 1, at = pretty(Xaxis), labels = pretty(Xaxis), cex.axis = cexaxis); } + if (BOOL_Dates) { + axis(side = 1, at = Seq1, labels = FALSE, cex.axis = cex.axis, ...); + axis(side = 1, at = Seq2, labels = Labels2, lwd.ticks = 1.5, cex.axis = cex.axis, ...); + } else { axis(side = 1, at = pretty(Xaxis), labels = pretty(Xaxis), cex.axis = cex.axis, ...); } } ##SnowPack - if(BOOLPLOT_SnowPack){ + if (BOOLPLOT_SnowPack) { kPlot <- kPlot+1; mar <- c(3, 5, 1, 5); par(new = FALSE, mar = mar, las = 0) ylim1 <- c(+99999, -99999) - for(iLayer in 1:NLayers){ + for(iLayer in 1:NLayers) { ylim1[1] <- min(ylim1[1], OutputsModel$CemaNeigeLayers[[iLayer]]$SnowPack); ylim1[2] <- max(ylim1[2], OutputsModel$CemaNeigeLayers[[iLayer]]$SnowPack); - if(iLayer == 1){ SnowPackLayerMean <- OutputsModel$CemaNeigeLayers[[iLayer]]$SnowPack/NLayers; + if (iLayer == 1) { SnowPackLayerMean <- OutputsModel$CemaNeigeLayers[[iLayer]]$SnowPack/NLayers; } else { SnowPackLayerMean <- SnowPackLayerMean + OutputsModel$CemaNeigeLayers[[iLayer]]$SnowPack/NLayers; } } - plot(SnowPackLayerMean[IndPeriod_Plot], type = "l", ylim = ylim1, lwd = lwdLine*1.2, col = "royalblue", xlab = "", ylab = "", xaxt = "n", yaxt = "n", ...) - for(iLayer in 1:NLayers){ lines(OutputsModel$CemaNeigeLayers[[iLayer]]$SnowPack[IndPeriod_Plot], lty = 3, col = "royalblue", lwd = lwdLine*0.8); } - axis(side = 2, at = pretty(ylim1), labels = pretty(ylim1), cex.axis = cexaxis) - par(las = 0); mtext(side = 2, paste("snow pack", "[mm]", sep = " "), line = lineY, cex = cexlab); par(las = 0); - legend("topright", c("mean", "layers"), col = c("royalblue", "royalblue"), lty = c(1, 3), lwd = c(lwdLine*1.2, lwdLine*0.8), bty = "o", bg = bgleg, box.col = bgleg, cex = cexleg) + plot(SnowPackLayerMean[IndPeriod_Plot], type = "l", ylim = ylim1, lwd = lwd * lwdk *1.2, col = "royalblue", xlab = "", ylab = "", xaxt = "n", yaxt = "n", ...) + for(iLayer in 1:NLayers) { lines(OutputsModel$CemaNeigeLayers[[iLayer]]$SnowPack[IndPeriod_Plot], lty = 3, col = "royalblue", lwd = lwd * lwdk *0.8); } + axis(side = 2, at = pretty(ylim1), labels = pretty(ylim1), cex.axis = cex.axis, ...) + par(las = 0); mtext(side = 2, paste0("snow pack ", "[mm]"), line = line, cex = cex.lab); par(las = 0); + legend("topright", legend = c("mean", "layers"), col = c("royalblue", "royalblue"), lty = c(1, 3), lwd = c(lwd*1.2, lwd*0.8), bty = "o", bg = bg, box.col = bg, cex = cex.leg) box() - if(BOOL_Dates){ - axis(side = 1, at = Seq1, labels = FALSE, cex.axis = cexaxis); - axis(side = 1, at = Seq2, labels = Labels2, lwd.ticks = 1.5, cex.axis = cexaxis); - } else { axis(side = 1, at = pretty(Xaxis), labels = pretty(Xaxis), cex.axis = cexaxis); } + if (BOOL_Dates) { + axis(side = 1, at = Seq1, labels = FALSE, cex.axis = cex.axis, ...); + axis(side = 1, at = Seq2, labels = Labels2, lwd.ticks = 1.5, cex.axis = cex.axis, ...); + } else { axis(side = 1, at = pretty(Xaxis), labels = pretty(Xaxis), cex.axis = cex.axis, ...); } } ##Flows - if(BOOLPLOT_Flows & log_scale) { + if (BOOLPLOT_Flows & log_scale) { kPlot <- kPlot+1; mar <- c(3, 5, 1, 5); par(new = FALSE, mar = mar, las = 0) @@ -239,75 +260,75 @@ plot.OutputsModel <- function(x, Qobs = NULL, IndPeriod_Plot = NULL, BasinArea = DATA3 <- log(DATA3) ylim1 <- range(DATA3[IndPeriod_Plot], na.rm = TRUE); - if(BOOL_Qobs){ ylim1 <- range(c(ylim1, DATA2[IndPeriod_Plot]), na.rm = TRUE); } + if (BOOL_Qobs) { ylim1 <- range(c(ylim1, DATA2[IndPeriod_Plot]), na.rm = TRUE); } ylim2 <- c(ylim1[1], 1.2*ylim1[2]); plot(Xaxis, rep(NA, length(Xaxis)), type = "n", ylim = ylim2, xlab = "", ylab = "", xaxt = "n", yaxt = "n", ...); txtleg <- NULL; colleg <- NULL; - if(BOOL_Qobs){ lines(Xaxis, DATA2[IndPeriod_Plot], lwd = lwdLine, lty = 1, col = par("fg")); txtleg <- c(txtleg, "observed"); colleg <- c(colleg, par("fg")); } - if(BOOL_Qsim){ lines(Xaxis, DATA3[IndPeriod_Plot], lwd = lwdLine, lty = 1, col = "orangered"); txtleg <- c(txtleg, "simulated"); colleg <- c(colleg, "orangered"); } - axis(side = 2, at = seqDATA1, labels = seqDATA2, cex.axis = cexaxis) - par(las = 0); mtext(side = 2, paste("flow", plotunit, sep = " "), line = lineY, cex = cexlab); par(las = 0); - if(!is.null(BasinArea)){ + if (BOOL_Qobs) { lines(Xaxis, DATA2[IndPeriod_Plot], lwd = lwd * lwdk , lty = 1, col = par("fg")); txtleg <- c(txtleg, "observed"); colleg <- c(colleg, par("fg")); } + if (BOOL_Qsim) { lines(Xaxis, DATA3[IndPeriod_Plot], lwd = lwd * lwdk , lty = 1, col = "orangered"); txtleg <- c(txtleg, "simulated"); colleg <- c(colleg, "orangered"); } + axis(side = 2, at = seqDATA1, labels = seqDATA2, cex.axis = cex.axis, ...) + par(las = 0); mtext(side = 2, paste0("flow ", plotunit), line = line, cex = cex.lab); par(las = 0); + if (!is.null(BasinArea)) { Factor <- Factor_UNIT_M3S; - axis(side = 4, at = pretty(ylim1*Factor)/Factor, labels = pretty(ylim1*Factor), cex.axis = cexaxis); - par(las = 0); mtext(side = 4, paste("flow", "m3/s", sep = " "), line = lineY, cex = cexlab); par(las = 0); } - if(BOOL_Dates){ - axis(side = 1, at = Seq1, labels = FALSE, cex.axis = cexaxis); - axis(side = 1, at = Seq2, labels = Labels2, lwd.ticks = 1.5, cex.axis = cexaxis); - } else { axis(side = 1, at = pretty(Xaxis), labels = pretty(Xaxis), cex.axis = cexaxis); } - legend("topright", txtleg, col = colleg, lty = 1, lwd = lwdLine, bty = "o", bg = bgleg, box.col = bgleg, cex = cexleg) - legend("bottomright", "log scale", lty = 1, col = NA, bty = "o", bg = bgleg, box.col = bgleg, cex = cexleg) + axis(side = 4, at = pretty(ylim1*Factor)/Factor, labels = pretty(ylim1*Factor), cex.axis = cex.axis, ...); + par(las = 0); mtext(side = 4, paste0("flow ", "m3/s"), line = line, cex = cex.lab); par(las = 0); } + if (BOOL_Dates) { + axis(side = 1, at = Seq1, labels = FALSE, cex.axis = cex.axis, ...); + axis(side = 1, at = Seq2, labels = Labels2, lwd.ticks = 1.5, cex.axis = cex.axis, ...); + } else { axis(side = 1, at = pretty(Xaxis), labels = pretty(Xaxis), cex.axis = cex.axis, ...); } + legend("topright", txtleg, col = colleg, lty = 1, lwd = lwd * lwdk , bty = "o", bg = bg, box.col = bg, cex = cex.leg) + legend("bottomright", "log scale", lty = 1, col = NA, bty = "o", bg = bg, box.col = bg, cex = cex.leg) box() } - if(BOOLPLOT_Flows & !log_scale){ + if (BOOLPLOT_Flows & !log_scale) { kPlot <- kPlot+1; mar <- c(3, 5, 1, 5); par(new = FALSE, mar = mar, las = 0) ylim1 <- range(OutputsModel$Qsim[IndPeriod_Plot], na.rm = TRUE); - if(BOOL_Qobs){ ylim1 <- range(c(ylim1, Qobs[IndPeriod_Plot]), na.rm = TRUE); } + if (BOOL_Qobs) { ylim1 <- range(c(ylim1, Qobs[IndPeriod_Plot]), na.rm = TRUE); } ylim2 <- c(ylim1[1], 1.2*ylim1[2]); plot(Xaxis, rep(NA, length(Xaxis)), type = "n", ylim = ylim2, xlab = "", ylab = "", xaxt = "n", yaxt = "n", ...); txtleg <- NULL; colleg <- NULL; - if(BOOL_Qobs){ lines(Xaxis, Qobs[IndPeriod_Plot], lwd = lwdLine, lty = 1, col = par("fg")); txtleg <- c(txtleg, "observed"); colleg <- c(colleg, par("fg")); } - if(BOOL_Qsim){ lines(Xaxis, OutputsModel$Qsim[IndPeriod_Plot], lwd = lwdLine, lty = 1, col = "orangered"); txtleg <- c(txtleg, "simulated"); colleg <- c(colleg, "orangered"); } - axis(side = 2, at = pretty(ylim1), labels = pretty(ylim1), cex.axis = cexaxis) - par(las = 0); mtext(side = 2, paste("flow", plotunit, sep = " "), line = lineY, cex = cexlab); par(las = 0); - if(!is.null(BasinArea)){ + if (BOOL_Qobs) { lines(Xaxis, Qobs[IndPeriod_Plot], lwd = lwd * lwdk , lty = 1, col = par("fg")); txtleg <- c(txtleg, "observed"); colleg <- c(colleg, par("fg")); } + if (BOOL_Qsim) { lines(Xaxis, OutputsModel$Qsim[IndPeriod_Plot], lwd = lwd * lwdk , lty = 1, col = "orangered"); txtleg <- c(txtleg, "simulated"); colleg <- c(colleg, "orangered"); } + axis(side = 2, at = pretty(ylim1), labels = pretty(ylim1), cex.axis = cex.axis, ...) + par(las = 0); mtext(side = 2, paste("flow", plotunit, sep = " "), line = line, cex = cex.lab); par(las = 0); + if (!is.null(BasinArea)) { Factor <- Factor_UNIT_M3S; - axis(side = 4, at = pretty(ylim1*Factor)/Factor, labels = pretty(ylim1*Factor), cex.axis = cexaxis); - par(las = 0); mtext(side = 4, paste("flow", "m3/s", sep = " "), line = lineY, cex = cexlab); par(las = 0); } - if(BOOL_Dates){ - axis(side = 1, at = Seq1, labels = FALSE, cex.axis = cexaxis); - axis(side = 1, at = Seq2, labels = Labels2, lwd.ticks = 1.5, cex.axis = cexaxis); - } else { axis(side = 1, at = pretty(Xaxis), labels = pretty(Xaxis), cex.axis = cexaxis); } - legend("topright", txtleg, col = colleg, lty = 1, lwd = lwdLine, bty = "o", bg = bgleg, box.col = bgleg, cex = cexleg) + axis(side = 4, at = pretty(ylim1*Factor)/Factor, labels = pretty(ylim1*Factor), cex.axis = cex.axis, ...); + par(las = 0); mtext(side = 4, paste("flow", "m3/s", sep = " "), line = line, cex = cex.lab); par(las = 0); } + if (BOOL_Dates) { + axis(side = 1, at = Seq1, labels = FALSE, cex.axis = cex.axis, ...); + axis(side = 1, at = Seq2, labels = Labels2, lwd.ticks = 1.5, cex.axis = cex.axis, ...); + } else { axis(side = 1, at = pretty(Xaxis), labels = pretty(Xaxis), cex.axis = cex.axis, ...); } + legend("topright", txtleg, col = colleg, lty = 1, lwd = lwd * lwdk , bty = "o", bg = bg, box.col = bg, cex = cex.leg) box() } ##Regime - if(BOOLPLOT_Regime){ + if (BOOLPLOT_Regime) { kPlot <- kPlot+1; mar <- c(6, 5, 1, 5); plotunitregime <- "[mm/month]"; par(new = FALSE, mar = mar, las = 0) ##Data_formating_as_table DataModel <- as.data.frame(matrix(as.numeric(NA), nrow = length(IndPeriod_Plot), ncol = 5)); DataModel[, 1] <- as.numeric(format(OutputsModel$DatesR[IndPeriod_Plot], format = "%Y%m%d%H")); - if(BOOL_Pobs){ DataModel[, 2] <- OutputsModel$Precip[IndPeriod_Plot]; } - if(BOOL_Psol){ DataModel[, 3] <- PsolLayerMean[IndPeriod_Plot]; } - if(BOOL_Qobs){ DataModel[, 4] <- Qobs[IndPeriod_Plot]; } - if(BOOL_Qsim){ DataModel[, 5] <- OutputsModel$Qsim[IndPeriod_Plot]; } + if (BOOL_Pobs) { DataModel[, 2] <- OutputsModel$Precip[IndPeriod_Plot]; } + if (BOOL_Psol) { DataModel[, 3] <- PsolLayerMean[IndPeriod_Plot]; } + if (BOOL_Qobs) { DataModel[, 4] <- Qobs[IndPeriod_Plot]; } + if (BOOL_Qsim) { DataModel[, 5] <- OutputsModel$Qsim[IndPeriod_Plot]; } colnames(DataModel) <- c("Dates", "Precip", "Psol", "Qobs", "Qsim"); TxtDatesDataModel <- formatC(DataModel$Dates, format = "d", width = 8, flag = "0"); ##Building_of_daily_time_series_if_needed - if(NameTS == "month"){ DataDaily <- NULL; } - if(NameTS == "day" ){ DataDaily <- DataModel; } - if(NameTS == "hour" ){ DataDaily <- as.data.frame(aggregate(DataModel[, 2:5], by = list(as.numeric(substr(TxtDatesDataModel, 1, 8))), FUN = sum, na.rm = T)); } - if(NameTS %in% c("hour", "day")){ + if (NameTS == "month") { DataDaily <- NULL; } + if (NameTS == "day" ) { DataDaily <- DataModel; } + if (NameTS == "hour" ) { DataDaily <- as.data.frame(aggregate(DataModel[, 2:5], by = list(as.numeric(substr(TxtDatesDataModel, 1, 8))), FUN = sum, na.rm = T)); } + if (NameTS %in% c("hour", "day")) { colnames(DataDaily) <- c("Dates", "Precip", "Psol", "Qobs", "Qsim"); TxtDatesDataDaily <- formatC(DataDaily$Dates, format = "d", width = 8, flag = "0"); } ##Building_of_monthly_time_series_if_needed - if(NameTS == "month"){ DataMonthly <- DataModel; } - if(NameTS == "day" ){ DataMonthly <- as.data.frame(aggregate(DataDaily[, 2:5], by = list(as.numeric(substr(TxtDatesDataDaily, 1, 6))), FUN = sum, na.rm = T)); } - if(NameTS == "hour" ){ DataMonthly <- as.data.frame(aggregate(DataDaily[, 2:5], by = list(as.numeric(substr(TxtDatesDataDaily, 1, 6))), FUN = sum, na.rm = T)); } + if (NameTS == "month") { DataMonthly <- DataModel; } + if (NameTS == "day" ) { DataMonthly <- as.data.frame(aggregate(DataDaily[, 2:5], by = list(as.numeric(substr(TxtDatesDataDaily, 1, 6))), FUN = sum, na.rm = T)); } + if (NameTS == "hour" ) { DataMonthly <- as.data.frame(aggregate(DataDaily[, 2:5], by = list(as.numeric(substr(TxtDatesDataDaily, 1, 6))), FUN = sum, na.rm = T)); } colnames(DataMonthly) <- c("Dates", "Precip", "Psol", "Qobs", "Qsim"); TxtDatesDataMonthly <- formatC(DataMonthly$Dates, format = "d", width = 6, flag = "0"); ##Computation_of_interannual_mean_series @@ -319,14 +340,14 @@ plot.OutputsModel <- function(x, Qobs = NULL, IndPeriod_Plot = NULL, BasinArea = colnames(DataDailyInterAn) <- c("Dates", "Precip", "Psol", "Qobs", "Qsim") DataDailyInterAn <- merge(SeqY, DataDailyInterAn, by = "Dates", all.x = TRUE, all.y = FALSE) } - if(!is.null(DataMonthly)){ + if (!is.null(DataMonthly)) { SeqM <- data.frame(Dates = 1:12) DataMonthlyInterAn <- as.data.frame(aggregate(DataMonthly[, 2:5], by = list(as.numeric(substr(TxtDatesDataMonthly, 5, 6))), FUN = mean, na.rm = T)); colnames(DataMonthlyInterAn) <- c("Dates", "Precip", "Psol", "Qobs", "Qsim") DataMonthlyInterAn <- merge(SeqM, DataMonthlyInterAn, by = "Dates", all.x = TRUE, all.y = FALSE) } ##Smoothing_of_daily_series_and_scale_conversion_to_make_them_become_a_monthly_regime - if(!is.null(DataDaily)){ + if (!is.null(DataDaily)) { ##Smoothing NDaysWindow <- 30; DataDailyInterAn <- as.data.frame(cbind(DataDailyInterAn$Dates, @@ -334,14 +355,14 @@ plot.OutputsModel <- function(x, Qobs = NULL, IndPeriod_Plot = NULL, BasinArea = MyRollMean3(DataDailyInterAn$Qobs , NDaysWindow), MyRollMean3(DataDailyInterAn$Qsim, NDaysWindow))); colnames(DataDailyInterAn) <- c("Dates", "Precip", "Psol", "Qobs", "Qsim"); ##Scale_conversion_to_make_them_become_a_monthly_regime - if(plotunitregime != "[mm/month]"){ stop(paste("incorrect unit for regime plot \n", sep = "")); return(NULL); } + if (plotunitregime != "[mm/month]") { stop(paste("incorrect unit for regime plot \n", sep = "")); return(NULL); } DataDailyInterAn <- as.data.frame(cbind(DataDailyInterAn[1], DataDailyInterAn[2:5]*30)); } ##Plot_preparation DataPlotP <- DataMonthlyInterAn; - if(!is.null(DataDaily)){ + if (!is.null(DataDaily)) { DataPlotQ <- DataDailyInterAn; - SeqX1 <- c( 1, 32, 61, 92, 122, 153, 183, 214, 245, 275, 306, 336, 366); + SeqX1 <- c( 1, 32, 61, 92, 122, 153, 183, 214, 245, 275, 306, 336, 366); SeqX2 <- c( 15, 46, 75, 106, 136, 167, 197, 228, 259, 289, 320, 350); labX <- "30-days rolling mean"; } else { @@ -353,103 +374,109 @@ plot.OutputsModel <- function(x, Qobs = NULL, IndPeriod_Plot = NULL, BasinArea = xLabels1 <- rep("", 13); xLabels2 <- month.abb ylimQ <- range(c(DataPlotQ$Qobs, DataPlotQ$Qsim), na.rm = TRUE); - if(BOOL_Pobs){ ylimP <- c(max(DataPlotP$Precip, na.rm = TRUE), 0); } + if (BOOL_Pobs) { ylimP <- c(max(DataPlotP$Precip, na.rm = TRUE), 0); } txtleg <- NULL; colleg <- NULL; lwdleg <- NULL; lwdP = 10; ##Plot_forcings - if(BOOL_Pobs){ - plot(SeqX2[DataMonthlyInterAn$Dates], DataPlotP$Precip, type = "h", xlim = range(SeqX1), ylim = c(3*ylimP[1], ylimP[2]), lwd = lwdP, lend = 1, lty = 1, col = "royalblue", xlab = "", ylab = "", xaxt = "n", yaxt = "n", yaxs = "i", bty = "n", ...) + if (BOOL_Pobs) { + plot(SeqX2[DataMonthlyInterAn$Dates], DataPlotP$Precip, type = "h", xlim = range(SeqX1), ylim = c(3*ylimP[1], ylimP[2]), lwd = lwdk * lwdP, lend = 1, lty = 1, col = "royalblue", xlab = "", ylab = "", xaxt = "n", yaxt = "n", yaxs = "i", bty = "n", ...) txtleg <- c(txtleg, "Ptot" ); colleg <- c(colleg, "royalblue"); lwdleg <- c(lwdleg, lwdP/3); - axis(side = 2, at = pretty(0.8*ylimP, n = 3), labels = pretty(0.8*ylimP, n = 3), cex.axis = cexaxis, col.axis = "royalblue", col.ticks = "royalblue"); + axis(side = 2, at = pretty(0.8*ylimP, n = 3), labels = pretty(0.8*ylimP, n = 3), col.axis = "royalblue", col.ticks = "royalblue", cex.axis = cex.axis, ...); par(new = TRUE); } - if(BOOL_Psol){ - plot(SeqX2, DataPlotP$Psol[DataMonthlyInterAn$Dates], type = "h", xlim = range(SeqX1), ylim = c(3*ylimP[1], ylimP[2]), lwd = lwdP, lend = 1, lty = 1, col = "lightblue", xlab = "", ylab = "", xaxt = "n", yaxt = "n", yaxs = "i", bty = "n", ...); + if (BOOL_Psol) { + plot(SeqX2, DataPlotP$Psol[DataMonthlyInterAn$Dates], type = "h", xlim = range(SeqX1), ylim = c(3*ylimP[1], ylimP[2]), lwd = lwdk * lwdP, lend = 1, lty = 1, col = "lightblue", xlab = "", ylab = "", xaxt = "n", yaxt = "n", yaxs = "i", bty = "n", ...); txtleg <- c(txtleg, "Psol" ); colleg <- c(colleg, "lightblue"); lwdleg <- c(lwdleg, lwdP/3); par(new = TRUE); } ##Plot_flows plot(NULL, type = "n", xlim = range(SeqX1), ylim = c(ylimQ[1], 2*ylimQ[2]), xlab = "", ylab = "", xaxt = "n", yaxt = "n", ...) - if(BOOL_Qobs){ lines(1:nrow(DataPlotQ), DataPlotQ$Qobs, lwd = lwdLine, lty = 1, col = par("fg") ); txtleg <- c(txtleg, "Qobs" ); colleg <- c(colleg, par("fg") ); lwdleg <- c(lwdleg, lwdLine); } - if(BOOL_Qsim){ lines(1:nrow(DataPlotQ), DataPlotQ$Qsim, lwd = lwdLine, lty = 1, col = "orangered"); txtleg <- c(txtleg, "Qsim"); colleg <- c(colleg, "orangered"); lwdleg <- c(lwdleg, lwdLine); } + if (BOOL_Qobs) { lines(1:nrow(DataPlotQ), DataPlotQ$Qobs, lwd = lwd * lwdk , lty = 1, col = par("fg") ); txtleg <- c(txtleg, "Qobs" ); colleg <- c(colleg, par("fg") ); lwdleg <- c(lwdleg, lwd); } + if (BOOL_Qsim) { lines(1:nrow(DataPlotQ), DataPlotQ$Qsim, lwd = lwd * lwdk , lty = 1, col = "orangered"); txtleg <- c(txtleg, "Qsim"); colleg <- c(colleg, "orangered"); lwdleg <- c(lwdleg, lwd); } ##Axis_and_legend - axis(side = 1, at = SeqX1, tick = TRUE , labels = xLabels1, cex.axis = cexaxis) - axis(side = 1, at = SeqX2, tick = FALSE, labels = xLabels2, cex.axis = cexaxis) - axis(side = 2, at = pretty(ylimQ), labels = pretty(ylimQ), cex.axis = cexaxis) - par(las = 0); mtext(side = 1, labX, line = lineX, cex = cexlab); par(las = 0); + axis(side = 1, at = SeqX1, tick = TRUE , labels = xLabels1, cex.axis = cex.axis, ...) + axis(side = 1, at = SeqX2, tick = FALSE, labels = xLabels2, cex.axis = cex.axis, ...) + axis(side = 2, at = pretty(ylimQ), labels = pretty(ylimQ), cex.axis = cex.axis, ...) + par(las = 0); mtext(side = 1, labX, line = line, cex = cex.lab); par(las = 0); posleg <- "topright"; txtlab <- "flow regime"; - if(BOOL_Pobs){ posleg <- "right"; txtlab <- "precip. & flow regime"; } - par(las = 0); mtext(side = 2, paste(txtlab, " ", plotunitregime, sep = ""), line = lineY, cex = cexlab); par(las = 0); - if(!is.null(BasinArea)){ + if (BOOL_Pobs) { posleg <- "right"; txtlab <- "precip. & flow regime"; } + par(las = 0); mtext(side = 2, paste0(txtlab, " ", plotunitregime), line = line, cex = cex.lab); par(las = 0); + if (!is.null(BasinArea)) { Factor <- Factor_MMM_M3S; - axis(side = 4, at = pretty(ylimQ*Factor)/Factor, labels = pretty(ylimQ*Factor), cex.axis = cexaxis); - par(las = 0); mtext(side = 4, paste("flow ", "m3/s", sep = ""), line = lineY, cex = cexlab); par(las = 0); } - ### posleg <- "topright"; if(BOOL_Pobs){ posleg <- "right"; } - ### legend(posleg, txtleg, col = colleg, lty = 1, lwd = lwdleg, bty = "o", bg = bgleg, box.col = bgleg, cex = cexleg) + axis(side = 4, at = pretty(ylimQ*Factor)/Factor, labels = pretty(ylimQ*Factor), cex.axis = cex.axis, ...); + par(las = 0); mtext(side = 4, paste0("flow ", "m3/s"), line = line, cex = cex.lab); par(las = 0); } + ### posleg <- "topright"; if (BOOL_Pobs) { posleg <- "right"; } + ### legend(posleg, txtleg, col = colleg, lty = 1, lwd = lwdleg, bty = "o", bg = bg, box.col = bg, cex = cex.leg) box() } ##Cumulative_frequency - if(BOOLPLOT_CumFreq){ + if (BOOLPLOT_CumFreq) { kPlot <- kPlot+1; mar <- c(6, 5, 1, 5); par(new = FALSE, mar = mar, las = 0) xlim <- c(0, 1); - if(BOOL_Qobs & !BOOL_Qsim){ SelectNotZero <- SelectQobsNotZero; + if (BOOL_Qobs & !BOOL_Qsim) { SelectNotZero <- SelectQobsNotZero; ylim <- range(log(Qobs[IndPeriod_Plot][SelectNotZero]), na.rm = TRUE); } - if(BOOL_Qsim & !BOOL_Qobs){ SelectNotZero <- SelectQsimNotZero; + if (BOOL_Qsim & !BOOL_Qobs) { SelectNotZero <- SelectQsimNotZero; ylim <- range(log(OutputsModel$Qsim[IndPeriod_Plot][SelectNotZero]), na.rm = TRUE); } - if(BOOL_Qobs & BOOL_Qsim){ SelectNotZero <- SelectQobsNotZero & SelectQsimNotZero; + if (BOOL_Qobs & BOOL_Qsim) { SelectNotZero <- SelectQobsNotZero & SelectQsimNotZero; ylim <- range(log(c(Qobs[IndPeriod_Plot][SelectNotZero], OutputsModel$Qsim[IndPeriod_Plot][SelectNotZero])), na.rm = TRUE); } - plot(0, 0, type = "n", xlim = xlim, ylim = ylim, xaxt = "n", yaxt = "n", xlab = "", ylab = "", main = "", ...); + plot(0, 0, type = "n", + xlim = xlim, ylim = ylim, + xaxt = "n", yaxt = "n", + xlab = "", ylab = "", main = "", ...); ### abline(h = 0, lty = 2, col = grey(0.5)); ### abline(h = 1, lty = 2, col = grey(0.5)); - axis(side = 1, at = pretty(xlim), labels = pretty(xlim), cex.axis = cexaxis); - par(las = 0); mtext(side = 1, text = "non-exceedance prob. [-]", line = lineY, cex = cexlab); par(las = 0); - axis(side = 2, at = seqDATA1, labels = seqDATA2, cex.axis = cexaxis) - par(las = 0); mtext(side = 2, text = paste("flow ", plotunit, "", sep = ""), line = lineY, cex = cexlab); par(las = 0); + axis(side = 1, at = pretty(xlim), labels = pretty(xlim), cex.axis = cex.axis, ...); + par(las = 0); mtext(side = 1, text = "non-exceedance prob. [-]", line = line, cex = cex.lab); par(las = 0); + axis(side = 2, at = seqDATA1, labels = seqDATA2, cex.axis = cex.axis, ...) + par(las = 0); mtext(side = 2, text = paste("flow ", plotunit, "", sep = ""), line = line, cex = cex.lab); par(las = 0); txtleg <- NULL; colleg <- NULL; - if(BOOL_Qobs){ + if (BOOL_Qobs) { DATA2 <- log(Qobs[IndPeriod_Plot][SelectNotZero]); SeqQuant <- seq(0, 1, by = 1/(length(DATA2))); Quant <- as.numeric(quantile(DATA2, SeqQuant, na.rm = TRUE)); Fn <- ecdf(DATA2); YY <- DATA2; YY <- YY[order( Fn(DATA2) )]; XX <- Fn(DATA2); XX <- XX[order( Fn(DATA2) )]; - lines(XX, YY, lwd = 1, col = par("fg")); + lines(XX, YY, lwd = lwd, col = par("fg")); txtleg <- c(txtleg, "observed"); colleg <- c(colleg, par("fg")); } - if(BOOL_Qsim){ + if (BOOL_Qsim) { DATA2 <- log(OutputsModel$Qsim[IndPeriod_Plot][SelectNotZero]); SeqQuant <- seq(0, 1, by = 1/(length(DATA2))); Quant <- as.numeric(quantile(DATA2, SeqQuant, na.rm = TRUE)); Fn <- ecdf(DATA2); YY <- DATA2; YY <- YY[order( Fn(DATA2) )]; XX <- Fn(DATA2); XX <- XX[order( Fn(DATA2) )]; - lines(XX, YY, lwd = 1, col = "orangered"); + lines(XX, YY, lwd = lwd, col = "orangered"); txtleg <- c(txtleg, "simulated"); colleg <- c(colleg, "orangered"); } - if(!is.null(BasinArea)){ + if (!is.null(BasinArea)) { Factor <- Factor_UNIT_M3S; - axis(side = 4, at = seqDATA1, labels = round(seqDATA2*Factor), cex.axis = cexaxis) - par(las = 0); mtext(side = 4, paste("flow ", "m3/s", sep = ""), line = lineY, cex = cexlab); par(las = 0); } - legend("topleft", txtleg, col = colleg, lty = 1, lwd = lwdLine, bty = "o", bg = bgleg, box.col = bgleg, cex = cexleg) - legend("bottomright", "log scale", lty = 1, col = NA, bty = "o", bg = bgleg, box.col = bgleg, cex = cexleg) + axis(side = 4, at = seqDATA1, labels = round(seqDATA2*Factor), cex.axis = cex.axis, ...) + par(las = 0); mtext(side = 4, paste0("flow ", "m3/s"), line = line, cex = cex.lab); par(las = 0); } + legend("topleft", txtleg, col = colleg, lty = 1, lwd = lwd, bty = "o", bg = bg, box.col = bg, cex = cex.leg) + legend("bottomright", "log scale", lty = 1, col = NA, bty = "o", bg = bg, box.col = bg, cex = cex.leg) box() } ##Correlation_QQ - if(BOOLPLOT_CorQQ){ + if (BOOLPLOT_CorQQ) { kPlot <- kPlot+1; mar <- c(6, 5, 1, 5); par(new = FALSE, mar = mar, las = 0) ylim <- log(range(c(Qobs[IndPeriod_Plot][SelectQobsNotZero & SelectQsimNotZero], OutputsModel$Qsim[IndPeriod_Plot][SelectQobsNotZero & SelectQsimNotZero]), na.rm = TRUE)); - plot(log(Qobs[IndPeriod_Plot][SelectQobsNotZero & SelectQsimNotZero]), log(OutputsModel$Qsim[IndPeriod_Plot][SelectQobsNotZero & SelectQsimNotZero]), type = "p", pch = 1, cex = 0.9, col = par("fg"), xlim = ylim, ylim = ylim, xaxt = "n", yaxt = "n", xlab = "", ylab = "", ...) - abline(a = 0, b = 1, col = "royalblue"); - axis(side = 1, at = seqDATA1, labels = seqDATA2, cex = cexaxis); - axis(side = 2, at = seqDATA1, labels = seqDATA2, cex = cexaxis); - par(las = 0); mtext(side = 1, paste("observed flow ", plotunit, "", sep = ""), line = lineX, cex = cexlab); par(las = 0); - par(las = 0); mtext(side = 2, paste("simulated flow ", plotunit, "", sep = ""), line = lineY, cex = cexlab); par(las = 0); - if(!is.null(BasinArea)){ + plot(log(Qobs[IndPeriod_Plot][SelectQobsNotZero & SelectQsimNotZero]), + log(OutputsModel$Qsim[IndPeriod_Plot][SelectQobsNotZero & SelectQsimNotZero]), + type = "p", pch = 1, cex = 0.9, col = par("fg"), lwd = lwd, + xlim = ylim, ylim = ylim, xaxt = "n", yaxt = "n", xlab = "", ylab = "", ...) + abline(a = 0, b = 1, col = "royalblue", lwd = lwd); + axis(side = 1, at = seqDATA1, labels = seqDATA2, cex = cex.leg, cex.axis = cex.axis, ...); + axis(side = 2, at = seqDATA1, labels = seqDATA2, cex = cex.leg, cex.axis = cex.axis, ...); + par(las = 0); mtext(side = 1, paste0("observed flow ", plotunit), line = line, cex = cex.lab); par(las = 0); + par(las = 0); mtext(side = 2, paste0("simulated flow ", plotunit), line = line, cex = cex.lab); par(las = 0); + if (!is.null(BasinArea)) { Factor <- Factor_UNIT_M3S; - axis(side = 4, at = seqDATA1, labels = round(seqDATA2*Factor), cex.axis = cexaxis); - par(las = 0); mtext(side = 4, paste("flow ", "m3/s", sep = ""), line = lineY, cex = cexlab); par(las = 0); } - legend("bottomright", "log scale", lty = 1, col = NA, bty = "o", bg = bgleg, box.col = bgleg, cex = cexleg) + axis(side = 4, at = seqDATA1, labels = round(seqDATA2*Factor), cex.axis = cex.axis, ...); + par(las = 0); mtext(side = 4, paste0("flow ", "m3/s"), line = line, cex = cex.lab); par(las = 0); } + legend("bottomright", "log scale", lty = 1, col = NA, bty = "o", bg = bg, box.col = bg, cex = cex.leg) box() } ##Empty_plots - while(kPlot < iPlotMax){ + while (kPlot < iPlotMax) { kPlot <- kPlot+1; par(new = FALSE) plot(0, 0, type = "n", xlab = "", ylab = "", axes = FALSE, ...) diff --git a/man/plot.OutputsModel.Rd b/man/plot.OutputsModel.Rd index 036b1119..f1a805f6 100644 --- a/man/plot.OutputsModel.Rd +++ b/man/plot.OutputsModel.Rd @@ -5,7 +5,8 @@ \title{Default preview of model outputs} \usage{ \method{plot}{OutputsModel}(x, Qobs = NULL, IndPeriod_Plot = NULL, - BasinArea = NULL, which = "all", log_scale = FALSE, verbose = TRUE, ...) + BasinArea = NULL, which = "all", log_scale = FALSE, + cex.axis = 1, cex.lab = 0.9, cex.leg = 0.9, lwd = 1, verbose = TRUE, ...) } \arguments{ \item{x}{[object of class \emph{OutputsModel}] list of model outputs (which must at least include DatesR, Precip and Qsim) [POSIXlt, mm, mm]} @@ -18,9 +19,17 @@ \item{which}{(optional) [character] choice of plots \cr (e.g. c(\code{"Precip"}, \code{"Temp"}, \code{"SnowPack"}, \code{"Flows"}, \code{"Regime"}, \code{"CumFreq"}, \code{"CorQQ"})), default = \code{"all"}} -\item{log_scale}{(optional) [boolean] boolean indicating if the flow axis is to be logarithmic, default = \code{FALSE}} +\item{log_scale}{(optional) [boolean] indicating if the flow axis is to be logarithmic, default = \code{FALSE}} -\item{verbose}{(optional) [boolean] boolean indicating if the function is run in verbose mode or not, default = \code{TRUE}} +\item{cex.axis}{(optional) [numeric] the magnification to be used for axis annotation relative to the current setting of \code{cex}} + +\item{cex.lab}{(optional) [numeric] the magnification to be used for x and y labels relative to the current setting of \code{cex}} + +\item{cex.leg}{(optional) [numeric] the magnification to be used for the legend labels relative to the current setting of \code{cex}} + +\item{lwd}{(optional) [numeric] the line width (a positive number)} + +\item{verbose}{(optional) [boolean] indicating if the function is run in verbose mode or not, default = \code{TRUE}} \item{...}{other parameters to be passed through to plotting functions} } -- GitLab