diff --git a/NAMESPACE b/NAMESPACE index b196f73d46805ea5225d923cd586a4865307a525..9357f631e08a05372bb0cc82671b401ca13a8397 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -1,5 +1,6 @@ # Generated by roxygen2: do not edit by hand +export(AMV) export(AnimateMap) export(Ano) export(Clim) @@ -19,6 +20,8 @@ export(ConfigShowSimilarEntries) export(ConfigShowTable) export(Corr) export(Eno) +export(GMST) +export(GSAT) export(InsertDim) export(LeapYear) export(Load) @@ -35,8 +38,10 @@ export(RMSSS) export(RandomWalkTest) export(Regression) export(Reorder) +export(SPOD) export(Season) export(Smoothing) +export(TPI) export(ToyModel) export(Trend) export(clim.colors) @@ -52,7 +57,9 @@ import(multiApply) import(ncdf4) import(parallel) import(plyr) +importFrom(ClimProjDiags,CombineIndices) importFrom(ClimProjDiags,Subset) +importFrom(ClimProjDiags,WeightedMean) importFrom(abind,abind) importFrom(abind,adrop) importFrom(grDevices,bmp) diff --git a/R/AMV.R b/R/AMV.R new file mode 100644 index 0000000000000000000000000000000000000000..6657f6fe409602a6a37bdb22a4a2fd6ceec58e45 --- /dev/null +++ b/R/AMV.R @@ -0,0 +1,242 @@ +#'Compute the Atlantic Multidecadal Variability (AMV) index +#' +#'The Atlantic Multidecadal Variability (AMV), also known as Atlantic +#'Multidecadal Oscillation (AMO), is a mode of natural variability of the sea +#'surface temperatures (SST) over the North Atlantic Ocean on multi-decadal +#'time scales. The AMV index is computed as the difference of weighted-averaged +#'SST anomalies over the North Atlantic region (0ºN-60ºN, 280ºE-360ºE) and the +#'weighted-averaged SST anomalies over 60ºS-60ºN, 0ºE-360ºE (Trenberth & +#'Dennis, 2005; Doblas-Reyes et al., 2013). +#' +#'@param data A numerical array to be used for the index computation with the +#' dimensions: 1) latitude, longitude, start date, forecast month, and member +#' (in case of decadal predictions), 2) latitude, longitude, year, month and +#' member (in case of historical simulations), or 3) latitude, longitude, year +#' and month (in case of observations or reanalyses). This data has to be +#' provided, at least, over the whole region needed to compute the index. +#'@param data_lats A numeric vector indicating the latitudes of the data. +#'@param data_lons A numeric vector indicating the longitudes of the data. +#'@param type A character string indicating the type of data ('dcpp' for +#' decadal predictions, 'hist' for historical simulations, or 'obs' for +#' observations or reanalyses). +#'@param lat_dim A character string of the name of the latitude dimension. The +#' default value is 'lat'. +#'@param lon_dim A character string of the name of the longitude dimension. The +#' default value is 'lon'. +#'@param mask An array of a mask (with 0's in the grid points that have to be +#' masked) or NULL (i.e., no mask is used). This parameter allows to remove +#' the values over land in case the dataset is a combination of surface air +#' temperature over land and sea surface temperature over the ocean. Also, it +#' can be used to mask those grid points that are missing in the observational +#' dataset for a fair comparison between the forecast system and the reference +#' dataset. The default value is NULL. +#'@param monini An integer indicating the month in which the forecast system is +#' initialized. Only used when parameter 'type' is 'dcpp'. The default value +#' is 11, i.e., initialized in November. +#'@param fmonth_dim A character string indicating the name of the forecast +#' month dimension. Only used if parameter 'type' is 'dcpp'. The default value +#' is 'fmonth'. +#'@param sdate_dim A character string indicating the name of the start date +#' dimension. Only used if parameter 'type' is 'dcpp'. The default value is +#' 'sdate'. +#'@param indices_for_clim A numeric vector of the indices of the years to +#' compute the climatology for calculating the anomalies, or NULL so the +#' climatology is calculated over the whole period. If the data are already +#' anomalies, set it to FALSE. The default value is NULL.\cr +#' In case of parameter 'type' is 'dcpp', 'indices_for_clim' must be relative +#' to the first forecast year, and the climatology is automatically computed +#' over the actual common period for the different forecast years. +#'@param year_dim A character string indicating the name of the year dimension +#' The default value is 'year'. Only used if parameter 'type' is 'hist' or +#' 'obs'. +#'@param month_dim A character string indicating the name of the month +#' dimension. The default value is 'month'. Only used if parameter 'type' is +#' 'hist' or 'obs'. +#'@param member_dim A character string indicating the name of the member +#' dimension. The default value is 'member'. Only used if parameter 'type' is +#' 'dcpp' or 'hist'. +#' +#'@return A numerical array of the AMV index with the dimensions of: +#' 1) sdate, forecast year, and member (in case of decadal predictions); +#' 2) year and member (in case of historical simulations); or +#' 3) year (in case of observations or reanalyses). +#' +#'@examples +#' ## Observations or reanalyses +#' obs <- array(1:100, dim = c(year = 5, lat = 19, lon = 37, month = 12)) +#' lat <- seq(-90, 90, 10) +#' lon <- seq(0, 360, 10) +#' index_obs <- AMV(data = obs, data_lats = lat, data_lons = lon, type = 'obs') +#' +#' ## Historical simulations +#' hist <- array(1:100, dim = c(year = 5, lat = 19, lon = 37, month = 12, member = 5)) +#' lat <- seq(-90, 90, 10) +#' lon <- seq(0, 360, 10) +#' index_hist <- AMV(data = hist, data_lats = lat, data_lons = lon, type = 'hist') +#' +#' ## Decadal predictions +#' dcpp <- array(1:100, dim = c(sdate = 5, lat = 19, lon = 37, fmonth = 24, member = 5)) +#' lat <- seq(-90, 90, 10) +#' lon <- seq(0, 360, 10) +#' index_dcpp <- AMV(data = dcpp, data_lats = lat, data_lons = lon, type = 'dcpp', monini = 1) +#' +#'@author Carlos Delgado-Torres, \email{carlos.delgado@bsc.es} +#'@author Roberto Bilbao, \email{roberto.bilbao@bsc.es} +#'@author Núria Pérez-Zanón, \email{nuria.perez@bsc.es} +#' +#'@importFrom ClimProjDiags WeightedMean +#'@import multiApply +#'@export +AMV <- function(data, data_lats, data_lons, type, lat_dim = 'lat', lon_dim = 'lon', + mask = NULL, monini = 11, fmonth_dim = 'fmonth', sdate_dim = 'sdate', + indices_for_clim = NULL, year_dim = 'year', month_dim = 'month', + member_dim = 'member') { + + ## Input Checks + # data + if (is.null(data)) { + stop("Parameter 'data' cannot be NULL.") + } + if (!is.numeric(data)) { + stop("Parameter 'data' must be a numeric array.") + } + # data_lats and data_lons part1 + if (!(class(data_lats) == 'numeric' | class(data_lats) == 'integer')) { + stop("Parameter 'data_lats' must be a numeric vector.") + } + if (!(class(data_lons) == 'numeric' | class(data_lons) == 'integer')) { + stop("Parameter 'data_lons' must be a numeric vector.") + } + # type + if (!type %in% c('dcpp', 'hist', 'obs')) { + stop("Parameter 'type' must be 'dcpp', 'hist', or 'obs'.") + } + # lat_dim + if (!(is.character(lat_dim) & length(lat_dim) == 1)) { + stop("Parameter 'lat_dim' must be a character string.") + } + if (!lat_dim %in% names(dim(data))) { + stop("Parameter 'lat_dim' is not found in 'data' dimension.") + } + # lon_dim + if (!(is.character(lon_dim) & length(lon_dim) == 1)) { + stop("Parameter 'lon_dim' must be a character string.") + } + if (!lon_dim %in% names(dim(data))) { + stop("Parameter 'lon_dim' is not found in 'data' dimension.") + } + # data_lats and data_lons part2 + if (dim(data)[lat_dim] != length(data_lats)){ + stop(paste0("The latitude dimension of parameter 'data' must be the same", + " length of parameter 'data_lats'.")) + } + if (dim(data)[lon_dim] != length(data_lons)){ + stop(paste0("The longitude dimension of parameter 'data' must be the same", + " length of parameter 'data_lons'.")) + } + # mask + if (!is.null(mask)) { + if (is.array(mask) & identical(names(dim(mask)), c(lat_dim,lon_dim)) & + identical(as.integer(dim(mask)), c(length(data_lats), length(data_lons)))) { + ## To mask those grid point that are missing in the observations + mask <- s2dv::Reorder(data = mask, order = c(lat_dim, lon_dim)) + fun_mask <- function(data, mask){ + data[mask == 0] <- NA + return(data) + } + data <- multiApply::Apply(data = data, target_dims = c(lat_dim, lon_dim), + fun = fun_mask, mask = mask)$output1 + } else { + stop(paste0("Parameter 'mask' must be NULL (no mask) or a numerical array ", + "with c(lat_dim, lon_dim) dimensions and 0 in those grid ", + "points that have to be masked.")) + } + } + # monini + if (type == 'dcpp') { + if (!is.numeric(monini) | monini %% 1 != 0 | monini < 1 | + monini > 12) { + stop("Parameter 'monini' must be an integer from 1 to 12.") + } + } + # fmonth_dim + if (type == 'dcpp') { + if (!(is.character(fmonth_dim) & length(fmonth_dim) == 1)) { + stop("Parameter 'fmonth_dim' must be a character string.") + } + if (!fmonth_dim %in% names(dim(data))) { + stop("Parameter 'fmonth_dim' is not found in 'data' dimension.") + } + } + # sdate_dim + if (type == 'dcpp') { + if (!(is.character(sdate_dim) & length(sdate_dim) == 1)) { + stop("Parameter 'sdate_dim' must be a character string.") + } + if (!sdate_dim %in% names(dim(data))) { + stop("Parameter 'sdate_dim' is not found in 'data' dimension.") + } + } + # indices_for_clim + if (!is.null(indices_for_clim)) { + if (!class(indices_for_clim) %in% c('numeric', 'integer') + & !(is.logical(indices_for_clim) & !any(indices_for_clim))) { + stop(paste0("The parameter 'indices_for_clim' must be a numeric vector ", + "or NULL to compute the anomalies based on the whole period, ", + "or FALSE if data are already anomalies")) + } + } + # year_dim + if (type == 'hist' | type == 'obs') { + if (!(is.character(year_dim) & length(year_dim) == 1)) { + stop("Parameter 'year_dim' must be a character string.") + } + if (!year_dim %in% names(dim(data))) { + stop("Parameter 'year_dim' is not found in 'data' dimension.") + } + } + # month_dim + if (type == 'hist' | type == 'obs') { + if (!(is.character(month_dim) & length(month_dim) == 1)) { + stop("Parameter 'month_dim' must be a character string.") + } + if (!month_dim %in% names(dim(data))) { + stop("Parameter 'month_dim' is not found in 'data' dimension.") + } + } + # member_dim + if (type == 'hist' | type == 'dcpp') { + if (!(is.character(member_dim) & length(member_dim) == 1)) { + stop("Parameter 'member_dim' must be a character string.") + } + if (!member_dim %in% names(dim(data))) { + stop("Parameter 'member_dim' is not found in 'data' dimension.") + } + } + + ## Regions for AMV (Doblas-Reyes et al., 2013) + lat_min_1 <- 0; lat_max_1 <- 60 + lon_min_1 <- 280; lon_max_1 <- 359.9 + lat_min_2 <- -60; lat_max_2 <- 60 + lon_min_2 <- 0; lon_max_2 <- 359.9 + regions <- NULL + regions$reg1 <- c(lon_min_1, lon_max_1, lat_min_1, lat_max_1) + regions$reg2 <- c(lon_min_2, lon_max_2, lat_min_2, lat_max_2) + + mean_1 <- ClimProjDiags::WeightedMean(data = data, lon = data_lons, lat = data_lats, + region = regions$reg1, + londim = which(names(dim(data)) == lon_dim), + latdim = which(names(dim(data)) == lat_dim)) + mean_2 <- ClimProjDiags::WeightedMean(data = data, lon = data_lons, lat = data_lats, + region = regions$reg2, + londim = which(names(dim(data)) == lon_dim), + latdim = which(names(dim(data)) == lat_dim)) + data <- ClimProjDiags::CombineIndices(indices = list(mean_1,mean_2), weights = NULL, + operation = 'subtract') # (mean_1 - mean_2) + + INDEX <- .Indices(data = data, type = type, monini = monini, + indices_for_clim = indices_for_clim, fmonth_dim = fmonth_dim, + sdate_dim = sdate_dim, year_dim = year_dim, + month_dim = month_dim, member_dim = member_dim) + return(INDEX) +} diff --git a/R/GMST.R b/R/GMST.R new file mode 100644 index 0000000000000000000000000000000000000000..5e6f4af1b83daca1d4f14dfb566e0af5f9b927bb --- /dev/null +++ b/R/GMST.R @@ -0,0 +1,280 @@ +#'Compute the Global Mean Surface Temperature (GMST) anomalies +#' +#'The Global Mean Surface Temperature (GMST) anomalies are computed as the +#'weighted-averaged surface air temperature anomalies over land and sea surface +#'temperature anomalies over the ocean. +#' +#'@param data_tas A numerical array indicating the surface air temperature data +#' to be used for the index computation with the dimensions: 1) latitude, +#' longitude, start date, forecast month, and member (in case of decadal +#' predictions), 2) latitude, longitude, year, month and member (in case of +#' historical simulations), or 3) latitude, longitude, year and month (in case +#' of observations or reanalyses). This data has to be provided, at least, +#' over the whole region needed to compute the index. The dimensions must be +#' identical to those of data_tos. +#'@param data_tos A numerical array indicating the sea surface temperature data +#' to be used for the index computation with the dimensions: 1) latitude, +#' longitude, start date, forecast month, and member (in case of decadal +#' predictions), 2) latitude, longitude, year, month and member (in case of +#' historical simulations), or 3) latitude, longitude, year and month (in case +#' of observations or reanalyses). This data has to be provided, at least, +#' over the whole region needed to compute the index. The dimensions must be +#' identical to those of data_tas. +#'@param data_lats A numeric vector indicating the latitudes of the data. +#'@param data_lons A numeric vector indicating the longitudes of the data. +#'@param mask_sea_land An array with dimensions [lat_dim = data_lats, lon_dim = +#' data_lons] for blending 'data_tas' and 'data_tos'. +#'@param sea_value A numeric value indicating the sea grid points in +#' 'mask_sea_land'. +#'@param type A character string indicating the type of data ('dcpp' for +#' decadal predictions, 'hist' for historical simulations, or 'obs' for +#' observations or reanalyses). +#'@param lat_dim A character string of the name of the latitude dimension. The +#' default value is 'lat'. +#'@param lon_dim A character string of the name of the longitude dimension. The +#' default value is 'lon'. +#'@param mask An array of a mask (with 0's in the grid points that have to be +#' masked) or NULL (i.e., no mask is used). This parameter allows to remove +#' the values over land in case the dataset is a combination of surface air +#' temperature over land and sea surface temperature over the ocean. Also, it +#' can be used to mask those grid points that are missing in the observational +#' dataset for a fair comparison between the forecast system and the reference +#' dataset. The default value is NULL. +#'@param monini An integer indicating the month in which the forecast system is +#' initialized. Only used when parameter 'type' is 'dcpp'. The default value +#' is 11, i.e., initialized in November. +#'@param fmonth_dim A character string indicating the name of the forecast +#' month dimension. Only used if parameter 'type' is 'dcpp'. The default value +#' is 'fmonth'. +#'@param sdate_dim A character string indicating the name of the start date +#' dimension. Only used if parameter 'type' is 'dcpp'. The default value is +#' 'sdate'. +#'@param indices_for_clim A numeric vector of the indices of the years to +#' compute the climatology for calculating the anomalies, or NULL so the +#' climatology is calculated over the whole period. If the data are already +#' anomalies, set it to FALSE. The default value is NULL.\cr +#' In case of parameter 'type' is 'dcpp', 'indices_for_clim' must be relative +#' to the first forecast year, and the climatology is automatically computed +#' over the actual common period for the different forecast years. +#'@param year_dim A character string indicating the name of the year dimension +#' The default value is 'year'. Only used if parameter 'type' is 'hist' or +#' 'obs'. +#'@param month_dim A character string indicating the name of the month +#' dimension. The default value is 'month'. Only used if parameter 'type' is +#' 'hist' or 'obs'. +#'@param member_dim A character string indicating the name of the member +#' dimension. The default value is 'member'. Only used if parameter 'type' is +#' 'dcpp' or 'hist'. +#' +#'@return A numerical array of the GMST anomalies with the dimensions of: +#' 1) sdate, forecast year, and member (in case of decadal predictions); +#' 2) year and member (in case of historical simulations); or +#' 3) year (in case of observations or reanalyses). +#' +#'@examples +#' ## Observations or reanalyses +#' obs_tas <- array(1:100, dim = c(year = 5, lat = 19, lon = 37, month = 12)) +#' obs_tos <- array(2:101, dim = c(year = 5, lat = 19, lon = 37, month = 12)) +#' mask_sea_land <- array(c(1,0,1), dim = c(lat = 19, lon = 37)) +#' sea_value <- 1 +#' lat <- seq(-90, 90, 10) +#' lon <- seq(0, 360, 10) +#' index_obs <- GMST(data_tas = obs_tas, data_tos = obs_tos, data_lats = lat, data_lons = lon, type = 'obs', mask_sea_land = mask_sea_land, sea_value = sea_value) +#' +#' ## Historical simulations +#' hist_tas <- array(1:100, dim = c(year = 5, lat = 19, lon = 37, month = 12, member = 5)) +#' hist_tos <- array(2:101, dim = c(year = 5, lat = 19, lon = 37, month = 12, member = 5)) +#' mask_sea_land <- array(c(1,0,1), dim = c(lat = 19, lon = 37)) +#' sea_value <- 1 +#' lat <- seq(-90, 90, 10) +#' lon <- seq(0, 360, 10) +#' index_hist <- GMST(data_tas = hist_tas, data_tos = hist_tos, data_lats = lat, data_lons = lon, type = 'hist', mask_sea_land = mask_sea_land, sea_value = sea_value) +#' +#' ## Decadal predictions +#' dcpp_tas <- array(1:100, dim = c(sdate = 5, lat = 19, lon = 37, fmonth = 24, member = 5)) +#' dcpp_tos <- array(2:101, dim = c(sdate = 5, lat = 19, lon = 37, fmonth = 24, member = 5)) +#' mask_sea_land <- array(c(1,0,1), dim = c(lat = 19, lon = 37)) +#' sea_value <- 1 +#' lat <- seq(-90, 90, 10) +#' lon <- seq(0, 360, 10) +#' index_dcpp <- GMST(data_tas = dcpp_tas, data_tos = dcpp_tos, data_lats = lat, data_lons = lon, type = 'dcpp', monini = 1, mask_sea_land = mask_sea_land, sea_value = sea_value) +#' +#'@author Carlos Delgado-Torres, \email{carlos.delgado@bsc.es} +#'@author Roberto Bilbao, \email{roberto.bilbao@bsc.es} +#'@author Núria Pérez-Zanón, \email{nuria.perez@bsc.es} +#' +#'@importFrom ClimProjDiags WeightedMean +#'@import multiApply +#'@export +GMST <- function(data_tas, data_tos, data_lats, data_lons, mask_sea_land, sea_value, + type, mask = NULL, lat_dim = 'lat', lon_dim = 'lon', monini = 11, + fmonth_dim = 'fmonth', sdate_dim = 'sdate', indices_for_clim = NULL, + year_dim = 'year', month_dim = 'month', member_dim = 'member') { + + ## Input Checks + # data_tas and data_tos + if (is.null(data_tas) | is.null(data_tos)) { + stop("Parameter 'data_tas' and 'data_tos' cannot be NULL.") + } + if (!is.numeric(data_tas) | !is.numeric(data_tos)) { + stop("Parameter 'data_tas' and 'data_tos' must be a numeric array.") + } + if (!identical(dim(data_tas), dim(data_tos))) { + stop("The dimension of data_tas and data_tos must be identical.") + } + # data_lats and data_lons part1 + if (!(class(data_lats) == 'numeric' | class(data_lats) == 'integer')) { + stop("Parameter 'data_lats' must be a numeric vector.") + } + if (!(class(data_lons) == 'numeric' | class(data_lons) == 'integer')) { + stop("Parameter 'data_lons' must be a numeric vector.") + } + # lat_dim + if (!(is.character(lat_dim) & length(lat_dim) == 1)) { + stop("Parameter 'lat_dim' must be a character string.") + } + if (!lat_dim %in% names(dim(data_tas)) | !lat_dim %in% names(dim(data_tos))) { + stop("Parameter 'lat_dim' is not found in 'data_tas' or 'data_tos' dimension.") + } + # lon_dim + if (!(is.character(lon_dim) & length(lon_dim) == 1)) { + stop("Parameter 'lon_dim' must be a character string.") + } + if (!lon_dim %in% names(dim(data_tas)) | !lon_dim %in% names(dim(data_tos))) { + stop("Parameter 'lon_dim' is not found in 'data_tas' or 'data_tos' dimension.") + } + # data_lats and data_lons part2 + if (dim(data_tas)[lat_dim] != length(data_lats) | + dim(data_tos)[lat_dim] != length(data_lats)) { + stop(paste0("The latitude dimension of parameter 'data_tas' and 'data_tos'", + " must be the same length of parameter 'data_lats'.")) + } + if (dim(data_tas)[lon_dim] != length(data_lons) | + dim(data_tos)[lon_dim] != length(data_lons)) { + stop(paste0("The longitude dimension of parameter 'data_tas' and 'data_tos'", + " must be the same length of parameter 'data_lons'.")) + } + # sea_value + if (!is.numeric(sea_value) | length(sea_value) != 1) { + stop("Parameter 'sea_value' must be a numeric value.") + } + # mask_sea_land + if (!is.array(mask_sea_land)) { + stop("Parameter 'mask_sea_land' must be an array with dimensions [lat_dim, lon_dim].") + } else if (!identical(names(dim(mask_sea_land)), c(lat_dim, lon_dim))) { + stop("Parameter 'mask_sea_land' must be an array with dimensions [lat_dim, lon_dim].") + } else if (!identical(as.integer(dim(mask_sea_land)), + c(length(data_lats), length(data_lons)))) { + stop("Parameter 'mask_sea_land' dimensions must be equal to the length of 'data_lats' and 'data_lons'.") + } + # type + if (!type %in% c('dcpp', 'hist', 'obs')) { + stop("Parameter 'type' must be 'dcpp', 'hist', or 'obs'.") + } + # mask + if (!is.null(mask)) { + if (!is.array(mask) | !identical(names(dim(mask)), c(lat_dim,lon_dim)) | + !identical(as.integer(dim(mask)), c(length(data_lats), length(data_lons)))) { + stop(paste0("Parameter 'mask' must be NULL (no mask) or a numerical array ", + "with c(lat_dim, lon_dim) dimensions and 0 in those grid ", + "points that have to be masked.")) + } + } + # monini + if (type == 'dcpp') { + if (!is.numeric(monini) | monini %% 1 != 0 | monini < 1 | + monini > 12) { + stop("Parameter 'monini' must be an integer from 1 to 12.") + } + } + # fmonth_dim + if (type == 'dcpp') { + if (!(is.character(fmonth_dim) & length(fmonth_dim) == 1)) { + stop("Parameter 'fmonth_dim' must be a character string.") + } + if (!fmonth_dim %in% names(dim(data_tas)) | !fmonth_dim %in% names(dim(data_tos))) { + stop("Parameter 'fmonth_dim' is not found in 'data_tas' or 'data_tos' dimension.") + } + } + # sdate_dim + if (type == 'dcpp') { + if (!(is.character(sdate_dim) & length(sdate_dim) == 1)) { + stop("Parameter 'sdate_dim' must be a character string.") + } + if (!sdate_dim %in% names(dim(data_tas)) | !sdate_dim %in% names(dim(data_tos))) { + stop("Parameter 'sdate_dim' is not found in 'data_tas' or 'data_tos' dimension.") + } + } + # indices_for_clim + if (!is.null(indices_for_clim)) { + if (!class(indices_for_clim) %in% c('numeric', 'integer') + & !(is.logical(indices_for_clim) & !any(indices_for_clim))) { + stop(paste0("The parameter 'indices_for_clim' must be a numeric vector ", + "or NULL to compute the anomalies based on the whole period, ", + "or FALSE if data are already anomalies")) + } + } + # year_dim + if (type == 'hist' | type == 'obs') { + if (!(is.character(year_dim) & length(year_dim) == 1)) { + stop("Parameter 'year_dim' must be a character string.") + } + if (!year_dim %in% names(dim(data_tas)) | !year_dim %in% names(dim(data_tos))) { + stop("Parameter 'year_dim' is not found in 'data_tas' or 'data_tos' dimension.") + } + } + # month_dim + if (type == 'hist' | type == 'obs') { + if (!(is.character(month_dim) & length(month_dim) == 1)) { + stop("Parameter 'month_dim' must be a character string.") + } + if (!month_dim %in% names(dim(data_tas)) | !month_dim %in% names(dim(data_tos))) { + stop("Parameter 'month_dim' is not found in 'data_tas' or 'data_tos' dimension.") + } + } + # member_dim + if (type == 'hist' | type == 'dcpp') { + if (!(is.character(member_dim) & length(member_dim) == 1)) { + stop("Parameter 'member_dim' must be a character string.") + } + if (!member_dim %in% names(dim(data_tas)) | !member_dim %in% names(dim(data_tos))) { + stop("Parameter 'member_dim' is not found in 'data_tas' or 'data_tos' dimension.") + } + } + + ## combination of tas and tos (data) + mask_tas_tos <- function(data_tas, data_tos, mask_sea_land, sea_value) { + data <- data_tas + data[mask_sea_land == sea_value] <- data_tos[mask_sea_land == sea_value] + return(data) + } + mask_sea_land <- s2dv::Reorder(data = mask_sea_land, order = c(lat_dim,lon_dim)) + data <- multiApply::Apply(data = list(data_tas, data_tos), + target_dims = c(lat_dim, lon_dim), + fun = mask_tas_tos, mask_sea_land = mask_sea_land, + sea_value = sea_value)$output1 + data <- drop(data) + rm(data_tas, data_tos) + + ## To mask those grid point that are missing in the observations + if (!is.null(mask)) { + mask <- s2dv::Reorder(data = mask, order = c(lat_dim, lon_dim)) + fun_mask <- function(data, mask) { + data[mask == 0] <- NA + return(data) + } + data <- multiApply::Apply(data = data, target_dims = c(lat_dim, lon_dim), + fun = fun_mask, mask = mask)$output1 + } + + data <- ClimProjDiags::WeightedMean(data = data, lon = data_lons, lat = data_lats, + region = NULL, + londim = which(names(dim(data)) == lon_dim), + latdim = which(names(dim(data)) == lat_dim)) + + INDEX <- .Indices(data = data, type = type, monini = monini, + indices_for_clim = indices_for_clim, fmonth_dim = fmonth_dim, + sdate_dim = sdate_dim, year_dim = year_dim, + month_dim = month_dim, member_dim = member_dim) + return(INDEX) +} diff --git a/R/GSAT.R b/R/GSAT.R new file mode 100644 index 0000000000000000000000000000000000000000..b9d0f6f3de7f280a574b6979afd02f1c1f55042e --- /dev/null +++ b/R/GSAT.R @@ -0,0 +1,222 @@ +#'Compute the Global Surface Air Temperature (GSAT) anomalies +#' +#'The Global Surface Air Temperature (GSAT) anomalies are computed as the +#'weighted-averaged surface air temperature anomalies over the global region. +#' +#'@param data A numerical array to be used for the index computation with the +#' dimensions: 1) latitude, longitude, start date, forecast month, and member +#' (in case of decadal predictions), 2) latitude, longitude, year, month and +#' member (in case of historical simulations), or 3) latitude, longitude, year +#' and month (in case of observations or reanalyses). This data has to be +#' provided, at least, over the whole region needed to compute the index. +#'@param data_lats A numeric vector indicating the latitudes of the data. +#'@param data_lons A numeric vector indicating the longitudes of the data. +#'@param type A character string indicating the type of data ('dcpp' for +#' decadal predictions, 'hist' for historical simulations, or 'obs' for +#' observations or reanalyses). +#'@param lat_dim A character string of the name of the latitude dimension. The +#' default value is 'lat'. +#'@param lon_dim A character string of the name of the longitude dimension. The +#' default value is 'lon'. +#'@param mask An array of a mask (with 0's in the grid points that have to be +#' masked) or NULL (i.e., no mask is used). This parameter allows to remove +#' the values over land in case the dataset is a combination of surface air +#' temperature over land and sea surface temperature over the ocean. Also, it +#' can be used to mask those grid points that are missing in the observational +#' dataset for a fair comparison between the forecast system and the reference +#' dataset. The default value is NULL. +#'@param monini An integer indicating the month in which the forecast system is +#' initialized. Only used when parameter 'type' is 'dcpp'. The default value +#' is 11, i.e., initialized in November. +#'@param fmonth_dim A character string indicating the name of the forecast +#' month dimension. Only used if parameter 'type' is 'dcpp'. The default value +#' is 'fmonth'. +#'@param sdate_dim A character string indicating the name of the start date +#' dimension. Only used if parameter 'type' is 'dcpp'. The default value is +#' 'sdate'. +#'@param indices_for_clim A numeric vector of the indices of the years to +#' compute the climatology for calculating the anomalies, or NULL so the +#' climatology is calculated over the whole period. If the data are already +#' anomalies, set it to FALSE. The default value is NULL.\cr +#' In case of parameter 'type' is 'dcpp', 'indices_for_clim' must be relative +#' to the first forecast year, and the climatology is automatically computed +#' over the actual common period for the different forecast years. +#'@param year_dim A character string indicating the name of the year dimension +#' The default value is 'year'. Only used if parameter 'type' is 'hist' or +#' 'obs'. +#'@param month_dim A character string indicating the name of the month +#' dimension. The default value is 'month'. Only used if parameter 'type' is +#' 'hist' or 'obs'. +#'@param member_dim A character string indicating the name of the member +#' dimension. The default value is 'member'. Only used if parameter 'type' is +#' 'dcpp' or 'hist'. +#' +#'@return A numerical array of the GSAT anomalies with the dimensions of: +#' 1) sdate, forecast year, and member (in case of decadal predictions); +#' 2) year and member (in case of historical simulations); or +#' 3) year (in case of observations or reanalyses). +#' +#'@examples +#' ## Observations or reanalyses +#' obs <- array(1:100, dim = c(year = 5, lat = 19, lon = 37, month = 12)) +#' lat <- seq(-90, 90, 10) +#' lon <- seq(0, 360, 10) +#' index_obs <- GSAT(data = obs, data_lats = lat, data_lons = lon, type = 'obs') +#' +#' ## Historical simulations +#' hist <- array(1:100, dim = c(year = 5, lat = 19, lon = 37, month = 12, member = 5)) +#' lat <- seq(-90, 90, 10) +#' lon <- seq(0, 360, 10) +#' index_hist <- GSAT(data = hist, data_lats = lat, data_lons = lon, type = 'hist') +#' +#' ## Decadal predictions +#' dcpp <- array(1:100, dim = c(sdate = 5, lat = 19, lon = 37, fmonth = 24, member = 5)) +#' lat <- seq(-90, 90, 10) +#' lon <- seq(0, 360, 10) +#' index_dcpp <- GSAT(data = dcpp, data_lats = lat, data_lons = lon, type = 'dcpp', monini = 1) +#' +#'@author Carlos Delgado-Torres, \email{carlos.delgado@bsc.es} +#'@author Núria Pérez-Zanón, \email{nuria.perez@bsc.es} +#'@author Roberto Bilbao, \email{roberto.bilbao@bsc.es} +#' +#'@importFrom ClimProjDiags WeightedMean +#'@import multiApply +#'@export +GSAT <- function(data, data_lats, data_lons, type, lat_dim = 'lat', lon_dim = 'lon', + mask = NULL, monini = 11, fmonth_dim = 'fmonth', sdate_dim = 'sdate', + indices_for_clim = NULL, year_dim = 'year', month_dim = 'month', + member_dim = 'member') { + + ## Input Checks + # data + if (is.null(data)) { + stop("Parameter 'data' cannot be NULL.") + } + if (!is.numeric(data)) { + stop("Parameter 'data' must be a numeric array.") + } + # data_lats and data_lons part1 + if (!(class(data_lats) == 'numeric' | class(data_lats) == 'integer')) { + stop("Parameter 'data_lats' must be a numeric vector.") + } + if (!(class(data_lons) == 'numeric' | class(data_lons) == 'integer')) { + stop("Parameter 'data_lons' must be a numeric vector.") + } + # type + if (!type %in% c('dcpp', 'hist', 'obs')) { + stop("Parameter 'type' must be 'dcpp', 'hist', or 'obs'.") + } + # lat_dim + if (!(is.character(lat_dim) & length(lat_dim) == 1)) { + stop("Parameter 'lat_dim' must be a character string.") + } + if (!lat_dim %in% names(dim(data))) { + stop("Parameter 'lat_dim' is not found in 'data' dimension.") + } + # lon_dim + if (!(is.character(lon_dim) & length(lon_dim) == 1)) { + stop("Parameter 'lon_dim' must be a character string.") + } + if (!lon_dim %in% names(dim(data))) { + stop("Parameter 'lon_dim' is not found in 'data' dimension.") + } + # data_lats and data_lons part2 + if (dim(data)[lat_dim] != length(data_lats)){ + stop(paste0("The latitude dimension of parameter 'data' must be the same", + " length of parameter 'data_lats'.")) + } + if (dim(data)[lon_dim] != length(data_lons)){ + stop(paste0("The longitude dimension of parameter 'data' must be the same", + " length of parameter 'data_lons'.")) + } + # mask + if (!is.null(mask)) { + if (is.array(mask) & identical(names(dim(mask)), c(lat_dim,lon_dim)) & + identical(as.integer(dim(mask)), c(length(data_lats), length(data_lons)))) { + ## To mask those grid point that are missing in the observations + mask <- s2dv::Reorder(data = mask, order = c(lat_dim, lon_dim)) + fun_mask <- function(data, mask){ + data[mask == 0] <- NA + return(data) + } + data <- multiApply::Apply(data = data, target_dims = c(lat_dim, lon_dim), + fun = fun_mask, mask = mask)$output1 + } else { + stop(paste0("Parameter 'mask' must be NULL (no mask) or a numerical array ", + "with c(lat_dim, lon_dim) dimensions and 0 in those grid ", + "points that have to be masked.")) + } + } + # monini + if (type == 'dcpp') { + if (!is.numeric(monini) | monini %% 1 != 0 | monini < 1 | + monini > 12) { + stop("Parameter 'monini' must be an integer from 1 to 12.") + } + } + # fmonth_dim + if (type == 'dcpp') { + if (!(is.character(fmonth_dim) & length(fmonth_dim) == 1)) { + stop("Parameter 'fmonth_dim' must be a character string.") + } + if (!fmonth_dim %in% names(dim(data))) { + stop("Parameter 'fmonth_dim' is not found in 'data' dimension.") + } + } + # sdate_dim + if (type == 'dcpp') { + if (!(is.character(sdate_dim) & length(sdate_dim) == 1)) { + stop("Parameter 'sdate_dim' must be a character string.") + } + if (!sdate_dim %in% names(dim(data))) { + stop("Parameter 'sdate_dim' is not found in 'data' dimension.") + } + } + # indices_for_clim + if (!is.null(indices_for_clim)) { + if (!class(indices_for_clim) %in% c('numeric', 'integer') + & !(is.logical(indices_for_clim) & !any(indices_for_clim))) { + stop(paste0("The parameter 'indices_for_clim' must be a numeric vector ", + "or NULL to compute the anomalies based on the whole period, ", + "or FALSE if data are already anomalies")) + } + } + # year_dim + if (type == 'hist' | type == 'obs') { + if (!(is.character(year_dim) & length(year_dim) == 1)) { + stop("Parameter 'year_dim' must be a character string.") + } + if (!year_dim %in% names(dim(data))) { + stop("Parameter 'year_dim' is not found in 'data' dimension.") + } + } + # month_dim + if (type == 'hist' | type == 'obs') { + if (!(is.character(month_dim) & length(month_dim) == 1)) { + stop("Parameter 'month_dim' must be a character string.") + } + if (!month_dim %in% names(dim(data))) { + stop("Parameter 'month_dim' is not found in 'data' dimension.") + } + } + # member_dim + if (type == 'hist' | type == 'dcpp') { + if (!(is.character(member_dim) & length(member_dim) == 1)) { + stop("Parameter 'member_dim' must be a character string.") + } + if (!member_dim %in% names(dim(data))) { + stop("Parameter 'member_dim' is not found in 'data' dimension.") + } + } + + data <- ClimProjDiags::WeightedMean(data = data, lon = data_lons, lat = data_lats, + region = NULL, + londim = which(names(dim(data)) == lon_dim), + latdim = which(names(dim(data)) == lat_dim)) + + INDEX <- .Indices(data = data, type = type, monini = monini, + indices_for_clim = indices_for_clim, fmonth_dim = fmonth_dim, + sdate_dim = sdate_dim, year_dim = year_dim, + month_dim = month_dim, member_dim = member_dim) + return(INDEX) +} diff --git a/R/SPOD.R b/R/SPOD.R new file mode 100644 index 0000000000000000000000000000000000000000..d939df79a8ecd67cb99fa450c68a18e8b77b5b0b --- /dev/null +++ b/R/SPOD.R @@ -0,0 +1,241 @@ +#'Compute the South Pacific Ocean Dipole (SPOD) index +#' +#'The South Pacific Ocean Dipole (SPOD) index is related to the El +#'Nino-Southern Oscillation (ENSO) and the Inderdecadal Pacific Oscillation +#'(IPO). The SPOD index is computed as the difference of weighted-averaged SST +#'anomalies over 20ºS-48ºS, 165ºE-190ºE (NW pole) and the weighted-averaged SST +#' anomalies over 44ºS-65ºS, 220ºE-260ºE (SE pole) (Saurral et al., 2020). +#' +#'@param data A numerical array to be used for the index computation with the +#' dimensions: 1) latitude, longitude, start date, forecast month, and member +#' (in case of decadal predictions), 2) latitude, longitude, year, month and +#' member (in case of historical simulations), or 3) latitude, longitude, year +#' and month (in case of observations or reanalyses). This data has to be +#' provided, at least, over the whole region needed to compute the index. +#'@param data_lats A numeric vector indicating the latitudes of the data. +#'@param data_lons A numeric vector indicating the longitudes of the data. +#'@param type A character string indicating the type of data ('dcpp' for +#' decadal predictions, 'hist' for historical simulations, or 'obs' for +#' observations or reanalyses). +#'@param lat_dim A character string of the name of the latitude dimension. The +#' default value is 'lat'. +#'@param lon_dim A character string of the name of the longitude dimension. The +#' default value is 'lon'. +#'@param mask An array of a mask (with 0's in the grid points that have to be +#' masked) or NULL (i.e., no mask is used). This parameter allows to remove +#' the values over land in case the dataset is a combination of surface air +#' temperature over land and sea surface temperature over the ocean. Also, it +#' can be used to mask those grid points that are missing in the observational +#' dataset for a fair comparison between the forecast system and the reference +#' dataset. The default value is NULL. +#'@param monini An integer indicating the month in which the forecast system is +#' initialized. Only used when parameter 'type' is 'dcpp'. The default value +#' is 11, i.e., initialized in November. +#'@param fmonth_dim A character string indicating the name of the forecast +#' month dimension. Only used if parameter 'type' is 'dcpp'. The default value +#' is 'fmonth'. +#'@param sdate_dim A character string indicating the name of the start date +#' dimension. Only used if parameter 'type' is 'dcpp'. The default value is +#' 'sdate'. +#'@param indices_for_clim A numeric vector of the indices of the years to +#' compute the climatology for calculating the anomalies, or NULL so the +#' climatology is calculated over the whole period. If the data are already +#' anomalies, set it to FALSE. The default value is NULL.\cr +#' In case of parameter 'type' is 'dcpp', 'indices_for_clim' must be relative +#' to the first forecast year, and the climatology is automatically computed +#' over the actual common period for the different forecast years. +#'@param year_dim A character string indicating the name of the year dimension +#' The default value is 'year'. Only used if parameter 'type' is 'hist' or +#' 'obs'. +#'@param month_dim A character string indicating the name of the month +#' dimension. The default value is 'month'. Only used if parameter 'type' is +#' 'hist' or 'obs'. +#'@param member_dim A character string indicating the name of the member +#' dimension. The default value is 'member'. Only used if parameter 'type' is +#' 'dcpp' or 'hist'. +#' +#'@return A numerical array of the SPOD index with the dimensions of: +#' 1) sdate, forecast year, and member (in case of decadal predictions); +#' 2) year and member (in case of historical simulations); or +#' 3) year (in case of observations or reanalyses). +#' +#'@examples +#' ## Observations or reanalyses +#' obs <- array(1:100, dim = c(year = 5, lat = 19, lon = 37, month = 12)) +#' lat <- seq(-90, 90, 10) +#' lon <- seq(0, 360, 10) +#' index_obs <- SPOD(data = obs, data_lats = lat, data_lons = lon, type = 'obs') +#' +#' ## Historical simulations +#' hist <- array(1:100, dim = c(year = 5, lat = 19, lon = 37, month = 12, member = 5)) +#' lat <- seq(-90, 90, 10) +#' lon <- seq(0, 360, 10) +#' index_hist <- SPOD(data = hist, data_lats = lat, data_lons = lon, type = 'hist') +#' +#' ## Decadal predictions +#' dcpp <- array(1:100, dim = c(sdate = 5, lat = 19, lon = 37, fmonth = 24, member = 5)) +#' lat <- seq(-90, 90, 10) +#' lon <- seq(0, 360, 10) +#' index_dcpp <- SPOD(data = dcpp, data_lats = lat, data_lons = lon, type = 'dcpp', monini = 1) +#' +#'@author Carlos Delgado-Torres, \email{carlos.delgado@bsc.es} +#'@author Núria Pérez-Zanón, \email{nuria.perez@bsc.es} +#'@author Roberto Bilbao, \email{roberto.bilbao@bsc.es} +#' +#'@importFrom ClimProjDiags WeightedMean CombineIndices +#'@import multiApply +#'@export +SPOD <- function(data, data_lats, data_lons, type, lat_dim = 'lat', lon_dim = 'lon', + mask = NULL, monini = 11, fmonth_dim = 'fmonth', sdate_dim = 'sdate', + indices_for_clim = NULL, year_dim = 'year', month_dim = 'month', + member_dim = 'member') { + + ## Input Checks + # data + if (is.null(data)) { + stop("Parameter 'data' cannot be NULL.") + } + if (!is.numeric(data)) { + stop("Parameter 'data' must be a numeric array.") + } + # data_lats and data_lons part1 + if (!(class(data_lats) == 'numeric' | class(data_lats) == 'integer')) { + stop("Parameter 'data_lats' must be a numeric vector.") + } + if (!(class(data_lons) == 'numeric' | class(data_lons) == 'integer')) { + stop("Parameter 'data_lons' must be a numeric vector.") + } + # type + if (!type %in% c('dcpp', 'hist', 'obs')) { + stop("Parameter 'type' must be 'dcpp', 'hist', or 'obs'.") + } + # lat_dim + if (!(is.character(lat_dim) & length(lat_dim) == 1)) { + stop("Parameter 'lat_dim' must be a character string.") + } + if (!lat_dim %in% names(dim(data))) { + stop("Parameter 'lat_dim' is not found in 'data' dimension.") + } + # lon_dim + if (!(is.character(lon_dim) & length(lon_dim) == 1)) { + stop("Parameter 'lon_dim' must be a character string.") + } + if (!lon_dim %in% names(dim(data))) { + stop("Parameter 'lon_dim' is not found in 'data' dimension.") + } + # data_lats and data_lons part2 + if (dim(data)[lat_dim] != length(data_lats)){ + stop(paste0("The latitude dimension of parameter 'data' must be the same", + " length of parameter 'data_lats'.")) + } + if (dim(data)[lon_dim] != length(data_lons)){ + stop(paste0("The longitude dimension of parameter 'data' must be the same", + " length of parameter 'data_lons'.")) + } + # mask + if (!is.null(mask)) { + if (is.array(mask) & identical(names(dim(mask)), c(lat_dim,lon_dim)) & + identical(as.integer(dim(mask)), c(length(data_lats), length(data_lons)))) { + ## To mask those grid point that are missing in the observations + mask <- s2dv::Reorder(data = mask, order = c(lat_dim, lon_dim)) + fun_mask <- function(data, mask){ + data[mask == 0] <- NA + return(data) + } + data <- multiApply::Apply(data = data, target_dims = c(lat_dim, lon_dim), + fun = fun_mask, mask = mask)$output1 + } else { + stop(paste0("Parameter 'mask' must be NULL (no mask) or a numerical array ", + "with c(lat_dim, lon_dim) dimensions and 0 in those grid ", + "points that have to be masked.")) + } + } + # monini + if (type == 'dcpp') { + if (!is.numeric(monini) | monini %% 1 != 0 | monini < 1 | + monini > 12) { + stop("Parameter 'monini' must be an integer from 1 to 12.") + } + } + # fmonth_dim + if (type == 'dcpp') { + if (!(is.character(fmonth_dim) & length(fmonth_dim) == 1)) { + stop("Parameter 'fmonth_dim' must be a character string.") + } + if (!fmonth_dim %in% names(dim(data))) { + stop("Parameter 'fmonth_dim' is not found in 'data' dimension.") + } + } + # sdate_dim + if (type == 'dcpp') { + if (!(is.character(sdate_dim) & length(sdate_dim) == 1)) { + stop("Parameter 'sdate_dim' must be a character string.") + } + if (!sdate_dim %in% names(dim(data))) { + stop("Parameter 'sdate_dim' is not found in 'data' dimension.") + } + } + # indices_for_clim + if (!is.null(indices_for_clim)) { + if (!class(indices_for_clim) %in% c('numeric', 'integer') + & !(is.logical(indices_for_clim) & !any(indices_for_clim))) { + stop(paste0("The parameter 'indices_for_clim' must be a numeric vector ", + "or NULL to compute the anomalies based on the whole period, ", + "or FALSE if data are already anomalies")) + } + } + # year_dim + if (type == 'hist' | type == 'obs') { + if (!(is.character(year_dim) & length(year_dim) == 1)) { + stop("Parameter 'year_dim' must be a character string.") + } + if (!year_dim %in% names(dim(data))) { + stop("Parameter 'year_dim' is not found in 'data' dimension.") + } + } + # month_dim + if (type == 'hist' | type == 'obs') { + if (!(is.character(month_dim) & length(month_dim) == 1)) { + stop("Parameter 'month_dim' must be a character string.") + } + if (!month_dim %in% names(dim(data))) { + stop("Parameter 'month_dim' is not found in 'data' dimension.") + } + } + # member_dim + if (type == 'hist' | type == 'dcpp') { + if (!(is.character(member_dim) & length(member_dim) == 1)) { + stop("Parameter 'member_dim' must be a character string.") + } + if (!member_dim %in% names(dim(data))) { + stop("Parameter 'member_dim' is not found in 'data' dimension.") + } + } + + ## Regions for IPO_SPOD (Saurral et al., 2020) + lat_min_1 <- -48; lat_max_1 <- -20 + lon_min_1 <- 165; lon_max_1 <- 190 + lat_min_2 <- -65; lat_max_2 <- -44 + lon_min_2 <- 220; lon_max_2 <- 260 + regions <- NULL + regions$reg1 <- c(lon_min_1, lon_max_1, lat_min_1, lat_max_1) + regions$reg2 <- c(lon_min_2, lon_max_2, lat_min_2, lat_max_2) + + mean_1 <- ClimProjDiags::WeightedMean(data = data, lon = data_lons, lat = data_lats, + region = regions$reg1, + londim = which(names(dim(data)) == lon_dim), + latdim = which(names(dim(data)) == lat_dim)) + mean_2 <- ClimProjDiags::WeightedMean(data = data, lon = data_lons, lat = data_lats, + region = regions$reg2, + londim = which(names(dim(data)) == lon_dim), + latdim = which(names(dim(data)) == lat_dim)) + + data <- ClimProjDiags::CombineIndices(indices = list(mean_1,mean_2), + weights = NULL, operation = 'subtract') # (mean_1 - mean_2) + + INDEX <- .Indices(data = data, type = type, monini = monini, + indices_for_clim = indices_for_clim, fmonth_dim = fmonth_dim, + sdate_dim = sdate_dim, year_dim = year_dim, + month_dim = month_dim, member_dim = member_dim) + return(INDEX) +} diff --git a/R/TPI.R b/R/TPI.R new file mode 100644 index 0000000000000000000000000000000000000000..6733bdd718432a59aec393fe129b5e87656785de --- /dev/null +++ b/R/TPI.R @@ -0,0 +1,250 @@ +#'Compute the Tripole Index (TPI) for the Interdecadal Pacific Oscillation (IPO) +#' +#'The Tripole Index (TPI) for the Interdecadal Pacific Oscillation (IPO) is +#'computed as the difference of weighted-averaged SST anomalies over 10ºS-10ºN, +#'170ºE-270ºE minus the mean of the weighted-averaged SST anomalies over +#'25ºN-45ºN, 140ºE-215ºE and 50ºS-15ºS, 150ºE-200ºE (Henley et al., 2015). +#' +#'@param data A numerical array to be used for the index computation with the +#' dimensions: 1) latitude, longitude, start date, forecast month, and member +#' (in case of decadal predictions), 2) latitude, longitude, year, month and +#' member (in case of historical simulations), or 3) latitude, longitude, year +#' and month (in case of observations or reanalyses). This data has to be +#' provided, at least, over the whole region needed to compute the index. +#'@param data_lats A numeric vector indicating the latitudes of the data. +#'@param data_lons A numeric vector indicating the longitudes of the data. +#'@param A character string indicating the type of data ('dcpp' for +#' decadal predictions, 'hist' for historical simulations, or 'obs' for +#' observations or reanalyses). +#'@param lat_dim A character string of the name of the latitude dimension. The +#' default value is 'lat'. +#'@param lon_dim A character string of the name of the longitude dimension. The +#' default value is 'lon'. +#'@param mask An array of a mask (with 0's in the grid points that have to be +#' masked) or NULL (i.e., no mask is used). This parameter allows to remove +#' the values over land in case the dataset is a combination of surface air +#' temperature over land and sea surface temperature over the ocean. Also, it +#' can be used to mask those grid points that are missing in the observational +#' dataset for a fair comparison between the forecast system and the reference +#' dataset. The default value is NULL. +#'@param monini An integer indicating the month in which the forecast system is +#' initialized. Only used when parameter 'type' is 'dcpp'. The default value +#' is 11, i.e., initialized in November. +#'@param fmonth_dim A character string indicating the name of the forecast +#' month dimension. Only used if parameter 'type' is 'dcpp'. The default value +#' is 'fmonth'. +#'@param sdate_dim A character string indicating the name of the start date +#' dimension. Only used if parameter 'type' is 'dcpp'. The default value is +#' 'sdate'. +#'@param indices_for_clim A numeric vector of the indices of the years to +#' compute the climatology for calculating the anomalies, or NULL so the +#' climatology is calculated over the whole period. If the data are already +#' anomalies, set it to FALSE. The default value is NULL.\cr +#' In case of parameter 'type' is 'dcpp', 'indices_for_clim' must be relative +#' to the first forecast year, and the climatology is automatically computed +#' over the actual common period for the different forecast years. +#'@param year_dim A character string indicating the name of the year dimension +#' The default value is 'year'. Only used if parameter 'type' is 'hist' or +#' 'obs'. +#'@param month_dim A character string indicating the name of the month +#' dimension. The default value is 'month'. Only used if parameter 'type' is +#' 'hist' or 'obs'. +#'@param member_dim A character string indicating the name of the member +#' dimension. The default value is 'member'. Only used if parameter 'type' is +#' 'dcpp' or 'hist'. +#' +#'@return A numerical array of the TPI index with the dimensions of: +#' 1) sdate, forecast year, and member (in case of decadal predictions); +#' 2) year and member (in case of historical simulations); or +#' 3) year (in case of observations or reanalyses). +#' +#'@examples +#' ## Observations or reanalyses +#' obs = array(1:100, dim = c(year = 5, lat = 19, lon = 37, month = 12)) +#' lat = seq(-90, 90, 10) +#' lon = seq(0, 360, 10) +#' index_obs = TPI(data = obs, data_lats = lat, data_lons = lon, type = 'obs') +#' +#' ## Historical simulations +#' hist = array(1:100, dim = c(year = 5, lat = 19, lon = 37, month = 12, member = 5)) +#' lat = seq(-90, 90, 10) +#' lon = seq(0, 360, 10) +#' index_hist = TPI(data = hist, data_lats = lat, data_lons = lon, type = 'hist') +#' +#' ## Decadal predictions +#' dcpp = array(1:100, dim = c(sdate = 5, lat = 19, lon = 37, fmonth = 24, member = 5)) +#' lat = seq(-90, 90, 10) +#' lon = seq(0, 360, 10) +#' index_dcpp = TPI(data = dcpp, data_lats = lat, data_lons = lon, type = 'dcpp', monini = 1) +#' +#'@author Carlos Delgado-Torres, \email{carlos.delgado@bsc.es} +#'@author Núria Pérez-Zanón, \email{nuria.perez@bsc.es} +#'@author Roberto Bilbao, \email{roberto.bilbao@bsc.es} +#' +#'@importFrom ClimProjDiags WeightedMean CombineIndices +#'@import multiApply +#'@export +TPI <- function(data, data_lats, data_lons, index, type, lat_dim = 'lat', lon_dim = 'lon', + mask = NULL, monini = 11, fmonth_dim = 'fmonth', sdate_dim = 'sdate', + indices_for_clim = NULL, year_dim = 'year', month_dim = 'month', + member_dim = 'member') { + + ## Input Checks + # data + if (is.null(data)) { + stop("Parameter 'data' cannot be NULL.") + } + if (!is.numeric(data)) { + stop("Parameter 'data' must be a numeric array.") + } + # data_lats and data_lons part1 + if (!(class(data_lats) == 'numeric' | class(data_lats) == 'integer')) { + stop("Parameter 'data_lats' must be a numeric vector.") + } + if (!(class(data_lons) == 'numeric' | class(data_lons) == 'integer')) { + stop("Parameter 'data_lons' must be a numeric vector.") + } + # type + if (!type %in% c('dcpp', 'hist', 'obs')) { + stop("Parameter 'type' must be 'dcpp', 'hist', or 'obs'.") + } + # lat_dim + if (!(is.character(lat_dim) & length(lat_dim) == 1)) { + stop("Parameter 'lat_dim' must be a character string.") + } + if (!lat_dim %in% names(dim(data))) { + stop("Parameter 'lat_dim' is not found in 'data' dimension.") + } + # lon_dim + if (!(is.character(lon_dim) & length(lon_dim) == 1)) { + stop("Parameter 'lon_dim' must be a character string.") + } + if (!lon_dim %in% names(dim(data))) { + stop("Parameter 'lon_dim' is not found in 'data' dimension.") + } + # data_lats and data_lons part2 + if (dim(data)[lat_dim] != length(data_lats)){ + stop(paste0("The latitude dimension of parameter 'data' must be the same", + " length of parameter 'data_lats'.")) + } + if (dim(data)[lon_dim] != length(data_lons)){ + stop(paste0("The longitude dimension of parameter 'data' must be the same", + " length of parameter 'data_lons'.")) + } + # mask + if (!is.null(mask)) { + if (is.array(mask) & identical(names(dim(mask)), c(lat_dim,lon_dim)) & + identical(as.integer(dim(mask)), c(length(data_lats), length(data_lons)))) { + ## To mask those grid point that are missing in the observations + mask <- s2dv::Reorder(data = mask, order = c(lat_dim, lon_dim)) + fun_mask <- function(data, mask){ + data[mask == 0] <- NA + return(data) + } + data <- multiApply::Apply(data = data, target_dims = c(lat_dim, lon_dim), + fun = fun_mask, mask = mask)$output1 + } else { + stop(paste0("Parameter 'mask' must be NULL (no mask) or a numerical array ", + "with c(lat_dim, lon_dim) dimensions and 0 in those grid ", + "points that have to be masked.")) + } + } + # monini + if (type == 'dcpp') { + if (!is.numeric(monini) | monini %% 1 != 0 | monini < 1 | + monini > 12) { + stop("Parameter 'monini' must be an integer from 1 to 12.") + } + } + # fmonth_dim + if (type == 'dcpp') { + if (!(is.character(fmonth_dim) & length(fmonth_dim) == 1)) { + stop("Parameter 'fmonth_dim' must be a character string.") + } + if (!fmonth_dim %in% names(dim(data))) { + stop("Parameter 'fmonth_dim' is not found in 'data' dimension.") + } + } + # sdate_dim + if (type == 'dcpp') { + if (!(is.character(sdate_dim) & length(sdate_dim) == 1)) { + stop("Parameter 'sdate_dim' must be a character string.") + } + if (!sdate_dim %in% names(dim(data))) { + stop("Parameter 'sdate_dim' is not found in 'data' dimension.") + } + } + # indices_for_clim + if (!is.null(indices_for_clim)) { + if (!class(indices_for_clim) %in% c('numeric', 'integer') + & !(is.logical(indices_for_clim) & !any(indices_for_clim))) { + stop(paste0("The parameter 'indices_for_clim' must be a numeric vector ", + "or NULL to compute the anomalies based on the whole period, ", + "or FALSE if data are already anomalies")) + } + } + # year_dim + if (type == 'hist' | type == 'obs') { + if (!(is.character(year_dim) & length(year_dim) == 1)) { + stop("Parameter 'year_dim' must be a character string.") + } + if (!year_dim %in% names(dim(data))) { + stop("Parameter 'year_dim' is not found in 'data' dimension.") + } + } + # month_dim + if (type == 'hist' | type == 'obs') { + if (!(is.character(month_dim) & length(month_dim) == 1)) { + stop("Parameter 'month_dim' must be a character string.") + } + if (!month_dim %in% names(dim(data))) { + stop("Parameter 'month_dim' is not found in 'data' dimension.") + } + } + # member_dim + if (type == 'hist' | type == 'dcpp') { + if (!(is.character(member_dim) & length(member_dim) == 1)) { + stop("Parameter 'member_dim' must be a character string.") + } + if (!member_dim %in% names(dim(data))) { + stop("Parameter 'member_dim' is not found in 'data' dimension.") + } + } + + # Regions for IPO_TPI (psl.noaa.gov/data/timeseries/IPOTPI) + lat_min_1 <- 25; lat_max_1 <- 45 + lon_min_1 <- 140; lon_max_1 <- 215 + lat_min_2 <- -10; lat_max_2 <- 10 + lon_min_2 <- 170; lon_max_2 <- 270 + lat_min_3 <- -50; lat_max_3 <- -15 + lon_min_3 <- 150; lon_max_3 <- 200 + regions <- NULL + regions$reg1 <- c(lon_min_1, lon_max_1, lat_min_1, lat_max_1) + regions$reg2 <- c(lon_min_2, lon_max_2, lat_min_2, lat_max_2) + regions$reg3 <- c(lon_min_3, lon_max_3, lat_min_3, lat_max_3) + + mean_1 <- ClimProjDiags::WeightedMean(data = data, lon = data_lons, lat = data_lats, + region = regions$reg1, + londim = which(names(dim(data)) == lon_dim), + latdim = which(names(dim(data)) == lat_dim)) + mean_2 <- ClimProjDiags::WeightedMean(data = data, lon = data_lons, lat = data_lats, + region = regions$reg2, + londim = which(names(dim(data)) == lon_dim), + latdim = which(names(dim(data)) == lat_dim)) + mean_3 <- ClimProjDiags::WeightedMean(data = data, lon = data_lons, lat = data_lats, + region = regions$reg3, + londim = which(names(dim(data)) == lon_dim), + latdim = which(names(dim(data)) == lat_dim)) + + mean_1_3 <- ClimProjDiags::CombineIndices(indices = list(mean_1, mean_3), + weights = NULL, operation = 'mean') + + data <- ClimProjDiags::CombineIndices(indices = list(mean_2, mean_1_3), + weights = NULL, operation = 'subtract') # mean_2 - ((mean_1 + mean_3)/2) + + INDEX <- .Indices(data = data, type = type, monini = monini, + indices_for_clim = indices_for_clim, fmonth_dim = fmonth_dim, + sdate_dim = sdate_dim, year_dim = year_dim, + month_dim = month_dim, member_dim = member_dim) + return(INDEX) +} diff --git a/R/Utils.R b/R/Utils.R index 6e8e7135c470b4e33548305bfc79cf664715203a..b815bc5525eb4f7fd816cddb37779f8f98454290 100644 --- a/R/Utils.R +++ b/R/Utils.R @@ -1669,3 +1669,108 @@ return(output) } +# to be used in AMV.R, TPI.R, SPOD.R, GSAT.R and GMST.R +.Indices <- function(data, type, monini, indices_for_clim, + fmonth_dim, sdate_dim, year_dim, month_dim, member_dim) { + + data = drop(data) + + if(member_dim %in% names(dim(data))){ + if (type == 'dcpp'){ + data = s2dv::Reorder(data = data, order = c(sdate_dim,fmonth_dim,member_dim)) + } else if (type %in% c('hist','obs')){ + data = s2dv::Reorder(data = data, order = c(year_dim,month_dim,member_dim)) + } + } + + if (type == 'dcpp'){ + + data = s2dv::Season(data = data, time_dim = fmonth_dim, + monini = monini, moninf = 1, monsup = 12, + method = mean, na.rm = FALSE) + names(dim(data))[which(names(dim(data))==fmonth_dim)] = 'fyear' + if (member_dim %in% names(dim(data))){ + data = s2dv::Reorder(data = data, order = c('fyear',sdate_dim,member_dim)) + } else { + data = s2dv::Reorder(data = data, order = c('fyear',sdate_dim)) + } + + if (is.logical(indices_for_clim)) { + if(!any(indices_for_clim)) { + # indices_for_clim == FALSE -> anomalies are directly given + anom = data + } + + } else { + + ## Different indices_for_clim for each forecast year (same actual years) + + n_fyears = as.numeric(dim(data)['fyear']) + n_sdates = as.numeric(dim(data)[sdate_dim]) + + if (is.null(indices_for_clim)){ + + # indices_for_clim == NULL -> anomalies based on the whole (common) period + first_years_for_clim = n_fyears : 1 + last_years_for_clim = n_sdates : (n_sdates - n_fyears + 1) + + } else { + + first_years_for_clim = seq(from = indices_for_clim[1], by = -1, length.out = n_fyears) + last_years_for_clim = seq(from = indices_for_clim[length(indices_for_clim)], by = -1, length.out = n_fyears) + + } + + anom = array(data = NA, dim = dim(data)) + if (member_dim %in% names(dim(data))){ + clim = array(data = NA, dim = c(dim(data)['fyear'],dim(data)[member_dim])) + } else { + clim = array(data = NA, dim = c(dim(data)['fyear'])) + } + for (i in 1:n_fyears){ + if (member_dim %in% names(dim(data))){ + for (m in 1:as.numeric(dim(data)[member_dim])){ + clim[i,m] = mean(data[i,first_years_for_clim[i]:last_years_for_clim[i],m]) + anom[i,,m] = data[i,,m] - clim[i,m] + } + } else { + clim = mean(data[i,first_years_for_clim[i]:last_years_for_clim[i]]) + anom[i,] = data[i,] - clim + } + } + } + + } else if (type %in% c('obs','hist')){ + + data = multiApply::Apply(data = data, target_dims = month_dim, fun = mean)$output1 + + if (is.logical(indices_for_clim)) { + if(!any(indices_for_clim)) { + anom = data + } + + } else { + + if (is.null(indices_for_clim)){ + + clim = multiApply::Apply(data = data, target_dims = year_dim, fun = mean)$output1 + + } else { + + if (member_dim %in% names(dim(data))){ + target_dims = c(year_dim,member_dim) + } else { + target_dims = year_dim + } + clim = multiApply::Apply(data = ClimProjDiags::Subset(x = data, along = year_dim, indices = indices_for_clim), + target_dims = target_dims, fun = mean)$output1 + } + anom = multiApply::Apply(data = data, target_dims = year_dim, + fun = function(data,clim){data-clim}, clim = clim)$output1 + } + + } else {stop('type must be dcpp, hist or obs')} + + return(anom) +} + diff --git a/man/AMV.Rd b/man/AMV.Rd new file mode 100644 index 0000000000000000000000000000000000000000..7cb7503a8d8f38c9d9738c9f368423218d3c29ca --- /dev/null +++ b/man/AMV.Rd @@ -0,0 +1,116 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/AMV.R +\name{AMV} +\alias{AMV} +\title{Compute the Atlantic Multidecadal Variability (AMV) index} +\usage{ +AMV(data, data_lats, data_lons, type, lat_dim = "lat", lon_dim = "lon", + mask = NULL, monini = 11, fmonth_dim = "fmonth", sdate_dim = "sdate", + indices_for_clim = NULL, year_dim = "year", month_dim = "month", + member_dim = "member") +} +\arguments{ +\item{data}{A numerical array to be used for the index computation with the +dimensions: 1) latitude, longitude, start date, forecast month, and member +(in case of decadal predictions), 2) latitude, longitude, year, month and +member (in case of historical simulations), or 3) latitude, longitude, year +and month (in case of observations or reanalyses). This data has to be +provided, at least, over the whole region needed to compute the index.} + +\item{data_lats}{A numeric vector indicating the latitudes of the data.} + +\item{data_lons}{A numeric vector indicating the longitudes of the data.} + +\item{type}{A character string indicating the type of data ('dcpp' for +decadal predictions, 'hist' for historical simulations, or 'obs' for +observations or reanalyses).} + +\item{lat_dim}{A character string of the name of the latitude dimension. The +default value is 'lat'.} + +\item{lon_dim}{A character string of the name of the longitude dimension. The +default value is 'lon'.} + +\item{mask}{An array of a mask (with 0's in the grid points that have to be +masked) or NULL (i.e., no mask is used). This parameter allows to remove +the values over land in case the dataset is a combination of surface air +temperature over land and sea surface temperature over the ocean. Also, it +can be used to mask those grid points that are missing in the observational +dataset for a fair comparison between the forecast system and the reference +dataset. The default value is NULL.} + +\item{monini}{An integer indicating the month in which the forecast system is +initialized. Only used when parameter 'type' is 'dcpp'. The default value +is 11, i.e., initialized in November.} + +\item{fmonth_dim}{A character string indicating the name of the forecast +month dimension. Only used if parameter 'type' is 'dcpp'. The default value +is 'fmonth'.} + +\item{sdate_dim}{A character string indicating the name of the start date +dimension. Only used if parameter 'type' is 'dcpp'. The default value is +'sdate'.} + +\item{indices_for_clim}{A numeric vector of the indices of the years to +compute the climatology for calculating the anomalies, or NULL so the +climatology is calculated over the whole period. If the data are already +anomalies, set it to FALSE. The default value is NULL.\cr +In case of parameter 'type' is 'dcpp', 'indices_for_clim' must be relative +to the first forecast year, and the climatology is automatically computed +over the actual common period for the different forecast years.} + +\item{year_dim}{A character string indicating the name of the year dimension +The default value is 'year'. Only used if parameter 'type' is 'hist' or +'obs'.} + +\item{month_dim}{A character string indicating the name of the month +dimension. The default value is 'month'. Only used if parameter 'type' is +'hist' or 'obs'.} + +\item{member_dim}{A character string indicating the name of the member +dimension. The default value is 'member'. Only used if parameter 'type' is +'dcpp' or 'hist'.} +} +\value{ +A numerical array of the AMV index with the dimensions of: + 1) sdate, forecast year, and member (in case of decadal predictions); + 2) year and member (in case of historical simulations); or + 3) year (in case of observations or reanalyses). +} +\description{ +The Atlantic Multidecadal Variability (AMV), also known as Atlantic +Multidecadal Oscillation (AMO), is a mode of natural variability of the sea +surface temperatures (SST) over the North Atlantic Ocean on multi-decadal +time scales. The AMV index is computed as the difference of weighted-averaged +SST anomalies over the North Atlantic region (0ºN-60ºN, 280ºE-360ºE) and the +weighted-averaged SST anomalies over 60ºS-60ºN, 0ºE-360ºE (Trenberth & +Dennis, 2005; Doblas-Reyes et al., 2013). +} +\examples{ +## Observations or reanalyses +obs <- array(1:100, dim = c(year = 5, lat = 19, lon = 37, month = 12)) +lat <- seq(-90, 90, 10) +lon <- seq(0, 360, 10) +index_obs <- AMV(data = obs, data_lats = lat, data_lons = lon, type = 'obs') + +## Historical simulations +hist <- array(1:100, dim = c(year = 5, lat = 19, lon = 37, month = 12, member = 5)) +lat <- seq(-90, 90, 10) +lon <- seq(0, 360, 10) +index_hist <- AMV(data = hist, data_lats = lat, data_lons = lon, type = 'hist') + +## Decadal predictions +dcpp <- array(1:100, dim = c(sdate = 5, lat = 19, lon = 37, fmonth = 24, member = 5)) +lat <- seq(-90, 90, 10) +lon <- seq(0, 360, 10) +index_dcpp <- AMV(data = dcpp, data_lats = lat, data_lons = lon, type = 'dcpp', monini = 1) + +} +\author{ +Carlos Delgado-Torres, \email{carlos.delgado@bsc.es} + +Roberto Bilbao, \email{roberto.bilbao@bsc.es} + +Núria Pérez-Zanón, \email{nuria.perez@bsc.es} +} + diff --git a/man/GMST.Rd b/man/GMST.Rd new file mode 100644 index 0000000000000000000000000000000000000000..3c9ace76f88fc47780cbf65b12c88ae5524134d2 --- /dev/null +++ b/man/GMST.Rd @@ -0,0 +1,138 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/GMST.R +\name{GMST} +\alias{GMST} +\title{Compute the Global Mean Surface Temperature (GMST) anomalies} +\usage{ +GMST(data_tas, data_tos, data_lats, data_lons, mask_sea_land, sea_value, type, + mask = NULL, lat_dim = "lat", lon_dim = "lon", monini = 11, + fmonth_dim = "fmonth", sdate_dim = "sdate", indices_for_clim = NULL, + year_dim = "year", month_dim = "month", member_dim = "member") +} +\arguments{ +\item{data_tas}{A numerical array indicating the surface air temperature data +to be used for the index computation with the dimensions: 1) latitude, +longitude, start date, forecast month, and member (in case of decadal +predictions), 2) latitude, longitude, year, month and member (in case of +historical simulations), or 3) latitude, longitude, year and month (in case +of observations or reanalyses). This data has to be provided, at least, +over the whole region needed to compute the index. The dimensions must be +identical to those of data_tos.} + +\item{data_tos}{A numerical array indicating the sea surface temperature data +to be used for the index computation with the dimensions: 1) latitude, +longitude, start date, forecast month, and member (in case of decadal +predictions), 2) latitude, longitude, year, month and member (in case of +historical simulations), or 3) latitude, longitude, year and month (in case +of observations or reanalyses). This data has to be provided, at least, +over the whole region needed to compute the index. The dimensions must be +identical to those of data_tas.} + +\item{data_lats}{A numeric vector indicating the latitudes of the data.} + +\item{data_lons}{A numeric vector indicating the longitudes of the data.} + +\item{mask_sea_land}{An array with dimensions [lat_dim = data_lats, lon_dim = +data_lons] for blending 'data_tas' and 'data_tos'.} + +\item{sea_value}{A numeric value indicating the sea grid points in +'mask_sea_land'.} + +\item{type}{A character string indicating the type of data ('dcpp' for +decadal predictions, 'hist' for historical simulations, or 'obs' for +observations or reanalyses).} + +\item{mask}{An array of a mask (with 0's in the grid points that have to be +masked) or NULL (i.e., no mask is used). This parameter allows to remove +the values over land in case the dataset is a combination of surface air +temperature over land and sea surface temperature over the ocean. Also, it +can be used to mask those grid points that are missing in the observational +dataset for a fair comparison between the forecast system and the reference +dataset. The default value is NULL.} + +\item{lat_dim}{A character string of the name of the latitude dimension. The +default value is 'lat'.} + +\item{lon_dim}{A character string of the name of the longitude dimension. The +default value is 'lon'.} + +\item{monini}{An integer indicating the month in which the forecast system is +initialized. Only used when parameter 'type' is 'dcpp'. The default value +is 11, i.e., initialized in November.} + +\item{fmonth_dim}{A character string indicating the name of the forecast +month dimension. Only used if parameter 'type' is 'dcpp'. The default value +is 'fmonth'.} + +\item{sdate_dim}{A character string indicating the name of the start date +dimension. Only used if parameter 'type' is 'dcpp'. The default value is +'sdate'.} + +\item{indices_for_clim}{A numeric vector of the indices of the years to +compute the climatology for calculating the anomalies, or NULL so the +climatology is calculated over the whole period. If the data are already +anomalies, set it to FALSE. The default value is NULL.\cr +In case of parameter 'type' is 'dcpp', 'indices_for_clim' must be relative +to the first forecast year, and the climatology is automatically computed +over the actual common period for the different forecast years.} + +\item{year_dim}{A character string indicating the name of the year dimension +The default value is 'year'. Only used if parameter 'type' is 'hist' or +'obs'.} + +\item{month_dim}{A character string indicating the name of the month +dimension. The default value is 'month'. Only used if parameter 'type' is +'hist' or 'obs'.} + +\item{member_dim}{A character string indicating the name of the member +dimension. The default value is 'member'. Only used if parameter 'type' is +'dcpp' or 'hist'.} +} +\value{ +A numerical array of the GMST anomalies with the dimensions of: + 1) sdate, forecast year, and member (in case of decadal predictions); + 2) year and member (in case of historical simulations); or + 3) year (in case of observations or reanalyses). +} +\description{ +The Global Mean Surface Temperature (GMST) anomalies are computed as the +weighted-averaged surface air temperature anomalies over land and sea surface +temperature anomalies over the ocean. +} +\examples{ +## Observations or reanalyses +obs_tas <- array(1:100, dim = c(year = 5, lat = 19, lon = 37, month = 12)) +obs_tos <- array(2:101, dim = c(year = 5, lat = 19, lon = 37, month = 12)) +mask_sea_land <- array(c(1,0,1), dim = c(lat = 19, lon = 37)) +sea_value <- 1 +lat <- seq(-90, 90, 10) +lon <- seq(0, 360, 10) +index_obs <- GMST(data_tas = obs_tas, data_tos = obs_tos, data_lats = lat, data_lons = lon, type = 'obs', mask_sea_land = mask_sea_land, sea_value = sea_value) + +## Historical simulations +hist_tas <- array(1:100, dim = c(year = 5, lat = 19, lon = 37, month = 12, member = 5)) +hist_tos <- array(2:101, dim = c(year = 5, lat = 19, lon = 37, month = 12, member = 5)) +mask_sea_land <- array(c(1,0,1), dim = c(lat = 19, lon = 37)) +sea_value <- 1 +lat <- seq(-90, 90, 10) +lon <- seq(0, 360, 10) +index_hist <- GMST(data_tas = hist_tas, data_tos = hist_tos, data_lats = lat, data_lons = lon, type = 'hist', mask_sea_land = mask_sea_land, sea_value = sea_value) + +## Decadal predictions +dcpp_tas <- array(1:100, dim = c(sdate = 5, lat = 19, lon = 37, fmonth = 24, member = 5)) +dcpp_tos <- array(2:101, dim = c(sdate = 5, lat = 19, lon = 37, fmonth = 24, member = 5)) +mask_sea_land <- array(c(1,0,1), dim = c(lat = 19, lon = 37)) +sea_value <- 1 +lat <- seq(-90, 90, 10) +lon <- seq(0, 360, 10) +index_dcpp <- GMST(data_tas = dcpp_tas, data_tos = dcpp_tos, data_lats = lat, data_lons = lon, type = 'dcpp', monini = 1, mask_sea_land = mask_sea_land, sea_value = sea_value) + +} +\author{ +Carlos Delgado-Torres, \email{carlos.delgado@bsc.es} + +Roberto Bilbao, \email{roberto.bilbao@bsc.es} + +Núria Pérez-Zanón, \email{nuria.perez@bsc.es} +} + diff --git a/man/GSAT.Rd b/man/GSAT.Rd new file mode 100644 index 0000000000000000000000000000000000000000..bb43cfb516cafe0064ee4f1d8a0cd38c0e388f42 --- /dev/null +++ b/man/GSAT.Rd @@ -0,0 +1,111 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/GSAT.R +\name{GSAT} +\alias{GSAT} +\title{Compute the Global Surface Air Temperature (GSAT) anomalies} +\usage{ +GSAT(data, data_lats, data_lons, type, lat_dim = "lat", lon_dim = "lon", + mask = NULL, monini = 11, fmonth_dim = "fmonth", sdate_dim = "sdate", + indices_for_clim = NULL, year_dim = "year", month_dim = "month", + member_dim = "member") +} +\arguments{ +\item{data}{A numerical array to be used for the index computation with the +dimensions: 1) latitude, longitude, start date, forecast month, and member +(in case of decadal predictions), 2) latitude, longitude, year, month and +member (in case of historical simulations), or 3) latitude, longitude, year +and month (in case of observations or reanalyses). This data has to be +provided, at least, over the whole region needed to compute the index.} + +\item{data_lats}{A numeric vector indicating the latitudes of the data.} + +\item{data_lons}{A numeric vector indicating the longitudes of the data.} + +\item{type}{A character string indicating the type of data ('dcpp' for +decadal predictions, 'hist' for historical simulations, or 'obs' for +observations or reanalyses).} + +\item{lat_dim}{A character string of the name of the latitude dimension. The +default value is 'lat'.} + +\item{lon_dim}{A character string of the name of the longitude dimension. The +default value is 'lon'.} + +\item{mask}{An array of a mask (with 0's in the grid points that have to be +masked) or NULL (i.e., no mask is used). This parameter allows to remove +the values over land in case the dataset is a combination of surface air +temperature over land and sea surface temperature over the ocean. Also, it +can be used to mask those grid points that are missing in the observational +dataset for a fair comparison between the forecast system and the reference +dataset. The default value is NULL.} + +\item{monini}{An integer indicating the month in which the forecast system is +initialized. Only used when parameter 'type' is 'dcpp'. The default value +is 11, i.e., initialized in November.} + +\item{fmonth_dim}{A character string indicating the name of the forecast +month dimension. Only used if parameter 'type' is 'dcpp'. The default value +is 'fmonth'.} + +\item{sdate_dim}{A character string indicating the name of the start date +dimension. Only used if parameter 'type' is 'dcpp'. The default value is +'sdate'.} + +\item{indices_for_clim}{A numeric vector of the indices of the years to +compute the climatology for calculating the anomalies, or NULL so the +climatology is calculated over the whole period. If the data are already +anomalies, set it to FALSE. The default value is NULL.\cr +In case of parameter 'type' is 'dcpp', 'indices_for_clim' must be relative +to the first forecast year, and the climatology is automatically computed +over the actual common period for the different forecast years.} + +\item{year_dim}{A character string indicating the name of the year dimension +The default value is 'year'. Only used if parameter 'type' is 'hist' or +'obs'.} + +\item{month_dim}{A character string indicating the name of the month +dimension. The default value is 'month'. Only used if parameter 'type' is +'hist' or 'obs'.} + +\item{member_dim}{A character string indicating the name of the member +dimension. The default value is 'member'. Only used if parameter 'type' is +'dcpp' or 'hist'.} +} +\value{ +A numerical array of the GSAT anomalies with the dimensions of: + 1) sdate, forecast year, and member (in case of decadal predictions); + 2) year and member (in case of historical simulations); or + 3) year (in case of observations or reanalyses). +} +\description{ +The Global Surface Air Temperature (GSAT) anomalies are computed as the +weighted-averaged surface air temperature anomalies over the global region. +} +\examples{ +## Observations or reanalyses +obs <- array(1:100, dim = c(year = 5, lat = 19, lon = 37, month = 12)) +lat <- seq(-90, 90, 10) +lon <- seq(0, 360, 10) +index_obs <- GSAT(data = obs, data_lats = lat, data_lons = lon, type = 'obs') + +## Historical simulations +hist <- array(1:100, dim = c(year = 5, lat = 19, lon = 37, month = 12, member = 5)) +lat <- seq(-90, 90, 10) +lon <- seq(0, 360, 10) +index_hist <- GSAT(data = hist, data_lats = lat, data_lons = lon, type = 'hist') + +## Decadal predictions +dcpp <- array(1:100, dim = c(sdate = 5, lat = 19, lon = 37, fmonth = 24, member = 5)) +lat <- seq(-90, 90, 10) +lon <- seq(0, 360, 10) +index_dcpp <- GSAT(data = dcpp, data_lats = lat, data_lons = lon, type = 'dcpp', monini = 1) + +} +\author{ +Carlos Delgado-Torres, \email{carlos.delgado@bsc.es} + +Núria Pérez-Zanón, \email{nuria.perez@bsc.es} + +Roberto Bilbao, \email{roberto.bilbao@bsc.es} +} + diff --git a/man/SPOD.Rd b/man/SPOD.Rd new file mode 100644 index 0000000000000000000000000000000000000000..09621c793845835585e5c3073dba550b5c395211 --- /dev/null +++ b/man/SPOD.Rd @@ -0,0 +1,114 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/SPOD.R +\name{SPOD} +\alias{SPOD} +\title{Compute the South Pacific Ocean Dipole (SPOD) index} +\usage{ +SPOD(data, data_lats, data_lons, type, lat_dim = "lat", lon_dim = "lon", + mask = NULL, monini = 11, fmonth_dim = "fmonth", sdate_dim = "sdate", + indices_for_clim = NULL, year_dim = "year", month_dim = "month", + member_dim = "member") +} +\arguments{ +\item{data}{A numerical array to be used for the index computation with the +dimensions: 1) latitude, longitude, start date, forecast month, and member +(in case of decadal predictions), 2) latitude, longitude, year, month and +member (in case of historical simulations), or 3) latitude, longitude, year +and month (in case of observations or reanalyses). This data has to be +provided, at least, over the whole region needed to compute the index.} + +\item{data_lats}{A numeric vector indicating the latitudes of the data.} + +\item{data_lons}{A numeric vector indicating the longitudes of the data.} + +\item{type}{A character string indicating the type of data ('dcpp' for +decadal predictions, 'hist' for historical simulations, or 'obs' for +observations or reanalyses).} + +\item{lat_dim}{A character string of the name of the latitude dimension. The +default value is 'lat'.} + +\item{lon_dim}{A character string of the name of the longitude dimension. The +default value is 'lon'.} + +\item{mask}{An array of a mask (with 0's in the grid points that have to be +masked) or NULL (i.e., no mask is used). This parameter allows to remove +the values over land in case the dataset is a combination of surface air +temperature over land and sea surface temperature over the ocean. Also, it +can be used to mask those grid points that are missing in the observational +dataset for a fair comparison between the forecast system and the reference +dataset. The default value is NULL.} + +\item{monini}{An integer indicating the month in which the forecast system is +initialized. Only used when parameter 'type' is 'dcpp'. The default value +is 11, i.e., initialized in November.} + +\item{fmonth_dim}{A character string indicating the name of the forecast +month dimension. Only used if parameter 'type' is 'dcpp'. The default value +is 'fmonth'.} + +\item{sdate_dim}{A character string indicating the name of the start date +dimension. Only used if parameter 'type' is 'dcpp'. The default value is +'sdate'.} + +\item{indices_for_clim}{A numeric vector of the indices of the years to +compute the climatology for calculating the anomalies, or NULL so the +climatology is calculated over the whole period. If the data are already +anomalies, set it to FALSE. The default value is NULL.\cr +In case of parameter 'type' is 'dcpp', 'indices_for_clim' must be relative +to the first forecast year, and the climatology is automatically computed +over the actual common period for the different forecast years.} + +\item{year_dim}{A character string indicating the name of the year dimension +The default value is 'year'. Only used if parameter 'type' is 'hist' or +'obs'.} + +\item{month_dim}{A character string indicating the name of the month +dimension. The default value is 'month'. Only used if parameter 'type' is +'hist' or 'obs'.} + +\item{member_dim}{A character string indicating the name of the member +dimension. The default value is 'member'. Only used if parameter 'type' is +'dcpp' or 'hist'.} +} +\value{ +A numerical array of the SPOD index with the dimensions of: + 1) sdate, forecast year, and member (in case of decadal predictions); + 2) year and member (in case of historical simulations); or + 3) year (in case of observations or reanalyses). +} +\description{ +The South Pacific Ocean Dipole (SPOD) index is related to the El +Nino-Southern Oscillation (ENSO) and the Inderdecadal Pacific Oscillation +(IPO). The SPOD index is computed as the difference of weighted-averaged SST +anomalies over 20ºS-48ºS, 165ºE-190ºE (NW pole) and the weighted-averaged SST +anomalies over 44ºS-65ºS, 220ºE-260ºE (SE pole) (Saurral et al., 2020). +} +\examples{ +## Observations or reanalyses +obs <- array(1:100, dim = c(year = 5, lat = 19, lon = 37, month = 12)) +lat <- seq(-90, 90, 10) +lon <- seq(0, 360, 10) +index_obs <- SPOD(data = obs, data_lats = lat, data_lons = lon, type = 'obs') + +## Historical simulations +hist <- array(1:100, dim = c(year = 5, lat = 19, lon = 37, month = 12, member = 5)) +lat <- seq(-90, 90, 10) +lon <- seq(0, 360, 10) +index_hist <- SPOD(data = hist, data_lats = lat, data_lons = lon, type = 'hist') + +## Decadal predictions +dcpp <- array(1:100, dim = c(sdate = 5, lat = 19, lon = 37, fmonth = 24, member = 5)) +lat <- seq(-90, 90, 10) +lon <- seq(0, 360, 10) +index_dcpp <- SPOD(data = dcpp, data_lats = lat, data_lons = lon, type = 'dcpp', monini = 1) + +} +\author{ +Carlos Delgado-Torres, \email{carlos.delgado@bsc.es} + +Núria Pérez-Zanón, \email{nuria.perez@bsc.es} + +Roberto Bilbao, \email{roberto.bilbao@bsc.es} +} + diff --git a/man/TPI.Rd b/man/TPI.Rd new file mode 100644 index 0000000000000000000000000000000000000000..d7938ffde4223c8231fddf5cc8f47f0f896e34b7 --- /dev/null +++ b/man/TPI.Rd @@ -0,0 +1,113 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/TPI.R +\name{TPI} +\alias{TPI} +\title{Compute the Tripole Index (TPI) for the Interdecadal Pacific Oscillation (IPO)} +\usage{ +TPI(data, data_lats, data_lons, index, type, lat_dim = "lat", + lon_dim = "lon", mask = NULL, monini = 11, fmonth_dim = "fmonth", + sdate_dim = "sdate", indices_for_clim = NULL, year_dim = "year", + month_dim = "month", member_dim = "member") +} +\arguments{ +\item{data}{A numerical array to be used for the index computation with the +dimensions: 1) latitude, longitude, start date, forecast month, and member +(in case of decadal predictions), 2) latitude, longitude, year, month and +member (in case of historical simulations), or 3) latitude, longitude, year +and month (in case of observations or reanalyses). This data has to be +provided, at least, over the whole region needed to compute the index.} + +\item{data_lats}{A numeric vector indicating the latitudes of the data.} + +\item{data_lons}{A numeric vector indicating the longitudes of the data.} + +\item{lat_dim}{A character string of the name of the latitude dimension. The +default value is 'lat'.} + +\item{lon_dim}{A character string of the name of the longitude dimension. The +default value is 'lon'.} + +\item{mask}{An array of a mask (with 0's in the grid points that have to be +masked) or NULL (i.e., no mask is used). This parameter allows to remove +the values over land in case the dataset is a combination of surface air +temperature over land and sea surface temperature over the ocean. Also, it +can be used to mask those grid points that are missing in the observational +dataset for a fair comparison between the forecast system and the reference +dataset. The default value is NULL.} + +\item{monini}{An integer indicating the month in which the forecast system is +initialized. Only used when parameter 'type' is 'dcpp'. The default value +is 11, i.e., initialized in November.} + +\item{fmonth_dim}{A character string indicating the name of the forecast +month dimension. Only used if parameter 'type' is 'dcpp'. The default value +is 'fmonth'.} + +\item{sdate_dim}{A character string indicating the name of the start date +dimension. Only used if parameter 'type' is 'dcpp'. The default value is +'sdate'.} + +\item{indices_for_clim}{A numeric vector of the indices of the years to +compute the climatology for calculating the anomalies, or NULL so the +climatology is calculated over the whole period. If the data are already +anomalies, set it to FALSE. The default value is NULL.\cr +In case of parameter 'type' is 'dcpp', 'indices_for_clim' must be relative +to the first forecast year, and the climatology is automatically computed +over the actual common period for the different forecast years.} + +\item{year_dim}{A character string indicating the name of the year dimension +The default value is 'year'. Only used if parameter 'type' is 'hist' or +'obs'.} + +\item{month_dim}{A character string indicating the name of the month +dimension. The default value is 'month'. Only used if parameter 'type' is +'hist' or 'obs'.} + +\item{member_dim}{A character string indicating the name of the member +dimension. The default value is 'member'. Only used if parameter 'type' is +'dcpp' or 'hist'.} + +\item{A}{character string indicating the type of data ('dcpp' for +decadal predictions, 'hist' for historical simulations, or 'obs' for +observations or reanalyses).} +} +\value{ +A numerical array of the TPI index with the dimensions of: + 1) sdate, forecast year, and member (in case of decadal predictions); + 2) year and member (in case of historical simulations); or + 3) year (in case of observations or reanalyses). +} +\description{ +The Tripole Index (TPI) for the Interdecadal Pacific Oscillation (IPO) is +computed as the difference of weighted-averaged SST anomalies over 10ºS-10ºN, +170ºE-270ºE minus the mean of the weighted-averaged SST anomalies over +25ºN-45ºN, 140ºE-215ºE and 50ºS-15ºS, 150ºE-200ºE (Henley et al., 2015). +} +\examples{ +## Observations or reanalyses +obs = array(1:100, dim = c(year = 5, lat = 19, lon = 37, month = 12)) +lat = seq(-90, 90, 10) +lon = seq(0, 360, 10) +index_obs = TPI(data = obs, data_lats = lat, data_lons = lon, type = 'obs') + +## Historical simulations +hist = array(1:100, dim = c(year = 5, lat = 19, lon = 37, month = 12, member = 5)) +lat = seq(-90, 90, 10) +lon = seq(0, 360, 10) +index_hist = TPI(data = hist, data_lats = lat, data_lons = lon, type = 'hist') + +## Decadal predictions +dcpp = array(1:100, dim = c(sdate = 5, lat = 19, lon = 37, fmonth = 24, member = 5)) +lat = seq(-90, 90, 10) +lon = seq(0, 360, 10) +index_dcpp = TPI(data = dcpp, data_lats = lat, data_lons = lon, type = 'dcpp', monini = 1) + +} +\author{ +Carlos Delgado-Torres, \email{carlos.delgado@bsc.es} + +Núria Pérez-Zanón, \email{nuria.perez@bsc.es} + +Roberto Bilbao, \email{roberto.bilbao@bsc.es} +} + diff --git a/tests/testthat/test-AMV.R b/tests/testthat/test-AMV.R new file mode 100644 index 0000000000000000000000000000000000000000..562605b6cf04efc348c27456fdbbc63250306b7b --- /dev/null +++ b/tests/testthat/test-AMV.R @@ -0,0 +1,240 @@ +context("s2dv::AMV tests") + +############################################## + lat <- seq(-90, 90, 10) + lon <- seq(0, 360, 10) + set.seed(1) + mask <- array(sample(c(0, 1), 19*37, replace = T), dim = c(lat = 19, lon = 37)) + +## Observations or reanalyses + dat1 <- array(1:100, dim = c(year = 5, lat = 19, lon = 37, month = 12)) + + ## Historical simulations + dat2 <- array(1:100, dim = c(year = 5, lat = 19, lon = 37, month = 12, member = 5)) + + ## Decadal predictions + dat3 <- array(1:100, dim = c(sdate = 5, lat = 19, lon = 37, fmonth = 24, member = 5)) + +############################################## +test_that("1. Input checks", { + + expect_error( + AMV(data = c()), + "Parameter 'data' cannot be NULL." + ) + expect_error( + AMV(data = 'a'), + "Parameter 'data' must be a numeric array." + ) + expect_error( + AMV(data = dat1, data_lats = c()), + "Parameter 'data_lats' must be a numeric vector." + ) + expect_error( + AMV(data = dat1, data_lats = lat, data_lons = c()), + "Parameter 'data_lons' must be a numeric vector." + ) + expect_error( + AMV(data = dat1, data_lats = lat, data_lons = lon, type = TRUE), + "Parameter 'type' must be 'dcpp', 'hist', or 'obs'." + ) + expect_error( + AMV(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + lat_dim = 2), + "Parameter 'lat_dim' must be a character string." + ) + expect_error( + AMV(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + lat_dim = 'lats'), + "Parameter 'lat_dim' is not found in 'data' dimension." + ) + expect_error( + AMV(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + lon_dim = 2), + "Parameter 'lon_dim' must be a character string." + ) + expect_error( + AMV(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + lon_dim = 'lons'), + "Parameter 'lon_dim' is not found in 'data' dimension." + ) + expect_error( + AMV(data = dat3, data_lats = seq(-90, 90, 11), data_lons = lon, type = 'dcpp'), + paste0("The latitude dimension of parameter 'data' must be the same", + " length of parameter 'data_lats'.") + ) + expect_error( + AMV(data = dat3, data_lats = lat, data_lons = c(1:10), type = 'dcpp'), + paste0("The longitude dimension of parameter 'data' must be the same", + " length of parameter 'data_lons'.") + ) + expect_error( + AMV(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + mask = c(1, 0, 1, 0, 0)), + paste0("Parameter 'mask' must be NULL \\(no mask\\) or a numerical array ", + "with c\\(lat_dim, lon_dim\\) dimensions and 0 in those grid ", + "points that have to be masked.") + ) + expect_error( + AMV(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + mask = array(rep(c(1, 0), 100), dim = c(lat = 20, lon = 5))), + paste0("Parameter 'mask' must be NULL \\(no mask\\) or a numerical array ", + "with c\\(lat_dim, lon_dim\\) dimensions and 0 in those grid ", + "points that have to be masked.") + ) + expect_error( + AMV(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + monini = 0), + "Parameter 'monini' must be an integer from 1 to 12." + ) + expect_error( + AMV(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + fmonth_dim = c('lat', 'lon')), + "Parameter 'fmonth_dim' must be a character string." + ) + expect_error( + AMV(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + fmonth_dim = 'asd'), + "Parameter 'fmonth_dim' is not found in 'data' dimension." + ) + expect_error( + AMV(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + sdate_dim = 2), + "Parameter 'sdate_dim' must be a character string." + ) + expect_error( + AMV(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + sdate_dim = 'date'), + "Parameter 'sdate_dim' is not found in 'data' dimension." + ) + expect_error( + AMV(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + indices_for_clim = TRUE), + paste0("The parameter 'indices_for_clim' must be a numeric vector ", + "or NULL to compute the anomalies based on the whole period, ", + "or FALSE if data are already anomalies") + ) + expect_error( + AMV(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + year_dim = 1), + "Parameter 'year_dim' must be a character string." + ) + expect_error( + AMV(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + year_dim = 'yr'), + "Parameter 'year_dim' is not found in 'data' dimension." + ) + expect_error( + AMV(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + month_dim = 0), + "Parameter 'month_dim' must be a character string." + ) + expect_error( + AMV(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + month_dim = 'mth'), + "Parameter 'month_dim' is not found in 'data' dimension." + ) + expect_error( + AMV(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + member_dim = 1), + "Parameter 'member_dim' must be a character string." + ) + expect_error( + AMV(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + member_dim = 'ens'), + "Parameter 'member_dim' is not found in 'data' dimension." + ) + +}) + +############################################## +test_that("2. Output checks: dat1", { + res1_1 <- AMV(data = dat1, data_lats = lat, data_lons = lon, type = 'obs') + + expect_equal( + dim(res1_1), + c(year = 5) + ) + expect_equal( + mean(res1_1)*10^17, + -4.440892, + tolerance = 0.00001 + ) + + res1_2 <- AMV(data = dat1, data_lats = lat, data_lons = lon, type = 'obs', + indices_for_clim = F) + expect_equal( + mean(res1_2), + -1.78994, + tolerance = 0.00001 + ) + +}) +############################################## +test_that("3. Output checks: dat2", { + res2_1 <- AMV(data = dat2, data_lats = lat, data_lons = lon, type = 'hist') + expect_equal( + dim(res2_1), + c(year = 5, member = 5) + ) + expect_equal( + mean(res2_1)*10^17, + -3.549383, + tolerance = 0.00001 + ) + + res2_2 <- AMV(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + indices_for_clim = 2:3) + expect_equal( + mean(res2_2)*10^16, + 8.616545, + tolerance = 0.00001 + ) + + res2_3 <- AMV(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + mask = mask) + expect_equal( + mean(res2_3)*10^17, + -4.437162, + tolerance = 0.00001 + ) + +}) +############################################## +test_that("4. Output checks: dat3", { + res3_1 <- AMV(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp') + expect_equal( + dim(res3_1), + c(fyear = 1, sdate = 5, member = 5) + ) + expect_equal( + mean(res3_1)*10^18, + -8.881784, + tolerance = 0.00001 + ) + + res3_2 <- AMV(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + indices_for_clim = 2:3) + expect_equal( + mean(res3_2)*10^16, + 8.459899, + tolerance = 0.00001 + ) + + res3_3 <- AMV(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + monini = 5) + expect_equal( + mean(res3_3)*10^17, + -4.884981, + tolerance = 0.00001 + ) + + res3_4 <- AMV(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + mask = mask) + expect_equal( + mean(res3_4)*10^17, + 2.220446, + tolerance = 0.00001 + ) + +}) diff --git a/tests/testthat/test-GMST.R b/tests/testthat/test-GMST.R new file mode 100644 index 0000000000000000000000000000000000000000..d6e3a4f882f02d38fa5d22164b7b7fc66398191b --- /dev/null +++ b/tests/testthat/test-GMST.R @@ -0,0 +1,288 @@ +context("s2dv::GMST tests") + +############################################## + lat <- seq(-90, 90, 10) + lon <- seq(0, 360, 10) + set.seed(1) + mask <- array(sample(c(0, 1), 19*37, replace = T), dim = c(lat = 19, lon = 37)) + mask_sea_land <- array(c(1, 0, 1), dim = c(lat = 19, lon = 37)) + sea_value <- 1 + +## Observations or reanalyses + dat_tas1 <- array(1:100, dim = c(year = 5, lat = 19, lon = 37, month = 12)) + dat_tos1 <- array(2:101, dim = c(year = 5, lat = 19, lon = 37, month = 12)) + + ## Historical simulations + dat_tas2 <- array(1:100, dim = c(year = 5, lat = 19, lon = 37, month = 12, member = 5)) + dat_tos2 <- array(2:101, dim = c(year = 5, lat = 19, lon = 37, month = 12, member = 5)) + + ## Decadal predictions + dat_tas3 <- array(1:100, dim = c(sdate = 5, lat = 19, lon = 37, fmonth = 24, member = 5)) + dat_tos3 <- array(2:101, dim = c(sdate = 5, lat = 19, lon = 37, fmonth = 24, member = 5)) + +############################################## +test_that("1. Input checks", { + + expect_error( + GMST(data_tas = c(), data_tos = c()), + "Parameter 'data_tas' and 'data_tos' cannot be NULL." + ) + expect_error( + GMST(data_tas = 'a', data_tos = dat_tos1), + "Parameter 'data_tas' and 'data_tos' must be a numeric array." + ) + expect_error( + GMST(data_tas = dat_tas1, data_tos = dat_tos2), + "The dimension of data_tas and data_tos must be identical." + ) + expect_error( + GMST(data_tas = dat_tas1, data_tos = dat_tos1, data_lats = c()), + "Parameter 'data_lats' must be a numeric vector." + ) + expect_error( + GMST(data_tas = dat_tas1, data_tos = dat_tos1, data_lats = lat, data_lons = c()), + "Parameter 'data_lons' must be a numeric vector." + ) + expect_error( + GMST(data_tas = dat_tas1, data_tos = dat_tos1, data_lats = lat, data_lons = lon, + sea_value = sea_value, mask_sea_land = TRUE), + "Parameter 'mask_sea_land' must be an array with dimensions \\[lat_dim, lon_dim\\]." + ) + expect_error( + GMST(data_tas = dat_tas1, data_tos = dat_tos1, data_lats = lat, data_lons = lon, + mask_sea_land = mask_sea_land, sea_value = TRUE), + "Parameter 'sea_value' must be a numeric value." + ) + expect_error( + GMST(data_tas = dat_tas1, data_tos = dat_tos1, data_lats = lat, data_lons = lon, + mask_sea_land = mask_sea_land, sea_value = sea_value, type = TRUE), + "Parameter 'type' must be 'dcpp', 'hist', or 'obs'." + ) + expect_error( + GMST(data_tas = dat_tas3, data_tos = dat_tos3, data_lats = lat, data_lons = lon, + mask_sea_land = mask_sea_land, sea_value = sea_value, type = 'dcpp', + lat_dim = 2), + "Parameter 'lat_dim' must be a character string." + ) + expect_error( + GMST(data_tas = dat_tas3, data_tos = dat_tos3, data_lats = lat, data_lons = lon, + mask_sea_land = mask_sea_land, sea_value = sea_value, type = 'dcpp', + lat_dim = 'lats'), + "Parameter 'lat_dim' is not found in 'data_tas' or 'data_tos' dimension." + ) + expect_error( + GMST(data_tas = dat_tas3, data_tos = dat_tos3, data_lats = lat, data_lons = lon, + mask_sea_land = mask_sea_land, sea_value = sea_value, type = 'dcpp', + lon_dim = 2), + "Parameter 'lon_dim' must be a character string." + ) + expect_error( + GMST(data_tas = dat_tas3, data_tos = dat_tos3, data_lats = lat, data_lons = lon, + mask_sea_land = mask_sea_land, sea_value = sea_value, type = 'dcpp', + lon_dim = 'lons'), + "Parameter 'lon_dim' is not found in 'data_tas' or 'data_tos' dimension." + ) + expect_error( + GMST(data_tas = dat_tas3, data_tos = dat_tos3, data_lats = seq(-90, 90, 11), data_lons = lon, mask_sea_land = mask_sea_land, sea_value = sea_value, type = 'dcpp'), + paste0("The latitude dimension of parameter 'data_tas' and 'data_tos'", + " must be the same length of parameter 'data_lats'.") + ) + expect_error( + GMST(data_tas = dat_tas3, data_tos = dat_tos3, data_lats = lat, data_lons = c(1:10), + mask_sea_land = mask_sea_land, sea_value = sea_value, type = 'dcpp'), + paste0("The longitude dimension of parameter 'data_tas' and 'data_tos'", + " must be the same length of parameter 'data_lons'.") + ) + expect_error( + GMST(data_tas = dat_tas3, data_tos = dat_tos3, data_lats = lat, data_lons = lon, + mask_sea_land = mask_sea_land, sea_value = sea_value, type = 'dcpp', + mask = c(1, 0, 1, 0, 0)), + paste0("Parameter 'mask' must be NULL \\(no mask\\) or a numerical array ", + "with c\\(lat_dim, lon_dim\\) dimensions and 0 in those grid ", + "points that have to be masked.") + ) + expect_error( + GMST(data_tas = dat_tas3, data_tos = dat_tos3, data_lats = lat, data_lons = lon, + mask_sea_land = mask_sea_land, sea_value = sea_value, type = 'dcpp', + mask = array(rep(c(1, 0), 100), dim = c(lat = 20, lon = 5))), + paste0("Parameter 'mask' must be NULL \\(no mask\\) or a numerical array ", + "with c\\(lat_dim, lon_dim\\) dimensions and 0 in those grid ", + "points that have to be masked.") + ) + expect_error( + GMST(data_tas = dat_tas3, data_tos = dat_tos3, data_lats = lat, data_lons = lon, + mask_sea_land = mask_sea_land, sea_value = sea_value, type = 'dcpp', + monini = 0), + "Parameter 'monini' must be an integer from 1 to 12." + ) + expect_error( + GMST(data_tas = dat_tas3, data_tos = dat_tos3, data_lats = lat, data_lons = lon, + mask_sea_land = mask_sea_land, sea_value = sea_value, type = 'dcpp', + fmonth_dim = c('lat', 'lon')), + "Parameter 'fmonth_dim' must be a character string." + ) + expect_error( + GMST(data_tas = dat_tas3, data_tos = dat_tos3, data_lats = lat, data_lons = lon, + mask_sea_land = mask_sea_land, sea_value = sea_value, type = 'dcpp', + fmonth_dim = 'asd'), + "Parameter 'fmonth_dim' is not found in 'data_tas' or 'data_tos' dimension." + ) + expect_error( + GMST(data_tas = dat_tas3, data_tos = dat_tos3, data_lats = lat, data_lons = lon, + mask_sea_land = mask_sea_land, sea_value = sea_value, type = 'dcpp', + sdate_dim = 2), + "Parameter 'sdate_dim' must be a character string." + ) + expect_error( + GMST(data_tas = dat_tas3, data_tos = dat_tos3, data_lats = lat, data_lons = lon, + mask_sea_land = mask_sea_land, sea_value = sea_value, type = 'dcpp', + sdate_dim = 'date'), + "Parameter 'sdate_dim' is not found in 'data_tas' or 'data_tos' dimension." + ) + expect_error( + GMST(data_tas = dat_tas3, data_tos = dat_tos3, data_lats = lat, data_lons = lon, + mask_sea_land = mask_sea_land, sea_value = sea_value, type = 'dcpp', + indices_for_clim = TRUE), + paste0("The parameter 'indices_for_clim' must be a numeric vector ", + "or NULL to compute the anomalies based on the whole period, ", + "or FALSE if data are already anomalies") + ) + expect_error( + GMST(data_tas = dat_tas2, data_tos = dat_tos2, data_lats = lat, data_lons = lon, + mask_sea_land = mask_sea_land, sea_value = sea_value, type = 'hist', + year_dim = 1), + "Parameter 'year_dim' must be a character string." + ) + expect_error( + GMST(data_tas = dat_tas2, data_tos = dat_tos2, data_lats = lat, data_lons = lon, + mask_sea_land = mask_sea_land, sea_value = sea_value, type = 'hist', + year_dim = 'yr'), + "Parameter 'year_dim' is not found in 'data_tas' or 'data_tos' dimension." + ) + expect_error( + GMST(data_tas = dat_tas2, data_tos = dat_tos2, data_lats = lat, data_lons = lon, + mask_sea_land = mask_sea_land, sea_value = sea_value, type = 'hist', + month_dim = 0), + "Parameter 'month_dim' must be a character string." + ) + expect_error( + GMST(data_tas = dat_tas2, data_tos = dat_tos2, data_lats = lat, data_lons = lon, + mask_sea_land = mask_sea_land, sea_value = sea_value, type = 'hist', + month_dim = 'mth'), + "Parameter 'month_dim' is not found in 'data_tas' or 'data_tos' dimension." + ) + expect_error( + GMST(data_tas = dat_tas2, data_tos = dat_tos2, data_lats = lat, data_lons = lon, + mask_sea_land = mask_sea_land, sea_value = sea_value, type = 'hist', + member_dim = 1), + "Parameter 'member_dim' must be a character string." + ) + expect_error( + GMST(data_tas = dat_tas2, data_tos = dat_tos2, data_lats = lat, data_lons = lon, + mask_sea_land = mask_sea_land, sea_value = sea_value, type = 'hist', + member_dim = 'ens'), + "Parameter 'member_dim' is not found in 'data_tas' or 'data_tos' dimension." + ) + +}) + +############################################## +test_that("2. Output checks: dat1", { + res1_1 <- GMST(data_tas = dat_tas1, data_tos = dat_tos1, data_lats = lat, data_lons = lon, + mask_sea_land = mask_sea_land, sea_value = sea_value, type = 'obs') + + expect_equal( + dim(res1_1), + c(year = 5) + ) + expect_equal( + res1_1, + array(seq(-2, 2, 1), dim = c(year = 5)) + ) + + res1_2 <- GMST(data_tas = dat_tas1, data_tos = dat_tos1, data_lats = lat, data_lons = lon, + mask_sea_land = mask_sea_land, sea_value = sea_value, type = 'obs', + indices_for_clim = F) + expect_equal( + mean(res1_2), + 51.09007, + tolerance = 0.00001 + ) + +}) +############################################## +test_that("3. Output checks: dat2", { + res2_1 <- GMST(data_tas = dat_tas2, data_tos = dat_tos2, data_lats = lat, data_lons = lon, + mask_sea_land = mask_sea_land, sea_value = sea_value, type = 'hist') + expect_equal( + dim(res2_1), + c(year = 5, member = 5) + ) + expect_equal( + mean(res2_1), + 0 + ) + expect_equal( + range(res2_1), + c(-2.028458, 2.028458), + tolerance = 0.00001 + ) + + res2_2 <- GMST(data_tas = dat_tas2, data_tos = dat_tos2, data_lats = lat, data_lons = lon, + mask_sea_land = mask_sea_land, sea_value = sea_value, type = 'hist', + indices_for_clim = 2:3) + expect_equal( + mean(res2_2), + 0.5 + ) + + res2_3 <- GMST(data_tas = dat_tas2, data_tos = dat_tos2, data_lats = lat, data_lons = lon, + mask_sea_land = mask_sea_land, sea_value = sea_value, type = 'hist', + mask = mask) + expect_equal( + mean(res2_3)*10^16, + 2.842084, + tolerance = 0.00001 + ) + +}) +############################################## +test_that("4. Output checks: dat3", { + res3_1 <- GMST(data_tas = dat_tas3, data_tos = dat_tos3, data_lats = lat, data_lons = lon, + mask_sea_land = mask_sea_land, sea_value = sea_value, type = 'dcpp') + expect_equal( + dim(res3_1), + c(fyear = 1, sdate = 5, member = 5) + ) + expect_equal( + mean(res3_1)*10^16, + -2.842084, + tolerance = 0.00001 + ) + + res3_2 <- GMST(data_tas = dat_tas3, data_tos = dat_tos3, data_lats = lat, data_lons = lon, + mask_sea_land = mask_sea_land, sea_value = sea_value, type = 'dcpp', + indices_for_clim = 2:3) + expect_equal( + mean(res3_2), + 0.5 + ) + + res3_3 <- GMST(data_tas = dat_tas3, data_tos = dat_tos3, data_lats = lat, data_lons = lon, + mask_sea_land = mask_sea_land, sea_value = sea_value, type = 'dcpp', + monini = 10) + expect_equal( + mean(res3_3), + 0 + ) + + res3_4 <- GMST(data_tas = dat_tas3, data_tos = dat_tos3, data_lats = lat, data_lons = lon, + mask_sea_land = mask_sea_land, sea_value = sea_value, type = 'dcpp', + mask = mask) + expect_equal( + mean(res3_4)*10^16, + -8.52686, + tolerance = 0.00001 + ) + +}) diff --git a/tests/testthat/test-GSAT.R b/tests/testthat/test-GSAT.R new file mode 100644 index 0000000000000000000000000000000000000000..bd53e8ab062e232098cdcdb0d0db4e33dd912919 --- /dev/null +++ b/tests/testthat/test-GSAT.R @@ -0,0 +1,240 @@ +context("s2dv::GSAT tests") + +############################################## + lat <- seq(-90, 90, 10) + lon <- seq(0, 360, 10) + set.seed(1) + mask <- array(sample(c(0, 1), 19*37, replace = T), dim = c(lat = 19, lon = 37)) + +## Observations or reanalyses + dat1 <- array(1:100, dim = c(year = 5, lat = 19, lon = 37, month = 12)) + + ## Historical simulations + dat2 <- array(1:100, dim = c(year = 5, lat = 19, lon = 37, month = 12, member = 5)) + + ## Decadal predictions + dat3 <- array(1:100, dim = c(sdate = 5, lat = 19, lon = 37, fmonth = 24, member = 5)) + +############################################## +test_that("1. Input checks", { + + expect_error( + GSAT(data = c()), + "Parameter 'data' cannot be NULL." + ) + expect_error( + GSAT(data = 'a'), + "Parameter 'data' must be a numeric array." + ) + expect_error( + GSAT(data = dat1, data_lats = c()), + "Parameter 'data_lats' must be a numeric vector." + ) + expect_error( + GSAT(data = dat1, data_lats = lat, data_lons = c()), + "Parameter 'data_lons' must be a numeric vector." + ) + expect_error( + GSAT(data = dat1, data_lats = lat, data_lons = lon, type = TRUE), + "Parameter 'type' must be 'dcpp', 'hist', or 'obs'." + ) + expect_error( + GSAT(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + lat_dim = 2), + "Parameter 'lat_dim' must be a character string." + ) + expect_error( + GSAT(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + lat_dim = 'lats'), + "Parameter 'lat_dim' is not found in 'data' dimension." + ) + expect_error( + GSAT(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + lon_dim = 2), + "Parameter 'lon_dim' must be a character string." + ) + expect_error( + GSAT(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + lon_dim = 'lons'), + "Parameter 'lon_dim' is not found in 'data' dimension." + ) + expect_error( + GSAT(data = dat3, data_lats = seq(-90, 90, 11), data_lons = lon, type = 'dcpp'), + paste0("The latitude dimension of parameter 'data' must be the same", + " length of parameter 'data_lats'.") + ) + expect_error( + GSAT(data = dat3, data_lats = lat, data_lons = c(1:10), type = 'dcpp'), + paste0("The longitude dimension of parameter 'data' must be the same", + " length of parameter 'data_lons'.") + ) + expect_error( + GSAT(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + mask = c(1, 0, 1, 0, 0)), + paste0("Parameter 'mask' must be NULL \\(no mask\\) or a numerical array ", + "with c\\(lat_dim, lon_dim\\) dimensions and 0 in those grid ", + "points that have to be masked.") + ) + expect_error( + GSAT(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + mask = array(rep(c(1, 0), 100), dim = c(lat = 20, lon = 5))), + paste0("Parameter 'mask' must be NULL \\(no mask\\) or a numerical array ", + "with c\\(lat_dim, lon_dim\\) dimensions and 0 in those grid ", + "points that have to be masked.") + ) + expect_error( + GSAT(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + monini = 0), + "Parameter 'monini' must be an integer from 1 to 12." + ) + expect_error( + GSAT(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + fmonth_dim = c('lat', 'lon')), + "Parameter 'fmonth_dim' must be a character string." + ) + expect_error( + GSAT(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + fmonth_dim = 'asd'), + "Parameter 'fmonth_dim' is not found in 'data' dimension." + ) + expect_error( + GSAT(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + sdate_dim = 2), + "Parameter 'sdate_dim' must be a character string." + ) + expect_error( + GSAT(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + sdate_dim = 'date'), + "Parameter 'sdate_dim' is not found in 'data' dimension." + ) + expect_error( + GSAT(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + indices_for_clim = TRUE), + paste0("The parameter 'indices_for_clim' must be a numeric vector ", + "or NULL to compute the anomalies based on the whole period, ", + "or FALSE if data are already anomalies") + ) + expect_error( + GSAT(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + year_dim = 1), + "Parameter 'year_dim' must be a character string." + ) + expect_error( + GSAT(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + year_dim = 'yr'), + "Parameter 'year_dim' is not found in 'data' dimension." + ) + expect_error( + GSAT(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + month_dim = 0), + "Parameter 'month_dim' must be a character string." + ) + expect_error( + GSAT(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + month_dim = 'mth'), + "Parameter 'month_dim' is not found in 'data' dimension." + ) + expect_error( + GSAT(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + member_dim = 1), + "Parameter 'member_dim' must be a character string." + ) + expect_error( + GSAT(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + member_dim = 'ens'), + "Parameter 'member_dim' is not found in 'data' dimension." + ) + +}) + +############################################## +test_that("2. Output checks: dat1", { + res1_1 <- GSAT(data = dat1, data_lats = lat, data_lons = lon, type = 'obs') + + expect_equal( + dim(res1_1), + c(year = 5) + ) + expect_equal( + mean(res1_1)*10^15, + 1.421042, + tolerance = 0.00001 + ) + + res1_2 <- GSAT(data = dat1, data_lats = lat, data_lons = lon, type = 'obs', + indices_for_clim = F) + expect_equal( + mean(res1_2), + 50.42349, + tolerance = 0.00001 + ) + +}) +############################################## +test_that("3. Output checks: dat2", { + res2_1 <- GSAT(data = dat2, data_lats = lat, data_lons = lon, type = 'hist') + expect_equal( + dim(res2_1), + c(year = 5, member = 5) + ) + expect_equal( + mean(res2_1)*10^16, + 2.841867, + tolerance = 0.00001 + ) + + res2_2 <- GSAT(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + indices_for_clim = 2:3) + expect_equal( + mean(res2_2), + 0.5, + tolerance = 0.00001 + ) + + res2_3 <- GSAT(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + mask = mask) + expect_equal( + mean(res2_3)*10^15, + 1.136842, + tolerance = 0.00001 + ) + +}) +############################################## +test_that("4. Output checks: dat3", { + res3_1 <- GSAT(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp') + expect_equal( + dim(res3_1), + c(fyear = 1, sdate = 5, member = 5) + ) + expect_equal( + mean(res3_1)*10^16, + -5.684776, + tolerance = 0.00001 + ) + + res3_2 <- GSAT(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + indices_for_clim = 2:3) + expect_equal( + mean(res3_2), + 0.5, + tolerance = 0.00001 + ) + + res3_3 <- GSAT(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + monini = 5) + expect_equal( + mean(res3_3)*10^16, + 2.841867, + tolerance = 0.00001 + ) + + res3_4 <- GSAT(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + mask = mask) + expect_equal( + mean(res3_4), + 0, + tolerance = 0.00001 + ) + +}) diff --git a/tests/testthat/test-SPOD.R b/tests/testthat/test-SPOD.R new file mode 100644 index 0000000000000000000000000000000000000000..6d23a59745158a7d3c36807244cf96ceb453608c --- /dev/null +++ b/tests/testthat/test-SPOD.R @@ -0,0 +1,240 @@ +context("s2dv::SPOD tests") + +############################################## + lat <- seq(-90, 90, 10) + lon <- seq(0, 360, 10) + set.seed(1) + mask <- array(sample(c(0, 1), 19*37, replace = T), dim = c(lat = 19, lon = 37)) + +## Observations or reanalyses + dat1 <- array(1:100, dim = c(year = 5, lat = 19, lon = 37, month = 12)) + + ## Historical simulations + dat2 <- array(1:100, dim = c(year = 5, lat = 19, lon = 37, month = 12, member = 5)) + + ## Decadal predictions + dat3 <- array(1:100, dim = c(sdate = 5, lat = 19, lon = 37, fmonth = 24, member = 5)) + +############################################## +test_that("1. Input checks", { + + expect_error( + SPOD(data = c()), + "Parameter 'data' cannot be NULL." + ) + expect_error( + SPOD(data = 'a'), + "Parameter 'data' must be a numeric array." + ) + expect_error( + SPOD(data = dat1, data_lats = c()), + "Parameter 'data_lats' must be a numeric vector." + ) + expect_error( + SPOD(data = dat1, data_lats = lat, data_lons = c()), + "Parameter 'data_lons' must be a numeric vector." + ) + expect_error( + SPOD(data = dat1, data_lats = lat, data_lons = lon, type = TRUE), + "Parameter 'type' must be 'dcpp', 'hist', or 'obs'." + ) + expect_error( + SPOD(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + lat_dim = 2), + "Parameter 'lat_dim' must be a character string." + ) + expect_error( + SPOD(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + lat_dim = 'lats'), + "Parameter 'lat_dim' is not found in 'data' dimension." + ) + expect_error( + SPOD(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + lon_dim = 2), + "Parameter 'lon_dim' must be a character string." + ) + expect_error( + SPOD(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + lon_dim = 'lons'), + "Parameter 'lon_dim' is not found in 'data' dimension." + ) + expect_error( + SPOD(data = dat3, data_lats = seq(-90, 90, 11), data_lons = lon, type = 'dcpp'), + paste0("The latitude dimension of parameter 'data' must be the same", + " length of parameter 'data_lats'.") + ) + expect_error( + SPOD(data = dat3, data_lats = lat, data_lons = c(1:10), type = 'dcpp'), + paste0("The longitude dimension of parameter 'data' must be the same", + " length of parameter 'data_lons'.") + ) + expect_error( + SPOD(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + mask = c(1, 0, 1, 0, 0)), + paste0("Parameter 'mask' must be NULL \\(no mask\\) or a numerical array ", + "with c\\(lat_dim, lon_dim\\) dimensions and 0 in those grid ", + "points that have to be masked.") + ) + expect_error( + SPOD(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + mask = array(rep(c(1, 0), 100), dim = c(lat = 20, lon = 5))), + paste0("Parameter 'mask' must be NULL \\(no mask\\) or a numerical array ", + "with c\\(lat_dim, lon_dim\\) dimensions and 0 in those grid ", + "points that have to be masked.") + ) + expect_error( + SPOD(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + monini = 0), + "Parameter 'monini' must be an integer from 1 to 12." + ) + expect_error( + SPOD(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + fmonth_dim = c('lat', 'lon')), + "Parameter 'fmonth_dim' must be a character string." + ) + expect_error( + SPOD(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + fmonth_dim = 'asd'), + "Parameter 'fmonth_dim' is not found in 'data' dimension." + ) + expect_error( + SPOD(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + sdate_dim = 2), + "Parameter 'sdate_dim' must be a character string." + ) + expect_error( + SPOD(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + sdate_dim = 'date'), + "Parameter 'sdate_dim' is not found in 'data' dimension." + ) + expect_error( + SPOD(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + indices_for_clim = TRUE), + paste0("The parameter 'indices_for_clim' must be a numeric vector ", + "or NULL to compute the anomalies based on the whole period, ", + "or FALSE if data are already anomalies") + ) + expect_error( + SPOD(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + year_dim = 1), + "Parameter 'year_dim' must be a character string." + ) + expect_error( + SPOD(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + year_dim = 'yr'), + "Parameter 'year_dim' is not found in 'data' dimension." + ) + expect_error( + SPOD(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + month_dim = 0), + "Parameter 'month_dim' must be a character string." + ) + expect_error( + SPOD(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + month_dim = 'mth'), + "Parameter 'month_dim' is not found in 'data' dimension." + ) + expect_error( + SPOD(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + member_dim = 1), + "Parameter 'member_dim' must be a character string." + ) + expect_error( + SPOD(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + member_dim = 'ens'), + "Parameter 'member_dim' is not found in 'data' dimension." + ) + +}) + +############################################## +test_that("2. Output checks: dat1", { + res1_1 <- SPOD(data = dat1, data_lats = lat, data_lons = lon, type = 'obs') + + expect_equal( + dim(res1_1), + c(year = 5) + ) + expect_equal( + mean(res1_1)*10^16, + 1.776357, + tolerance = 0.00001 + ) + + res1_2 <- SPOD(data = dat1, data_lats = lat, data_lons = lon, type = 'obs', + indices_for_clim = F) + expect_equal( + mean(res1_2), + 6.073599, + tolerance = 0.00001 + ) + +}) +############################################## +test_that("3. Output checks: dat2", { + res2_1 <- SPOD(data = dat2, data_lats = lat, data_lons = lon, type = 'hist') + expect_equal( + dim(res2_1), + c(year = 5, member = 5) + ) + expect_equal( + mean(res2_1)*10^16, + 2.754588, + tolerance = 0.00001 + ) + + res2_2 <- SPOD(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + indices_for_clim = 2:3) + expect_equal( + mean(res2_2)*10^16, + 2.645063, + tolerance = 0.00001 + ) + + res2_3 <- SPOD(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + mask = mask) + expect_equal( + mean(res2_3)*10^17, + 5.3553, + tolerance = 0.00001 + ) + +}) +############################################## +test_that("4. Output checks: dat3", { + res3_1 <- SPOD(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp') + expect_equal( + dim(res3_1), + c(fyear = 1, sdate = 5, member = 5) + ) + expect_equal( + mean(res3_1)*10^17, + -3.552714, + tolerance = 0.00001 + ) + + res3_2 <- SPOD(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + indices_for_clim = 2:3) + expect_equal( + mean(res3_2)*10^16, + 1.421085, + tolerance = 0.00001 + ) + + res3_3 <- SPOD(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + monini = 5) + expect_equal( + mean(res3_3)*10^16, + 2.087219, + tolerance = 0.00001 + ) + + res3_4 <- SPOD(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + mask = mask) + expect_equal( + mean(res3_4)*10^17, + 7.549517, + tolerance = 0.00001 + ) + +}) diff --git a/tests/testthat/test-TPI.R b/tests/testthat/test-TPI.R new file mode 100644 index 0000000000000000000000000000000000000000..3b584c6f5a3eaf5d038f9f20a04c00fcdfe2e39f --- /dev/null +++ b/tests/testthat/test-TPI.R @@ -0,0 +1,245 @@ +context("s2dv::TPI tests") + +############################################## + lat <- seq(-90, 90, 10) + lon <- seq(0, 360, 10) + set.seed(1) + mask <- array(sample(c(0, 1), 19*37, replace = T), dim = c(lat = 19, lon = 37)) + +## Observations or reanalyses + dat1 <- array(1:100, dim = c(year = 5, lat = 19, lon = 37, month = 12)) + + ## Historical simulations + dat2 <- array(1:100, dim = c(year = 5, lat = 19, lon = 37, month = 12, member = 5)) + + ## Decadal predictions + dat3 <- array(1:100, dim = c(sdate = 5, lat = 19, lon = 37, fmonth = 24, member = 5)) + +############################################## +test_that("1. Input checks", { + + expect_error( + TPI(data = c()), + "Parameter 'data' cannot be NULL." + ) + expect_error( + TPI(data = 'a'), + "Parameter 'data' must be a numeric array." + ) + expect_error( + TPI(data = dat1, data_lats = c()), + "Parameter 'data_lats' must be a numeric vector." + ) + expect_error( + TPI(data = dat1, data_lats = lat, data_lons = c()), + "Parameter 'data_lons' must be a numeric vector." + ) + expect_error( + TPI(data = dat1, data_lats = lat, data_lons = lon, type = TRUE), + "Parameter 'type' must be 'dcpp', 'hist', or 'obs'." + ) + expect_error( + TPI(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + lat_dim = 2), + "Parameter 'lat_dim' must be a character string." + ) + expect_error( + TPI(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + lat_dim = 'lats'), + "Parameter 'lat_dim' is not found in 'data' dimension." + ) + expect_error( + TPI(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + lon_dim = 2), + "Parameter 'lon_dim' must be a character string." + ) + expect_error( + TPI(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + lon_dim = 'lons'), + "Parameter 'lon_dim' is not found in 'data' dimension." + ) + expect_error( + TPI(data = dat3, data_lats = seq(-90, 90, 11), data_lons = lon, type = 'dcpp'), + paste0("The latitude dimension of parameter 'data' must be the same", + " length of parameter 'data_lats'.") + ) + expect_error( + TPI(data = dat3, data_lats = lat, data_lons = c(1:10), type = 'dcpp'), + paste0("The longitude dimension of parameter 'data' must be the same", + " length of parameter 'data_lons'.") + ) + expect_error( + TPI(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + mask = c(1, 0, 1, 0, 0)), + paste0("Parameter 'mask' must be NULL \\(no mask\\) or a numerical array ", + "with c\\(lat_dim, lon_dim\\) dimensions and 0 in those grid ", + "points that have to be masked.") + ) + expect_error( + TPI(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + mask = array(rep(c(1, 0), 100), dim = c(lat = 20, lon = 5))), + paste0("Parameter 'mask' must be NULL \\(no mask\\) or a numerical array ", + "with c\\(lat_dim, lon_dim\\) dimensions and 0 in those grid ", + "points that have to be masked.") + ) + expect_error( + TPI(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + monini = 0), + "Parameter 'monini' must be an integer from 1 to 12." + ) + expect_error( + TPI(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + fmonth_dim = c('lat', 'lon')), + "Parameter 'fmonth_dim' must be a character string." + ) + expect_error( + TPI(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + fmonth_dim = 'asd'), + "Parameter 'fmonth_dim' is not found in 'data' dimension." + ) + expect_error( + TPI(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + sdate_dim = 2), + "Parameter 'sdate_dim' must be a character string." + ) + expect_error( + TPI(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + sdate_dim = 'date'), + "Parameter 'sdate_dim' is not found in 'data' dimension." + ) + expect_error( + TPI(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + indices_for_clim = TRUE), + paste0("The parameter 'indices_for_clim' must be a numeric vector ", + "or NULL to compute the anomalies based on the whole period, ", + "or FALSE if data are already anomalies") + ) + expect_error( + TPI(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + year_dim = 1), + "Parameter 'year_dim' must be a character string." + ) + expect_error( + TPI(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + year_dim = 'yr'), + "Parameter 'year_dim' is not found in 'data' dimension." + ) + expect_error( + TPI(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + month_dim = 0), + "Parameter 'month_dim' must be a character string." + ) + expect_error( + TPI(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + month_dim = 'mth'), + "Parameter 'month_dim' is not found in 'data' dimension." + ) + expect_error( + TPI(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + member_dim = 1), + "Parameter 'member_dim' must be a character string." + ) + expect_error( + TPI(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + member_dim = 'ens'), + "Parameter 'member_dim' is not found in 'data' dimension." + ) + +}) + +############################################## +test_that("2. Output checks: dat1", { + res1_1 <- TPI(data = dat1, data_lats = lat, data_lons = lon, type = 'obs') + + expect_equal( + dim(res1_1), + c(year = 5) + ) + expect_equal( + mean(res1_1), + 0, + tolerance = 0.00001 + ) + expect_equal( + res1_1*10^16, + array(c(0.000000, -12.21245, -5.551115, 12.21245, 5.551115), dim = c(year = 5)), + tolerance = 0.00001 + ) + + res1_2 <- TPI(data = dat1, data_lats = lat, data_lons = lon, type = 'obs', + indices_for_clim = F) + expect_equal( + mean(res1_2), + 0.9376566, + tolerance = 0.00001 + ) + +}) +############################################## +test_that("3. Output checks: dat2", { + res2_1 <- TPI(data = dat2, data_lats = lat, data_lons = lon, type = 'hist') + expect_equal( + dim(res2_1), + c(year = 5, member = 5) + ) + expect_equal( + mean(res2_1)*10^17, + -2.220099, + tolerance = 0.00001 + ) + + res2_2 <- TPI(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + indices_for_clim = 2:3) + expect_equal( + mean(res2_2)*10^16, + 3.840678, + tolerance = 0.00001 + ) + + res2_3 <- TPI(data = dat2, data_lats = lat, data_lons = lon, type = 'hist', + mask = mask) + expect_equal( + mean(res2_3)*10^17, + -2.288135, + tolerance = 0.00001 + ) + +}) +############################################## +test_that("4. Output checks: dat3", { + res3_1 <- TPI(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp') + expect_equal( + dim(res3_1), + c(fyear = 1, sdate = 5, member = 5) + ) + expect_equal( + mean(res3_1)*10^17, + 1.776357, + tolerance = 0.00001 + ) + + res3_2 <- TPI(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + indices_for_clim = 2:3) + expect_equal( + mean(res3_2)*10^16, + 3.064216, + tolerance = 0.00001 + ) + + res3_3 <- TPI(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + monini = 5) + expect_equal( + mean(res3_3)*10^17, + -2.664535, + tolerance = 0.00001 + ) + + res3_4 <- TPI(data = dat3, data_lats = lat, data_lons = lon, type = 'dcpp', + mask = mask) + expect_equal( + mean(res3_4)*10^17, + -1.360023, + tolerance = 0.00001 + ) + +})