diff --git a/NEWS.md b/NEWS.md
index ae7955dd26f648927f46437124eddcebf62384f6..9a4c78e066f915b46bc8b7ee9db664f512ca455d 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -9,6 +9,7 @@
     + Calibration retains correlation absolute value 
     + Calibration fixed when cal.methodi == rpc-based, apply_to == sign, 
         eval.method == 'leave-one-out' and the correlation is not significant
+    + PlotMostLikelyQuantileMap reoder latitudes of an array provided in 'dots' parameter.
 
 ### CSTools 4.0.0 
 **Submission date to CRAN: 23-02-2021**
diff --git a/R/PlotCombinedMap.R b/R/PlotCombinedMap.R
index 2aee71b8ef7a3861a3834e05f0110ea10bc1751d..7169e4c43a8bab1e9a30accca7433b0702760d12 100644
--- a/R/PlotCombinedMap.R
+++ b/R/PlotCombinedMap.R
@@ -14,6 +14,12 @@
 #'@param col_unknown_map Colour to use to paint the grid cells for which a map is not possible to be chosen according to 'map_select_fun' or for those values that go beyond 'display_range'. Takes the value 'white' by default.
 #'@param mask Optional numeric array with dimensions (latitude, longitude), with values in the range [0, 1], indicating the opacity of the mask over each grid point. Cells with a 0 will result in no mask, whereas cells with a 1 will result in a totally opaque superimposed pixel coloured in 'col_mask'.
 #'@param col_mask Colour to be used for the superimposed mask (if specified in 'mask'). Takes the value 'grey' by default.
+#'@param dots Array of same dimensions as 'var' or with dimensions 
+#'  c(n, dim(var)), where n is the number of dot/symbol layers to add to the 
+#'  plot. A value of TRUE at a grid cell will draw a dot/symbol on the 
+#'  corresponding square of the plot. By default all layers provided in 'dots' 
+#'  are plotted with dots, but a symbol can be specified for each of the 
+#'  layers via the parameter 'dot_symbol'.
 #'@param bar_titles Optional vector of character strings providing the titles to be shown on top of each of the colour bars.
 #'@param legend_scale Scale factor for the size of the colour bar labels. Takes 1 by default.
 #'@param fileout File where to save the plot. If not specified (default) a graphics device will pop up. Extensions allowed: eps/ps, jpeg, png, pdf, bmp and tiff
@@ -61,6 +67,7 @@ PlotCombinedMap <- function(maps, lon, lat,
                             brks = NULL, cols = NULL,  
                             col_unknown_map = 'white',
                             mask = NULL, col_mask = 'grey',
+                            dots = NULL,
                             bar_titles = NULL, legend_scale = 1,
                             fileout = NULL, width = 8, height = 5, 
                             size_units = 'in', res = 100, 
@@ -280,7 +287,17 @@ PlotCombinedMap <- function(maps, lon, lat,
       stop("Parameter 'mask' must have dimensions c(lat, lon).")
     }
   }
-  
+  # Check dots
+  if (!is.null(dots)) {
+    if (length(dim(dots)) != 2) {
+      stop("Parameter 'mask' must have two dimensions.")
+    }
+    if ((dim(dots)[1] != dim(maps)[lat_dim]) ||
+        (dim(dots)[2] != dim(maps)[lon_dim])) {
+      stop("Parameter 'mask' must have dimensions c(lat, lon).")
+    }
+  }
+ 
   #----------------------
   # Identify the most likely map
   #----------------------
@@ -327,6 +344,9 @@ PlotCombinedMap <- function(maps, lon, lat,
     if (!is.null(mask)){
       mask <- mask[nlat:1, ]
     }
+    if (!is.null(dots)){
+      dots <- dots[nlat:1,]
+    }
   }
   
   #----------------------
@@ -353,7 +373,7 @@ PlotCombinedMap <- function(maps, lon, lat,
   tbrks <- c(-1, brks_norm + rep(1:nmap, each = length(brks)))
   PlotEquiMap(var = ml_map, lon = lon, lat = lat, 
               brks = tbrks, cols = tcols, drawleg = FALSE, 
-              filled.continents = FALSE, ...)
+              filled.continents = FALSE, dots = dots, ...)
   
   #----------------------
   # Add overplot on top
diff --git a/man/Calibration.Rd b/man/Calibration.Rd
index f61a3cd3537844714533a899fa196125e9917319..7ac9cc2d11f272b7aa2bd05a07cca7c7c4843a42 100644
--- a/man/Calibration.Rd
+++ b/man/Calibration.Rd
@@ -48,7 +48,7 @@ Calibration(
 an array containing the calibrated forecasts with the same dimensions as the \code{exp} array.
 }
 \description{
-Four types of member-by-member bias correction can be performed. The \code{"bias"} method corrects the bias only, the \code{"evmos"} method applies a variance inflation technique to ensure the correction of the bias and the correspondence of variance between forecast and observation (Van Schaeybroeck and Vannitsem, 2011). The ensemble calibration methods \code{"mse_min"} and \code{"crps_min"} correct the bias, the overall forecast variance and the ensemble spread as described in Doblas-Reyes et al. (2005) and Van Schaeybroeck and Vannitsem (2015), respectively. While the \code{"mse_min"} method minimizes a constrained mean-squared error using three parameters, the \code{"crps_min"} method features four parameters and minimizes the Continuous Ranked Probability Score (CRPS). The \code{"rpc-based"} method adjusts the forecast variance ensuring that the ratio of predictable components (RPC) is equal to one, as in Eade et al. (2014).
+Five types of member-by-member bias correction can be performed. The \code{"bias"} method corrects the bias only, the \code{"evmos"} method applies a variance inflation technique to ensure the correction of the bias and the correspondence of variance between forecast and observation (Van Schaeybroeck and Vannitsem, 2011). The ensemble calibration methods \code{"mse_min"} and \code{"crps_min"} correct the bias, the overall forecast variance and the ensemble spread as described in Doblas-Reyes et al. (2005) and Van Schaeybroeck and Vannitsem (2015), respectively. While the \code{"mse_min"} method minimizes a constrained mean-squared error using three parameters, the \code{"crps_min"} method features four parameters and minimizes the Continuous Ranked Probability Score (CRPS). The \code{"rpc-based"} method adjusts the forecast variance ensuring that the ratio of predictable components (RPC) is equal to one, as in Eade et al. (2014).
 
 Both in-sample or our out-of-sample (leave-one-out cross validation) calibration are possible.
 }
diff --git a/man/PlotCombinedMap.Rd b/man/PlotCombinedMap.Rd
index c45d1afb594c362a6ecf8e4e641cf7223f8a0e1c..e631761e8188e591b045a905852a0be76c0904fa 100644
--- a/man/PlotCombinedMap.Rd
+++ b/man/PlotCombinedMap.Rd
@@ -16,6 +16,7 @@ PlotCombinedMap(
   col_unknown_map = "white",
   mask = NULL,
   col_mask = "grey",
+  dots = NULL,
   bar_titles = NULL,
   legend_scale = 1,
   fileout = NULL,
@@ -49,6 +50,13 @@ PlotCombinedMap(
 
 \item{col_mask}{Colour to be used for the superimposed mask (if specified in 'mask'). Takes the value 'grey' by default.}
 
+\item{dots}{Array of same dimensions as 'var' or with dimensions 
+c(n, dim(var)), where n is the number of dot/symbol layers to add to the 
+plot. A value of TRUE at a grid cell will draw a dot/symbol on the 
+corresponding square of the plot. By default all layers provided in 'dots' 
+are plotted with dots, but a symbol can be specified for each of the 
+layers via the parameter 'dot_symbol'.}
+
 \item{bar_titles}{Optional vector of character strings providing the titles to be shown on top of each of the colour bars.}
 
 \item{legend_scale}{Scale factor for the size of the colour bar labels. Takes 1 by default.}