diff --git a/R/CST_MergeDims.R b/R/CST_MergeDims.R index 4b66629ef724460287e01e26a22b098ae705ab05..dabdc57f3ae0aafd543f45d9f2594acdd5738c88 100644 --- a/R/CST_MergeDims.R +++ b/R/CST_MergeDims.R @@ -15,8 +15,6 @@ #' \code{merge_dims} will be used. #'@param na.rm A logical indicating if the NA values should be removed or not. #' -#'@import abind -#'@importFrom ClimProjDiags Subset #'@examples #'data <- 1 : c(2 * 3 * 4 * 5 * 6 * 7) #'dim(data) <- c(time = 7, lat = 2, lon = 3, monthly = 4, member = 6, @@ -35,8 +33,32 @@ CST_MergeDims <- function(data, merge_dims = c('ftime', 'monthly'), if (!inherits(data, 's2dv_cube')) { stop("Parameter 'data' must be of the class 's2dv_cube'.") } + if (is.null(rename_dim)) { + rename_dim <- merge_dims[1] + } + # data data$data <- MergeDims(data$data, merge_dims = merge_dims, rename_dim = rename_dim, na.rm = na.rm) + # dims + data$dims <- dim(data$data) + + # rename_dim + if (length(rename_dim) > 1) { + rename_dim <- as.character(rename_dim[1]) + } + # coords + data$coords[merge_dims] <- NULL + data$coords[[rename_dim]] <- 1:dim(data$data)[rename_dim] + attr(data$coords[[rename_dim]], 'indices') <- TRUE + + # attrs + if (all(merge_dims %in% names(dim(data$attrs$Dates)))) { + dim(data$attrs$Dates) <- dim(data$data)[rename_dim] + } else if (any(merge_dims %in% names(dim(data$attrs$Dates)))) { + warning("The dimensions of 'Dates' array will be different from ", + "the temporal dimensions in 'data'. Parameter 'merge_dims' ", + "only includes one temporal dimension of 'Dates'.") + } return(data) } #'Function to Split Dimension @@ -55,12 +77,12 @@ CST_MergeDims <- function(data, merge_dims = c('ftime', 'monthly'), #' \code{merge_dims} will be used. #'@param na.rm A logical indicating if the NA values should be removed or not. #' -#'@import abind -#'@importFrom ClimProjDiags Subset #'@examples #'data <- 1 : 20 #'dim(data) <- c(time = 10, lat = 2) #'new_data <- MergeDims(data, merge_dims = c('time', 'lat')) +#'@import abind +#'@importFrom ClimProjDiags Subset #'@export MergeDims <- function(data, merge_dims = c('time', 'monthly'), rename_dim = NULL, na.rm = FALSE) { diff --git a/R/CST_Subset.R b/R/CST_Subset.R index 28c8498cbcb52cc984a37a21f189700fa7484fc8..098df3dff90af695fe097e6aca83cdb06b26d418 100644 --- a/R/CST_Subset.R +++ b/R/CST_Subset.R @@ -79,9 +79,9 @@ CST_Subset <- function(x, along, indices, drop = FALSE, var_dim = NULL, } # Subset data - x$data <- ClimProjDiags::Subset(x$data, along = along, - indices = indices, - drop = drop) + x$data <- Subset(x$data, along = along, + indices = indices, + drop = drop) # Adjust dimensions x$dims <- dim(x$data) # Adjust coordinates @@ -113,10 +113,10 @@ CST_Subset <- function(x, along, indices, drop = FALSE, var_dim = NULL, } if ((!is.null(x$attrs$source_files)) && (dim_name %in% names(dim(x$attrs$source_files)))) { - x$attrs$source_files <- ClimProjDiags::Subset(x$attrs$source_files, - along = dim_name, - indices = index, - drop = drop) + x$attrs$source_files <- Subset(x$attrs$source_files, + along = dim_name, + indices = index, + drop = drop) } } # Remove metadata from variables that were dropped @@ -128,10 +128,10 @@ CST_Subset <- function(x, along, indices, drop = FALSE, var_dim = NULL, if (!(length(time_along) == 0)) { time_indices <- indices[match(time_along, along)] original_dates <- x$attrs$Dates - x$attrs$Dates <- ClimProjDiags::Subset(x$attrs$Dates, - along = time_along, - indices = time_indices, - drop = drop) + x$attrs$Dates <- Subset(x$attrs$Dates, + along = time_along, + indices = time_indices, + drop = drop) } # Subset metadata for (variable in 1:length(names(x$attrs$Variable$metadata))) { @@ -153,9 +153,9 @@ CST_Subset <- function(x, along, indices, drop = FALSE, var_dim = NULL, if (any(is.null(dim(x)), length(dim(x)) == 1)) { l <- x[args_subset[['indices']][[1]]] } else { - l <- ClimProjDiags::Subset(x, along = args_subset[['along']], - indices = args_subset[['indices']], - drop = args_subset[['drop']]) + l <- Subset(x, along = args_subset[['along']], + indices = args_subset[['indices']], + drop = args_subset[['drop']]) } attr.names <- names(attributes(x)) attr.names <- attr.names[attr.names != 'names'] diff --git a/README.md b/README.md index 2dfbb5b7e20d691c468c889909a7c50180ba9297..76c82d59739eedf9e15fb3cf81c996efa759c585 100644 --- a/README.md +++ b/README.md @@ -19,6 +19,7 @@ A part from this GitLab project, that allows you to monitor CSTools progress, to - The CRAN repository [https://CRAN.R-project.org/package=CSTools](https://CRAN.R-project.org/package=CSTools) which includes the user manual and vignettes. - Video tutorials [https://www.medscope-project.eu/products/tool-box/cstools-video-tutorials/](https://www.medscope-project.eu/products/tool-box/cstools-video-tutorials/). - Other resources are under-development such [training material](https://earth.bsc.es/gitlab/external/cstools/-/tree/MEDCOF2022/inst/doc/MEDCOF2022) and a [full reproducible use case for forecast calibration](https://earth.bsc.es/gitlab/external/cstools/-/tree/develop-CalibrationVignette/FOCUS_7_2). +- See and run package [**use cases**](inst/doc/usecase.md) Installation ------------ @@ -46,43 +47,55 @@ Overview The CSTools package functions can be distributed in the following methods: -- **Data retrieval and formatting:** CST_Load, CST_Anomaly, CST_MergeDims, CST_SplitDims, CST_Subset, as.s2dv_cube, s2dv_cube, CST_SaveExp. -- **Classification:** CST_MultiEOF, CST_WeatherRegimes, CST_RegimsAssign, CST_CategoricalEnsCombination, CST_EnsClustering. -- **Downscaling:** CST_Analogs, CST_RainFARM, CST_RFTemp, CST_AdamontAnalog, CST_AnalogsPredictors. -- **Correction:** CST_BEI_Weighting, CST_BiasCorrection, CST_Calibration, CST_QuantileMapping, CST_DynBiasCorrection. -- **Assessment:** CST_MultiMetric, CST_MultivarRMSE -- **Visualization:** PlotCombinedMap, PlotForecastPDF, PlotMostLikelyQuantileMap, PlotPDFsOLE, PlotTriangles4Categories, PlotWeeklyClim. +- **Data retrieval and formatting:** [CST_Start](R/CST_Start.R), [CST_SaveExp](R/CST_SaveExp.R), [CST_MergeDims](R/CST_MergeDims.R), [CST_SplitDim](R/CST_SplitDim.R), [CST_Subset](R/CST_Subset), [CST_InsertDim](R/CST_InsertDim.R), [CST_ChangeDimNames](R/CST_ChangeDimNames.R), [as.s2dv_cube](R/as.s2dv_cube.R) and [s2dv_cube](R/s2dv_cube.R). +- **Classification:** [CST_MultiEOF](R/CST_MultiEOF.R), [CST_WeatherRegimes](R/CST_WeatherRegimes.R), [CST_RegimsAssign](R/CST_RegimesAssign.R), [CST_CategoricalEnsCombination](R/CST_CategoricalEnsCombination.R), [CST_EnsClustering](R/CST_EnsClustering.R). +- **Downscaling:** [CST_Analogs](R/CST_Analogs.R), [CST_RainFARM](R/CST_RainFARM.R), [CST_RFTemp](R/CST_RFTemp.R), [CST_AdamontAnalog](R/CST_AdamontAnalog.R), [CST_AnalogsPredictors](R/CST_AnalogsPredictors.R). +- **Correction and transformation:** [CST_BiasCorrection](R/CST_BiasCorrection.R), [CST_Calibration](R/CST_Calibration.R), [CST_QuantileMapping](R/CST_QuantileMapping.R), [CST_Anomaly](R/CST_Anomaly.R), [CST_BEI_Weighting](R/CST_BEI_Weighting.R), [CST_DynBiasCorrection](R/CST_DynBiasCorrection.R). +- **Assessment:** [CST_MultiMetric](R/CST_MultiMetric.R), [CST_MultivarRMSE](R/CST_MultivarRMSE.R) +- **Visualization:** [PlotCombinedMap](R/PlotCombinedMap.R), [PlotForecastPDF](R/PlotForecastPDF.R), [PlotMostLikelyQuantileMap](R/PlotMostLikelyQuantileMap.R), [PlotPDFsOLE](R/PlotPDFsOLE.R), [PlotTriangles4Categories](R/PlotTriangles4Categories.R), [PlotWeeklyClim](R/PlotWeeklyClim.R). -This package is designed to be compatible with other R packages such as [s2dv](https://CRAN.R-project.org/package=s2dv), [startR](https://CRAN.R-project.org/package=startR), [CSIndicators](https://CRAN.R-project.org/package=CSIndicators), [CSDownscale](https://earth.bsc.es/gitlab/es/csdownscale). Functions with the prefix **CST_** deal with a common object called `s2dv_cube` as inputs. Also, this object can be created from Load (s2dv) and from Start (startR) directly. Multiple functions from different packages can operate on this common data structure to easily define a complete post-processing workflow. - -The class `s2dv_cube` is mainly a list of named elements to keep data and metadata in a single object. Basic structure of the object: +An `s2dv_cube` is an object to store ordered multidimensional array with named dimensions, specific coordinates and stored metadata (in-memory representation of a NetCDF file). Its “methods” are the **CST** prefix functions. The basic structure of the class `s2dv_cube` is a list of lists. The first level elements are: `data`, `dims`, `coords` and `attrs`. To access any specific element it will be done using the `$` operator. +As an example, this is how an `s2dv_cube` looks like (see `lonlat_temp_st$exp`): ```r -$ data: [data array] -$ dims: [dimensions vector] -$ coords: [List of coordinates vectors] - $ sdate - $ time - $ lon - [...] -$ attrs: [List of the attributes] - $ Variable: - $ varName - $ metadata - $ Datasets - $ Dates - $ source_files - $ when - $ load_parameters +'s2dv_cube' +Data [ 279.99, 280.34, 279.45, 281.99, 280.92, ... ] +Dimensions ( dataset = 1, var = 1, member = 15, sdate = 6, ftime = 3, lat = 22, lon = 53 ) +Coordinates + * dataset : dat1 + * var : tas + member : 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15 + * sdate : 20001101, 20011101, 20021101, 20031101, 20041101, 20051101 + ftime : 1, 2, 3 + * lat : 48, 47, 46, 45, 44, 43, 42, 41, 40, 39, 38, 37, 36, 35, 34, ... + * lon : 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, ... +Attributes + Dates : 2000-11-01 2001-11-01 2002-11-01 2003-11-01 2004-11-01 ... + varName : tas + metadata : + lat + units : degrees_north + long name : latitude + lon + units : degrees_east + long name : longitude + ftime + units : hours since 2000-11-01 00:00:00 + tas + units : K + long name : 2 metre temperature + Datasets : dat1 + when : 2023-10-02 10:11:06 + source_files : "/ecmwf/system5c3s/monthly_mean/tas_f6h/tas_20001101.nc" ... + load_parameters : + ( dat1 ) : dataset = dat1, var = tas, sdate = 20001101 ... ``` -More information about the `s2dv_cube` object class can be found here: [description of the s2dv_cube object structure document](https://docs.google.com/document/d/1ko37JFl_h6mOjDKM5QSQGikfLBKZq1naL11RkJIwtMM/edit?usp=sharing). +This package is designed to be compatible with other R packages such as [s2dv](https://CRAN.R-project.org/package=s2dv), [startR](https://CRAN.R-project.org/package=startR), [CSIndicators](https://CRAN.R-project.org/package=CSIndicators), [CSDownscale](https://earth.bsc.es/gitlab/es/csdownscale). -The current `s2dv_cube` object (CSTools 5.0.0) differs from the original object used in the previous versions of the packages. If you have **questions** on this change you can follow some of the points below: +More information about the `s2dv_cube` object class can be found here: [description of the s2dv_cube object structure document](https://docs.google.com/document/d/1ko37JFl_h6mOjDKM5QSQGikfLBKZq1naL11RkJIwtMM/edit?usp=sharing). -- [New s2dv_cube object discussion](https://earth.bsc.es/gitlab/external/cstools/-/issues/94) -- [How to deal with the compatibility break](https://earth.bsc.es/gitlab/external/cstools/-/issues/112) -- [Testing issue and specifications](https://earth.bsc.es/gitlab/external/cstools/-/issues/110) +**Note:** The current `s2dv_cube` object (CSTools version > 5.0.0) differs from the original object used in the previous versions of the packages. If you have doubts on this change you can follow some of the issues: [New s2dv_cube object discussion](https://earth.bsc.es/gitlab/external/cstools/-/issues/94), [How to deal with the compatibility break](https://earth.bsc.es/gitlab/external/cstools/-/issues/112) and [Testing issue and specifications](https://earth.bsc.es/gitlab/external/cstools/-/issues/110) Contribute ---------- diff --git a/inst/doc/usecase.md b/inst/doc/usecase.md index 8ebf047840f0f886738e09d1d2e01dc0ff0d6153..272964abe77352b93a6691b8b4b544705cdb561c 100644 --- a/inst/doc/usecase.md +++ b/inst/doc/usecase.md @@ -1,12 +1,15 @@ -# Usecase scripts +# Use case and example scripts -In this document, you can link to the example scripts for different usage of the function: +In this document, you will find example scripts of the package. The first ones are use cases of cliimate data assessment. The second ones are example scripts on the use of the 's2dv_cube' object. -1. **Climate data assesment and downscaling** +1. **Use cases of climate data assesment and downscaling** 1. [Bias adjustment for assessment of an extreme event](inst/doc/usecase/UseCase1_WindEvent_March2018.R) 2. [Precipitation Downscaling with RainFARM RF 4](inst/doc/usecase/UseCase2_PrecipitationDownscaling_RainFARM_RF4.R) 3. [Precipitation Downscaling with RainFARM RF 100](inst/doc/usecase/UseCase2_PrecipitationDownscaling_RainFARM_RF100.R) 4. [Seasonal forecasts for a river flow](inst/doc/usecase/UseCase3_data_preparation_SCHEME_model.R) -2. **Examples using 's2dv_cube'** - 1. [Save 's2dv_cube'](inst/doc/usecase/UseCase4_SaveExp.R) \ No newline at end of file +2. **Examples on how to use 's2dv_cube'** + 1. [Create an 's2dv_cube'](inst/doc/usecase/ex1_create.R) + 2. [Save 's2dv_cube'](inst/doc/usecase/ex2_save.R) + 3. [Modify any 's2dv_cube' dimension](inst/doc/usecase/ex3_modify_dims.R) + 4. [Subset any 's2dv_cube' dimension](inst/doc/usecase/ex4_subset.R) diff --git a/inst/doc/usecase/ex1_create.R b/inst/doc/usecase/ex1_create.R new file mode 100644 index 0000000000000000000000000000000000000000..3c07db4118167144b1bf633c28c334915838c21f --- /dev/null +++ b/inst/doc/usecase/ex1_create.R @@ -0,0 +1,212 @@ +#******************************************************************************* +# Title: Example script to create 's2dv_cube' objects +# Author: Eva Rifà Rovira +# Date: 16/01/2024 +#******************************************************************************* +# This example shows how to create an 's2dv_cube' object. +# There are two ways of creating an 's2dv_cube' object. +# (1) With the function s2dv_cube(): create it from scratch with any data. +# (2) With the function CST_Start(). This function returns an 's2dv_cube' +# from an 'startR_array'. + +# Needed packages +library(CSTools) +library(startR) +################################################################################ +#----------------------------------------------------- +# Example 1: Function s2dv_cube() from defined data +#----------------------------------------------------- +# Minimal use case, with s2dv_cube function. + +# In this example we use the function s2dv_cube() to create an object of class +# 's2dv_cube' with the correct structure. + +# (1) We define the array with named dimensions: +dat <- array(1:100, dim = c(time = 10, lat = 4, lon = 10)) +# (2) We define the coordinates as a list of vectors: +coords <- list(time = 1:10, lat = 43:40, lon = 0:9) +# (3) The metadata: +metadata <- list(tas = list(level = '2m'), + lon = list(cdo_grid_name = 'r360x181'), + lat = list(cdo_grid_name = 'r360x181')) +# (4) The creation of Dates array. +# First the initial date: +ini_date <- as.POSIXct('2010-01-01', format = '%Y-%m-%d') +# The sequence of dates +dates <- seq(ini_date, by = 'days', length.out = 10) +# We define the dates dimensions +dim(dates) <- c(time = 10) +# (5) We call the function s2dv_cube() +dat_cube <- s2dv_cube(data = dat, coords = coords, + varName = 'tas', metadata = metadata, + Dates = dates, + when = "2019-10-23 19:15:29 CET", + source_files = c("/path/to/file1.nc", "/path/to/file2.nc"), + Datasets = 'test_dataset') + +# We print the result to see the 's2dv_cube' structure: +# > dat_cube +# 's2dv_cube' +# Data [ 1, 2, 3, 4, 5, 6, 7, 8 ... ] +# Dimensions ( time = 10, lat = 4, lon = 10 ) +# Coordinates +# * time : 1, 2, 3, 4, 5, 6, 7, 8, 9, 10 +# * lat : 43, 42, 41, 40 +# * lon : 0, 1, 2, 3, 4, 5, 6, 7, 8, 9 +# Attributes +# Dates : 2010-01-01 2010-01-02 2010-01-03 2010-01-04 2010-01-05 ... +# varName : tas +# metadata : +# tas +# other : level +# lon +# other : cdo_grid_name +# lat +# other : cdo_grid_name +# Datasets : test_dataset +# when : 2019-10-23 19:15:29 CET +# source_files : /path/to/file1.nc ... + +#----------------------------------------------------- +# Example 2: Function as.s2dv_cube() +#----------------------------------------------------- +# (1) Example using CST_Start + +# NOTE 1: CST_Start() is just a wrapper of function Start() with the +# transformation to 's2dv_cube' object. +# NOTE 2: In order that the input argument auxiliary functions from startR +# work, we need to call them explicitly the startR namespace. +# (e.g. startR::indices()) + +# We just need to define a CST_Start call with all the information: + +repos1 <- "/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" + +res <- CST_Start(dat = list(list(name = 'system4_m1', path = repos2), + list(name = 'system5_m1', path = repos1)), + var = c('tas', 'sfcWind'), + sdate = c('20160101', '20170101'), + ensemble = startR::indices(1:2), + time = startR::indices(1:2), + lat = startR::indices(1:10), + lon = startR::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 = TRUE) + + +# Now we can explore the object: + +# 1st level +names(res) +# "data" "dims" "coords" "attrs" + +dim(res$data) +# dat var sdate ensemble time lat lon +# 2 2 2 2 2 10 10 + +res$coords$lon +# [1] 0.000000 0.703125 1.406250 2.109375 2.812500 3.515625 4.218750 4.921875 +# [9] 5.625000 6.328125 +attr(res$coords$lon, 'indices') +# [1] FALSE +# NOTE: The attribute 'indices' is FALSE, it means that the longitude elements +# are the actual values of longitude coordinate. + +res$coords$ensemble +# [1] 1 2 +# attr(,"indices") +# [1] TRUE + +# Now we take a look into the Dates array. It must have the time dimensions +# of the data. +dim(res$attrs$Dates) +# sdate time +# 2 2 + +# To see the nested list structure of the object, we just need to use the +# function str(): +str(res) + +#----------------------------------------------------- + +# (2) Example using as.s2dv_cube() function + +# We'll load the data with Start and then we'll transform the 'startR_array' +# to 's2dv_cube' object with the function as.s2dv_cube(). We are going +# to load the same data as before, with the same call: + +repos1 <- "/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" + +res <- Start(dat = list(list(name = 'system4_m1', path = repos2), + list(name = 'system5_m1', path = repos1)), + var = c('tas', 'sfcWind'), + sdate = c('20160101', '20170101'), + ensemble = startR::indices(1:2), + time = startR::indices(1:2), + lat = startR::indices(1:10), + lon = startR::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 = TRUE) + +# Now, we use the function as.s2dv_cube() to transform the 'startR_array' +# into an 's2dv_cube': +res_cube <- as.s2dv_cube(res) + +# If we call directly the object directly into the terminal, we can see +# all the elements nicely: + +# > res_cube +# 's2dv_cube' +# Data [ 248.241973876953, 247.365753173828, 6.80753087997437, 5.46453714370728, 247.256896972656, 248.500869750977, 6.25862503051758, 5.76889991760254 ... ] +# Dimensions ( dat = 2, var = 2, sdate = 2, ensemble = 2, time = 2, lat = 10, lon = 10 ) +# Coordinates +# * dat : system4_m1, system5_m1 +# * var : tas, sfcWind +# * sdate : 20160101, 20170101 +# ensemble : 1, 2 +# time : 1, 2 +# * lat : 89.4628215685774, 88.7669513528422, 88.0669716474306, 87.366063433082, 86.6648030134408, 85.9633721608804, 85.2618460607126, 84.5602613830534, 83.8586381286076, 83.1569881285417 +# * lon : 0, 0.703125, 1.40625, 2.109375, 2.8125, 3.515625, 4.21875, 4.921875, 5.625, 6.328125 +# Attributes +# Dates : 2016-02-01 2017-02-01 2016-03-01 2017-03-01 +# varName : tas sfcWind +# metadata : +# time +# units : hours since 2016-01-01 00:00:00 +# other : ndims, size, standard_name, calendar +# lon +# units : degrees_east +# long name : longitude +# other : ndims, size, standard_name, axis +# lat +# units : degrees_north +# long name : latitude +# other : ndims, size, standard_name, axis +# tas +# units : K +# long name : 2 metre temperature +# other : prec, dim, unlim, make_missing_value, missval, hasAddOffset, hasScaleFact, code, table, grid_type +# sfcWind +# units : m s**-1 +# long name : 10 meter windspeed +# other : prec, dim, unlim, make_missing_value, missval, hasAddOffset, hasScaleFact, code, table, grid_type +# Datasets : system4_m1 ... +# when : 2024-01-17 11:38:27 +# source_files : /esarchive/exp/ecmwf/system4_m1/monthly_mean/tas_f6h/tas_20160101.nc ... +# load_parameters : +# ( system4_m1 ) : dat = system4_m1, var = tas ..., sdate = 20160101 ... +# ... + +################################################################################ \ No newline at end of file diff --git a/inst/doc/usecase/UseCase4_CST_SaveExp.R b/inst/doc/usecase/ex2_save.R similarity index 96% rename from inst/doc/usecase/UseCase4_CST_SaveExp.R rename to inst/doc/usecase/ex2_save.R index 926754f3c1eb4f66a0bb81c5e325e4d8e60afb11..7fedd3492a6cde143b595995173406d0d1bb8c9e 100644 --- a/inst/doc/usecase/UseCase4_CST_SaveExp.R +++ b/inst/doc/usecase/ex2_save.R @@ -1,25 +1,20 @@ #******************************************************************************* -# Script to test examples of CST_SaveExp -# Eva Rifà Rovira -# 29/11/2024 +# Title: Example script to save 's2dv_cube' to NetCDF using CST_SaveExp +# Author: Eva Rifà Rovira +# Date: 29/11/2024 #******************************************************************************* +# In this script, we'll see multiple ways to store the 's2dv_cube' (CST_SaveExp) +# or the multidimensional array (SaveExp) to NetCDF. -#------------------------------------------------------------------------------- -# Needed packages before a new version is installed +# Needed packages: +library(CSTools) library(CSIndicators) -library(multiApply) -library(easyNCDF) library(s2dv) -library(ClimProjDiags) -library(CSTools) library(startR) -source("https://earth.bsc.es/gitlab/external/cstools/-/raw/develop-SaveCube/R/CST_SaveExp.R") -source("https://earth.bsc.es/gitlab/external/cstools/-/raw/develop-SaveCube/R/zzz.R") ################################################################################ -# Tests: #----------------------------------------------------- -# Tests 1: Multidimensional array and Dates, without metadata and coordinates +# Example 1: Multidimensional array and Dates, without metadata and coordinates #----------------------------------------------------- # (1.1) Minimal use case, without Dates data <- array(1:5, dim = c(sdate = 5, lon = 4, lat = 4)) @@ -70,7 +65,7 @@ SaveExp(data, ftime_dim = 'ftime', memb_dim = NULL, dat_dim = NULL, extra_string = 'test', global_attrs = list(system = 'tes1', reference = 'test2')) #----------------------------------------------------- -# Tests 2: Test sample data from Start and from Load +# Example 2: Test sample data from Start and from Load #----------------------------------------------------- # (2.1) Test SaveExp @@ -267,7 +262,7 @@ CST_SaveExp(data = data, ftime_dim = NULL, single_file = FALSE, units_hours_since = FALSE) #----------------------------------------------------- -# Test 3: Special cases +# Example 3: Special cases #----------------------------------------------------- # (3.1) Two variables and two datasets in separated files @@ -377,7 +372,7 @@ CST_SaveExp(data = obscube, ftime_dim = 'time', var_dim = 'var', single_file = FALSE, extra_string = 'obs_tas') #----------------------------------------------------- -# Test 4: Time bounds: +# Example 4: Time bounds: #----------------------------------------------------- # example: /esarchive/exp/ncep/cfs-v2/weekly_mean/s2s/tas_f24h/tas_20231128.nc @@ -440,7 +435,7 @@ CST_SaveExp(data = res, ftime_dim = 'time', var_dim = 'var', units_hours_since = FALSE) #----------------------------------------------------- -# Test 5: Read data with Load +# Example 5: Read data with Load #----------------------------------------------------- data <- lonlat_temp$exp @@ -462,5 +457,7 @@ CST_SaveExp(data = data, ftime_dim = 'ftime', # latmin = 27, latmax = 48, # lonmin = -12, lonmax = 40, # output = "lonlat") -# Error + +# NOTE: This case hasn't been developed since the function to load data +# that will be maintianed will be CST_Start. ################################################################################ \ No newline at end of file diff --git a/inst/doc/usecase/ex3_modify_dims.R b/inst/doc/usecase/ex3_modify_dims.R new file mode 100644 index 0000000000000000000000000000000000000000..6d1fa5b87b2300231bb954cf43a242b9362aeeea --- /dev/null +++ b/inst/doc/usecase/ex3_modify_dims.R @@ -0,0 +1,165 @@ +#******************************************************************************* +# Title: Script to modify the dimensions of the 's2dv_cube' +# Author: Eva Rifà Rovira +# Date: 18/01/2024 +#******************************************************************************* +# In this example, we will explore different methods to modify the dimensions +# of the 's2dv_cube': +# (1) Changing dimension names +# (2) Adding new dimensions +# (3) Merge 2 dimensions +# (4) Split a dimension + +# Needed packages: +library(CSTools) + +################################################################################ +#----------------------------------------------------- +# Example 1: Change dimension names with CST_ChangeDimNames +#----------------------------------------------------- +# With using this function, we can change the dimension names in all elements +# of the 's2dv_cube' object: + +# (1) Check original dimensions and coordinates +lonlat_temp$exp$dims +names(lonlat_temp$exp$coords) +dim(lonlat_temp$exp$attrs$Dates) +# (2) Change 'dataset' to 'dat' and 'ftime' to 'time' +exp <- CST_ChangeDimNames(lonlat_temp$exp, + original_names = c("dataset", "ftime", "lon", "lat"), + new_names = c("dat", "time", "longitude", "latitude")) +# (3) Check new dimensions and coordinates +exp$dims +names(exp$coords) +dim(exp$attrs$Dates) + +#----------------------------------------------------- +# Example 2: Insert a new dimension with CST_InsertDim +#----------------------------------------------------- +# With this function, we can add a dimension into the 's2dv_cube'. +# NOTE: When the dimension that we want to add has length greater than 1, the +# values of the data are repeated for that new dimension. + +# (1) Check original dimensions and coordinates +lonlat_temp$exp$dims +names(lonlat_temp$exp$coords) +# (2) Add 'variable' dimension +exp <- CST_InsertDim(lonlat_temp$exp, + posdim = 2, + lendim = 2, + name = "variable", + values = c("tas", "tos")) +# (3) Check new dimensions and coordinates +exp$dims +exp$coords$variable +# We see that the values will be repeated along the new dimension: +exp$data[, , 1, 1, 1, 1, 1] + +#----------------------------------------------------- +# Example 3: Merge two dimensions with CST_MergeDims +#----------------------------------------------------- +# In this example, we will merge the dimensions corresponding to the latitude +# and the longitude of the data. The new dimension will be named 'grid'. + +# (1) Call the function: +new_data <- CST_MergeDims(lonlat_temp$exp, merge_dims = c('lat', 'lon'), + rename_dim = 'grid') + +# (2) Check the dimensions of the data: +dim(new_data$data) +# dataset member sdate ftime grid +# 1 15 6 3 1166 + +# (3) Check the names of the coordinates: +names(new_data$coords) +# [1] "dataset" "member" "sdate" "ftime" "grid" + +# (4) Explore the object by printing it in the terminal: +new_data +# NOTE: Be aware that when we print the object, we see that its name in +# "Coordinates" field appears without the asterisk (*) at its left. This means +# that the values of that coordinate, are indices, not the actual values. We +# can also find this information with the attribute "indices": +attributes(new_data$coords$grid) +# $indices +# [1] TRUE + +# (5) Now, we want to merge time dimensions start date and forecast time: +new_data <- CST_MergeDims(data = lonlat_temp_st$exp, merge_dims = c('sdate', 'ftime')) +# In this case, the Dates dimensions will be merged too. +# (6) Check the dimensions of Dates: +dim(new_data$attrs$Dates) +# sdate +# 18 + +# NOTE: When we want to merge temporal and other dimensions nature, +# the Dates dimensions are kept as the original. In this case, the function +# returns a Warning Message, we must pay attention! +new_data <- CST_MergeDims(data = lonlat_temp$exp, + merge_dims = c('lat', 'ftime'), + rename_dim = 'newdim') + +#----------------------------------------------------- +# Example 4: Split two dimensions with SplitDim and CST_SplitDim +#----------------------------------------------------- +# In this example, we will start working with the function SplitDim, +# that it can be used to split dimensions of an array. + +# NOTE: Take into account that time dimensions will be treated differently than +# other dimensions: + +# (1) Decadal example: We define an array of consecutive days of different years: +dates <- seq(as.Date("01-01-2000", "%d-%m-%Y", tz = 'UTC'), + as.Date("31-12-2005","%d-%m-%Y", tz = 'UTC'), "day") +dim(dates) <- c(time = 2192) + +# (2) Now, we will split the array in a new 'year' dimension: +dates_year <- SplitDim(dates, indices = dates, + split_dim = 'time', freq = 'year') +# time year +# 366 6 + +# (3) Now, we can try: freq = 'month' and 'day' +dates_month <- SplitDim(dates, indices = dates, + split_dim = 'time', freq = 'month') + +dates_day <- SplitDim(dates, indices = dates, + split_dim = 'time', freq = 'day') + +# (4) Finnally, we need to convert them again from numeric to 'POSIXct': +dates_year <- as.POSIXct(dates_year * 24 * 3600, origin = '1970-01-01', tz = 'UTC') +dates_month <- as.POSIXct(dates_month * 24 * 3600, origin = '1970-01-01', tz = 'UTC') +dates_day <- as.POSIXct(dates_day * 24 * 3600, origin = '1970-01-01', tz = 'UTC') + +#----------------------------------------------------- + +# In the following example, we will use the sample data of the package. We +# will use lonlat_prec_st because it is daily data: + +# NOTE: By Jan 2024, a development is needed regarding updates in other fields +# of the 's2dv_cube' + +# (1) Call the function CST_SplitDim with adding 'day' dimension: +data_day <- CST_SplitDim(lonlat_prec_st, indices = lonlat_prec_st$attrs$Dates[1,], + split_dim = 'ftime', freq = 'day') +# (2) Explore the dimensions of the data array +dim(data_day$data) +# dataset var member sdate ftime lat lon day +# 1 1 6 3 1 4 4 31 + +# (3) Call the function CST_SplitDim with adding 'month' dimension: +data_month <- CST_SplitDim(lonlat_prec_st, indices = lonlat_prec_st$attrs$Dates[1,], + split_dim = 'ftime', freq = 'month') + +dim(data_month$data) +# dataset var member sdate ftime lat lon month +# 1 1 6 3 31 4 4 1 + +# (4) Call the function CST_SplitDim with adding 'year' dimension: +data_year <- CST_SplitDim(lonlat_prec_st, indices = lonlat_prec_st$attrs$Dates[,1], + split_dim = 'sdate', freq = 'year') +dim(data_year$data) +# dataset var member sdate ftime lat lon year +# 1 1 6 1 31 4 4 3 + +################################################################################ \ No newline at end of file diff --git a/inst/doc/usecase/ex4_subset.R b/inst/doc/usecase/ex4_subset.R new file mode 100644 index 0000000000000000000000000000000000000000..360ca08fc479f42b12d421ce8ffa5d82b8a12e47 --- /dev/null +++ b/inst/doc/usecase/ex4_subset.R @@ -0,0 +1,140 @@ +#******************************************************************************* +# Title: Example script to subset any dimension of an 's2dv_cube' +# Author: Eva Rifà Rovira +# Date: 16/11/2024 +#******************************************************************************* +# This example shows how to subset any dimension of an 's2dv_cube'. To do it, +# we will use the function CST_Subset. This function is the 's2dv_cube' method +# version of Subset from the package ClimProjDiags. +# (1) First we will see how Subset works. +# (2) Then, we will use CST_Subset with an 's2dv_cube' + +# Needed packages: +library(CSTools) +library(ClimProjDiags) + +################################################################################ +#----------------------------------------------------- +# Example 1: Subset an example array +#----------------------------------------------------- +# This is a minimal use case about spatial coordinates subset. + +# (1) We create the array amd we print it: +dat <- array(1:100, dim = c(lat = 10, lon = 10)) +dat +# [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] +# [1,] 1 11 21 31 41 51 61 71 81 91 +# [2,] 2 12 22 32 42 52 62 72 82 92 +# [3,] 3 13 23 33 43 53 63 73 83 93 +# [4,] 4 14 24 34 44 54 64 74 84 94 +# [5,] 5 15 25 35 45 55 65 75 85 95 +# [6,] 6 16 26 36 46 56 66 76 86 96 +# [7,] 7 17 27 37 47 57 67 77 87 97 +# [8,] 8 18 28 38 48 58 68 78 88 98 +# [9,] 9 19 29 39 49 59 69 79 89 99 +# [10,] 10 20 30 40 50 60 70 80 90 100 + +# (2) We call the function Subset from ClimProjDiags and we see the result: +dat_subset <- Subset(x = dat, along = c('lat', 'lon'), indices = list(1:5, 1:7), + drop = 'all') +dat_subset +# [,1] [,2] [,3] [,4] [,5] [,6] [,7] +# [1,] 1 11 21 31 41 51 61 +# [2,] 2 12 22 32 42 52 62 +# [3,] 3 13 23 33 43 53 63 +# [4,] 4 14 24 34 44 54 64 +# [5,] 5 15 25 35 45 55 65 + +#----------------------------------------------------- +# Example 2: Subset an 's2dv_cube' using sample data +#----------------------------------------------------- +# In this example we will not drop any dimension, we will select only the first +# member, the first and the second start dates, and also subset the longitude and +# keep only the values from [0, 21]: + +# (1) Explore the sample data: +dat <- lonlat_temp_st$exp + +dat$dims +# dataset var member sdate ftime lat lon +# 1 1 15 6 3 22 53 + +dat +# 's2dv_cube' +# Data [ 279.994110107422, 280.337463378906, 279.450866699219, ... ] +# Dimensions ( dataset = 1, var = 1, member = 15, sdate = 6, ftime = 3, +# lat = 22, lon = 53 ) +# Coordinates +# * dataset : dat1 +# * var : tas +# member : 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15 +# * sdate : 20001101, 20011101, 20021101, 20031101, 20041101, 20051101 +# ftime : 1, 2, 3 +# * lat : 48, 47, 46, 45, 44, 43, 42, 41, 40, 39, 38, 37, 36, 35, 34, ... +# * lon : 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, ... +# Attributes +# Dates : 2000-11-01 2001-11-01 2002-11-01 2003-11-01 2004-11-01 ... +# varName : tas +# metadata : +# lat +# units : degrees_north +# long name : latitude +# lon +# units : degrees_east +# long name : longitude +# ftime +# units : hours since 2000-11-01 00:00:00 +# tas +# units : K +# long name : 2 metre temperature +# Datasets : dat1 +# when : 2023-10-02 10:11:06 +# source_files : /monthly_mean/tas_f6h/tas_20001101.nc ... +# load_parameters : +# ( dat1 ) : dataset = dat1, var = tas, sdate = 20001101 ... +# ... + +# (2) Call the function CST_Subset: +dat_subset <- CST_Subset(x = dat, along = c('member', 'sdate', 'lon'), + indices = list(1, 1:2, 1:22), drop = 'none') + +# (3) Explore the 's2dv_cube' +dat_subset +# 's2dv_cube' +# Data [ 279.994110107422, 277.161102294922, 278.825836181641, 276.8271484375, 276.052703857422, 276.950805664062, 280.677215576172, 277.285247802734 ... ] +# Dimensions ( dataset = 1, var = 1, member = 1, sdate = 2, ftime = 3, lat = 22, lon = 22 ) +# Coordinates +# * dataset : dat1 +# * var : tas +# member : 1 +# * sdate : 20001101, 20011101 +# ftime : 1, 2, 3 +# * lat : 48, 47, 46, 45, 44, 43, 42, 41, 40, 39, 38, 37, 36, 35, 34, 33, 32, 31, 30, 29, 28, 27 +# * lon : 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21 +# Attributes +# Dates : 2000-11-01 2001-11-01 2000-12-01 2001-12-01 2001-01-01 ... +# varName : tas +# metadata : +# ftime +# units : hours since 2000-11-01 00:00:00 +# other : ndims, size, standard_name, calendar +# lat +# units : degrees_north +# long name : latitude +# other : ndims, size, standard_name, axis +# lon +# units : degrees_east +# long name : longitude +# other : ndims, size, standard_name, axis +# tas +# units : K +# long name : 2 metre temperature +# other : prec, dim, unlim, make_missing_value, missval, hasAddOffset, hasScaleFact, code, table +# Datasets : dat1 +# when : 2023-10-02 10:11:06 +# source_files : /esarchive/exp/ecmwf/system5c3s/monthly_mean/tas_f6h/tas_20001101.nc ... +# load_parameters : +# ( dat1 ) : dataset = dat1, var = tas, sdate = 20001101 ... +# ... + +################################################################################ \ No newline at end of file diff --git a/tests/testthat/test-CST_MergeDims.R b/tests/testthat/test-CST_MergeDims.R index f7eac6acf5369988d1bdff0f0a38921810b8f215..a99d5eba65e56aec072727dd582230f5685d34f3 100644 --- a/tests/testthat/test-CST_MergeDims.R +++ b/tests/testthat/test-CST_MergeDims.R @@ -1,65 +1,96 @@ ############################################## -test_that("Sanity checks", { +# data1 +data1 <- list(data = 1:10) +class(data1) <- 's2dv_cube' + +# data2 +data <- 1 : 20 +dim(data) <- c(time = 20) +data2 <- list(data = data) +data2$dims <- dim(data) +data2$coords <- list(time = 1:20) +attr(data2$coords$time, 'indices') <- TRUE +class(data2) <- 's2dv_cube' + +# exp +exp <- 1 : 20 +dim(exp) <- c(time = 10, lat = 2) +exp <- list(data = exp) +class(exp) <- 's2dv_cube' + +# data3 +data3 <- data2 +names(dim(data3$data)) <- 'Dim1' +data3$dims <- dim(data3$data) +names(data3$coords) <- 'Dim1' + +############################################## + +test_that("1. Sanity checks", { expect_error( CST_MergeDims(data = 1), - paste0("Parameter 'data' must be of the class 's2dv_cube'.")) -data <- list(data = 1:10) -class(data) <- 's2dv_cube' + paste0("Parameter 'data' must be of the class 's2dv_cube'.") + ) expect_error( - CST_MergeDims(data = data), - paste0("Parameter 'data' must have dimensions.")) + CST_MergeDims(data = data1), + paste0("Parameter 'data' must have dimensions.") + ) - data <- 1 : 20 - dim(data) <- c(time = 20) - data <- list(data = data) - class(data) <- 's2dv_cube' expect_error( - CST_MergeDims(data = data), - "Parameter 'merge_dims' must match with dimension names in parameter 'data'.") + CST_MergeDims(data = data2), + "Parameter 'merge_dims' must match with dimension names in parameter 'data'." + ) expect_error( - CST_MergeDims(data = data, merge_dims = 1), - paste0("Parameter 'merge_dims' must be a character vector indicating the names", - " of the dimensions to be merged.")) + CST_MergeDims(data = data2, merge_dims = 1), + paste0("Parameter 'merge_dims' must be a character vector indicating the ", + "names of the dimensions to be merged.") + ) expect_error( - CST_MergeDims(data = data, merge_dims = 'time'), - "Parameter 'merge_dims' must be of length two.") + CST_MergeDims(data = data2, merge_dims = 'time'), + "Parameter 'merge_dims' must be of length two." + ) expect_error( - CST_MergeDims(data = data, merge_dims = c('time', 'sdates')), + CST_MergeDims(data = data2, merge_dims = c('time', 'sdates')), paste0("Parameter 'merge_dims' must match with dimension ", - "names in parameter 'data'.")) + "names in parameter 'data'.") + ) +}) + +############################################## - exp <- 1 : 20 - dim(exp) <- c(time = 10, lat = 2) - exp <- list(data = exp) - class(exp) <- 's2dv_cube' +test_that("2. Output checks", { expect_equal( - CST_MergeDims(data = exp, merge_dims = c('time', 'lat')), data) - + CST_MergeDims(data = exp, merge_dims = c('time', 'lat')), + data2 + ) expect_warning( CST_MergeDims(data = exp, merge_dims = c('time', 'lat', 'lon')), paste0("Only two dimensions can be merge, only the first two dimension", " will be used. To merge further dimensions consider to use this ", - "function multiple times.")) + "function multiple times.") + ) expect_warning( - CST_MergeDims(data = exp, merge_dims = c('time', 'lat'), rename_dim = c('lat', 'lon')), + CST_MergeDims(data = exp, merge_dims = c('time', 'lat'), + rename_dim = c('lat', 'lon')), paste0("Parameter 'rename_dim' has length greater than 1 and only the ", - "first element will be used.")) - names(dim(data$data)) <- 'Dim1' + "first element will be used.") + ) expect_equal( - CST_MergeDims(data = exp, merge_dims = c('time', 'lat'), rename_dim = 'Dim1'), - data) - + CST_MergeDims(data = exp, merge_dims = c('time', 'lat'), + rename_dim = 'Dim1'), + data3 + ) expect_equal( CST_MergeDims(data = exp, merge_dims = c('time', 'lat'), - rename_dim = 'Dim1', na.rm = TRUE), data) - - exp$data[1,] <- NA - data <- c(2 : 10, 12 : 20) - dim(data) <- c(Dim1 = 18) - data <- list(data = data) - class(data) <- 's2dv_cube' + rename_dim = 'Dim1', na.rm = TRUE), + data3 + ) expect_equal( CST_MergeDims(data = exp, merge_dims = c('time', 'lat'), - rename_dim = 'Dim1', na.rm = TRUE), data) + rename_dim = 'Dim1', na.rm = TRUE), + data3 + ) }) + +############################################## \ No newline at end of file