forked from ARCCSS-extremes/climpact2
-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathclimpact2.batch.stations.r
159 lines (141 loc) · 6.89 KB
/
climpact2.batch.stations.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
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
# ------------------------------------------------
# This file contains functionality to batch process station files stored in "RClimdex format" (i.e. 6 column format, see climpact2.sampledata.1d.time-series.txt for an example)
# University of New South Wales
# nherold, May 2016.
# ------------------------------------------------
#
# CALLING THIS FILE:
# Rscript climpact2.batch.stations.r /full/path/to/station/files/ /full/path/to/metadata.txt base_period_begin base_period_end cores_to_use
# e.g. Rscript climpact2.batch.stations.r ./sample_data/Renovados_hasta_2010 ./sample_data/climpact2.sample.batch.metadata.txt 1971 2000 4
#
# NOTE: This file expects that all of your individual station files are kept in /directory/path/to/station/files/ and that each file name (excluding path directory) is found in
# column one of metadata.txt, with corresponding metadata in subsequent columns.
#
# COMMAND LINE FLAGS:
# - /full/path/to/station/files/ : directory path to where files listed in column one of metadata.txt are kept.
# - /full/path/to/metadata.txt : text file containing 12 columns; station file, latitude, longitude, wsdin, csdin, Tb_HDD, Tb_CDD, Tb_GDD, rx_ud, rnnmm_ud, txtn_ud, SPEI
# - base_period_begin : beginning year of base period
# - base_period_end : end year of base period
# - cores_to_use : number of cores to use in parallel
# load and source and specify cores
library(foreach)
library(climdex.pcic)
library(doParallel)
source("climpact2.GUI.r")
source("ancillary/climpact2.etsci-functions.r")
# return a nice list of station metadata
read.file.list.metadata <- function(file.list.metadata)
{
file.list.metadata <- read.table(file.list.metadata,header=T,col.names=c("station","latitude","longitude","wsdin","csdin","Tb_HDD","Tb_CDD","Tb_GDD","rxnday","rnnmm","txtn","SPEI"),
colClasses=c("character","real","real","integer","integer","real","real","integer","real","real","integer"))
return(file.list.metadata)
}
# call QC and index calculation functionality for each file specified in metadata.txt
batch <- function(input.directory,output.directory,base.start,base.end) {
metadata <- read.file.list.metadata(file.list.metadata)
# foreach does not support 'next'. This code is removed from the dopar loop in order to provide 'next' like functionality, in the form of return(NA) calls.
func = function(file.number) {
file=paste(input.directory,"/",metadata$station[file.number],sep="")
print(file)
user.data <- read.user.file(file,graphics=FALSE)
user.data <- check.and.create.dates(user.data)
get.file.path(file)
create.dir(file)
# define variables for indices
lat <- as.numeric(metadata$latitude[file.number])
lon <- as.numeric(metadata$longitude[file.number])
wsdi_ud <<- metadata$wsdin[file.number]
csdi_ud <<- metadata$csdin[file.number]
Tb_HDD <<- metadata$Tb_HDD[file.number]
Tb_CDD <<- metadata$Tb_CDD[file.number]
Tb_GDD <<- metadata$Tb_GDD[file.number]
rx_ud <<- metadata$rxnday[file.number]
rnnmm_ud <<- metadata$rnnmm[file.number]
txtn_ud <<- metadata$txtn[file.number]
custom_SPEI <<- metadata$SPEI[file.number]
# global variables needed for calling climpact2.GUI.r functionality
station.metadata <- create.metadata(lat,lon,base.start,base.end,user.data$dates,"ofile_filler")
assign("metadata",station.metadata,envir=.GlobalEnv)
version.climpact <<- software_id
quantiles <<- NULL
if(lat<0) lat_text = "°S" else lat_text = "°N"
if(lon<0) lon_text = "°W" else lon_text = "°E"
Encoding(lon_text) <- "UTF-8" # to ensure proper plotting of degree symbol in Windows (which uses Latin encoding by default)
Encoding(lat_text) <- "UTF-8"
title.station <- paste(ofilename, " [", lat,lat_text, ", ", lon,lon_text, "]", sep = "")
assign("title.station", title.station, envir = .GlobalEnv)
plot.title<-gsub('\\#',title.station,"Station: #"); assign('plot.title',plot.title,envir=.GlobalEnv)
barplot_flag <<- TRUE
min_trend <<- 10
temp.quantiles <<- c(0.05,0.1,0.5,0.9,0.95)
prec.quantiles <<- c(0.05,0.1,0.5,0.9,0.95,0.99)
op.choice <<- NULL
skip <<- FALSE
if(file_test("-f",paste(file,".error.txt",sep=""))) { file.remove(paste(file,".error.txt",sep="")) }
# run quality control and create climdex input object
catch1 <- tryCatch(QC.wrapper(station.metadata,user.data,file,graphics=FALSE),
error=function(msg) {
fileConn<-file(paste(file,".error.txt",sep=""))
writeLines(toString(msg), fileConn)
close(fileConn)
if(file_test("-f",paste0(file,".temporary"))) { file.remove(paste0(file,".temporary")) }
})
if(skip) { return(NA) }
# calculate indices
catch2 <- tryCatch(index.calc(station.metadata,graphics=FALSE),
error=function(msg) {
fileConn<-file(paste(file,".error.txt",sep=""))
writeLines(toString(msg), fileConn)
close(fileConn)
if(file_test("-f",paste0(file,".temporary"))) { file.remove(paste0(file,".temporary")) }
})
if(skip) { return(NA) }
# RJHD - NH addition for pdf error 2-aug-17
graphics.off()
print(paste(file," done",sep=""))
rm(list=ls())
}
# capture.output(stdout <- foreach(file.number=1:length(metadata$station)) %dopar%
# for (file.number in 1:length(metadata$station))
foreach(file.number=1:length(metadata$station)) %dopar%
{
func(file.number)
}
print("",quote=FALSE)
print("",quote=FALSE)
print("",quote=FALSE)
print("",quote=FALSE)
print("",quote=FALSE)
print("",quote=FALSE)
print("*********************************************************************************************",quote=FALSE)
print("*********************************************************************************************",quote=FALSE)
print("*********************************************************************************************",quote=FALSE)
print("PROCESSING COMPLETE.",quote=FALSE)
print("",quote=FALSE)
print("",quote=FALSE)
print("",quote=FALSE)
print("Any errors encountered during processing are listed below by input file. Assess these files carefully and correct any errors.",quote=FALSE)
print("",quote=FALSE)
error.files <- suppressWarnings(list.files(path=input.directory,pattern=paste("*error.txt",sep="")))
if(length(error.files)==0) { print("... no errors detected in processing your files. That doesn't mean there aren't any!",quote=FALSE) }
else {
for (i in 1:length(error.files)) { #system(paste("ls ",input.directory,"*error.txt | wc -l",sep=""))) {
print(error.files[i],quote=FALSE)
#system(paste("cat ",input.directory,"/",error.files[i],sep=""))
}
}
}
# Enable reading of command line arguments
args<-commandArgs(TRUE)
# where one or more station files are kept
input.directory = toString(args[1])
# metadata text file
file.list.metadata = toString(args[2])
# begin base period
base.start = as.numeric(args[3])
# end base period
base.end = as.numeric(args[4])
# set up multiple cores
registerDoParallel(cores=as.numeric(args[5]))
# call main function
batch(input.directory,output.directory,base.start,base.end)