diff --git a/DESCRIPTION b/DESCRIPTION
index 9123c589b1411638ea21e364a1592dc7be137f30..3d62d705d17f0ad119154b9faff8ee026a1cb241 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -76,7 +76,9 @@ Imports:
     grDevices,
     stats,
     utils,
-    verification
+    verification,
+    lubridate,
+    scales
 Suggests:
     zeallot,
     testthat,
diff --git a/NAMESPACE b/NAMESPACE
index 5884b2e411691d6b4ac1984dc618f7349ad10e0e..3fa054b2a747109be678c4898880e4e410ecc874 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -43,6 +43,7 @@ export(PlotForecastPDF)
 export(PlotMostLikelyQuantileMap)
 export(PlotPDFsOLE)
 export(PlotTriangles4Categories)
+export(PlotWeeklyClim)
 export(Predictability)
 export(ProxiesAttractor)
 export(QuantileMapping)
@@ -57,13 +58,16 @@ export(WeatherRegime)
 export(as.s2dv_cube)
 export(s2dv_cube)
 export(training_analogs)
+import(RColorBrewer)
 import(abind)
 import(ggplot2)
+import(lubridate)
 import(multiApply)
 import(ncdf4)
 import(qmap)
 import(rainfarmr)
 import(s2dv)
+import(scales)
 import(stats)
 importFrom(ClimProjDiags,SelBox)
 importFrom(ClimProjDiags,Subset)
diff --git a/R/PlotWeeklyClim.R b/R/PlotWeeklyClim.R
new file mode 100644
index 0000000000000000000000000000000000000000..6f103d9641ca7f02d19f853dfabd4bfbc0587408
--- /dev/null
+++ b/R/PlotWeeklyClim.R
@@ -0,0 +1,231 @@
+#'Plots the observed weekly means and climatology of a timeseries data
+#' 
+#'@description This function plots the observed weekly means and climatology of 
+#'a timeseries data using ggplot package. It compares the weekly climatology in 
+#'a specified period (reference period) to the observed conditions during the 
+#'target period analyzed in the case study (included in the reference period).
+#' 
+#'@param data A multidimensional array with named dimensions with at least sdate 
+#'  and time dimensions containing observed daily data. It can also be a 
+#'  dataframe with computed percentiles as input for ggplot. The target year 
+#'  must be included in the input data.
+#'@param first_date The first date of the target period of timeseries. It can be 
+#'  of class 'Date', 'POSIXct' or a character string in the format 'yyyy-mm-dd'. 
+#'  It must be a date included in the reference period.
+#'@param ref_period_ini A numeric value indicating the first year of the 
+#'  reference period.
+#'@param ref_period_end A numeric value indicating the last year of the 
+#'  reference period.
+#'@param time_dim A character string indicating the daily time dimension name. 
+#'  The default value is 'time'.
+#'@param sdate_dim A character string indicating the start year dimension name. 
+#'  The default value is 'sdate'.
+#'@param title The text for the top title of the plot. 
+#'@param palette A palette name from the R Color Brewer’s package. The default 
+#'  value is 'Blues'.
+#'@param fileout A character string indicating the file name where to save the 
+#'  plot. If not specified (default) a graphics device will pop up.
+#'@param device A character string indicating the device to use. Can either be 
+#'  a device function (e.g. png), or one of "eps", "ps", "tex" (pictex), 
+#'  "pdf", "jpeg", "tiff", "png", "bmp", "svg" or "wmf" (windows only).
+#'@param width A numeric value of the plot width in units ("in", "cm", "mm", or 
+#'  "px"). It is set to 8 by default.
+#'@param height A numeric value of the plot height in units ("in", "cm", "mm", 
+#'  or "px"). It is set to 6 by default.
+#'@param units Units of the size of the device (file or window) to plot in. 
+#'  Inches (’in’) by default.
+#'@param dpi A numeric value of the plot resolution. It is set to 300 by 
+#'  default.
+#' 
+#'@return A ggplot object containing the plot.
+#' 
+#'@examples
+#'data <- array(rnorm(49*20, 274, 7), dim = c(time = 49, sdate = 20))
+#'PlotWeeklyClim(data = data, first_date = '2010-08-09', 
+#'               ref_period_ini = 1998,
+#'               ref_period_end = 2020)
+#' 
+#'@import multiApply
+#'@import lubridate
+#'@import ggplot2
+#'@import RColorBrewer
+#'@import scales
+#'@importFrom ClimProjDiags Subset
+#'@importFrom s2dv MeanDims
+#'@export
+PlotWeeklyClim <- function(data, first_date, ref_period_ini, ref_period_end, 
+                           time_dim = 'time', sdate_dim = 'sdate',
+                           title = "Observed weekly means and climatology",
+                           palette = "Blues", fileout = NULL, 
+                           device = NULL, width = 8, height = 6,
+                           units = 'in', dpi = 300) {
+  # Check input arguments
+  # data
+  if (is.array(data)) {
+    if (is.null(names(dim(data)))) {
+      stop("Parameter 'data' must have named dimensions.")
+    }
+    is_array <- TRUE
+  } else if (is.data.frame(data)) {
+    col_names <- c("week", "clim", "p10", "p90", "p33", "p66", 
+                   "week_mean", "day", "data")
+    if (!all(col_names %in% names(data))) {
+      stop(paste0("If parameter 'data' is a data frame, it must contain the ",
+                  "following column names: 'week', 'clim', 'p10', 'p90', 'p33', ", 
+                  "'p66', 'week_mean', 'day' and 'data'."))
+    }
+    is_array <- FALSE
+  } else {
+    stop("Parameter 'data' must be an array or a data frame.")
+  }
+  if (is_array) {
+    # time_dim
+    if (!is.character(time_dim)) {
+      stop("Parameter 'time_dim' must be a character string.")
+    }
+    if (!all(time_dim %in% names(dim(data)))) {
+      stop("Parameter 'time_dim' is not found in 'data' dimension.")
+    }
+    if (dim(data)[time_dim] < 7) {
+      stop(paste0("Parameter 'data' must have the dimension 'time_dim' of ",
+                  "length equal or grater than 7 to compute the weekly means."))
+    }
+    # sdate_dim
+    if (!is.character(sdate_dim)) {
+      stop("Parameter 'sdate_dim' must be a character string.")
+    }
+    if (!sdate_dim %in% names(dim(data))) {
+      warning(paste0("Parameter 'sdate_dim' is not found in 'data' dimension. ",
+                     "A dimension of length 1 has been added."))
+      data <- InsertDim(data, 1, lendim = 1, name = sdate_dim)
+    }
+    # ref_period_ini and ref_period_end
+    if (!is.numeric(ref_period_ini) | !is.numeric(ref_period_end)) {
+      stop("Parameters 'ref_period_ini' and 'ref_period_end' must be numeric.")
+    }
+    # first_date
+    if ((!inherits(first_date, "POSIXct") & !inherits(first_date, "Date")) &&
+        (!is.character(first_date) | nchar(first_date) != 10)) {
+      stop(paste0("Parameter 'first_date' must be a character string ", 
+                  "indicating the date in the format 'yyyy-mm-dd', 'POSIXct' ", 
+                  "or 'Dates' class."))
+    }
+    first_date <- ymd(first_date)
+    target_year <- year(first_date)
+    if (target_year < ref_period_ini | target_year > ref_period_end) {
+      stop("Parameter 'first_date' must be a date included in the reference period.")
+    }
+
+    # Dates creation
+    dates <- seq(first_date, first_date + days(dim(data)[time_dim]-1), by = "1 day")
+    index_first_date <- which(dates == first_date)
+    index_last_date <- length(dates) - (length(dates) %% 7)
+    last_date <- dates[index_last_date]
+
+    # Weekly aggregations
+    data_subset <- Subset(data, along = time_dim, 
+                          indices = index_first_date:index_last_date)
+    weekly_aggre <- SplitDim(data_subset, split_dim = time_dim, 
+                             indices = sort(rep(1:(length(index_first_date:index_last_date)/7), 7)),
+                             new_dim_name = 'week')
+    weekly_means <- MeanDims(weekly_aggre, time_dim)
+    weekly_clim <- MeanDims(weekly_means, sdate_dim)
+
+    weekly_p10 <- Apply(weekly_means, target_dims = sdate_dim,
+                        fun = function(x) {quantile(x, 0.10)})$output1
+    weekly_p90 <- Apply(weekly_means, target_dims = sdate_dim,
+                        fun = function(x) {quantile(x, 0.90)})$output1
+    weekly_p33 <- Apply(weekly_means, target_dims = sdate_dim,
+                        fun = function(x) {quantile(x, 0.33)})$output1
+    weekly_p66 <- Apply(weekly_means, target_dims = sdate_dim,
+                        fun = function(x) {quantile(x, 0.66)})$output1
+
+    clim <- p10 <- p90 <- p33 <- p66 <- NULL
+    weekly_data <- data.frame(clim = as.vector(weekly_clim),
+                              p10 = as.vector(weekly_p10),
+                              p90 = as.vector(weekly_p90), 
+                              p33 = as.vector(weekly_p33),
+                              p66 = as.vector(weekly_p66),
+                              week = 1:(length(index_first_date:index_last_date)/7))
+
+    daily <- Subset(data_subset, along = sdate_dim, 
+                    indices = which(ref_period_ini:ref_period_end == target_year), 
+                    drop = TRUE)
+
+    dims_subset <- names(dim(daily))[which(!names(dim(daily)) %in% c(time_dim, sdate_dim))]
+
+    if (!identical(dims_subset, character(0))) {
+      daily <- Subset(daily, dims_subset, as.list(rep(1, length(dims_subset))), drop = TRUE)
+    }
+    
+    daily_data <- data.frame(day = seq(first_date, last_date, by = "1 day"),
+                             data = daily, 
+                             week = sort(rep(1:(length(index_first_date:index_last_date)/7), 7)))
+    week_mean <- aggregate(data ~ week, data = daily_data, mean)
+    weekly_data <- cbind(weekly_data, week_mean$data)
+    colnames(weekly_data)[7] <- 'week_mean'
+    all <- merge(weekly_data, daily_data, by = 'week')
+  } else {
+    all <- data
+  }
+
+  # Create a ggplot object
+  cols <- colorRampPalette(brewer.pal(9, palette))(6)
+
+  p <- ggplot(all, aes(x = day)) +
+       geom_ribbon(aes(ymin = p10, ymax = p90, group = week, fill = "p10-p90"), 
+                   alpha = 0.7) + # extremes clim
+       geom_ribbon(aes(ymin = p33, ymax = p66, group = week, fill = "p33-p66"),
+                   alpha = 0.7) +  # terciles clim
+       geom_line(aes(y = clim, group = week, color = "climatological mean",
+                 linetype = "climatological mean"),
+                 alpha = 1.0, size = 0.7) + # mean clim
+       geom_line(aes(y = data, color = "observed daily mean",
+                 linetype = "observed daily mean"),
+                 alpha = 1, size = 0.2) + # daily evolution
+       geom_line(aes(y = week_mean, group = week, color = "observed weekly mean",
+                 linetype = "observed weekly mean"),
+                 alpha = 1, size = 0.7) + 		 # weekly evolution
+       theme_bw() + ylab(paste0('tas', " (", "deg.C", ")")) + xlab(NULL) + 
+       ggtitle(title) +
+       scale_x_date(breaks = seq(min(all$day), max(all$day), by = "7 days"),
+                    minor_breaks = NULL, expand = c(0.03, 0.03),
+                    labels = date_format("%d %b %Y")) +
+       theme(axis.text.x = element_text(angle = 45, hjust = 1),
+             panel.grid.major = element_line(size = 0.5, linetype = 'solid',
+                                             colour = "gray92"),
+             panel.grid.minor = element_line(size = 0.25, linetype = 'solid',
+                                             colour = "gray92"),
+             legend.spacing = unit(-0.2, "cm")) +
+       scale_fill_manual(name = NULL,
+                         values = c("p10-p90" = cols[3], "p33-p66" = cols[4])) +
+       scale_color_manual(name = NULL, values = c("climatological mean" = cols[5],
+                                                  "observed daily mean" = "grey20",
+                                                  "observed weekly mean" = "black")) +
+       scale_linetype_manual(name = NULL, values = c("climatological mean" = "solid",
+                                                     "observed daily mean" = "dashed",
+                                                     "observed weekly mean" = "solid"),
+       guide = guide_legend(override.aes = list(lwd = c(0.7, 0.2, 0.7)))) +
+       guides(fill = guide_legend(order = 1))
+
+  # Return the ggplot object
+  if (is.null(fileout)) {
+    return(p)
+  } else {
+    ggsave(filename = fileout, plot = p, device = device, height = height, 
+           width = width, units = units, dpi = dpi) 
+  }
+}
+
+
+
+
+
+
+
+
+
+
+
+
+  
\ No newline at end of file
diff --git a/man/PlotWeeklyClim.Rd b/man/PlotWeeklyClim.Rd
new file mode 100644
index 0000000000000000000000000000000000000000..746c641ea4f411cfc36d5d52249f87b45ebbb5d2
--- /dev/null
+++ b/man/PlotWeeklyClim.Rd
@@ -0,0 +1,85 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/PlotWeeklyClim.R
+\name{PlotWeeklyClim}
+\alias{PlotWeeklyClim}
+\title{Plots the observed weekly means and climatology of a timeseries data}
+\usage{
+PlotWeeklyClim(
+  data,
+  first_date,
+  ref_period_ini,
+  ref_period_end,
+  time_dim = "time",
+  sdate_dim = "sdate",
+  title = "Observed weekly means and climatology",
+  palette = "Blues",
+  fileout = NULL,
+  device = NULL,
+  width = 8,
+  height = 6,
+  units = "in",
+  dpi = 300
+)
+}
+\arguments{
+\item{data}{A multidimensional array with named dimensions with at least sdate 
+and time dimensions containing observed daily data. It can also be a 
+dataframe with computed percentiles as input for ggplot. The target year 
+must be included in the input data.}
+
+\item{first_date}{The first date of the target period of timeseries. It can be 
+of class 'Date', 'POSIXct' or a character string in the format 'yyyy-mm-dd'. 
+It must be a date included in the reference period.}
+
+\item{ref_period_ini}{A numeric value indicating the first year of the 
+reference period.}
+
+\item{ref_period_end}{A numeric value indicating the last year of the 
+reference period.}
+
+\item{time_dim}{A character string indicating the daily time dimension name. 
+The default value is 'time'.}
+
+\item{sdate_dim}{A character string indicating the start year dimension name. 
+The default value is 'sdate'.}
+
+\item{title}{The text for the top title of the plot.}
+
+\item{palette}{A palette name from the R Color Brewer’s package. The default 
+value is 'Blues'.}
+
+\item{fileout}{A character string indicating the file name where to save the 
+plot. If not specified (default) a graphics device will pop up.}
+
+\item{device}{A character string indicating the device to use. Can either be 
+a device function (e.g. png), or one of "eps", "ps", "tex" (pictex), 
+"pdf", "jpeg", "tiff", "png", "bmp", "svg" or "wmf" (windows only).}
+
+\item{width}{A numeric value of the plot width in units ("in", "cm", "mm", or 
+"px"). It is set to 8 by default.}
+
+\item{height}{A numeric value of the plot height in units ("in", "cm", "mm", 
+or "px"). It is set to 6 by default.}
+
+\item{units}{Units of the size of the device (file or window) to plot in. 
+Inches (’in’) by default.}
+
+\item{dpi}{A numeric value of the plot resolution. It is set to 300 by 
+default.}
+}
+\value{
+A ggplot object containing the plot.
+}
+\description{
+This function plots the observed weekly means and climatology of 
+a timeseries data using ggplot package. It compares the weekly climatology in 
+a specified period (reference period) to the observed conditions during the 
+target period analyzed in the case study (included in the reference period).
+}
+\examples{
+data <- array(rnorm(49*20, 274, 7), dim = c(time = 49, sdate = 20))
+PlotWeeklyClim(data = data, first_date = '2010-08-09', 
+              ref_period_ini = 1998,
+              ref_period_end = 2020)
+
+}
diff --git a/tests/testthat/test-PlotWeeklyClim.R b/tests/testthat/test-PlotWeeklyClim.R
new file mode 100644
index 0000000000000000000000000000000000000000..e7886fbc6836b809ca039dc577a37a65fd42037b
--- /dev/null
+++ b/tests/testthat/test-PlotWeeklyClim.R
@@ -0,0 +1,90 @@
+context("CSTools::PlotWeeklyClim tests")
+
+##############################################
+
+# dat1
+dat1 <- array(rnorm(1*7), dim = c(dat = 1, var = 1, sdate = 1, time = 7))
+
+##############################################
+
+test_that("1. Input checks", {
+  # data
+  expect_error(
+    PlotWeeklyClim(data = array(1:92), first_date = '2020-03-01', 
+                       ref_period_ini = 1993, ref_period_end = 2021),
+    "Parameter 'data' must have named dimensions."
+  )
+  expect_error(
+    PlotWeeklyClim(data = data.frame(week = 1:92), first_date = '2020-03-01', 
+                       ref_period_ini = 1993, ref_period_end = 2021),
+    paste0("If parameter 'data' is a data frame, it must contain the ",
+           "following column names: 'week', 'clim', 'p10', 'p90', 'p33', ", 
+           "'p66', 'week_mean', 'day' and 'data'.")
+  )
+  expect_error(
+    PlotWeeklyClim(data = 1:92, first_date = '2020-03-01', 
+                       ref_period_ini = 1993, ref_period_end = 2021),
+    "Parameter 'data' must be an array or a data frame."
+  )
+  # time_dim
+  expect_error(
+    PlotWeeklyClim(data = dat1, first_date = '2020-03-01', 
+                       ref_period_ini = 2020, ref_period_end = 2020, time_dim = 1),
+    "Parameter 'time_dim' must be a character string."
+  )
+  expect_error(
+    PlotWeeklyClim(data = array(rnorm(1), dim = c(dat = 1)), 
+                       first_date = '2020-03-01', ref_period_ini = 2020, 
+                       ref_period_end = 2020),
+    "Parameter 'time_dim' is not found in 'data' dimension."
+  )
+  expect_error(
+    PlotWeeklyClim(data = array(rnorm(1*7), dim = c(time = 6)), 
+                       first_date = '2020-03-01', ref_period_ini = 2020, 
+                       ref_period_end = 2020),
+    paste0("Parameter 'data' must have the dimension 'time_dim' of length ",
+           "equal or grater than 7 to compute the weekly means.")
+  )
+  # sdate_dim
+  expect_error(
+    PlotWeeklyClim(data = dat1, first_date = '2020-03-01', 
+                       ref_period_ini = 2020, ref_period_end = 2020, 
+                       sdate_dim = 1),
+    "Parameter 'sdate_dim' must be a character string."
+  )
+  expect_warning(
+    PlotWeeklyClim(data = array(rnorm(7), dim = c(time = 7)), 
+                       first_date = '2020-03-01', ref_period_ini = 2020, 
+                       ref_period_end = 2020),
+    paste0("Parameter 'sdate_dim' is not found in 'data' dimension. ",
+           "A dimension of length 1 has been added.")
+  )
+  # ref_period_ini
+  expect_error(
+    PlotWeeklyClim(data = dat1, first_date = '2020-03-01', 
+                       ref_period_ini = "2020", ref_period_end = 2020),
+    "Parameters 'ref_period_ini' and 'ref_period_end' must be numeric."
+  )
+  # first_date
+  expect_error(
+    PlotWeeklyClim(data = dat1, first_date = 2020-03-01, 
+                       ref_period_ini = 2020, ref_period_end = 2020),
+    paste0("Parameter 'first_date' must be a character string ", 
+           "indicating the date in the format 'yyyy-mm-dd', 'POSIXct' ", 
+           "or 'Dates' class.")
+  )
+  expect_error(
+    PlotWeeklyClim(data = dat1, first_date = 'a', 
+                       ref_period_ini = 2020, ref_period_end = 2020),
+    paste0("Parameter 'first_date' must be a character string ", 
+           "indicating the date in the format 'yyyy-mm-dd', 'POSIXct' ", 
+           "or 'Dates' class.")
+  )
+  expect_error(
+    PlotWeeklyClim(data = dat1, first_date = '2020-03-01', ref_period_ini = 2021,
+                       ref_period_end = 2022),
+    "Parameter 'first_date' must be a date included in the reference period."
+  )
+})
+
+##############################################
\ No newline at end of file