#****************************************************************************
# Script to test examples of CST_SaveExp
# Eva Rifà Rovira 
# 29/11/2024
#****************************************************************************

# Tests 1: Multidimensional array and Dates, without metadata and coordinates
# (1) Minimal use case, without Dates
data <- array(1:5, dim = c(sdate = 5, lon = 4, lat = 4))
SaveExp(data, ftime_dim = NULL, memb_dim = NULL, dat_dim = NULL,
        var_dim = NULL, single_file = TRUE)
SaveExp(data, ftime_dim = NULL, memb_dim = NULL, dat_dim = NULL,
        var_dim = NULL, sdate_dim = NULL, single_file = TRUE) # same result

# (2) Forecast time dimension, without Dates
data <- array(1:5, dim = c(ftime = 5, lon = 4, lat = 4))
SaveExp(data, ftime_dim = 'ftime', memb_dim = NULL, dat_dim = NULL,
        var_dim = NULL, sdate_dim = NULL, single_file = TRUE)

# (2) Start date dimension, without Dates
data <- array(1:5, dim = c(sdate = 5, lon = 4, lat = 4))
SaveExp(data, ftime_dim = NULL, memb_dim = NULL, dat_dim = NULL,
        var_dim = NULL, sdate_dim = 'sdate', single_file = TRUE)

# (3) Only forecast time dimension (no sdate), with Dates
data <- array(1:5, dim = c(ftime = 5, lon = 4, lat = 4))
dates <- c('20000101', '20010102', '20020103', '20030104', '20040105')
dates <- as.Date(dates, format = "%Y%m%d", tz = "UTC")
dim(dates) <- c(ftime = 5)
SaveExp(data, ftime_dim = 'ftime', memb_dim = NULL, dat_dim = NULL,
        var_dim = NULL, sdate_dim = NULL, Dates = dates, single_file = FALSE)
# For this case we have the same result using: single_file = FALSE /TRUE.

# (4) Forecast time and 1 sdate, with Dates
data <- array(1:5, dim = c(sdate = 1, ftime = 5, lon = 4, lat = 4))
dates <- c('20000101', '20010102', '20020103', '20030104', '20040105')
dates <- as.Date(dates, format = "%Y%m%d", tz = "UTC")
dim(dates) <- c(sdate = 1, ftime = 5)
SaveExp(data, ftime_dim = 'ftime', memb_dim = NULL, dat_dim = NULL,
        var_dim = NULL, sdate_dim = 'sdate', Dates = dates, single_file = TRUE)
SaveExp(data, ftime_dim = 'ftime', memb_dim = NULL, dat_dim = NULL,
        var_dim = NULL, sdate_dim = 'sdate', Dates = dates, single_file = FALSE)

################################################################################

# Tests 2: Test sample data from CSTools

# (a) Data loaded with Start
# (1) lonlat_temp_st$exp in a single file with units 'hours since'

data <- lonlat_temp_st$exp
CST_SaveExp(data = data, ftime_dim = 'ftime', 
            var_dim = 'var', dat_dim = 'dataset', sdate_dim = 'sdate', 
            single_file = TRUE)

# Now we read the output with Start:
sdate <- as.vector(lonlat_temp_st$exp$coords$sdate)
path <- paste0(getwd(),"/$var$_", sdate[1], "_", sdate[length(sdate)], ".nc")
out <- Start(dat = path, 
             var = 'tas', 
             lon = 'all', 
             lat = 'all',
             ftime = 'all', 
             sdate = 'all', 
             member = 'all',
             return_vars = list(lon = 'dat',
                                lat = 'dat',
                                ftime = NULL,
                                sdate = NULL),
             retrieve = TRUE)

attributes(out)$Variables$common$ftime

# (2) lonlat_temp_st$exp in a single file with units of time frequency
data <- lonlat_temp_st$exp
CST_SaveExp(data = data, ftime_dim = 'ftime', 
            var_dim = 'var', dat_dim = 'dataset', sdate_dim = 'sdate', 
            single_file = TRUE, units_hours_since = FALSE)

# Now we read the output with Start:
sdate <- as.vector(lonlat_temp_st$exp$coords$sdate)
path <- paste0(getwd(),"/$var$_", sdate[1], "_", sdate[length(sdate)], ".nc")
out <- Start(dat = path,
             var = 'tas',
             lon = 'all',
             lat = 'all',
             ftime = 'all',
             sdate = 'all',
             member = 'all',
             return_vars = list(
                lon = 'dat',
                lat = 'dat',
                ftime = NULL,
                sdate = NULL),
             retrieve = TRUE)

attributes(out)$Variables$common$ftime
# [1] "1 months" "2 months" "3 months"

# (3) lonlat_temp_st$exp in separated files with units of hours since
data <- lonlat_temp_st$exp
CST_SaveExp(data = data, ftime_dim = 'ftime', 
            var_dim = 'var', dat_dim = 'dataset', sdate_dim = 'sdate', 
            single_file = FALSE)
sdate <- as.vector(lonlat_temp_st$exp$coords$sdate)
path <- paste0(getwd(),"/dat1/$var$/$var$_$sdate$.nc")

out <- Start(dat = path, var = 'tas',
             sdate = sdate,
             lon = 'all',
             lat = 'all',
             ftime = 'all',
             member = 'all',
             return_vars = list(lon = 'dat',
                                lat = 'dat',
                                ftime = NULL),
             retrieve = TRUE)

# (4) lonlat_prec_st$exp in a single file with units of time frequency
data <- lonlat_prec_st
CST_SaveExp(data = data, ftime_dim = 'ftime', 
            var_dim = 'var', dat_dim = 'dataset', sdate_dim = 'sdate', 
            single_file = TRUE, units_hours_since = FALSE)

# Now we read the output with Start:
sdate <- as.vector(data$coords$sdate)
path <- paste0(getwd(),"/$var$_", sdate[1], "_", sdate[length(sdate)], ".nc")
out <- Start(dat = path,
             var = 'prlr',
             lon = 'all',
             lat = 'all',
             ftime = 'all',
             sdate = 'all',
             member = 'all',
             return_vars = list(
                lon = 'dat',
                lat = 'dat',
                ftime = NULL,
                sdate = NULL),
             retrieve = TRUE)

attributes(out)$Variables$common$ftime
#  [1] "1 days"  "2 days"  "3 days"  "4 days"  "5 days"  "6 days"  "7 days" 
#  [8] "8 days"  "9 days"  "10 days" "11 days" "12 days" "13 days" "14 days"
# [15] "15 days" "16 days" "17 days" "18 days" "19 days" "20 days" "21 days"
# [22] "22 days" "23 days" "24 days" "25 days" "26 days" "27 days" "28 days"
# [29] "29 days" "30 days" "31 days"

# (5) Test observations
data <- lonlat_temp$obs
CST_SaveExp(data = data, ftime_dim = 'ftime', memb_dim = NULL, 
            var_dim = NULL, dat_dim = 'dataset', sdate_dim = 'sdate', 
            single_file = TRUE, units_hours_since = FALSE)
# Now we read the output with Start:
sdate <- c('20001101', '20051101')
path <- paste0(getwd(),"/$var$_", sdate[1], "_", sdate[length(sdate)], ".nc")
out <- Start(dat = path,
             var = 'tas', # tas
             lon = 'all',
             lat = 'all',
             ftime = 'all',
             member = 1,
             sdate = 'all',
             return_vars = list(
                lon = 'dat',
                lat = 'dat',
                ftime = NULL,
                sdate = NULL),
             retrieve = TRUE)
dim(out)
attributes(out)$Variables$common$ftime

# (b) Data loaded with Load
data <- lonlat_temp$exp
data <- lonlat_temp$obs
data <- lonlat_prec
CST_SaveExp(data = data, ftime_dim = 'ftime', 
            var_dim = NULL, dat_dim = 'dataset', sdate_dim = 'sdate', 
            single_file = TRUE, units_hours_since = FALSE)
# Now we read the output with Start:
sdate <- as.vector(data$coords$sdate)
path <- paste0(getwd(),"/$var$_", sdate[1], "_", sdate[length(sdate)], ".nc")
out <- Start(dat = path,
             var = 'prlr', # tas
             lon = 'all',
             lat = 'all',
             ftime = 'all',
             sdate = 'all',
             member = 'all',
             return_vars = list(
                lon = 'dat',
                lat = 'dat',
                ftime = NULL,
                sdate = NULL),
             retrieve = TRUE)
dim(out)
lonlat_prec$dims

# Test with ftime_dim NULL
data <- lonlat_temp$exp
data <- CST_Subset(data, along = 'ftime', indices = 1, drop = 'selected')

CST_SaveExp(data = data, ftime_dim = NULL, 
            var_dim = NULL, dat_dim = 'dataset', sdate_dim = 'sdate', 
            single_file = FALSE, units_hours_since = FALSE)

################################################################################
# Test 3: Special cases

# (1) two variables and two datasets in separated files

repos <- "/esarchive/exp/ecmwf/system5_m1/monthly_mean/$var$_f6h/$var$_$sdate$.nc"
repos2 <- "/esarchive/exp/ecmwf/system4_m1/monthly_mean/$var$_f6h/$var$_$sdate$.nc"

data3 <- Start(dat = list(list(name = 'system4_m1', path = repos2),
                          list(name = 'system5_m1', path = repos)),
              var = c('tas', 'sfcWind'),
              sdate = c('20160101', '20170101'),
              ensemble = indices(1),
              time = indices(1:2),
              lat = indices(1:10),
              lon = indices(1:10),
              synonims = list(lat = c('lat', 'latitude'),
                              lon = c('lon', 'longitude')),
              return_vars =  list(time = 'sdate',
                                  longitude = 'dat',
                                  latitude = 'dat'),
              metadata_dims = c('dat', 'var'),
              retrieve = T
              )
cube3 <- as.s2dv_cube(data3)

CST_SaveExp(data = cube3, ftime_dim = 'time', var_dim = 'var', 
             memb_dim = 'ensemble', dat_dim = 'dat')

# We read again the data with start
repos <- paste0(getwd(), "/system4_m1/$var$/$var$_$sdate$.nc")
repos2 <- paste0(getwd(), "/system5_m1/$var$/$var$_$sdate$.nc")

data3out <- Start(dat = list(list(name = 'system4_m1', path = repos2),
                          list(name = 'system5_m1', path = repos)),
              var = c('tas', 'sfcWind'),
              sdate = c('20160101', '20170101'),
              ensemble = indices(1),
              time = indices(1:2),
              lat = indices(1:10),
              lon = indices(1:10),
              synonims = list(lat = c('lat', 'latitude'),
                              lon = c('lon', 'longitude')),
              return_vars =  list(time = 'sdate',
                                  longitude = 'dat',
                                  latitude = 'dat'),
              metadata_dims = c('dat', 'var'),
              retrieve = T)

summary(data3out)
summary(data3)

dim(data3)
dim(data3out)

# (2) two variables and two datasets in the same file

CST_SaveExp(data = cube3, ftime_dim = 'time', var_dim = 'var', 
             memb_dim = 'ensemble', dat_dim = 'dat', 
             single_file = TRUE)

# (3) Observations
repos_exp <- paste0('/esarchive/exp/ecearth/a1tr/cmorfiles/CMIP/EC-Earth-Consortium/',
                    'EC-Earth3/historical/r24i1p1f1/Amon/$var$/gr/v20190312/',
                    '$var$_Amon_EC-Earth3_historical_r24i1p1f1_gr_$sdate$01-$sdate$12.nc')

exp <- Start(dat = repos_exp,
             var = 'tas',
             sdate = as.character(c(2005:2008)),
             time = indices(1:3),
             lat = 1:10,
             lat_reorder = Sort(),
             lon = 1:10,
             lon_reorder = CircularSort(0, 360),
             synonims = list(lat = c('lat', 'latitude'),
                             lon = c('lon', 'longitude')),
             return_vars = list(lon = NULL,
                                lat = NULL,
                                time = 'sdate'),
             retrieve = FALSE)
lats <- attr(exp, 'Variables')$common$lat
lons <- attr(exp, 'Variables')$common$lon
## The 'time' attribute is a two-dim array
dates <- attr(exp, 'Variables')$common$time
dim(dates)
repos_obs <- '/esarchive/recon/ecmwf/erainterim/monthly_mean/$var$_f6h/$var$_$date$.nc'

obs <- Start(dat = repos_obs,
             var = 'tas',
             date = unique(format(dates, '%Y%m')),
             time = values(dates),  #dim: [sdate = 4, time = 3]
             lat = 1:10,
             lat_reorder = Sort(),
             lon = 1:10,
             lon_reorder = CircularSort(0, 360),
             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 = TRUE)
dim(obs)
attributes(obs)$Variables$common$time
obscube <- as.s2dv_cube(obs)
CST_SaveExp(data = obscube, ftime_dim = 'time', var_dim = 'var', 
             memb_dim = 'ensemble', dat_dim = 'dat', 
             single_file = TRUE)



################################################################################