diff --git a/inst/doc/usecase.md b/inst/doc/usecase.md index 06f2c044eb64673cedaf1141a00c0cc69bc2b5b5..c5f8fe439311fb1df4f9a05a8a9a4f6e8903726b 100644 --- a/inst/doc/usecase.md +++ b/inst/doc/usecase.md @@ -12,13 +12,14 @@ In this document, you can link to the example scripts for various demands. For t data is one file per year, each file contains 12 months (time = 12). However, observational data is one file per month, each file contains only one time step. You can learn how to select all the required year and month for observation, and - tweak the dimension to make it consistent with experiment. + twist the dimensions to make them consistent with experiment. The highlight paramters used in this usecase are: **'*_across'**, **'merge_across_dims'**, and **'split_multiselected_dims'**. 3. [Use experimental data attribute to load in oberservational data](inst/doc/usecase/ex1_3_attr_loadin.R) Like ex1_2, it shows how to retrieve the experimental data and observational data in a comparable structure. It also shows how to use parameters `xxx_tolerance`, `xxx_across`, `merge_across_dims`, `merge_across_dims_narm`, and `split_multiselected_dims`. +It is recommended reading ex1_2 first since there is more explanation. 4. [Checking impact of start date order in the number of members](inst/doc/usecase/ex1_4_variable_nmember.R) Mixing start dates of different months can lead to load different number of members, check the code provided and the [FAQ 10](/inst/doc/faq.md). diff --git a/inst/doc/usecase/ex1_2_exp_obs_attr.R b/inst/doc/usecase/ex1_2_exp_obs_attr.R index 2e4b5da965ba642081bf0e362c7b873004f9d9e7..fa2323d6aaa5ace57195ea52e127a8ad1fa1b3c7 100644 --- a/inst/doc/usecase/ex1_2_exp_obs_attr.R +++ b/inst/doc/usecase/ex1_2_exp_obs_attr.R @@ -1,14 +1,24 @@ #--------------------------------------------------------------------- # This script tells you how to load experimental and observational data in a -# consistent way, facilating the following comparison. - -# First, we load the experimental data. Because the latitude order of observation -# is opposite with experiment, and the sdate/time dimension is also different, we -# use the attributes (sdate and latitude) of experimental data to define the -# selectors for observation. - -# You can see how to use parameter '*_across', 'merge_across_dims', and -# 'split_multiselected_dims' to create the consistent dimension as experiment. +# consistent way, facilating the following comparison. We use the attributes of +# the experimental data to define the selectors of obs Start() call, so they +# can have the same dimension structure. + +# Spatial dimensions: +# The exp and obs data happen to have the same spatial resolution (256x512) and +# the grids are not shifted, so we don't need to regrid them. However, their latitude +# orders are opposite. exp has ascending order while obs has descending order. +# To make them consistent, we cannot simply use 'all' as the selector of obs because +# for now, the reordering parameter '*_Sort' is only functional when the dimension is +# defined by values(). We can use either `indices(256:1)` or the exp attributes (`values()`) +# to define the latitude of obs. + +# Temporal dimensions: +# The exp and obs files have different date/time structure. exp has one file per year and +# each file has 12 time steps. obs has one file per month and each file has 1 time step. +# To shape the obs array as exp, we need to use the time attribute of exp to define +# the date/time selector of obs. You can see how to use parameter '*_across', 'merge_across_dims', +# and 'split_multiselected_dims' to achieve the goal. #--------------------------------------------------------------------- library(startR) @@ -31,30 +41,41 @@ exp <- Start(dat = repos_exp, time = 'sdate'), retrieve = FALSE) -# Retrieve attributes for the following observation. -# Because latitude order in experiment is [-90, 90] but in observation is [90, -90], -# latitude values need to be retrieved and used below. +attr(exp, 'Dimensions') +# dat var sdate time lat lon +# 1 1 4 3 256 512 + +# Retrieve attributes for observational data retrieval. +## Because latitude order in experiment is [-90, 90] but in observation is [90, -90], +## latitude values need to be retrieved and used below. lats <- attr(exp, 'Variables')$common$lat -# The 'time' attribute is dependent on 'sdate'. You can see the dimension below. +lons <- attr(exp, 'Variables')$common$lon +## The 'time' attribute is a two-dim array dates <- attr(exp, 'Variables')$common$time -# dim(dates) +dim(dates) #sdate time # 4 3 +dates +# [1] "2005-01-16 12:00:00 UTC" "2006-01-16 12:00:00 UTC" +# [3] "2007-01-16 12:00:00 UTC" "2008-01-16 12:00:00 UTC" +# [5] "2005-02-15 00:00:00 UTC" "2006-02-15 00:00:00 UTC" +# [7] "2007-02-15 00:00:00 UTC" "2008-02-15 12:00:00 UTC" +# [9] "2005-03-16 12:00:00 UTC" "2006-03-16 12:00:00 UTC" +#[11] "2007-03-16 12:00:00 UTC" "2008-03-16 12:00:00 UTC" #------------------------------------------- # obs -# 1. For lat, use experiment attribute. For lon, it is not necessary because they have -# same values. -# 2. For dimension 'date', it is a vector involving the first 3 months (ftime) of the four years (sdate). -# 3. Dimension 'time' is assigned by the matrix, so we can seperate 'sdate' and 'time' -# using 'split_multiselected_dims' later. -# 4. Because the 'time' is actually across all the files, so we need to specify -# 'time_across'. Then, use 'merge_across_dims' to make dimension 'date' disappears. +# 1. For lat, use the experiment attribute or reversed indices. For lon, it is not necessary +# because their lons are identical, but either way works. +# 2. For dimension 'date', it is a vector involving the 3 months (ftime) of the four years (sdate). +# 3. Dimension 'time' is assigned by the matrix, so we can split it into 'sdate' and 'time' +# by 'split_multiselected_dims'. +# 4. Because 'time' is actually across all the files, so we need to specify 'time_across'. +# Then, use 'merge_across_dims' to make dimension 'date' disappears. # At this moment, the dimension is 'time = 12'. # 5. However, we want to seperate year and month (which are 'sdate' and 'ftime' in -# experimental data). So we use 'split_multiselected_dims' to split the two dimensions -# of dimension 'time'. +# experimental data). So we use 'split_multiselected_dims' to split 'time' into the two dimensions. repos_obs <- '/esarchive/recon/ecmwf/erainterim/monthly_mean/$var$_f6h/$var$_$date$.nc' @@ -62,8 +83,8 @@ obs <- Start(dat = repos_obs, var = 'tas', date = unique(format(dates, '%Y%m')), time = values(dates), #dim: [sdate = 4, time = 3] - lat = values(lats), - lon = 'all', + lat = values(lats), # indices(256:1), + lon = values(lons), # 'all', time_across = 'date', merge_across_dims = TRUE, split_multiselected_dims = TRUE, @@ -74,9 +95,13 @@ obs <- Start(dat = repos_obs, time = 'date'), retrieve = FALSE) -#========================== -# Check attributes -#========================== +attr(obs, 'Dimensions') +# dat var sdate time lat lon +# 1 1 4 3 256 512 + +#==================================================== +# Check the attributes. They should be all identical. +#==================================================== ##-----dimension----- print(attr(exp, 'Dimensions')) @@ -87,7 +112,8 @@ print(attr(obs, 'Dimensions')) # dat var sdate time lat lon # 1 1 4 3 256 512 -##-----time----- +##-----time----- +## They're not identical but the years and months are. See below for more details. print(attr(exp, 'Variables')$common$time) # [1] "2005-01-16 12:00:00 UTC" "2006-01-16 12:00:00 UTC" # [3] "2007-01-16 12:00:00 UTC" "2008-01-16 12:00:00 UTC" @@ -114,3 +140,92 @@ print(attr(obs, 'Variables')$common$lat[1:3]) print(attr(obs, 'Variables')$common$lat[256]) #[1] 89.46282 +##-----lon----- +print(attr(exp, 'Variables')$common$lon[1:3]) +#[1] 0.000000 0.703125 1.406250 +print(attr(obs, 'Variables')$common$lon[1:3]) +#[1] 0.000000 0.703125 1.406250 + + +#======================= +# About time attributes +#======================= +# You may notice that the date and time of exp and obs are not the same. +# In this case, the data are monthly data, so only the years and months matter. +# The thing worth noticing is that the actual time values of obs are half month different +# from the values we assigned. For example, the first time from exp is "2005-01-16 12:00:00 UTC", +# and the obs time we get is "2005-01-31 18:00:00 UTC". +# If the provided selector is value, Start() looks for the closest value in the data. +# So for "2005-01-16 12:00:00 UTC", the two closest obs values are "2004-12-31 18:00:00 UTC" and +# "2005-01-31 18:00:00 UTC", and the later one is the closest and happen to be the desired one. +# It's fortunate that in this case, all the provided values are closer to the values we want. + +#----- 1. Manually adjust the values ----- +# It is always necessary to check the data attributes before and after data retrieval. +# If the provided exp values are quite in the middle of two values in obs, we can adjust a bit to +# make exp values closer to the desired obs values. +dates_adjust <- dates + 86400*15 +dates_adjust +# [1] "2005-01-31 12:00:00 UTC" "2006-01-31 12:00:00 UTC" +# [3] "2007-01-31 12:00:00 UTC" "2008-01-31 12:00:00 UTC" +# [5] "2005-03-02 00:00:00 UTC" "2006-03-02 00:00:00 UTC" +# [7] "2007-03-02 00:00:00 UTC" "2008-03-01 12:00:00 UTC" +# [9] "2005-03-31 12:00:00 UTC" "2006-03-31 12:00:00 UTC" +#[11] "2007-03-31 12:00:00 UTC" "2008-03-31 12:00:00 UTC" + +obs2 <- Start(dat = repos_obs, + var = 'tas', + date = unique(format(dates, '%Y%m')), + time = values(dates_adjust), # use the adjust ones + lat = values(lats), + lon = values(lons), + time_across = 'date', + merge_across_dims = TRUE, + split_multiselected_dims = TRUE, + synonims = list(lat = c('lat', 'latitude'), + lon = c('lon', 'longitude')), + return_vars = list(lon = NULL, + lat = NULL, + time = 'date'), + retrieve = FALSE) + +# The time should be the same as obs above. +print(attr(obs2, 'Variables')$common$time) +# [1] "2005-01-31 18:00:00 UTC" "2006-01-31 18:00:00 UTC" +# [3] "2007-01-31 18:00:00 UTC" "2008-01-31 18:00:00 UTC" +# [5] "2005-02-28 18:00:00 UTC" "2006-02-28 18:00:00 UTC" +# [7] "2007-02-28 18:00:00 UTC" "2008-02-29 18:00:00 UTC" +# [9] "2005-03-31 18:00:00 UTC" "2006-03-31 18:00:00 UTC" +#[11] "2007-03-31 18:00:00 UTC" "2008-03-31 18:00:00 UTC" + +#----- 2. Set the tolerance ----- +# Sometimes, it may be useful to set the tolerance. If the provided values are too much different +# from the values in obs, Start() returns an error directly (if none of the data found) or returns +# incorrect time attributes. + +obs3 <- Start(dat = repos_obs, + var = 'tas', + date = unique(format(dates, '%Y%m')), + time = values(dates), + lat = values(lats), + lon = values(lons), + time_across = 'date', + time_tolerance = as.difftime(15, units = 'days'), + merge_across_dims = TRUE, + split_multiselected_dims = TRUE, + synonims = list(lat = c('lat', 'latitude'), + lon = c('lon', 'longitude')), + return_vars = list(lon = NULL, + lat = NULL, + time = 'date'), + retrieve = FALSE) + +# We lose many data because there are no data within 15 days from the provided time values. +print(attr(obs3, 'Variables')$common$time) +[1] "2005-02-28 18:00:00 UTC" "2006-02-28 18:00:00 UTC" +[3] "2007-02-28 18:00:00 UTC" "2008-02-29 18:00:00 UTC" + +# If 'time_tolerance' is changed to "as.difftime(1, units = 'days')", an error shows: +# Selectors do not match any of the possible values for the dimension 'time'. + + diff --git a/inst/doc/usecase/ex1_3_attr_loadin.R b/inst/doc/usecase/ex1_3_attr_loadin.R index e2c821117a82612e842d29d3a90cdf3869f5c2cd..80690eecbdb84734f825af380c8ece97b62e5b9b 100644 --- a/inst/doc/usecase/ex1_3_attr_loadin.R +++ b/inst/doc/usecase/ex1_3_attr_loadin.R @@ -2,10 +2,13 @@ # This usecase shows you how to load experimental and observational data in a # consistent way. -# First, load the experimental data. Then, use the time attributes of experimental data to define the selectors for observational data. +# First, load the experimental data. Then, use the time attributes of experimental data to +# define the selectors for observational data. # You can see how to use parameter '*_across', 'merge_across_dims', 'merge_across_dims_narm', # and 'split_multiselected_dims' to create the same dimension structure. + +# If you haven't read ex1_2, it is recommended reading it first since there is more explanation. #--------------------------------------------------------------------- # experimental data diff --git a/tests/testthat/test-Compute-CDORemap.R b/tests/testthat/test-Compute-CDORemap.R new file mode 100644 index 0000000000000000000000000000000000000000..67f2a1005654ab7c0592cd88d6d808bf484a4ff5 --- /dev/null +++ b/tests/testthat/test-Compute-CDORemap.R @@ -0,0 +1,53 @@ +context("Compute use CDORemap") + +test_that("ex2_3", { + + repos <- '/esarchive/exp/ecmwf/system5_m1/monthly_mean/$var$_f6h/$var$_$sdate$.nc' + data <- Start(dat = repos, + var = 'tas', + sdate = c('20170101'), + ensemble = indices(1), + time = indices(1), + latitude = 'all', + longitude = 'all', + return_vars = list(latitude = 'dat', longitude = 'dat', time = 'sdate'), + retrieve = FALSE) + + + fun <- function(x) { + lons_data <- as.vector(attr(x, 'Variables')$dat1$longitude) + lats_data <- as.vector(attr(x, 'Variables')$dat1$latitude) + r <- s2dverification::CDORemap(x, lons_data, lats_data, "r360x181", + 'bil', crop = FALSE, force_remap = TRUE)[[1]] + return(r) + } + + step3 <- Step(fun = fun, + target_dims = c('latitude','longitude'), + output_dims = c('latitude', 'longitude'), + use_attributes = list(data = "Variables")) + wf3 <- AddStep(list(data = data), step3) + + res3 <- Compute(workflow = wf3, + chunks = list(ensemble = 1)) + +expect_equal( +attr(data, 'Dimensions'), +c(dat = 1, var = 1, sdate = 1, ensemble = 1, time = 1, latitude = 640, longitude = 1296) +) +expect_equal( +dim(res3$output), +c(latitude = 181, longitude = 360, dat = 1, var = 1, sdate = 1, ensemble = 1, time = 1) +) +expect_equal( +mean(res3$output), +277.0346, +tolerance = 0.0001 +) +expect_equal( +res3$output[20,11,1,,1,1,1], +c(265.5362), +tolerance = 0.0001 +) + +}) diff --git a/tests/testthat/test-Compute-extra_params.R b/tests/testthat/test-Compute-extra_params.R new file mode 100644 index 0000000000000000000000000000000000000000..9dae43c1be0614e09b112ea9896b49773b59c453 --- /dev/null +++ b/tests/testthat/test-Compute-extra_params.R @@ -0,0 +1,122 @@ +context("Compute, extra function arguments") + +test_that("ex2_6", { + + +#========================= +# Prepare sdates and paths +#========================= + dataset <- "/esarchive/recon/ecmwf/era5/daily_mean/$var$_f1h/$var$_$sdate$.nc" + sdates <- paste0(1981:1982, rep(10:12, 2)) +#=================== +# Get daily winds +#=================== + wind <- Start(dataset = dataset, + var = "sfcWind", + sdate = sdates, + time = 'all', + longitude = indices(1:3), + latitude = indices(1:2), + return_vars = list(time = NULL, latitude = NULL, longitude = NULL), + retrieve = FALSE, + synonims = list(longitude = c('lon', 'longitude'), + latitude = c('lat', 'latitude'))) + + # Synthetic MJO for 'season = OND': + set.seed(1) + MJO <- data.frame(vdate = 1:(30 * 2 + 31 * 4), + phase = c(rep(1:8, 23)), + amplitude = 10 * rnorm(31 * 4 + 30 * 2)) + + stratify_atomic <- function(field, MJO, season = c("JFM", "OND"), + lag = 0, ampl = 2, relative = TRUE, signif = 0.05) { + # Arrange wind in form (days) to match MJO + nmonths <- dim(field)[3] + field <- aperm(field, c(1, 2, 4, 3)) + dim(field) <- c(31 * nmonths) + if(season == "JFM") { + daysok <- rep(c(rep(TRUE, 31), rep(TRUE, 28), + rep(FALSE, 3), rep(TRUE, 31)), nmonths / 3) + } else if (season == "OND") { + daysok <- rep(c(rep(TRUE, 31), rep(TRUE, 30), + rep(FALSE, 1), rep(TRUE, 31)), nmonths / 3) + } + field <- field[daysok] + dim(field) <- c(days = length(field)) + + if(dim(field)[1] != dim(MJO)[1]) { + stop("MJO indices and wind data have different number of days") + } + + idx <- function(MJO, phase, ampl, lag){ + if(lag == 0) { + return(MJO$phase == phase & MJO$amplitude > ampl) + } + if(lag > 0) { + return(dplyr::lag(MJO$phase == phase & MJO$amplitude > ampl, + lag, default = FALSE)) + } + if(lag < 0) { + return(dplyr::lead(MJO$phase == phase & MJO$amplitude > ampl, + - 1 * lag, default = FALSE)) + } + } + strat <- plyr::laply(1:8, function(i) { + idx2 <- idx(MJO, i, ampl, lag) + if (relative) { + return(mean(field[idx2]) / mean(field) - 1) + } else { + return(mean(field[idx2]) - mean(field)) + }}) + strat.t.test <- plyr::laply(1:8, function(i) { + idx2 <- idx(MJO, i, ampl, lag) + return(t.test(field[idx2], field)$p.value)}) + return(list(strat = strat, t_test = strat.t.test)) + } + + step <- Step(stratify_atomic, + target_dims = list(field = c('dataset', 'var', 'sdate', 'time')), + output_dims = list(strat = c('phase'), t_test = c('phase'))) + workflow <- AddStep(wind, step, MJO = MJO, season = "OND", lag = "0", amp = 0) + + res <- Compute(workflow$strat, + chunks = list(latitude = 2)) + +expect_equal( +attr(wind, 'Dimensions'), +c(dataset = 1, var = 1, sdate = 6, time = 31, longitude = 3, latitude = 2) +) +expect_equal( +names(res), +c('strat', 't_test') +) +expect_equal( +dim(res$strat), +c(phase = 8, longitude = 3, latitude = 2) +) +expect_equal( +dim(res$t_test), +c(phase = 8, longitude = 3, latitude = 2) +) +expect_equal( +mean(res$strat), +-0.01373227, +tolerance = 0.0001 +) +expect_equal( +res$strat[1:6,2,1], +c(-0.002499522, 0.125437301, -0.044554040, -0.034862961, 0.019349007, -0.143963809), +tolerance = 0.0001 +) +expect_equal( +res$t_test[1:6,2,1], +c(0.9808923, 0.3378701, 0.6251017, 0.7305827, 0.8573760, 0.2473257), +tolerance = 0.0001 +) +expect_equal( +mean(res$t_test), +0.6419336, +tolerance = 0.0001 +) + +}) diff --git a/tests/testthat/test-Compute-timedim.R b/tests/testthat/test-Compute-timedim.R new file mode 100644 index 0000000000000000000000000000000000000000..4c17806c6adacbffef08020d5de243e311f9021d --- /dev/null +++ b/tests/testthat/test-Compute-timedim.R @@ -0,0 +1,55 @@ +context("Compute on time dimension") + +test_that("ex2_1", { + + repos <- '/esarchive/exp/ecmwf/system5_m1/monthly_mean/$var$_f6h/$var$_$sdate$.nc' + +suppressWarnings( + data <- Start(dat = repos, + var = 'tas', + sdate = c('20170101', '20180101'), + ensemble = indices(1:2), + time = 'all', + latitude = indices(1:10), + longitude = indices(1:15), + return_vars = list(latitude = 'dat', longitude = 'dat', time = 'sdate'), + retrieve = FALSE) +) + + fun_spring <- function(x) { + y <- s2dv::Season(x, time_dim = 'time', monini = 1, moninf = 3, monsup = 5) + return(y) + } + + step1 <- Step(fun = fun_spring, + target_dims = c('var', 'time'), + output_dims = c('var', 'time')) + + wf1 <- AddStep(data, step1) + + res1 <- Compute(wf1, + chunks = list(ensemble = 2)) + +expect_equal( +attr(data, 'Dimensions'), +c(dat = 1, var = 1, sdate = 2, ensemble = 2, time = 7, latitude = 10, longitude = 15) +) +expect_equal( +dim(res1$output), +c(var = 1, time = 1, dat = 1, sdate = 2, ensemble = 2, latitude = 10, longitude = 15) +) +expect_equal( +mean(res1$output), +258.3792, +tolerance = 0.0001 +) +expect_equal( +res1$output[1,1,1,,2,10,2], +c(256.4469, 260.3636), +tolerance = 0.0001 +) + + + + +}) diff --git a/tests/testthat/test-Compute-two_data.R b/tests/testthat/test-Compute-two_data.R new file mode 100644 index 0000000000000000000000000000000000000000..e55a1538141711197c27f3b35ff1538deb1d0d93 --- /dev/null +++ b/tests/testthat/test-Compute-two_data.R @@ -0,0 +1,78 @@ +context("Compute with two datasets") + +test_that("ex2_7", { + +# exp data + repos <- paste0('/esarchive/exp/ecmwf/system4_m1/monthly_mean/', + '$var$_f6h/$var$_$sdate$.nc') + sdates <- sapply(2013:2014, function(x) paste0(x, sprintf('%02d', 1:12), '01')) + + exp <- Start(dat = repos, + var = 'tas', + sdate = sdates, + time = indices(1), + ensemble = indices(1:2), + latitude = values(list(10, 12)), + latitude_reorder = Sort(), + longitude = values(list(0, 2)), + longitude_reorder = CircularSort(0, 360), + synonims = list(longitude = c('lon', 'longitude'), + latitude = c('lat', 'latitude')), + retrieve = F) + +# obs data + repos_obs <- paste0('/esarchive/recon/ecmwf/erainterim/monthly_mean/', + '$var$_f6h/$var$_$sdate$.nc') + sdates_obs <- (sapply(2012:2013, function(x) paste0(x, sprintf('%02d', 1:12)))) + + obs <- Start(dat = repos_obs, + var = 'tas', + sdate = sdates_obs, + time = indices(1), + latitude = values(list(10, 12)), + latitude_reorder = Sort(), + longitude = values(list(0, 2)), + longitude_reorder = CircularSort(0, 360), + synonims = list(longitude = c('lon', 'longitude'), + latitude = c('lat', 'latitude')), + retrieve = F) + + + func <- function(x, y) { + crps <- mean(SpecsVerification::EnsCrps(x, y, R.new = Inf)) + return(crps) + } + step <- Step(func, target_dims = list(c('sdate', 'ensemble'), c('sdate')), + output_dims = NULL) + wf <- AddStep(list(exp, obs), step) + + +# Compute() on fatnodes + res <- Compute(wf, + chunks = list(latitude = 2)) + +expect_equal( +attr(exp, 'Dimensions'), +c(dat = 1, var = 1, sdate = 24, time = 1, ensemble = 2, latitude = 3, longitude = 3) +) +expect_equal( +attr(obs, 'Dimensions'), +c(dat = 1, var = 1, sdate = 24, time = 1, latitude = 3, longitude = 3) +) +expect_equal( +dim(res$output), +c(dat = 1, var = 1, time = 1, latitude = 3, longitude = 3) +) +expect_equal( +mean(res$output), +0.8646249, +tolerance = 0.0001 +) +expect_equal( +res$output[1,1,1,2,1], +0.7980703, +tolerance = 0.0001 +) + + +}) diff --git a/tests/testthat/test-Compute-use_attribute.R b/tests/testthat/test-Compute-use_attribute.R new file mode 100644 index 0000000000000000000000000000000000000000..2ca73a770f040f19c650960e1be8e8c97b6e294e --- /dev/null +++ b/tests/testthat/test-Compute-use_attribute.R @@ -0,0 +1,56 @@ +context("Compute use attributes") + +test_that("ex2_2", { + + repos <- '/esarchive/exp/ecmwf/system5_m1/monthly_mean/$var$_f6h/$var$_$sdate$.nc' + +suppressWarnings( + data <- Start(dat = repos, + var = 'tas', + sdate = c('20170101', '20180101'), + ensemble = indices(1:2), + time = 'all', + latitude = indices(1:10), + longitude = indices(1:15), + return_vars = list(latitude = 'dat', longitude = 'dat', time = 'sdate'), + retrieve = FALSE) +) + funp <- function(x) { + lat <- attributes(x)$Variables$dat1$latitude + weight <- sqrt(cos(lat * pi / 180)) + corrected <- Apply(list(x), target_dims = "latitude", + fun = function(x) {x * weight}) + } + + + step2 <- Step(fun = funp, + target_dims = 'latitude', + output_dims = 'latitude', + use_attributes = list(data = "Variables")) + wf2 <- AddStep(list(data = data), step2) + +suppressWarnings( + res2 <- Compute(workflow = wf2, + chunks = list(sdate = 2)) +) + +expect_equal( +attr(data, 'Dimensions'), +c(dat = 1, var = 1, sdate = 2, ensemble = 2, time = 7, latitude = 10, longitude = 15) +) +expect_equal( +dim(res2$output), +c(latitude = 10, dat = 1, var = 1, sdate = 2, ensemble = 2, time = 7, longitude = 15) +) +expect_equal( +mean(res2$output), +39.84091, +tolerance = 0.0001 +) +expect_equal( +res2$output[2,1,1,,1,7,2], +c(25.40159, 25.40265), +tolerance = 0.0001 +) + +})