# invoke like this:
#
#  R --no-save --silent < cimis_graphs.R
#
# add to crontab:
# http://www.debian-administration.org/articles/56
#

# establish URLs
sfrec <- url('ftp://ftpcimis.water.ca.gov//pub/hourlyMetric/hlymet084.csv')
madera <- url('ftp://ftpcimis.water.ca.gov//pub/hourlyMetric/hlymet080.csv')
davis <- url('ftp://ftpcimis.water.ca.gov//pub/hourlyMetric/hlymet006.csv')

# read in relevent columns
x.sfrec <- read.csv(sfrec, header=FALSE, na.strings='       --')[, c(2,3,4,6,8,10,12,14,16,18,20,22,24)]
x.madera <- read.csv(madera, header=FALSE, na.strings='       --')[, c(2,3,4,6,8,10,12,14,16,18,20,22,24)]
x.davis <- read.csv(davis, header=FALSE, na.strings='       --')[, c(2,3,4,6,8,10,12,14,16,18,20,22,24)]


source('cimis_functions.R')

cimis_plot(x.sfrec, plot.title='SFREC', filename='sfrec-cimis.png')
cimis_plot(x.madera, plot.title='Station Fresno/Madera', filename='fno-madera-cimis.png')
cimis_plot(x.davis, plot.title='Davis', filename='davis-cimis.png')

cimis_plot <- function(cimis.data, plot.title='', filename=plot.title)
{

# we need a cairo device for this to work
require(Cairo)

# fix column names
names(cimis.data) <- c('date', 'time', 'julian_day','eto', 'precip', 'solar_rad', 'vapor_press', 'air_temp', 'rel_humidity', 'dew_pt', 'wind_speed', 'wind_dir', 'soil_temp')

# re-make date/time column
# convert time to better format by adding a leading 0 for hours < 10:00
# reduce time to hour 1300 -> 13 and slide back to 0-23 style notation
cimis.data$hour <- as.numeric(sprintf("%02i", (cimis.data$time/100) - 1))

# convert to datetime
cimis.data$datetime <- as.POSIXct(strptime(paste(cimis.data$date, cimis.data$hour), format='%m/%d/%Y %H'))

# get datetime range and make a list of tick marks
d.range <- range(cimis.data$datetime, na.rm=TRUE)
d.list <- seq(d.range[1], d.range[2], by='6 hours')

# get a range in eto and preci-- so that the graph can be scaled accordingly
# fetch the upper limit of rain
if(range(cimis.data$precip, na.rm=TRUE)[2] > range(cimis.data$eto, na.rm=TRUE)[2])
        ymax <- max(cimis.data$precip, na.rm=TRUE)
else
        ymax <- max(cimis.data$eto, na.rm=TRUE)


# start output file
Cairo(file=filename, type='png', bg='white', width=600, height=350)
# setup margins etc.
par(mar = c(0.5, 4, 0, 0.5), oma = c(4, 0, 1.5, 0), mfcol = c(2,1))
plot(eto ~ datetime, data=cimis.data, type='l', ylim=c(0, ymax), ylab='ETo / Rainfall (mm)', col='gray', xaxt='n', las=2, cex.axis=0.75)
lines(precip ~ datetime, data=cimis.data, type='h', col='blue')

plot(air_temp ~ datetime, data=cimis.data, type='l', ylab='Temp (C)', col='RoyalBlue', xaxt='n', las=2, cex.axis=0.75)
lines(soil_temp ~ datetime, data=cimis.data, type='l', col='orange')

# add shared axis
axis.POSIXct(at=d.list, side=1, format="%b-%d %H", cex.axis=0.75, las=2)

# add a title
mtext(plot.title, outer=TRUE, line=0.5, font=2)

# close file
dev.off()
}