r - avgokmts returning incorrect maximum rain from ok mesonet data -
i'm using okmesonet
package data on rainfall. i've tried using avgokmts
package calculate rainfall each day, i'm getting non-sensical values.
get rain data norman, ok (cumulative rain in mm on day @ 5 min intervals)
library(okmesonet) raindat <- okmts(begintime="2016-06-21 00:00:00", endtime="2016-07-04 00:00:00", station="nrmn", variables="rain", localtime=true)
calculate max rain per day
avgokmts(raindat, by="day", metric="max")
which returns these values
stid stnm day month year rain time date 1 nrmn 121 21 06 2016 0.00 23:55:00 2016-06-22 2 nrmn 121 22 06 2016 0.25 23:55:00 2016-06-23 3 nrmn 121 23 06 2016 59.70 23:55:00 2016-06-24 4 nrmn 121 24 06 2016 0.00 23:55:00 2016-06-25 5 nrmn 121 25 06 2016 0.00 23:55:00 2016-06-26 6 nrmn 121 26 06 2016 0.00 23:55:00 2016-06-27 7 nrmn 121 27 06 2016 0.00 23:55:00 2016-06-28 8 nrmn 121 28 06 2016 0.00 23:55:00 2016-06-29 9 nrmn 121 29 06 2016 0.00 23:55:00 2016-06-30 10 nrmn 121 30 06 2016 28.19 23:55:00 2016-07-01 11 nrmn 121 01 07 2016 0.00 23:55:00 2016-07-02 12 nrmn 121 02 07 2016 0.51 23:55:00 2016-07-03 13 nrmn 121 03 07 2016 0.00 23:55:00 2016-07-04 14 nrmn 121 04 07 2016 0.00 00:00:00 2016-07-04
but these rainfall values don't match rainfall graphed below (peak rainfall occurs on june 27th , july 3rd).
plot(raindat$time, raindat$rain, xlab="date", ylab="cumulative daily rain (mm)")
why isn't avgokmts
working in case? there error in how i'm calling function? there alternative way calculate daily rainfall using dataset?
i'm pretty sure pkg author did not deal utc<->cdt conversions precip readings. here's fragile way max precip per day if using single station. expansion of procedure handle multiple stations should adding 1 more group_by()
variable.
library(okmesonet) library(dplyr) library(ggplot2) library(gridextra) raindat <- okmts(begintime="2016-06-21 00:00:00", endtime="2016-07-04 00:00:00", station="nrmn", variables="rain", localtime=true) # use pkg calculation ------------------------------------------------- pkg_calc <- avgokmts(raindat, by="day", metric="max") # begin our own calculations ---------------------------------------------- raindat <- mutate(raindat, day=format(time, "%y-%m-%d")) day_precip_max <- function(x) { prev_day_last_reading_time <- as.posixct(sprintf("%s 23:55:00", x$day[1]), tz="america/chicago") - as.difftime(1, unit="days") prev_day_last_reading <- raindat[raindat$time==prev_day_last_reading_time, "rain"] if (length(prev_day_last_reading) == 0) prev_day_last_reading <- 0 x <- mutate(x, rain=rain - prev_day_last_reading) data_frame( stid=x$stid[1], stnm=x$stnm[1], day=substr(x$day[1], 9, 10), month=substr(x$day[1], 6, 7), year=substr(x$day[1], 1, 4), rain=max(x$rain) ) } new_calc <- group_by(raindat, day) %>% do(day_precip_max(.)) %>% ungroup() # convert posixct common plotting axis ------------------------------ pkg_calc <- mutate(pkg_calc, day=as.posixct(sprintf("%s-%s-%s 23:55:00", year, month, day), tz="america/chicago")) new_calc <- mutate(new_calc, day=as.posixct(sprintf("%s-%s-%s 23:55:00", year, month, day), tz="america/chicago")) grid.arrange( ggplot(raindat, aes(x=time, y=rain)) + geom_point() + scale_x_datetime(date_breaks="1 day", date_labels="%d") + labs(x=null, y="rain", title="raw readings") , ggplot(pkg_calc, aes(x=day, y=rain)) + geom_point() + scale_x_datetime(date_breaks="1 day", date_labels="%d", limits=range(raindat$time)) + labs(x=null, y="rain", title="package aggregation (max)") , ggplot(new_calc, aes(x=day, y=rain)) + geom_point() + scale_x_datetime(date_breaks="1 day", date_labels="%d", limits=range(raindat$time)) + labs(x=null, y="rain", title="manual aggregation (max)") , ncol=1 )
i have plot displaying max reading @ 23:55:00.
Comments
Post a Comment