2016-07-04 37 views
1

我正在使用okmesonet包獲取有關降雨量的數據。我試過使用這個軟件包中的avgokmts來計算每天的降雨量,但我得到的是非感性的值。avgokmts從ok mesonet數據返回不正確的最大雨量

獲取降雨數據諾曼,OK(以5分鐘間隔毫米累計降雨超過一天)

library(okmesonet) 
rainDat <- okmts(begintime="2016-06-21 00:00:00", endtime="2016-07-04 00:00:00", 
      station="NRMN", variables="RAIN", localtime=TRUE) 

計算每它返回這些值一天

avgokmts(rainDat, by="day", metric="max") 

最大雨

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 

但是,這些降雨量值非常明顯地與降雨量不匹配,如下圖所示(高峯期出現在6月27日和7月3日)。

plot(rainDat$TIME, rainDat$RAIN, xlab="Date", ylab="Cumulative Daily Rain (mm)") 

cumulative rain fall plot

爲什麼不avgokmts在這種情況下工作?我是如何調用函數的錯誤?有沒有其他方法可以使用這個數據集來計算每日降雨量?

回答

0

我敢肯定,pkg作者沒有正確處理UTC < - > CDT轉換正確的降水量讀數。如果您使用的是單站,那麼每天都可以獲得最大降水量。處理多個電臺的程序的擴展應該是增加一個group_by()變量。

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 the 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 to POSIXct for 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 
) 

enter image description here

我有顯示最大讀數23:55:00劇情。