From 50db921a6b4c73589c2087492c073ac869f55f30 Mon Sep 17 00:00:00 2001 From: aho Date: Mon, 4 Nov 2019 09:57:19 +0100 Subject: [PATCH] Create usecase folder and usecase.md --- inst/doc/usecase.md | 19 +++ inst/doc/usecase/ex1_1_tranform.R | 102 ++++++++++++++ inst/doc/usecase/ex2_1_timedim.R | 69 ++++++++++ inst/doc/usecase/ex2_2_attr.R | 71 ++++++++++ inst/doc/usecase/ex2_3_cdo.R | 76 +++++++++++ inst/doc/usecase/ex2_4_two_func.R | 77 +++++++++++ inst/doc/usecase/ex2_5_exp_and_obs.R | 50 +++++++ inst/doc/usecase/ex2_6_ext_param_func.R | 124 ++++++++++++++++++ .../ex2_7_seasonal_forecast_verification.R | 45 +++++++ 9 files changed, 633 insertions(+) create mode 100644 inst/doc/usecase.md create mode 100644 inst/doc/usecase/ex1_1_tranform.R create mode 100644 inst/doc/usecase/ex2_1_timedim.R create mode 100644 inst/doc/usecase/ex2_2_attr.R create mode 100644 inst/doc/usecase/ex2_3_cdo.R create mode 100644 inst/doc/usecase/ex2_4_two_func.R create mode 100644 inst/doc/usecase/ex2_5_exp_and_obs.R create mode 100644 inst/doc/usecase/ex2_6_ext_param_func.R create mode 100644 inst/doc/usecase/ex2_7_seasonal_forecast_verification.R diff --git a/inst/doc/usecase.md b/inst/doc/usecase.md new file mode 100644 index 0000000..a584f1e --- /dev/null +++ b/inst/doc/usecase.md @@ -0,0 +1,19 @@ +# Usecase scripts + +In this document, you can link to the example scripts for various demands. For the beginners, it is highly recommended to read the [practical guide](inst/doc/practical_guide.md) carefully first. You can find basic scripts and the configuration for different machines there. + +1. **Retrieve data (use `Start()` only)** + 1. [Interpolation in Start()](inst/doc/usecase/ex1_1_tranform.R) + Do the interpolation within Start(), and compare with Load() result. When the Start() parameter `transform_extra_cells = 2`, the two results will be the same. + 2. + +2. **Execute computation (use `Compute()`)** + 1. [Function working on time dimension](inst/doc/usecase/ex2_1_timedim.R) + 2. [Function using attributes of the data](inst/doc/usecase/ex2_2_attr.R) + Using attributes is only available in startR_v0.1.3 or above. + 3. [Use function CDORemap for interpolation](inst/doc/usecase/ex2_3_cdo.R) + Using parameter `CDO_module` is only available in startR_v0.1.3 or above. Interpolate data by using `s2dverification::CDORemap` in the workflow. + 4. [Use two functions in workflow](inst/doc/usecase/ex2_4_two_func.R) + 5. [Using experimental and (date-corresponding) observational data](inst/doc/usecase/ex2_5_exp_and_obs.R) + 6. [Use external parameters in atomic function](inst/doc/usecase/ex2_6_ext_param_func.R) + 7. [Seasonal forecast verification on cca](inst/doc/usecase/ex2_7_seasonal_forecast_verification.R) diff --git a/inst/doc/usecase/ex1_1_tranform.R b/inst/doc/usecase/ex1_1_tranform.R new file mode 100644 index 0000000..5fbeb6f --- /dev/null +++ b/inst/doc/usecase/ex1_1_tranform.R @@ -0,0 +1,102 @@ +# ------------------------------------------------------------------ +# Do the interpolation within Start(). Compare with Load(). +# ------------------------------------------------------------------ + +library(startR) +library(s2dverification) + +obs_path <- '/esarchive/recon/ecmwf/era5/monthly_mean/$var$_f1h/$var$_$sdate$.nc' +pobs <- paste0('/esarchive/recon/ecmwf/era5/monthly_mean/', + '$VAR_NAME$_f1h/$VAR_NAME$_$YEAR$$MONTH$.nc') +var_name<-'sfcWind' + +lons.min <- 0 +lons.max <- 10 +lats.min <- 0 +lats.max <- 10 + +#------------------------- +# Start() +#------------------------- +obs <- Start(dat = obs_path, + var = var_name, + sdate = '201811', + time = 'all', + latitude = values(list(lats.min, lats.max)), + latitude_reorder = Sort(), + longitude = values(list(lons.min, lons.max)), + synonims = list(longitude = c('lon', 'longitude'), + latitude = c('lat', 'latitude')), + transform = CDORemapper, + transform_extra_cells = 2, + transform_params = list(grid = 'r360x181', + method = 'conservative', + crop = c(lons.min, lons.max, lats.min, lats.max)), + transform_vars = c('latitude', 'longitude'), + return_vars = list(time = NULL, + latitude = 'dat', + longitude = 'dat'), + split_multiselected_dims = TRUE, + merge_across_dims = TRUE, + #num_procs = 4, + retrieve = T) + +# -- Result -- +dim(obs) +# dat var sdate time latitude longitude +# 1 1 1 1 11 11 +range(obs) +#[1] 1.090472 4.987585 +str(attr(obs, 'Variables')$common[[var_name]]$dim[[1]]) +# $ name : chr "longitude" +# $ vals : num [1:11(1d)] 0 1 2 3 4 5 6 7 8 9 ... +str(attr(obs, 'Variables')$common[[var_name]]$dim[[2]]) +# $ name : chr "latitude" +# $ vals : num [1:11(1d)] 0 1 2 3 4 5 6 7 8 9 ... +obs[1,1,1,1,,1] +# [1] 4.987585 4.299388 4.007287 3.758497 3.673881 3.908989 2.224407 1.774090 +# [9] 1.800998 1.718733 1.816074 + +#------------------------- +# s2dverification::Load() +#------------------------- +obs1 <- Load(var = var_name, + obs = list(list(path = pobs)), + sdates = '20181101', + leadtimemin = 1, + leadtimemax = 1, + output = 'lonlat', + grid = 'r360x181', + method = 'conservative', + storefreq = 'monthly', + latmin = lats.min, + latmax = lats.max, + lonmin = lons.min, + lonmax = lons.max) + +# -- Result -- +dim(obs1$obs) +#dataset member sdate ftime lat lon +# 1 1 1 1 11 11 +range(obs1$obs) +#[1] 1.0905 4.9876 +str(obs1) +# $ lon : num [1:11(1d)] 0 1 2 3 4 5 6 7 8 9 ... +# $ lat : num [1:11(1d)] 10 9 8 7 6 5 4 3 2 1 ... +obs1$obs[1,1,1,1,,1] #when lon = 0 +# [1] 1.8161 1.7187 1.8010 1.7741 2.2244 3.9090 3.6739 3.7585 4.0073 4.2994 +#[11] 4.9876 + + +#----------------------- +# Comparison +#----------------------- +diff <- drop(obs)[11:1, ] - drop(obs1$obs) #lat order is different +range(diff) +#[1] -4.918404e-05 4.693756e-05 + + + + + + diff --git a/inst/doc/usecase/ex2_1_timedim.R b/inst/doc/usecase/ex2_1_timedim.R new file mode 100644 index 0000000..2365e00 --- /dev/null +++ b/inst/doc/usecase/ex2_1_timedim.R @@ -0,0 +1,69 @@ +# ----------------------------------------------------------------- +# Function working on time dimension e.g.: Season +# ------------------------------------------------------------------ + + repos <- '/esarchive/exp/ecmwf/system5_m1/monthly_mean/$var$_f6h/$var$_$sdate$.nc' + data <- Start(dat = repos, + var = 'tas', + sdate = c('20170101', '20180101'), + ensemble = indices(1:20), + time = 'all', + latitude = 'all', + longitude = indices(1:40), + return_vars = list(latitude = 'dat', longitude = 'dat', time = 'sdate'), + retrieve = FALSE) + + library(multiApply) + fun_spring <- function(x) { + source("/esarchive/scratch/nperez/Season_v2.R") + y <- Season_v2(x, monini = 1, moninf = 3, monsup = 5) + return(y) + } + + step1 <- Step(fun = fun_spring, + target_dims = c('time'), + output_dims = c('time')) + + wf1 <- AddStep(data, step1) + +## locally + res1 <- Compute(wf1, + chunks = list(ensemble = 2, + sdate = 2)) + + dim(res1$output1) + str(res1$output1) + summary(res1$output1) +# ----------------------------------------------------------------- +#> dim(res1$output1) +# time dat var sdate ensemble latitude longitude +# 1 1 1 2 20 640 40 +#> str(res1$output1) +# num [1, 1, 1, 1:2, 1:20, 1:640, 1:40] 260 257 257 260 257 ... +#> summary(res1$output1) +# Min. 1st Qu. Median Mean 3rd Qu. Max. +# 220.5 273.3 282.7 279.9 295.6 306.8 +# ---------------------------------------------------------------- + +## on Power9 +#-----------modify according to your personal info--------- + queue_host = 'p1' #your own host name for power9 + temp_dir = '/gpfs/scratch/bsc32/bsc32734/startR_hpc/' + ecflow_suite_dir = '/home/Earth/aho/startR_local/' #your own local directory +#------------------------------------------------------------ + res <- Compute(wf1, + chunks = list(ensemble = 20, + sdate = 2), + threads_load = 2, + threads_compute = 4, + cluster = list(queue_host = queue_host, + queue_type = 'slurm', + cores_per_job = 2, + temp_dir = temp_dir, + r_module = 'R/3.5.0-foss-2018b', + polling_period = 10, + job_wallclock = '01:00:00', + max_jobs = 40, + bidirectional = FALSE), + ecflow_suite_dir = ecflow_suite_dir, + wait = TRUE) diff --git a/inst/doc/usecase/ex2_2_attr.R b/inst/doc/usecase/ex2_2_attr.R new file mode 100644 index 0000000..4a73ae2 --- /dev/null +++ b/inst/doc/usecase/ex2_2_attr.R @@ -0,0 +1,71 @@ +#----------------------------------------------------------------- +# Function using attributes of the data e.g.: latitudinal correction +# ------------------------------------------------------------------ + + repos <- '/esarchive/exp/ecmwf/system5_m1/monthly_mean/$var$_f6h/$var$_$sdate$.nc' + data <- Start(dat = repos, + var = 'tas', + sdate = c('20170101', '20180101'), + ensemble = indices(1:20), + time = 'all', + latitude = 'all', + longitude = indices(1:40), + 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(data, step2) + +## locally + res2 <- Compute(workflow = wf2, + chunks = list(ensemble = 2, + sdate = 2)) + + dim(res2$output1) + head(res2$output1) + summary(res2$output1) +# ------------------------------------------------------------------ +# Output: +#> dim(res2$output1) +# latitude dat var sdate ensemble time longitude +# 640 1 1 2 20 7 40 +#> head(res2$output1) +#[1] 15.22683 23.09543 28.94978 33.82982 38.11695 41.99680 +#> summary(res2$output1) +# Min. 1st Qu. Median Mean 3rd Qu. Max. +# 13.19 169.60 237.60 217.90 284.10 305.10 +# ----------------------------------------------------------------- + +## on Power9 +#-----------modify according to your personal info--------- + queue_host = 'p1' #your own host name for power9 + temp_dir = '/gpfs/scratch/bsc32/bsc32734/startR_hpc/' + ecflow_suite_dir = '/home/Earth/aho/startR_local/' #your own local directory +#------------------------------------------------------------ + res2_P <- Compute(wf2, + chunks = list(ensemble = 20, + sdate = 2), + threads_load = 2, + threads_compute = 4, + cluster = list(queue_host = queue_host, #your own host name for power9 + queue_type = 'slurm', + cores_per_job = 2, + temp_dir = temp_dir, + r_module = 'R/3.5.0-foss-2018b', + #extra_queue_params = list('#SBATCH --mem-per-cpu=3000'), + polling_period = 10, + job_wallclock = '01:00:00', + max_jobs = 40, + bidirectional = FALSE), + ecflow_suite_dir = ecflow_suite_dir, #your own local directory + wait = TRUE) diff --git a/inst/doc/usecase/ex2_3_cdo.R b/inst/doc/usecase/ex2_3_cdo.R new file mode 100644 index 0000000..cd7d4a3 --- /dev/null +++ b/inst/doc/usecase/ex2_3_cdo.R @@ -0,0 +1,76 @@ +# -------------------------------------------------------------- +# Function doing regridding CDO (e.g.: s2dverification::CDORemap) +#--------------------------------------------------------------- + + repos <- '/esarchive/exp/ecmwf/system5_m1/monthly_mean/$var$_f6h/$var$_$sdate$.nc' + data <- Start(dat = repos, + var = 'tas', + sdate = c('20170101', '20180101'), + ensemble = indices(1:2), + time = 'all', + latitude = 'all', + longitude = 'all', + return_vars = list(latitude = 'dat', longitude = 'dat', time = 'sdate'), + retrieve = FALSE) + + + fun_deb2 <- function(x) { + lons_data = as.vector(attr(x, 'Variables')$dat1$longitude) + lats_data = as.vector(attr(x, 'Variables')$dat1$latitude) + resgrid = "r360x180" # prlr + r <- s2dverification::CDORemap(x, lons_data, lats_data, resgrid, + 'bil', crop = FALSE, force_remap = TRUE)[[1]] + return(r) + } + + step3 <- Step(fun = fun_deb2, + target_dims = c('latitude','longitude'), + output_dims = c('latitude', 'longitude'), + use_attributes = list(data = "Variables")) + wf3 <- AddStep(data, step3) + +## locally + res3 <- Compute(workflow = wf3, + chunks = list(ensemble = 2, + sdate = 2)) + dim(res3$output1) + head(res3$output1) + summary(res3$output1) + +# -------------------------------------------------------------------- +# Output: +#> dim(res3$output1) +# latitude longitude dat var sdate ensemble time +# 180 360 1 1 2 2 7 +#> head(res3$output1) +#[1] 245.2156 244.6700 245.1513 244.7413 244.4198 245.0041 +#> summary(res3$output1) +# Min. 1st Qu. Median Mean 3rd Qu. Max. +# 206.3 269.2 283.3 278.8 296.4 314.2 +# -------------------------------------------------------------------- + +## on Power9 +#-----------modify according to your personal info--------- + queue_host = 'p1' #your own host name for power9 + temp_dir = '/gpfs/scratch/bsc32/bsc32734/startR_hpc/' + ecflow_suite_dir = '/home/Earth/aho/startR_local/' #your own local directory +#------------------------------------------------------------ + res3_P <- Compute(wf3, + chunks = list(ensemble = 2, + sdate = 2), + threads_load = 2, + threads_compute = 4, + cluster = list(queue_host = queue_host, #your own host name for power9 + queue_type = 'slurm', + cores_per_job = 1, + temp_dir = temp_dir, + r_module = 'R/3.5.0-foss-2018b', + CDO_module = 'CDO/1.9.5-foss-2018b', + extra_queue_params = list('#SBATCH --mem-per-cpu=3000'), + polling_period = 100, + job_wallclock = '01:00:00', + max_jobs = 4, + bidirectional = FALSE), + ecflow_suite_dir = ecflow_suite_dir, #your own local directory + wait = TRUE) + diff --git a/inst/doc/usecase/ex2_4_two_func.R b/inst/doc/usecase/ex2_4_two_func.R new file mode 100644 index 0000000..799bd3f --- /dev/null +++ b/inst/doc/usecase/ex2_4_two_func.R @@ -0,0 +1,77 @@ +# -------------------------------------------------------------- +# Two functions (e.g.: s2dverification::CDORemap and Season) +#--------------------------------------------------------------- + + repos <- '/esarchive/exp/ecmwf/system5_m1/monthly_mean/$var$_f6h/$var$_$sdate$.nc' + data <- Start(dat = repos, + var = 'tas', + sdate = c('20170101', '20180101'), + ensemble = indices(1:2), + time = 'all', + latitude = 'all', + longitude = 'all', + return_vars = list(latitude = 'dat', longitude = 'dat', time = 'sdate'), + retrieve = FALSE) + + fun_deb3 <- function(x) { + source("/esarchive/scratch/nperez/Season_v2.R") + lons_data = as.vector(attr(x, 'Variables')$dat1$longitude) + lats_data = as.vector(attr(x, 'Variables')$dat1$latitude) + resgrid = "r360x180" # prlr + y = Season_v2(x, posdim = 'time', monini = 1, moninf = 1, monsup = 3) + r <- s2dverification::CDORemap(y, lons_data, lats_data, resgrid, + 'bil', crop = FALSE, force_remap = TRUE)[[1]] + return(r) + } + + step4 <- Step(fun = fun_deb3, + target_dims = c('latitude','longitude', 'time'), + output_dims = c('latitude', 'longitude', 'time'), + use_attributes = list(data = "Variables")) + wf4 <- AddStep(data, step4) + +## locally + res4 <- Compute(workflow = wf4, + chunks = list(ensemble = 2, + sdate = 2)) + dim(res4$output1) + head(res4$output1) + summary(res4$output1) + +# ------------------------------------------------------------------ +# Output: +#> dim(res4$output1) +# latitude longitude time dat var sdate ensemble +# 180 360 1 1 1 2 2 +#> head(res4$output1) +#[1] 237.1389 237.2601 238.0882 238.0312 237.7883 238.4835 +#> summary(res4$output1) +# Min. 1st Qu. Median Mean 3rd Qu. Max. +# 227.3 259.6 280.8 277.1 296.2 306.7 +# ------------------------------------------------------------------ + +## on Power9 +#-----------modify according to your personal info----------- + queue_host = 'p1' #your own host name for power9 + temp_dir = '/gpfs/scratch/bsc32/bsc32734/startR_hpc/' + ecflow_suite_dir = '/home/Earth/aho/startR_local/' #your own local directory +#------------------------------------------------------------ + res4_P <- Compute(wf4, + chunks = list(ensemble = 2, + sdate = 2), + threads_load = 2, + threads_compute = 4, + cluster = list(queue_host = queue_host, #your own host name for power9 + queue_type = 'slurm', + cores_per_job = 1, + temp_dir = temp_dir, + r_module = 'R/3.5.0-foss-2018b', + CDO_module = 'CDO/1.9.5-foss-2018b', + extra_queue_params = list('#SBATCH --mem-per-cpu=3000'), + polling_period = 10, + job_wallclock = '01:00:00', + max_jobs = 6, + bidirectional = FALSE), + ecflow_suite_dir = ecflow_suite_dir, #your own local directory + wait = TRUE) + diff --git a/inst/doc/usecase/ex2_5_exp_and_obs.R b/inst/doc/usecase/ex2_5_exp_and_obs.R new file mode 100644 index 0000000..71874e6 --- /dev/null +++ b/inst/doc/usecase/ex2_5_exp_and_obs.R @@ -0,0 +1,50 @@ + repos <- paste0('/esnas/exp/ecmwf/system4_m1/6hourly/', + '$var$/$var$_$sdate$.nc') + + system4 <- Start(dat = repos, + var = 'sfcWind', + sdate = paste0(1981:1984, '1101'), + time = indices((30*4+1):(30*4+4)), + ensemble = 'all', + latitude = indices(1:10), + longitude = indices(1:10), + return_vars = list(latitude = NULL, + longitude = NULL, + time = c('sdate'))) + + repos <- paste0('/esnas/recon/ecmwf/erainterim/6hourly/', + '$var$/$var$_$file_date$.nc') + + dates <- attr(system4, 'Variables')$common$time + dates_file <- sort(unique(gsub('-', '', sapply(as.character(dates), + substr, 1, 7)))) + + erai <- Start(dat = repos, + var = 'sfcWind', + file_date = dates_file, + time = values(dates), + latitude = indices(1:10), + longitude = indices(1:10), + time_var = 'time', + time_tolerance = as.difftime(1, units = 'hours'), + time_across = 'file_date', + return_vars = list(latitude = NULL, + longitude = NULL, + time = 'file_date'), + merge_across_dims = TRUE, + split_multiselected_dims = TRUE) + + step <- Step(eqmcv_atomic, + list(a = c('ensemble', 'sdate'), + b = c('sdate')), + list(c = c('ensemble', 'sdate'))) + + res <- Compute(step, list(system4, erai), + chunks = list(latitude = 5, + longitude = 5, + time = 2), + cluster = list(queue_host = 'bsceslogin01.bsc.es', + max_jobs = 4, + cores_per_job = 2), + shared_dir = '/esnas/scratch/nmanuben/test_bychunk', + wait = FALSE) diff --git a/inst/doc/usecase/ex2_6_ext_param_func.R b/inst/doc/usecase/ex2_6_ext_param_func.R new file mode 100644 index 0000000..1bd74a8 --- /dev/null +++ b/inst/doc/usecase/ex2_6_ext_param_func.R @@ -0,0 +1,124 @@ + library(lubridate) + library(startR) + +#========================= +# Prepare sdates and paths +#========================= + dataset <- "/esarchive/recon/ecmwf/era5/daily_mean/$var$_f1h/$var$_$sdate$.nc" + firsty <- "1981" + lasty <- "1987" + months <- 10:12 + firstd <- paste0(firsty,"1001") + sdates <- ymd(firstd) + months(0:2) + rep(years(0:(1987-1981)), each = 3) + +#=================== +# Get daily winds +#=================== + wind <- Start(dataset = dataset, + var = "sfcWind", + sdate = format(sdates, "%Y%m"), + time = 'all', + longitude = indices(1:30), + latitude = indices(1:20), + return_vars = list(time = NULL, latitude = NULL, longitude = NULL), + num_procs = 4, retrieve = FALSE, + synonims = list(longitude = c('lon', 'longitude'), + latitude = c('lat', 'latitude'))) + + # Synthetic MJO for 'season = OND': + MJO <- data.frame(vdate = 1:(30 * 7 + 31 * 14), + phase = c(rep(1:8, 80), 1:4), + amplitude = 10 * rnorm(31 * 14 + 30 * 7)) + + 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, longitude = 2), + threads_load = 2, + threads_compute = 4) + +# ----------------------------------------------------------- +#names(res) +#[1] "strat" "t_test" +#> dim(res$strat) +# phase longitude latitude +# 8 30 20 +#> summary(res$strat) +# Min. 1st Qu. Median Mean 3rd Qu. Max. +#-0.133300 -0.032530 -0.001822 -0.005715 0.031700 0.094220 +#> res$strat[1:5, 1:2, 1] +# [,1] [,2] +#[1,] -0.04661354 -0.04661539 +#[2,] -0.01058483 -0.01053589 +#[3,] 0.06123498 0.06125660 +#[4,] 0.06371607 0.06374801 +#[5,] 0.07085255 0.07078778 + +# ------------------------------------------------------------- + +# ---- To be tried for comparison of results: + res <- Compute(workflow$strat, + chunks = list(latitude = 2, longitude = 2), + threads_load = 2, + threads_compute = 2, + cluster = list(queue_host = 'cte-power', #user-specific + queue_type = 'slurm', + temp_dir = '/gpfs/scratch/bsc32/bsc32339/', #user-specific + r_module = 'R/3.5.0-foss-2018b', + cores_per_job = 2, + job_wallclock = '10:00:00', + max_jobs = 4, + bidirectional = FALSE, + polling_period = 10) + ecflow_suite_dir = '/esarchive/scratch/nperez/ecFlow', #user-specific + wait = FALSE) + diff --git a/inst/doc/usecase/ex2_7_seasonal_forecast_verification.R b/inst/doc/usecase/ex2_7_seasonal_forecast_verification.R new file mode 100644 index 0000000..bcd11a5 --- /dev/null +++ b/inst/doc/usecase/ex2_7_seasonal_forecast_verification.R @@ -0,0 +1,45 @@ + library(startR) + + repos <- '/perm/ms/spesiccf/c3ah/qa4seas/data/seasonal/g1x1/ecmf-system4/msmm/atmos/seas/tprate/12/ecmf-system4_msmm_atmos_seas_sfc_$date$_tprate_g1x1_init12.nc' + + data <- Start(dat = repos, + var = 'tprate', + date = 'all', + time = 'all', + number = 'all', + latitude = 'all', + longitude = 'all', + return_vars = list(time = 'date')) + + dates <- attr(data, 'Variables')$common$time + + repos <- '/perm/ms/spesiccf/c3ah/qa4seas/data/ecmf-ei_msmm_atmos_seas_sfc_19910101-20161201_t2m_g1x1_init02.nc' + + obs <- Start(dat = repos, + var = 't2m', + time = values(dates), + latitude = 'all', + longitude = 'all', + split_multiselected_dims = TRUE) + + crps <- function(x, y) { + mean(SpecsVerification::EnsCrps(x, y, R.new = Inf)) + } + + s <- Step(crps, target_dims = list(c('date', 'number'), c('date')), + output_dims = NULL) + wf <- AddStep(list(data, obs), s) + + r <- Compute(wf, + chunks = list(latitude = 10, + longitude = 3), + cluster = list(queue_host = 'cca', + queue_type = 'pbs', + max_jobs = 10, + init_commands = list('module load ecflow'), + r_module = 'R/3.3.1', + extra_queue_params = list('#PBS -l EC_billing_account=spesiccf')), + ecflow_output_dir = '/perm/ms/spesiccf/c3ah/startR_test/', + is_ecflow_output_dir_shared = FALSE + ) + -- GitLab