diff --git a/DESCRIPTION b/DESCRIPTION
index c57b529b27e1f0ea5927345d2a8fa2e9ddd7df96..e8b597bde1f51637cddc6700e7b2fec18650e0de 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,7 +1,7 @@
 Package: airGR
 Type: Package
 Title: Suite of GR Hydrological Models for Precipitation-Runoff Modelling
-Version: 1.2.14.19
+Version: 1.2.15.0
 Date: 2019-05-02
 Authors@R: c(
   person("Laurent", "Coron", role = c("aut", "trl"), comment = c(ORCID = "0000-0002-1503-6204")),
diff --git a/NEWS.rmd b/NEWS.rmd
index 2b820d03f31db7de336752736cebe3d6f9b1102c..10c94f899a08b41aa085defcf7ac99c088bf9f8c 100644
--- a/NEWS.rmd
+++ b/NEWS.rmd
@@ -14,7 +14,7 @@ output:
 
 
 
-### 1.2.14.19 Release Notes (2019-05-02)
+### 1.2.15.0 Release Notes (2019-05-02)
 
 
 #### New features
@@ -22,6 +22,11 @@ output:
 - <code>CreateInputsCrit()</code> now allows power (as numeric or character values) and the Box-Cox transformations in the <code>transfo</code> argument.
 
 
+#### Major user-visible changes
+
+- <code>plot.OutputsModel()</code> can no drawn PET time serie if <code>which = "all"</code> or <code>"PotEvap"</code>.
+
+
 #### Minor user-visible changes
 
 - <code>.ErrorCrit()</code> private function added to check inputs into <code>ErrorCrit_&#42;()</code> functions. The <code>ErrorCrit_&#42;()</code> functions were simplified accordingly.
diff --git a/R/plot.OutputsModel.R b/R/plot.OutputsModel.R
index 0ea0e1a1514e282b1ba1440eb6f1762ee35c1290..05258173ed7b5d4da9e459673b8c8fd8417999b7 100644
--- a/R/plot.OutputsModel.R
+++ b/R/plot.OutputsModel.R
@@ -18,8 +18,13 @@ plot.OutputsModel <- function(x, Qobs = NULL, IndPeriod_Plot = NULL, BasinArea =
   if ("Precip" %in% names(OutputsModel)) {
     BOOL_Pobs <- TRUE
   }
-  BOOL_Qsim <- FALSE
   
+  BOOL_Eobs <- FALSE
+  if ("PotEvap" %in% names(OutputsModel)) {
+    BOOL_Eobs <- TRUE
+  }
+  
+  BOOL_Qsim <- FALSE
   if ("Qsim"   %in% names(OutputsModel)) {
     BOOL_Qsim <- TRUE
   }
@@ -33,6 +38,7 @@ plot.OutputsModel <- function(x, Qobs = NULL, IndPeriod_Plot = NULL, BasinArea =
   }
   BOOL_Snow <- FALSE
   
+  
   if ("CemaNeigeLayers" %in% names(OutputsModel)) {
     if ("SnowPack" %in% names(OutputsModel$CemaNeigeLayers[[1]])) {
       BOOL_Snow <- TRUE
@@ -45,6 +51,7 @@ plot.OutputsModel <- function(x, Qobs = NULL, IndPeriod_Plot = NULL, BasinArea =
     }
   }
   
+  
   if (is.null(which)) {
     stop("'which' must be a vector of character")
   }
@@ -54,18 +61,18 @@ plot.OutputsModel <- function(x, Qobs = NULL, IndPeriod_Plot = NULL, BasinArea =
   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 (any(!which %in% c("all", "PotEvap", "Precip", 'Temp', "SnowPack", "Flows", "Regime", "CumFreq", "CorQQ"))) {
+    stop("incorrect element found in argument 'which':\nit can only contain 'all', 'Precip', 'PotEvap', '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'")
+    stop("Incorrect element found in argument 'which':\nwithout CemaNeige it can only contain 'all', 'Precip', 'PotEvap', '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")
+    warning("Incorrect element found in argument 'which':\nit can only contain 'all', 'Precip', 'PotEvap', '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")
+    which <- c("Precip", "PotEvap", "Temp", "SnowPack", "Flows", "Regime", "CumFreq", "CorQQ")
   }
   
   
@@ -159,6 +166,7 @@ plot.OutputsModel <- function(x, Qobs = NULL, IndPeriod_Plot = NULL, BasinArea =
   
   ## Plots_choices
   BOOLPLOT_Precip   <- ( "Precip"   %in% which & BOOL_Pobs )
+  BOOLPLOT_PotEvap  <- ( "PotEvap"  %in% which & BOOL_Eobs )
   BOOLPLOT_Temp     <- ( "Temp"     %in% which & BOOL_Snow )
   BOOLPLOT_SnowPack <- ( "SnowPack" %in% which & BOOL_Snow )
   BOOLPLOT_Flows    <- ( "Flows"    %in% which & (BOOL_Qsim | BOOL_Qobs) )
@@ -183,6 +191,10 @@ plot.OutputsModel <- function(x, Qobs = NULL, IndPeriod_Plot = NULL, BasinArea =
       matlayout <- rbind(matlayout, c(iPlot+1, iPlot+1, iPlot+1))
       iPlot <- iPlot + 1
     }
+    if (BOOLPLOT_PotEvap) {
+      matlayout <- rbind(matlayout, c(iPlot + 1, iPlot + 1, iPlot + 1), c(iPlot + 1, iPlot + 1, iPlot + 1))
+      iPlot <- iPlot + 1
+    } 
     if (BOOLPLOT_Temp) {
       matlayout <- rbind(matlayout, c(iPlot+1, iPlot+1, iPlot+1), c(iPlot+1, iPlot+1, iPlot+1))
       iPlot <- iPlot+1
@@ -347,6 +359,32 @@ plot.OutputsModel <- function(x, Qobs = NULL, IndPeriod_Plot = NULL, BasinArea =
   }
   
   
+  ## PotEvap
+  if (BOOLPLOT_PotEvap) {
+    kPlot <- kPlot + 1
+    mar <- c(3, 5, 1, 5)
+    
+    par(new = FALSE, mar = mar, las = 0)
+    ylim1 <- range(OutputsModel$PotEvap[IndPeriod_Plot], na.rm = TRUE)
+    ylim2 <- ylim1 #* c(1.0, 1.1)
+    
+    plot(Xaxis, OutputsModel$PotEvap[IndPeriod_Plot],
+         type = "l", xaxt = "n", yaxt = "n", ylim = ylim2,
+         col = "green3", lwd = lwd * lwdk,
+         xlab = "", ylab = "", ...)
+    axis(side = 2, at = pretty(ylim1), labels = pretty(ylim1), cex.axis = cex.axis, ...)
+    par(las = 0)
+    mtext(side = 2, paste("pot. evap.", plotunit), cex = cex.lab, line = line)
+    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, ...)
+    }
+  }
+  
+  
   ## SnowPack
   if (BOOLPLOT_SnowPack) {
     kPlot <- kPlot + 1