Tortured R code to convert irregular timeseries of tipping bucket rain gage data into hourly rainfall accumulations.

tipping_bucket_to_hourly.r

rm(list=ls())

library(xts)

IN_TO_M = 0.0254

setwd('path/to/my/tipping/bucket/data')
gage1 = read.table("gage1_tipping.txt", sep="|", header=T)
gage2 = read.table("gage2_tipping.txt", sep="|", header=T)
data = list(gage1, gage2)

filenames = c("gage1_hourly.txt", "gage2_hourly.txt")
for (filename in filenames) {
	unlink(filename)
}

i = 0
for (datum in data) {
	i = i + 1
	outfile = file(filenames[i], "w")
	writeLines(paste("datetime", "rainfall_m", sep=","), outfile)
	
	time = as.POSIXct(datum$TimeStamp, format="%Y-%m-%d %H:%M:%S")
	rain = as.matrix(datum$Tips..in.)
	timeseries = xts(rain, order.by=time)
	timeseries_hourly = period.apply(timeseries, endpoints(timeseries, "hours"), sum)
	## Generate empty regular time series for water years of the entire record
	start = trunc(start(timeseries_hourly), "hour")
	fin = trunc(end(timeseries_hourly), "hour")
	unlist(unclass(start))
	unlist(unclass(fin))
	# Make sure we get the correct start and end of water years
	if (start$mon < 9) { # 0-indexed, so 9 = October
		year = start$year + 1900 - 1
	} else {
		year = start$year + 1900
	}
	if(fin$mon > 8) { # 0-indexed, so 8 = Sept.
		endYear = fin$year + 1900 + 1
	} else {
		endYear = fin$year + 1900
	}
	begin = as.POSIXct(sprintf("%d-%s-%s %s:%s", year, '10', '01', '00', '00'), format="%Y-%m-%d %H:%M")
	end = as.POSIXct(sprintf("%d-%s-%s %s:%s", endYear, '09', '30', '23', '00'), format="%Y-%m-%d %H:%M")
	interval = as.difftime(1, units="hours")
	regular_timeseries = seq(begin, end, by=interval)
	regular_timeseries_hourly = xts(rep(0, length(regular_timeseries)), order.by=regular_timeseries)
	## Merge irregular and regular timeseries
	# Drop minutes and seconds from input data
	index(timeseries_hourly) = as.POSIXct( gsub("(\\s+\\d{1,2}):.*","\\1:00", index(timeseries_hourly)), format="%Y-%m-%d %H:%M" )
	merged = merge(regular_timeseries_hourly, timeseries_hourly)
	merged[is.na(merged)] = 0
	final = merged[,2]
	# Convert from inches to m
	final$timeseries_hourly = final$timeseries_hourly * IN_TO_M
	## Output
	write.zoo(final, outfile, sep=",", col.names=F)
	
	close(outfile)
}

Input data are of format:

TimeStamp|Tips (in)|Quality code|interpolation|Remarks
2009-04-11 03:56:58|0.00|200|1| 
2009-04-11 03:58:38|0.01|200|1| 
2009-04-11 04:10:24|0.01|200|1|

Where Tips (in) is a tip equal to 0.01 inches.

Output data are of format:

datetime,rainfall_m
2009-10-07 00:00:00,0
2009-10-07 01:00:00,0
2009-10-07 02:00:00,0
2009-10-07 03:00:00,0
2009-10-07 04:00:00,0.000508

Where rainfall_m is rainfall accumulation for the hour, in meters.

There is probably a better way to do this (I'm not really an R programmer/user), but it works for me.