Commit d4a711c7 authored by aho's avatar aho
Browse files

Merge branch 'develop-transform_chunk_sri' into 'master'

Use destination grid instead of values to count sri

See merge request !190
parents edab7160 e715b34b
Pipeline #7694 passed with stage
in 52 minutes and 16 seconds
......@@ -3046,97 +3046,18 @@ Start <- function(..., # dim = indices/selectors,
sub_array_of_sri <- sub_array_of_sri[[1]]:sub_array_of_sri[[2]]
}
# Chunk sub_array_of_sri if this inner_dim needs to be chunked
#TODO: Potential problem: the transformed_subset_var value falls between
# the end of sub_sub_array_of_values of the 1st chunk and the beginning
# of sub_sub_array_of_values of the 2nd chunk. Then, one sub_array_of_sri
# will miss. 'previous_sri' is checked and will be included if this
# situation happens, but don't know if the transformed result is
# correct or not.
# NOTE: The chunking criteria may not be 100% correct. The current way
# is to pick the sri that larger than the minimal sub_sub_array_of_values
# and smaller than the maximal sub_sub_array_of_values; if it's
# the first chunk, make sure the 1st sri is included; if it's the
# last chunk, make sure the last sri is included.
if (chunks[[inner_dim]]["n_chunks"] > 1) {
sub_array_of_sri_complete <- sub_array_of_sri
if (is.list(sub_sub_array_of_values)) { # list
sub_array_of_sri <-
which(transformed_subset_var >= min(unlist(sub_sub_array_of_values)) &
transformed_subset_var <= max(unlist(sub_sub_array_of_values)))
# if it's 1st chunk & the first sri is not included, include it.
if (chunks[[inner_dim]]["chunk"] == 1 &
!(sub_array_of_sri_complete[1] %in% sub_array_of_sri)) {
sub_array_of_sri <- c(sub_array_of_sri_complete[1], sub_array_of_sri)
}
# if it's last chunk & the last sri is not included, include it.
if (chunks[[inner_dim]]["chunk"] == chunks[[inner_dim]]["n_chunks"] &
!(tail(sub_array_of_sri_complete, 1) %in% sub_array_of_sri)) {
sub_array_of_sri <- c(sub_array_of_sri, tail(sub_array_of_sri_complete, 1))
}
#========================================================
# Check if sub_array_of_sri perfectly connects to the previous sri.
# If not, inlclude the previous sri.
#NOTE 1: don't know if the transform for the previous sri is
# correct or not.
#NOTE 2: If crop = T, sub_array_of_sri always starts from 1.
# Don't know if the cropping will miss some sri or not.
if (sub_array_of_sri[1] != 1) {
if (!is.null(previous_sub_sub_array_of_values)) {
# if decreasing = F
if (transformed_subset_var[1] < transformed_subset_var[2]) {
previous_sri <- max(which(transformed_subset_var <= previous_sub_sub_array_of_values))
} else {
# if decreasing = T
previous_sri <- max(which(transformed_subset_var >= previous_sub_sub_array_of_values))
}
if (previous_sri + 1 != sub_array_of_sri[1]) {
sub_array_of_sri <- (previous_sri + 1):sub_array_of_sri[length(sub_array_of_sri)]
}
}
}
} else { # is vector
tmp <- which(transformed_subset_var >= min(sub_sub_array_of_values) &
transformed_subset_var <= max(sub_sub_array_of_values))
# Ensure tmp and sub_array_of_sri are both ascending or descending
if (is.unsorted(tmp) != is.unsorted(sub_array_of_sri)) {
tmp <- rev(tmp)
}
# Include first or last sri if tmp doesn't have. It's only for
# ""vectors"" because vectors look for the closest value.
#NOTE: The condition here is not correct. The criteria should be
# 'vector' instead of indices.
if (chunks[[inner_dim]]["chunk"] == 1) {
sub_array_of_sri <- unique(c(sub_array_of_sri[1], tmp))
} else if (chunks[[inner_dim]]["chunk"] ==
chunks[[inner_dim]]["n_chunks"]) { # last chunk
sub_array_of_sri <- unique(c(tmp, sub_array_of_sri[length(sub_array_of_sri)]))
} else {
sub_array_of_sri <- tmp
}
# Check if sub_array_of_sri perfectly connects to the previous sri.
# If not, inlclude the previous sri.
#NOTE 1: don't know if the transform for the previous sri is
# correct or not.
#NOTE 2: If crop = T, sub_array_of_sri always starts from 1.
# Don't know if the cropping will miss some sri or not.
if (sub_array_of_sri[1] != 1) {
if (!is.null(previous_sub_sub_array_of_values)) {
# if decreasing = F
if (transformed_subset_var[1] < transformed_subset_var[2]) {
previous_sri <- max(which(transformed_subset_var <= previous_sub_sub_array_of_values))
} else {
# if decreasing = T
previous_sri <- max(which(transformed_subset_var >= previous_sub_sub_array_of_values))
}
if (previous_sri + 1 != which(sub_array_of_sri[1] == sub_array_of_sri_complete)) {
sub_array_of_sri <- (previous_sri + 1):sub_array_of_sri[length(sub_array_of_sri)]
}
}
}
}
# Instead of using values to find sri, directly use the destination grid to count.
#NOTE: sub_array_of_sri seems to start at 1 always (because crop = c(lonmin, lonmax, latmin, latmax) already?)
if (chunks[[inner_dim]]["n_chunks"] > 1) {
sub_array_of_sri <- sub_array_of_sri[get_chunk_indices(
length(sub_array_of_sri),
chunks[[inner_dim]]["chunk"],
chunks[[inner_dim]]["n_chunks"],
inner_dim)]
}
#========================================================
ordered_sri <- sub_array_of_sri
sub_array_of_sri <- transformed_subset_var_unorder[sub_array_of_sri]
......@@ -3333,6 +3254,11 @@ Start <- function(..., # dim = indices/selectors,
selector_store_position <- chunk
}
sub_array_of_indices <- transformed_indices[which(indices_chunk == chunk)]
#NOTE: This 'with_transform' part is probably not tested because
# here is for the inner dim that goes across a file dim, which
# is normally not lat and lon dimension. If in the future, we
# can interpolate time, this part needs to be examined.
if (with_transform) {
# If the provided selectors are expressed in the world
# before transformation
......
......@@ -15,33 +15,36 @@ data <- Start(dat = path,
member = indices(1:2),
transform = CDORemapper,
transform_extra_cells = 2,
transform_params = list(grid = 'r100x50', method = 'conservative', crop = FALSE),
transform_vars = c('lat','lon'),
transform_params = list(grid = 'r100x50', method = 'conservative'),
transform_vars = c('lat', 'lon'),
synonims = list(lon = c('lon', 'longitude'), lat = c('lat', 'latitude')),
return_vars = list(lat = 'dat', lon = 'dat', time = 'sdate'),
retrieve = FALSE)
)
func <- function(x) {
a <- mean(x, na.rm = TRUE)
return(a)
}
step <- Step(func, target_dims = c('time'),
output_dims = NULL)
wf <- AddStep(data, step)
func <- function(x) {
return(x)
}
step <- Step(func, target_dims = c('time'),
output_dims = 'time')
wf <- AddStep(data, step)
#---
suppressWarnings(
res <- Compute(wf,
chunks = list(member = 2))
resT <- eval(data)
)
suppressWarnings(
res1 <- Compute(wf, chunks = list(member = 2))
)
expect_equal(
dim(res$output1),
c(dat = 1, var = 1, sdate = 1, lat = 50, lon = 100, fyear = 1, member = 2)
dim(res1$output1),
c(time = 2, dat = 1, var = 1, sdate = 1, lat = 50, lon = 100, fyear = 1, member = 2)
)
expect_equal(
res$output1[1, 1, 1, 10:12, 20, 1, 1],
c(274.2808, 275.8509, 277.7623),
tolerance = 0.0001
as.vector(resT),
as.vector(res1$output1)
)
})
......@@ -49,8 +52,7 @@ tolerance = 0.0001
test_that("2. chunk along lon", {
#skip_on_cran()
#!!!!!!!!!!!!!!!!!!!NOTE: the results are not identical when exp has extra cells = 2!!!!!!!!!!!!!!!!!!
# But exp2 (retrieve = T) has the same results with extra_cells = 2 and 8.
#NOTE: the results are not identical when exp has extra cells = 2
suppressWarnings(
exp <- Start(dat = '/esarchive/exp/ecmwf/system5_m1/monthly_mean/$var$_f6h/$var$_$sdate$.nc',
......@@ -64,8 +66,8 @@ exp <- Start(dat = '/esarchive/exp/ecmwf/system5_m1/monthly_mean/$var$_f6h/$var$
longitude_reorder = CircularSort(0, 360),
transform = CDORemapper,
transform_params = list(grid = 'r100x50',
method = 'con',
crop = FALSE),
method = 'con'),
# crop = FALSE),
transform_vars = c('latitude', 'longitude'),
transform_extra_cells = 8,
synonims = list(latitude = c('lat', 'latitude'),
......@@ -81,46 +83,35 @@ exp <- Start(dat = '/esarchive/exp/ecmwf/system5_m1/monthly_mean/$var$_f6h/$var$
}
step <- Step(func, target_dims = 'time', output_dims = 'time')
wf <- AddStep(exp, step)
#---
suppressWarnings(
res <- Compute(wf, chunks = list(longitude = 2))
resT <- eval(exp)
)
suppressWarnings(
res2 <- Compute(wf, chunks = list(ensemble = 1))
res1 <- Compute(wf, chunks = list(longitude = 2))$output1
)
expect_equal(
res$output1,
res2$output1
#suppressWarnings(
# res2 <- Compute(wf, chunks = list(ensemble = 1))
#)
suppressWarnings(
res3 <- Compute(wf, chunks = list(longitude = 2, latitude = 2))$output1
)
# Check with retrieve = TRUE
suppressWarnings(
exp2 <- Start(dat = '/esarchive/exp/ecmwf/system5_m1/monthly_mean/$var$_f6h/$var$_$sdate$.nc',
var = 'tas',
sdate = '20000101',
ensemble = indices(1),
time = indices(1),
latitude = 'all',
latitude_reorder = Sort(),
longitude = 'all',
longitude_reorder = CircularSort(0, 360),
transform = CDORemapper,
transform_params = list(grid = 'r100x50',
method = 'con',
crop = FALSE),
transform_vars = c('latitude', 'longitude'),
transform_extra_cells = 2,
synonims = list(latitude = c('lat', 'latitude'),
longitude = c('longitude', 'lon')),
return_vars = list(latitude = 'dat',#NULL,
longitude = 'dat',#NULL,
time = 'sdate'),
retrieve = T)
res4 <- Compute(wf, chunks = list(longitude = 3))$output1
)
expect_equal(
as.vector(res$output1),
as.vector(exp2)
as.vector(resT),
as.vector(res1)
)
expect_equal(
res1,
res3
)
expect_equal(
res1,
res4
)
})
......@@ -40,8 +40,8 @@ exp <- Start(dat = path,
lon_reorder = CircularSort(0, 360),
transform = CDORemapper,
transform_extra_cells = 8,
transform_params = list(grid = 'r100x50', method = 'conservative',
crop = c(0, 360, -90, 90)),
transform_params = list(grid = 'r100x50', method = 'conservative'),
# crop = c(0, 360, -90, 90)),
transform_vars = c('lat','lon'),
synonims = list(lon = c('lon', 'longitude'), lat = c('lat', 'latitude')),
return_vars = list(lat = 'dat', lon = 'dat', time = 'sdate'),
......@@ -53,20 +53,24 @@ func <- function(x) {
step <- Step(func, target_dims = 'time', output_dims = 'time')
wf <- AddStep(exp, step)
#---
suppressWarnings(
res1 <- Compute(wf, chunks = list(lon = 2))
resT <- eval(exp)
)
suppressWarnings(
res2 <- Compute(wf, chunks = list(ensemble = 1))
res1 <- Compute(wf, chunks = list(lon = 2))
)
#suppressWarnings(
#res2 <- Compute(wf, chunks = list(ensemble = 1))
#)
suppressWarnings(
res3 <- Compute(wf, chunks = list(lon = 3))
)
expect_equal(
res1$output1,
res2$output1
as.vector(res1$output1),
as.vector(resT)
)
expect_equal(
res1$output1,
......@@ -88,8 +92,8 @@ exp <- Start(dat = path,
lon_reorder = CircularSort(0, 360),
transform = CDORemapper,
transform_extra_cells = 8,
transform_params = list(grid = 'r100x50', method = 'conservative',
crop = c(0, 360, -90, 90)),
transform_params = list(grid = 'r100x50', method = 'conservative'),
# crop = c(0, 360, -90, 90)),
transform_vars = c('lat','lon'),
synonims = list(lon = c('lon', 'longitude'), lat = c('lat', 'latitude')),
return_vars = list(lat = 'dat', lon = 'dat', time = 'sdate'),
......@@ -101,6 +105,7 @@ func <- function(x) {
step <- Step(func, target_dims = 'time', output_dims = 'time')
wf <- AddStep(exp, step)
#---
suppressWarnings(
res1_list <- Compute(wf, chunks = list(lon = 2))
)
......@@ -109,116 +114,6 @@ res1$output1,
res1_list$output1
)
#-----------------------------------
# crop = FALSE
suppressWarnings(
exp <- Start(dat = path,
var = 'tas',
sdate = '20000101',
ensemble = indices(1),
time = indices(1),
lat = indices(1:640),
lon = indices(1:1296),
lat_reorder = Sort(),
lon_reorder = CircularSort(0, 360),
transform = CDORemapper,
transform_extra_cells = 8,
transform_params = list(grid = 'r100x50', method = 'conservative',
crop = FALSE),
transform_vars = c('lat','lon'),
synonims = list(lon = c('lon', 'longitude'), lat = c('lat', 'latitude')),
return_vars = list(lat = 'dat', lon = 'dat', time = 'sdate'),
retrieve = F)
)
func <- function(x) {
return(x)
}
step <- Step(func, target_dims = 'time', output_dims = 'time')
wf <- AddStep(exp, step)
suppressWarnings(
res_crop_F_1 <- Compute(wf, chunks = list(lon = 2))
)
suppressWarnings(
res_crop_F_2 <- Compute(wf, chunks = list(ensemble = 1))
)
suppressWarnings(
res_crop_F_3 <- Compute(wf, chunks = list(lon = 3))
)
expect_equal(
res1$output1,
res_crop_F_1$output1
)
expect_equal(
res_crop_F_1$output1,
res_crop_F_2$output1
)
expect_equal(
res_crop_F_1$output1,
res_crop_F_3$output1
)
#---------------------------------------------
#!!!!!!!!!!!!!!!!!!!!Problem when global + crop = T + chunk along lon!!!!!!!!!!!!!!!!
# crop = TRUE
suppressWarnings(
exp <- Start(dat = path,
var = 'tas',
sdate = '20000101',
ensemble = indices(1),
time = indices(1),
lat = indices(1:640),
lon = indices(1:1296),
lat_reorder = Sort(),
lon_reorder = CircularSort(0, 360),
transform = CDORemapper,
transform_extra_cells = 8,
transform_params = list(grid = 'r100x50', method = 'conservative',
crop = TRUE),
transform_vars = c('lat','lon'),
synonims = list(lon = c('lon', 'longitude'), lat = c('lat', 'latitude')),
return_vars = list(lat = 'dat', lon = 'dat', time = 'sdate'),
retrieve = F)
)
func <- function(x) {
return(x)
}
step <- Step(func, target_dims = 'time', output_dims = 'time')
wf <- AddStep(exp, step)
#WRONG??????
suppressWarnings(
res_crop_T_1 <- Compute(wf, chunks = list(lon = 2))
)
suppressWarnings(
res_crop_T_2 <- Compute(wf, chunks = list(ensemble = 1))
)
#WRONG?????
suppressWarnings(
res_crop_T_3 <- Compute(wf, chunks = list(lon = 3))
)
suppressWarnings(
res_crop_T_4 <- Compute(wf, chunks = list(lat = 2))
)
expect_equal(
res1$output1,
res_crop_T_1$output1
)
expect_equal(
res1$output1,
res_crop_T_2$output1
)
expect_equal(
res1$output1,
res_crop_T_3$output1
)
expect_equal(
res1$output1,
res_crop_T_4$output1
)
})
......@@ -281,149 +176,45 @@ func <- function(x) {
step <- Step(func, target_dims = 'time', output_dims = 'time')
wf <- AddStep(exp, step)
#---
suppressWarnings(
res <- Compute(wf, chunks = list(lon = 2))
resT <- eval(exp)
)
suppressWarnings(
res2 <- Compute(wf, chunks = list(ensemble = 1))
res1 <- Compute(wf, chunks = list(lon = 2))
)
#suppressWarnings(
#res2 <- Compute(wf, chunks = list(ensemble = 1))
#)
suppressWarnings(
res3 <- Compute(wf, chunks = list(lon = 3))
)
expect_equal(
res$output1,
res2$output1
as.vector(res1$output1),
as.vector(resT)
)
expect_equal(
res$output1,
res1$output1,
res3$output1
)
expect_equal(
drop(res$output1)[, 1],
drop(res1$output1)[, 1],
c(241.4042, 242.5804, 246.8507, 245.8008, 246.4318, 267.0983),
tolerance = 0.001
)
expect_equal(
drop(res$output1)[, 2],
drop(res1$output1)[, 2],
c(241.2223, 242.2564, 245.9863, 244.5377, 244.8937, 266.5749),
tolerance = 0.001
)
expect_equal(
drop(res$output1)[, 3],
drop(res1$output1)[, 3],
c(241.0895, 242.1896, 245.3183, 243.1169, 243.9446, 266.4386),
tolerance = 0.001
)
#------------------------------------------------------
# crop = FALSE
suppressWarnings(
exp <- Start(dat = path,
var = 'tas',
sdate = '20000101',
ensemble = indices(1),
time = indices(1),
lat = indices(1:80), # 1:80 = -89.78488:-67.58778
lon = indices(19:65), # 19:65 = 5.00000:17.7777778
lat_reorder = Sort(),
lon_reorder = CircularSort(0, 360),
transform = CDORemapper,
transform_extra_cells = 8,
transform_params = list(grid = 'r100x50',
method = 'conservative',
crop = FALSE),
transform_vars = c('lat','lon'),
synonims = list(lon = c('lon', 'longitude'), lat = c('lat', 'latitude')),
return_vars = list(lat = 'dat', lon = 'dat', time = 'sdate'),
retrieve = F)
)
func <- function(x) {
return(x)
}
step <- Step(func, target_dims = 'time', output_dims = 'time')
wf <- AddStep(exp, step)
suppressWarnings(
res_crop_F_1 <- Compute(wf, chunks = list(lon = 2))
)
suppressWarnings(
res_crop_F_2 <- Compute(wf, chunks = list(ensemble = 1))
)
suppressWarnings(
res_crop_F_3 <- Compute(wf, chunks = list(lon = 3))
)
expect_equal(
res$output1,
res_crop_F_1$output1
)
expect_equal(
res$output1,
res_crop_F_2$output1
)
expect_equal(
res$output1,
res_crop_F_3$output1
)
#------------------------------------------------------
# crop = TRUE
suppressWarnings(
exp <- Start(dat = path,
var = 'tas',
sdate = '20000101',
ensemble = indices(1),
time = indices(1),
lat = indices(1:80), # 1:80 = -89.78488:-67.58778
lon = indices(19:65), # 19:65 = 5.00000:17.7777778
lat_reorder = Sort(),
lon_reorder = CircularSort(0, 360),
transform = CDORemapper,
transform_extra_cells = 8,
transform_params = list(grid = 'r100x50',
method = 'conservative',
crop = T),
transform_vars = c('lat','lon'),
synonims = list(lon = c('lon', 'longitude'), lat = c('lat', 'latitude')),
return_vars = list(lat = 'dat', lon = 'dat', time = 'sdate'),
retrieve = F)
)
func <- function(x) {
return(x)
}
step <- Step(func, target_dims = 'time', output_dims = 'time')
wf <- AddStep(exp, step)
suppressWarnings(
res_crop_T_1 <- Compute(wf, chunks = list(lat = 2, lon = 2))
)
suppressWarnings(
res_crop_T_2 <- Compute(wf, chunks = list(ensemble = 1))
)