Newer
Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
#'CST_Analogs
#'
#'@author Carmen Alvarez-Castro, \email{carmen.alvarez-castro@cmcc.it}
#' adapted version of the method of Yiou et al 2013
#'
#'@references Yiou, P., T. Salameh, P. Drobinski, L. Menut, R. Vautard, and M. Vrac, 2013 :
#' Ensemble reconstruction of the atmospheric column from surface pressure using analogues.
#' Clim. Dyn., 41, 1419-1437. \email{pascal.yiou@lsce.ipsl.fr}
#'
#'@description search for days with similar atmospheric conditions based on the large scale slp (or geopotential height)
# and the local scale (precipitation or Temperature)
#'
#'@param month month of the analog day
#' day day of the analog day
#' yr year of the analog day
#' mAnalog month or list of months to search for analogs
#'@import
#'
#'@return list best.corr.time.ana,best.dist.time.ana,selec.dist.time,selec.corr.time
#' list file.dat with a list of days with the format yyyymmdd ordered by best analog with the dist (minima) and corr (maxima)
#' plot preliminary plot of the best analog selected
#' yr1 first year of the total period of study
#' yr2 last year of the total period of study
#' ical number of days per year in the calendar (360,365,366)
#'@example
# Analogs <- function(exp, obs) {
# if (!all(c('member', 'sdate') %in% names(dim(exp)))) {
# stop("Parameter 'exp' must have the dimensions 'member' and 'sdate'.")
# }
#
# if (!all(c('sdate') %in% names(dim(obs)))) {
# stop("Parameter 'obs' must have the dimension 'sdate'.")
# }
#
# if (any(is.na(exp))) {
# warning("Parameter 'exp' contains NA values.")
# }
#
# if (any(is.na(obs))) {
# warning("Parameter 'obs' contains NA values.")
# }
#
# target_dims_obs <- 'sdate'
# if ('member' %in% names(dim(obs))) {
# target_dims_obs <- c('member', target_dims_obs)
# }
#
# Analogs <- Apply(data = list(var_obs = obs, var_exp = exp),
# target_dims = list(target_dims_obs, c('member', 'sdate')),
# fun = .select)$output1
#
# return(Analogs)
# }
time_obsL <- as.Date(c("2005-01-01", "2005-02-01", "2005-03-01",
"2005-04-01", "2005-05-01"))
Analogs <- function(expL, obsL, time_obsL, expVar = NULL, obsVar = NULL,
criteria = "Large_dist",
lon_local = NULL, lat_local = NULL, region = NULL,
nAnalogs = 1, return_list = FALSE) {
# checks
metric <- Select(expL = expL, obsL = obsL, expVar = expVar, obsVar = obsVar,
criteria = criteria, lon_local = lon_local, lat_local = lat_local,
region = region)
best <- Apply(list(metric), target_dims = 'time', fun = BestAnalog,
criteria = criteria, return = return_list)
}
#'@example
#'expL <- 1 : (8*10*2*6*7)
dim(expL) <- c(lat = 8, lon = 10, time = 2, member = 6, sdate = 7)
obsL <- 1 : (8*10*5)
dim(obsL) <- c(lat = 8, lon = 10, time = 5)
#'lat_local <- lat <- seq(0, 19, 2.5)
lon_local <- lon <- seq(0, 23, 2.5)
met = Select(expL, obsL, criteria = "Local_dist", lon_local = lon,
lat_local = lat,
region = c(lonmin = 0, lonmax = 5, latmin = 0, latmax = 5))
#'dim(met$metric1)
#'str(met)
expL <- 1 : c(8*10*2)
dim(expL) <- c(lat = 8, lon = 10, time = 2)
obsL <- 1 : (8*10*5)
dim(obsL) <- c(lat = 8, lon = 10, time = 5)
met = Select(expL, obsL, criteria = "Local_dist", lon_local = lon,
lat_local = lat,
region = c(lonmin = 0, lonmax = 5, latmin = 0, latmax = 5))
#'pos <- BestAnalog(met, return_list = TRUE, nAnalogs = 2)
#'pos
res <- Apply(met, target_dims = 'time', fun = BestAnalog, criteria = 'Local_dist')$output1
BestAnalog <- function(metric, criteria = 'Large_dist', return_list = FALSE,
nAnalogs = 1) {
print(class(metric))
print(length(metric))
if (criteria == 'Large_dist') {
pos1 <- metric$pos1
if (return_list == FALSE) {
pos <- pos1[1]
} else {
pos <- pos1[1 : nAnalogs]
}
} else if (criteria== 'Local_dist') {
# pos1 <- c(7, 13, 5, 3, 6, 12, 10, 1, 8, 9, 11, 4, 2, 14)
# pos2 <- c(4, 8, 13, 6, 3, 1, 12, 5, 9, 7, 10, 2, 11, 14)
pos1 <- metric$pos1[1 : nAnalogs]
pos2 <- metric$pos2[1 : nAnalogs]
best <- match(pos1, pos2)
pos <- pos1[as.logical(best)]
pos <- pos[which(!is.na(pos))]
if (return_list == FALSE) {
pos <- pos[1]
}
} else if (criteria == 'Local_cor') {
pos1 <- metric$pos1[1 : nAnalogs]
pos2 <- metric$pos2[1 : nAnalogs]
best <- match(pos1, pos2)
pos <- pos1[as.logical(best)]
pos <- pos[which(!is.na(pos))]
# pos3 <- c(6, 11, 14, 3, 13, 7, 2, 5, 1, 12, 10, 9, 8, 4)
best <- match(pos, pos3)
pos <- pos[order(best, decreasing = F)]
pos <- pos[which(!is.na(pos))]
if (return_list == FALSE) {
pos[1]
}
return(pos)
}
expL <- (1 + 2): (4 * 3 * 2 + 2)
dim(expL) <- c(lat = 4, lon = 3, time = 2)
obsL <- 1 : c(4 * 3 * 5)
dim(obsL) <- c(lat = 4, lon = 3, time = 5)
res = Select(expL, obsL)
expL <- (1 + 2): (8 * 10 * 2 + 2)
dim(expL) <- c(lat = 8, lon = 10, time = 2)
obsL <- 1 : c(8 * 10 * 5)
dim(obsL) <- c(lat = 8, lon = 10, time = 5)
lat_local <- lat <- seq(0, 19, 2.5)
lon_local <- lon <- seq(0, 23, 2.5)
res = Select(expL, obsL, criteria = "Local_dist", lon_local = lon,
lat_local = lat,
region = c(lonmin = 0, lonmax = 5, latmin = 0, latmax = 5))
# probar mas ejemplos con diferentes criterios, latitudes, longitudes
Select <- function(expL, obsL, expVar = NULL, obsVar = NULL, criteria = "Large_dist",
lon_local = NULL, lat_local = NULL, region = NULL) {
#check expL
#check obsL
#check obsVar
if (any(names(dim(expL)) == 'time')) {
names(dim(expL))[which(names(dim(expL)) == 'time')] <- 'time_exp'
}
metric1 <- Apply(list(obsL), target_dims = list(c('lat', 'lon')),
fun = .select, expL, metric = "dist")$output1
dim_time_obs <- which(names(dim(metric1)) == 'time')
margins <- c(1 : length(dim(metric1)))[-dim_time_obs]
pos1 <- apply(metric1, margins, order)
names(dim(pos1))[1] <- 'time'
metric1 <- apply(metric1, margins, sort)
names(dim(metric1))[1] <- 'time'
if (criteria == "Large_dist") {
return(list(metric1 = metric1, pos1 = pos1))
}
if (criteria == "Local_dist" | criteria == "Local_cor") {
obs <- SelBox(obsL, lon = lon_local, lat = lat_local, region = region)$data
exp <- SelBox(expL, lon = lon_local, lat = lat_local, region = region)$data
metric2 <- Apply(list(obs), target_dims = list(c('lat', 'lon')),
fun = .select, exp, metric = "dist")$output1
margins <- c(1 : length(dim(metric2)))[-dim_time_obs]
pos2 <- apply(metric2, margins, order)
names(dim(pos2))[1] <- 'time'
metric2 <- apply(metric2, margins, sort)
names(dim(metric2))[1] <- 'time'
if (criteria == "Local_dist") {
return(list(metric1 = metric1, metric2 = metric2,
pos1 = pos1, pos2 = pos2))
}
}
if (criteria == "Local_cor") {
obs <- SelBox(obsVar, lon = lon_local, lat = lat_local, region = region)$data
exp <- SelBox(expVar, lon = lon_local, lat = lat_local, region = region)$data
metric3 <- Apply(list(obs), target_dims = list(c('lat', 'lon')),
fun = .select, exp, metric = "cor")$output1
margins <- c(1 : length(dim(metric3)))[-dim_time_obs]
pos3 <- apply(metric3, margins, order, decreasing = TRUE)
names(dim(pos3))[1] <- 'time'
metric3 <- apply(metric3, margins, sort)
names(dim(metric3))[1] <- 'time'
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
return(list(metric1 = metric1, metric2 = metric2, metric3 = metric3,
pos1 = pos1, pos2 = pos2, pos3 = pos3))
}
else {
stop("Parameter 'criteria' must to be one of the: 'Large_dist', ",
"'Local_dist','Local_cor'.")
}
}
# data <- 1:(20 * 3 * 2 * 4)
# dim(data) <- c(lon = 20, lat = 3, time = 2, model = 4)
# lon <- seq(2, 40, 2)
# lat <- c(1, 5, 10)
# a <- SelBox(data = data, lon = lon, lat = lat, region = c(2, 20, 1, 5),
# londim = 1, latdim = 2, mask = NULL)
# str(a)
#'@example
exp <- (1 + 2): (4 * 3 + 2)
dim(exp) <- c(lat = 4, lon = 3)
obs <- 1 : c(5 * 4 * 3)
dim(obs) <- c(time = 5, lat = 4, lon = 3)
res <- .select(exp, obs)
res
res <- .select(exp, obs, metric = 'cor')
dim(res)
.select <- function(exp, obs, metric = "dist") {
if (metric == "dist") {
#metric <- sum((obs - exp) ^ 2)
#metric <- apply(obs, "time", function(x) {sum((x - exp) ^ 2)})
result <- Apply(list(obs), target_dims = list(c('lat', 'lon')),
fun = function(x) {sum((x - exp) ^ 2)})$output1
} else if (metric == "cor") {
result <- Apply(list(obs), target_dims = list(c('lat', 'lon')),
fun = function(x) {cor(as.vector(x), as.vector(exp))})$output1
}
result
}
# Add '_exp' label to experiment dimmension if they are not lat and lon but in common with obs:
expL <- 1 : (8*10*2*6*7)
dim(expL) <- c(lat = 8, lon = 10, time = 2, member = 6, sdate = 7)
obsL <- 1 : (8*10*5*3)
dim(obsL) <- c(lat = 8, lon = 10, time = 5, member = 3)
dimnames_exp <- names(dim(expL))
#dimnames_exp <- dimnames_exp[-which(dimnames_exp == 'lat' | dimnames_exp == 'lon')]
dimnames_obs <- names(dim(obsL))
#dimnames_obs <- dimnames_obs[-which(dimnames_obs == 'lat' | dimnames_obs == 'lon')]
which(dimnames_exp == dimnames_obs & dimnames_exp)
dimnames_exp[which(dimnames_exp == dimnames_obs)] <- paste0(dimnames_exp[which(dimnames_exp == dimnames_obs)], "_exp")
names(dim(expL)) <- dimnames_exp
repeat_names <- function(names_exp, names_obs) {
latlon_dim_exp <- which(names_exp == 'lat' | names_exp == 'lon')
latlon_dim_obs <- which(names_obs == 'lat' | names_obs == 'lon')
dimnames_obs <- dimnames_obs[-which(dimnames_obs == 'lat' | dimnames_obs == 'lon')]
}