diff --git a/NAMESPACE b/NAMESPACE index efc5c8353c2d2e8a9f20f8e3980faf5d1d2027b1..af1faf4be6a91ba005ad601f00317ddea291ead4 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -35,6 +35,7 @@ export(PlotPDFsOLE) export(PlotTriangles4Categories) export(RFSlope) export(RFTemp) +export(RF_Weights) export(RainFARM) export(RegimesAssign) export(SaveExp) diff --git a/NEWS.md b/NEWS.md index c8dfc43f5e9d30c271aa6d240e53977051796fc7..cf9eae22058c01a73c340eec2157ec536967da04 100644 --- a/NEWS.md +++ b/NEWS.md @@ -2,7 +2,9 @@ **Submission date to CRAN: XX-12-2020** - New features: + PlotPDFsOLE includes parameters to modify legend style - + CST_RFSlope handless missing values in the temporal dimension and new ncores parameter for parallel computation + + CST_RFSlope handless missing values in the temporal dimension and new 'ncores' parameter allows parallel computation + + CST_RFWeights accepts s2dv_cube objects as input and new 'ncores' paramenter allows parallel computation + + RFWeights is exposed to users - Fixes: + PlotForecastPDF correctly displays terciles labels diff --git a/R/CST_RFWeights.R b/R/CST_RFWeights.R index 16415a572c4c84d3220fd740bae003f2c672f7fa..56d943da4f8bf1a63801c68d1d4a31339fa1e54e 100644 --- a/R/CST_RFWeights.R +++ b/R/CST_RFWeights.R @@ -17,16 +17,22 @@ #' high-resolution gridded climatology from observations, or a reconstruction such as those which #' can be downloaded from the WORLDCLIM (http://www.worldclim.org) or CHELSA (http://chelsa-climate.org) #' websites. The latter data will need to be converted to NetCDF format before being used (see for example the GDAL tools (https://www.gdal.org). +#' It could also be a 's2dv_cube' object. #' @param nf Refinement factor for downscaling (the output resolution is increased by this factor). -#' @param lon Vector of longitudes. +#' @param lon Vector of longitudes. #' @param lat Vector of latitudes. #' The number of longitudes and latitudes is expected to be even and the same. If not #' the function will perform a subsetting to ensure this condition. #' @param varname Name of the variable to be read from \code{climfile}. -#' @param fsmooth Logical to use smooth conservation (default) or large-scale box-average conservation. -#' @return A matrix containing the weights with dimensions (lon, lat). +#' @param fsmooth Logical to use smooth conservation (default) or large-scale box-average conservation. +#' @param lonname a character string indicating the name of the longitudinal dimension set as 'lon' by default. +#' @param latname a character string indicating the name of the latitudinal dimension set as 'lat' by default. +#' @param ncores an integer that indicates the number of cores for parallel computations using multiApply function. The default value is one. +#' +#' @return An object of class 's2dv_cube' containing in matrix \code{data} the weights with dimensions (lon, lat). #' @import ncdf4 #' @import rainfarmr +#' @import multiApply #' @importFrom utils tail #' @importFrom utils head #' @examples @@ -42,8 +48,15 @@ #' } #' @export -CST_RFWeights <- function(climfile, nf, lon, lat, varname = "", fsmooth=TRUE) { - +CST_RFWeights <- function(climfile, nf, lon, lat, varname = NULL, + fsmooth = TRUE, + lonname = 'lon', latname = 'lat', ncores = NULL) { + if (!inherits(climfile, "s2dv_cube")) { + if (!is.null(varname) & !is.character(varname)) { + stop("Parameter 'varname' must be a character string indicating the name", + " of the variable to be read from the file.") + } + } # Ensure input grid is square and with even dimensions if ((length(lat) != length(lon)) | (length(lon) %% 2 == 1)) { warning("Input data are expected to be on a square grid", @@ -57,17 +70,87 @@ CST_RFWeights <- function(climfile, nf, lon, lat, varname = "", fsmooth=TRUE) { " lat: [", lat[1], ", ", lat[length(lat)], "]")) } - ncin <- nc_open(climfile) - latin <- ncvar_get(ncin, grep("lat", attributes(ncin$dim)$names, - value = TRUE)) - lonin <- ncvar_get(ncin, grep("lon", attributes(ncin$dim)$names, - value = TRUE)) - if (varname == "") { - varname <- grep("bnds", attributes(ncin$var)$names, - invert = TRUE, value = TRUE)[1] + if (is.character(climfile)) { + ncin <- nc_open(climfile) + latin <- ncvar_get(ncin, grep(latname, attributes(ncin$dim)$names, + value = TRUE)) + lonin <- ncvar_get(ncin, grep(lonname, attributes(ncin$dim)$names, + value = TRUE)) + if (varname == "") { + varname <- grep("bnds", attributes(ncin$var)$names, + invert = TRUE, value = TRUE)[1] + } + zclim <- ncvar_get(ncin, varname) + nc_close(ncin) + } else if (inherits(climfile, "s2dv_cube")) { + zclim <- climfile$data + latin <- climfile$lat + lonin <- climfile$lon + } else { + stop("Parameter 'climfile' is expected to be a character string indicating", + " the path to the files or an object of class 's2dv_cube'.") + } + # Check dim names and order + if (length(names(dim(zclim))) < 1) { + stop("The dataset provided in 'climfile' requires dimension names.") + } + + result <- RF_Weights(zclim, latin, lonin, nf, lat, lon, fsmooth = fsmooth, + lonname = lonname, latname = latname, ncores = ncores) + if (inherits(climfile, "s2dv_cube")) { + climfile$data <- result$data + climfile$lon <- result$lon + climfile$lat <- result$lat + } else { + climfile <- s2dv_cube(data = result, lon = result$lon, lat = result$lat) } - zclim <- ncvar_get(ncin, varname) + return(climfile) +} +#' Compute climatological weights for RainFARM stochastic precipitation downscaling +#' +#' @author Jost von Hardenberg - ISAC-CNR, \email{j.vonhardenberg@isac.cnr.it} +#' +#' @description Compute climatological ("orographic") weights from a fine-scale precipitation climatology file. +#' @references Terzago, S., Palazzi, E., & von Hardenberg, J. (2018). +#' Stochastic downscaling of precipitation in complex orography: +#' A simple method to reproduce a realistic fine-scale climatology. +#' Natural Hazards and Earth System Sciences, 18(11), +#' 2825-2840. http://doi.org/10.5194/nhess-18-2825-2018 . +#' @param zclim a multi-dimensional array with named dimension containing at least one precipiation field with spatial dimensions. +#' @param lonin a vector indicating the longitudinal coordinates corresponding to the \code{zclim} parameter. +#' @param latin a vector indicating the latitudinal coordinates corresponding to the \code{zclim} parameter. +#' @param nf Refinement factor for downscaling (the output resolution is increased by this factor). +#' @param lon Vector of longitudes. +#' @param lat Vector of latitudes. +#' The number of longitudes and latitudes is expected to be even and the same. If not +#' the function will perform a subsetting to ensure this condition. +#' @param varname Name of the variable to be read from \code{climfile}. +#' @param fsmooth Logical to use smooth conservation (default) or large-scale box-average conservation. +#' @param lonname a character string indicating the name of the longitudinal dimension set as 'lon' by default. +#' @param latname a character string indicating the name of the latitudinal dimension set as 'lat' by default. +#' @param ncores an integer that indicates the number of cores for parallel computations using multiApply function. The default value is one. +#' +#' @return An object of class 's2dv_cube' containing in matrix \code{data} the weights with dimensions (lon, lat). +#' @import ncdf4 +#' @import rainfarmr +#' @import multiApply +#' @importFrom utils tail +#' @importFrom utils head +#' @examples +#' a <- array(1:2500, c(lat = 50, lon = 50)) +#' res <- RF_Weights(a, seq(0.1 ,5, 0.1), seq(0.1 ,5, 0.1), +#' nf = 5, lat = 1:5, lon = 1:5) +#' @export +RF_Weights <- function(zclim, latin, lonin, nf, lat, lon, fsmooth = TRUE, + lonname = 'lon', latname = 'lat', ncores = NULL) { + x <- Apply(list(zclim), target = c(lonname, latname), fun = rf_weights, + latin = latin, lonin = lonin, nf = nf, lat = lat, lon = lon, + fsmooth = fsmooth, ncores = ncores)$output1 + grid <- lon_lat_fine(lon, lat, nf) + return(list(data = x, lon = grid$lon, lat = grid$lat)) +} +rf_weights <- function(zclim, latin, lonin, nf, lat, lon, fsmooth = TRUE) { # Check if lon and lat need to be reversed if (lat[1] > lat[2]) { lat <- rev(lat) diff --git a/man/CST_RFWeights.Rd b/man/CST_RFWeights.Rd index ef5ebe4d5dd8585143c3fb78d2854e0783dbdd1d..acae8c6a512d775a3f0b0e892f4c121d12f907e0 100644 --- a/man/CST_RFWeights.Rd +++ b/man/CST_RFWeights.Rd @@ -4,7 +4,17 @@ \alias{CST_RFWeights} \title{Compute climatological weights for RainFARM stochastic precipitation downscaling} \usage{ -CST_RFWeights(climfile, nf, lon, lat, varname = "", fsmooth = TRUE) +CST_RFWeights( + climfile, + nf, + lon, + lat, + varname = NULL, + fsmooth = TRUE, + lonname = "lon", + latname = "lat", + ncores = NULL +) } \arguments{ \item{climfile}{Filename of a fine-scale precipitation climatology. @@ -15,7 +25,8 @@ Suitable climatology files could be for example a fine-scale precipitation clima from a high-resolution regional climate model (see e.g. Terzago et al. 2018), a local high-resolution gridded climatology from observations, or a reconstruction such as those which can be downloaded from the WORLDCLIM (http://www.worldclim.org) or CHELSA (http://chelsa-climate.org) -websites. The latter data will need to be converted to NetCDF format before being used (see for example the GDAL tools (https://www.gdal.org).} +websites. The latter data will need to be converted to NetCDF format before being used (see for example the GDAL tools (https://www.gdal.org). +It could also be a 's2dv_cube' object.} \item{nf}{Refinement factor for downscaling (the output resolution is increased by this factor).} @@ -28,9 +39,15 @@ the function will perform a subsetting to ensure this condition.} \item{varname}{Name of the variable to be read from \code{climfile}.} \item{fsmooth}{Logical to use smooth conservation (default) or large-scale box-average conservation.} + +\item{lonname}{a character string indicating the name of the longitudinal dimension set as 'lon' by default.} + +\item{latname}{a character string indicating the name of the latitudinal dimension set as 'lat' by default.} + +\item{ncores}{an integer that indicates the number of cores for parallel computations using multiApply function. The default value is one.} } \value{ -A matrix containing the weights with dimensions (lon, lat). +An object of class 's2dv_cube' containing in matrix \code{data} the weights with dimensions (lon, lat). } \description{ Compute climatological ("orographic") weights from a fine-scale precipitation climatology file. diff --git a/man/RF_Weights.Rd b/man/RF_Weights.Rd new file mode 100644 index 0000000000000000000000000000000000000000..7e8d735dedd35bc36e3d16a6eb0f3466099c9a04 --- /dev/null +++ b/man/RF_Weights.Rd @@ -0,0 +1,65 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/CST_RFWeights.R +\name{RF_Weights} +\alias{RF_Weights} +\title{Compute climatological weights for RainFARM stochastic precipitation downscaling} +\usage{ +RF_Weights( + zclim, + latin, + lonin, + nf, + lat, + lon, + fsmooth = TRUE, + lonname = "lon", + latname = "lat", + ncores = NULL +) +} +\arguments{ +\item{zclim}{a multi-dimensional array with named dimension containing at least one precipiation field with spatial dimensions.} + +\item{latin}{a vector indicating the latitudinal coordinates corresponding to the \code{zclim} parameter.} + +\item{lonin}{a vector indicating the longitudinal coordinates corresponding to the \code{zclim} parameter.} + +\item{nf}{Refinement factor for downscaling (the output resolution is increased by this factor).} + +\item{lat}{Vector of latitudes. +The number of longitudes and latitudes is expected to be even and the same. If not +the function will perform a subsetting to ensure this condition.} + +\item{lon}{Vector of longitudes.} + +\item{fsmooth}{Logical to use smooth conservation (default) or large-scale box-average conservation.} + +\item{lonname}{a character string indicating the name of the longitudinal dimension set as 'lon' by default.} + +\item{latname}{a character string indicating the name of the latitudinal dimension set as 'lat' by default.} + +\item{ncores}{an integer that indicates the number of cores for parallel computations using multiApply function. The default value is one.} + +\item{varname}{Name of the variable to be read from \code{climfile}.} +} +\value{ +An object of class 's2dv_cube' containing in matrix \code{data} the weights with dimensions (lon, lat). +} +\description{ +Compute climatological ("orographic") weights from a fine-scale precipitation climatology file. +} +\examples{ +a <- array(1:2500, c(lat = 50, lon = 50)) +res <- RF_Weights(a, seq(0.1 ,5, 0.1), seq(0.1 ,5, 0.1), + nf = 5, lat = 1:5, lon = 1:5) +} +\references{ +Terzago, S., Palazzi, E., & von Hardenberg, J. (2018). +Stochastic downscaling of precipitation in complex orography: +A simple method to reproduce a realistic fine-scale climatology. +Natural Hazards and Earth System Sciences, 18(11), +2825-2840. http://doi.org/10.5194/nhess-18-2825-2018 . +} +\author{ +Jost von Hardenberg - ISAC-CNR, \email{j.vonhardenberg@isac.cnr.it} +}