Newer
Older
#'Compute forecast or observation anomalies
#'
#'This function computes anomalies from a multidimensional data array and a
#'climatology array.
#'
#'@param data A numeric array with named dimensions, representing the model or
#' observational data to be calculated the anomalies. It should involve all
#' the dimensions in parameter 'clim', and it can have more other dimensions.
#'@param clim A numeric array with named dimensions, representing the
#' climatologies to be deducted from parameter 'data'. It can be generated by
#' Clim(). The dimensions should all be involved in parameter 'data' with the
#' same length.
#'@param ncores An integer indicating the number of cores to use for parallel
#' computation. The default value is NULL.
#'
aho
committed
#'@return An array with same dimensions as parameter 'data'.
#'
#'@examples
#'# Load sample data as in Load() example:
#'example(Load)
#'clim <- Clim(sampleData$mod, sampleData$obs)
#'ano_exp <- Ano(sampleData$mod, clim$clim_exp)
#'ano_obs <- Ano(sampleData$obs, clim$clim_obs)
#'\donttest{
#'PlotAno(ano_exp, ano_obs, startDates,
#' toptitle = 'Anomaly', ytitle = c('K', 'K', 'K'),
#' legends = 'ERSST', biglab = FALSE, fileout = 'tos_ano.png')
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
#'}
#'@import multiApply
#'@export
Ano <- function(data, clim, ncores = NULL) {
# Check inputs
## data
if (is.null(data)) {
stop("Parameter 'data' cannot be NULL.")
}
if (!is.numeric(data)) {
stop("Parameter 'data' must be a numeric array.")
}
if (is.null(dim(data))) { #is vector
dim(data) <- c(length(data))
names(dim(data)) <- 'tmp_name'
}
if(any(is.null(names(dim(data))))| any(nchar(names(dim(data))) == 0)) {
stop("Parameter 'data' must have dimension names.")
}
## clim
if (is.null(clim)) {
stop("Parameter 'clim' cannot be NULL.")
}
if (!is.numeric(clim)) {
stop("Parameter 'clim' must be a numeric array.")
}
if (is.null(dim(clim))) { #is vector
dim(clim) <- c(length(clim))
names(dim(clim)) <- 'tmp_name'
}
if (any(is.null(names(dim(clim))))| any(nchar(names(dim(clim))) == 0)) {
stop("Parameter 'clim' must have dimension names.")
}
aho
committed
## data and clim
if (!all(names(dim(clim)) %in% names(dim(data)))) {
stop("Parameter 'data' must have all the dimensions of parameter 'clim'.")
} else {
pos <- names(dim(data))[match(names(dim(clim)), names(dim(data)))]
if (any(dim(clim) != dim(data)[pos])) {
stop("Some dimensions of parameter 'clim' have different length from parameter 'data'.")
}
}
## ncores
if (!is.null(ncores)) {
aho
committed
if (!is.numeric(ncores)) {
stop("Parameter 'ncores' must be a positive integer.")
} else if (ncores %% 1 != 0 | ncores <= 0 | length(ncores) > 1) {
stop("Parameter 'ncores' must be a positive integer.")
}
}
###############################
# Calculate Ano
aho
committed
parallel_compute <- TRUE
if (is.null(ncores)) {
parallel_compute <- FALSE
} else if (ncores == 1) {
parallel_compute <- FALSE
}
if (!parallel_compute) {
target_dims_ind <- match(names(dim(clim)), names(dim(data)))
if (any(target_dims_ind != sort(target_dims_ind))) {
clim <- Reorder(clim, match(sort(target_dims_ind), target_dims_ind))
}
if (length(dim(data)) == length(dim(clim))) {
res <- data - clim
} else {
target_dims_ind <- match(dim(clim), dim(data))
margin_dims_ind <- c(1:length(dim(data)))[-target_dims_ind]
res <- apply(data, margin_dims_ind, .Ano, clim)
res <- array(res, dim = dim(data)[c(target_dims_ind, margin_dims_ind)])
}
} else {
res <- Apply(list(data),
target_dims = names(dim(clim)),
output_dims = names(dim(clim)),
fun = .Ano,
clim = clim,
ncores = ncores)$output1
}
aho
committed
# Reorder dim back to data
if (any(dim(res) != dim(data))) {
res <- Reorder(res, names(dim(data)))
}
aho
committed
return(invisible(res))
}
.Ano <- function(data, clim) {
aho
committed
data - clim