diff --git a/DESCRIPTION b/DESCRIPTION index af6d097..ec9b6dd 100755 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -28,8 +28,8 @@ Suggests: knitr Type: Package Title: Visualizing climate4R data and Communicating Uncertainty in Seasonal Climate Prediction -Version: 1.6.3 -Date: 2023-06-22 +Version: 1.6.4 +Date: 2023-10-26 Author: as.person(c( "Santander Meteorology Group [cph]", "María Dolores Frías [aut, cre]", diff --git a/NAMESPACE b/NAMESPACE index 3790e8a..27f8a18 100755 --- a/NAMESPACE +++ b/NAMESPACE @@ -19,6 +19,7 @@ export(skillMap) export(spatialMean) export(spatialPlot) export(spreadPlot) +export(stripePlot) export(temporalPlot) export(tercileBarplot) export(tercileMap) diff --git a/NEWS b/NEWS index 325d5df..f71be63 100755 --- a/NEWS +++ b/NEWS @@ -66,3 +66,5 @@ ## v1.6.3 (22 Jun 2023) * Update of the DESCRIPTION file +## v1.6.4 (26 Oct 2023) +* New function fo plot stripes diff --git a/R/map.hatching.R b/R/map.hatching.R index fcde3db..c43ce4a 100644 --- a/R/map.hatching.R +++ b/R/map.hatching.R @@ -36,6 +36,7 @@ #' # as example: #' #' clim <- climatology(tas.ncep) +#' attr(clim$xyCoords, "projection") <- NA #' spatialPlot(clim, backdrop.theme = "coastline", rev.colors = TRUE) #' #' # Basic usage diff --git a/R/stripePlot.R b/R/stripePlot.R new file mode 100644 index 0000000..30c759d --- /dev/null +++ b/R/stripePlot.R @@ -0,0 +1,105 @@ +# stripePlot.R Lattice plot methods for grids +# +# Copyright (C) 2017 Santander Meteorology Group (http://www.meteo.unican.es) +# +# This program is free software: you can redistribute it and/or modify +# it under the terms of the GNU General Public License as published by +# the Free Software Foundation, either version 3 of the License, or +# (at your option) any later version. +# +# This program is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. +# +# You should have received a copy of the GNU General Public License +# along with this program. If not, see . + + + +#' @title Lattice plot methods for producing stripe plots +#' @description A wrapper for the lattice (trellis) plot methods for grid and station data. +#' @param grid Input grid (or station data), ideally multimember and with time dimension > 1. +#' @param aggr.spatial list containing the function and the related arguments to perform spatial +#' aggregation. Default is \code{list(FUN = mean, na.rm = TRUE)}. +#' @param color.theme A character string indicating the color theme to use in the map. +#' Valid values are those available in the \code{\link{RColorBrewer}} themes. Additionally, +#' the \code{"jet.colors"} palette can be used (the rainbow colors, in general not advised, though), +#' for backwards compatibility. Default to the diverging, colorblind-friendly \code{"RdYlBu"} palette. +#' NOTE: the \code{color.theme} argument will be overriden if the \code{col.regions} option from \code{levelplot} is used. +#' @param rev.colors Should the chosen color theme be reversed? (default to \code{FALSE}, +#' leaving the palette \dQuote{as is}). +#' @param xlabels Character string of the labels to be shown in the x axis (Default to NULL) +#' @param ylabels Character string of the labels to be shown in the y axis (Default to NULL) +#' @param ... Additional arguments as passed to function \code{levelplot}. +#' @details The function applies the \code{\link[lattice]{levelplot}} method after computing spatial aggregation +#' (parameter \code{aggr.spatial}) to the imput grids (or station data). +#' +#' This function is appropriate for multimember grids, each row represents a member. +#' +#' \strong{Controlling graphical parameters} +#' +#' Many different aspects of the plot can be controlled passing the relevant arguments to +#' \code{levelplot}. Fine control of graphical parameters for the trellis display can +#' be also controlled using \code{\link[lattice]{trellis.par.set}}. +#' +#' +# +#' @return A lattice plot of class \dQuote{trellis}. +#' +#' +#'@author M. Iturbide +#'@export +#'@import lattice +#'@importFrom grDevices col2rgb rgb colors +#'@import latticeExtra +#'@importFrom transformeR aggregateGrid +#'@importFrom RColorBrewer brewer.pal +#' @examples \donttest{ +#' require(climate4R.datasets) +#' data("CFS_Iberia_pr") +#' stripePlot(CFS_Iberia_pr) +#' } + + +stripePlot <- function(grid, + aggr.spatial = list(FUN = "mean", na.rm = TRUE), + color.theme = "RdBu", + rev.colors = FALSE, + xlabels = NULL, + ylabels = NULL, + ...){ + lp <- list(...) + reg.mean <- aggregateGrid(grid, aggr.spatial = aggr.spatial) + mat <- t(reg.mean$Data) + if (is.null(xlabels)) { + rownames(mat) <- substring(grid$Dates$start, first = 1, last = 4) + } else { + if(length(xlabels) != nrow(mat)) stop("xlabels must have the same length as the time dimension in grid") + rownames(mat) <- xlabels + } + if (is.null(ylabels)) { + colnames(mat) <- grid$Members + } else { + if(length(ylabels) != ncol(mat)) stop("ylabels must have the same length as the member dimension in grid") + colnames(mat) <- ylabels + } + # PLOT + col <- brewer.pal(n = 9, name = color.theme) + if (rev.colors) col <- rev(col) + lp[["x"]] <- mat + if (is.null(lp[["col.regions"]])) lp[["col.regions"]] <- colorRampPalette(col) + if (is.null(lp[["xlab"]])) lp[["xlab"]] <- "" + if (is.null(lp[["ylab"]])) lp[["ylab"]] <- "" + if(nrow(mat) > 100) { xscales <- seq(1, nrow(mat), 100) + } else { + xscales <- seq(1, nrow(mat)) + } + if (is.null(lp[["scales"]])) lp[["scales"]] <- list(x = list(at = xscales, rot = 45, cex = 0.7), y = list(rot = 0, cex = 0.7)) + #if (is.null(lp[["scale"]])) lp[["scale"]] <- list(alternating = 1) + if (is.null(lp[["colorkey"]])) lp[["colorkey"]] <- list(width = 0.9, labels = list(cex = 0.7)) + if (is.null(lp[["aspect"]])) lp[["aspect"]] <- "fill" + do.call("levelplot", lp) +} + +#end diff --git a/R/temporalPlot.R b/R/temporalPlot.R index 3638f47..4402019 100644 --- a/R/temporalPlot.R +++ b/R/temporalPlot.R @@ -73,10 +73,9 @@ #' @examples \donttest{ #' require(climate4R.datasets) #' require(transformeR) -#' data("CFS_Iberia_pr") -#' data("EOBS_Iberia_pr") +#' data("CFS_Iberia_tas") #' data("EOBS_Iberia_tas") -#' data("VALUE_Iberia_pr") +#' data("VALUE_Iberia_tas") #' # Combine grids with members (CFS) and without members (EOBS) #' a <- subsetGrid(CFS_Iberia_tas, years = 1985:1992) #' b <- subsetGrid(EOBS_Iberia_tas, years = 1985:1992) diff --git a/man/map.hatching.Rd b/man/map.hatching.Rd index a8db280..4c6872f 100644 --- a/man/map.hatching.Rd +++ b/man/map.hatching.Rd @@ -69,6 +69,7 @@ data(tas.ncep) # as example: clim <- climatology(tas.ncep) +attr(clim$xyCoords, "projection") <- NA spatialPlot(clim, backdrop.theme = "coastline", rev.colors = TRUE) # Basic usage diff --git a/man/stripePlot.Rd b/man/stripePlot.Rd new file mode 100644 index 0000000..1e328ea --- /dev/null +++ b/man/stripePlot.Rd @@ -0,0 +1,65 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/stripePlot.R +\name{stripePlot} +\alias{stripePlot} +\title{Lattice plot methods for producing stripe plots} +\usage{ +stripePlot( + grid, + aggr.spatial = list(FUN = "mean", na.rm = TRUE), + color.theme = "RdBu", + rev.colors = FALSE, + xlabels = NULL, + ylabels = NULL, + ... +) +} +\arguments{ +\item{grid}{Input grid (or station data), ideally multimember and with time dimension > 1.} + +\item{aggr.spatial}{list containing the function and the related arguments to perform spatial +aggregation. Default is \code{list(FUN = mean, na.rm = TRUE)}.} + +\item{color.theme}{A character string indicating the color theme to use in the map. +Valid values are those available in the \code{\link{RColorBrewer}} themes. Additionally, +the \code{"jet.colors"} palette can be used (the rainbow colors, in general not advised, though), +for backwards compatibility. Default to the diverging, colorblind-friendly \code{"RdYlBu"} palette. +NOTE: the \code{color.theme} argument will be overriden if the \code{col.regions} option from \code{levelplot} is used.} + +\item{rev.colors}{Should the chosen color theme be reversed? (default to \code{FALSE}, +leaving the palette \dQuote{as is}).} + +\item{xlabels}{Character string of the labels to be shown in the x axis (Default to NULL)} + +\item{ylabels}{Character string of the labels to be shown in the y axis (Default to NULL)} + +\item{...}{Additional arguments as passed to function \code{levelplot}.} +} +\value{ +A lattice plot of class \dQuote{trellis}. +} +\description{ +A wrapper for the lattice (trellis) plot methods for grid and station data. +} +\details{ +The function applies the \code{\link[lattice]{levelplot}} method after computing spatial aggregation +(parameter \code{aggr.spatial}) to the imput grids (or station data). + +This function is appropriate for multimember grids, each row represents a member. + + \strong{Controlling graphical parameters} + + Many different aspects of the plot can be controlled passing the relevant arguments to + \code{levelplot}. Fine control of graphical parameters for the trellis display can + be also controlled using \code{\link[lattice]{trellis.par.set}}. +} +\examples{ +\donttest{ +require(climate4R.datasets) +data("CFS_Iberia_pr") +stripePlot(CFS_Iberia_pr) +} +} +\author{ +M. Iturbide +} diff --git a/man/temporalPlot.Rd b/man/temporalPlot.Rd index 3754a49..b1e30fa 100644 --- a/man/temporalPlot.Rd +++ b/man/temporalPlot.Rd @@ -79,10 +79,9 @@ and the range for plotting. The range is used to plot the shadow of the multimem \donttest{ require(climate4R.datasets) require(transformeR) -data("CFS_Iberia_pr") -data("EOBS_Iberia_pr") +data("CFS_Iberia_tas") data("EOBS_Iberia_tas") -data("VALUE_Iberia_pr") +data("VALUE_Iberia_tas") # Combine grids with members (CFS) and without members (EOBS) a <- subsetGrid(CFS_Iberia_tas, years = 1985:1992) b <- subsetGrid(EOBS_Iberia_tas, years = 1985:1992)