-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathsource_me.R
141 lines (107 loc) · 2.93 KB
/
source_me.R
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
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
rm(list = ls())
# library path
.libPaths('C:\\Users\\mbeck\\R\\library')
# startup message
cat('Case the joint...\n')
# packages to use
library(knitr)
library(reshape2)
library(plyr)
library(ggplot2)
library(scales)
library(doParallel)
library(foreach)
library(Metrics)
library(GGally)
library(gridExtra)
library(ggmap)
library(data.table)
setwd('M:/docs/SWMP/detiding_cases/')
# functions to use
source('case_funs.r')
cases <- c('ELKVM', 'PDBBY', 'RKBMB', 'SAPDC')
# setup parallel backend
cl <- makeCluster(8)
registerDoParallel(cl)
# window widths grid
dy_wins <- c(1, 3, 6, 9, 12)
hr_wins <- c(1, 3, 6, 9, 12)
td_wins <- c(0.2, 0.4, 0.6, 0.8, 1)
case_grds <- expand.grid(dy_wins, hr_wins, td_wins)
names(case_grds) <- c('dec_time', 'hour', 'Tide')
save(case_grds, file = 'case_grds.RData')
load('case_grds.RData')
# setup parallel backend
cl <- makeCluster(8)
registerDoParallel(cl)
# start time
strt <- Sys.time()
# do w/ tide
for(case in cases){
to_proc <- prep_wtreg(case)
subs <- format(to_proc$DateTimeStamp, '%Y') %in% '2012'
to_proc <- to_proc[subs, ]
foreach(i = 1:nrow(case_grds)) %dopar% {
# progress
sink('log.txt')
cat('Log entry time', as.character(Sys.time()), '\n')
cat(case, '\n')
cat(i, ' of ', nrow(case_grds), '\n')
print(Sys.time() - strt)
sink()
load('case_grds.RData')
# create wt reg contour surface
wtreg <- wtreg_fun(to_proc, wins = c(case_grds[i,]),
parallel = F,
progress = F)
# save results for each window
wtreg_nm <-paste0(case, '_wtreg_', i)
assign(wtreg_nm, wtreg)
save(
list = wtreg_nm,
file=paste0('wtreg/', wtreg_nm, '.RData')
)
# clear RAM
rm(list = c(wtreg_nm))
}
}
stopCluster(cl)
#####
# get metab ests before and after detiding
# one element per site, contains both metab ests in the same data frame
# setup parallel backend
cl <- makeCluster(8)
registerDoParallel(cl)
# start time
strt <- Sys.time()
cases <- list.files(path = paste0(getwd(), '/wtreg/'), pattern = '_wtreg_')
# metab ests as list
met_ls <- foreach(case = cases) %dopar% {
# progress
sink('log.txt')
cat('Log entry time', as.character(Sys.time()), '\n')
cat(which(case == cases), ' of ', length(cases), '\n')
print(Sys.time() - strt)
sink()
# get data for eval
load(paste0('wtreg/', case))
nm <- gsub('.RData', '', case)
stat <- gsub('_wtreg_[0-9]+$', '', nm)
dat_in <- get(nm)
# get metab for obs DO
met_obs <- nem.fun(dat_in, stat = stat,
DO_var = 'DO_obs')
met_dtd <- nem.fun(dat_in, stat = stat,
DO_var = 'DO_nrm')
# combine results
col_sel <- c('Pg', 'Rt', 'NEM')
met_obs <- met_obs[, c('Station', 'Date', 'Tide', col_sel)]
met_dtd <- met_dtd[, col_sel]
names(met_dtd) <- c('Pg_dtd', 'Rt_dtd', 'NEM_dtd')
met_out <- cbind(met_obs, met_dtd)
# return results
met_out
}
stopCluster(cl)
names(met_ls) <- cases
save(met_ls, file = 'met_ls.RData')