A function to read daily rainfall data from the GHCN data set available at http://www1.ncdc.noaa.gov/pub/data/ghcn/daily/all/

READ_GHCN_dly.R

READ_GHCN_dly <- function(filename){
  require(data.table)
  require(readr)
  require(magrittr)
  require(stringr)
  require(lubridate)
  # filename implies that the location of the file is already known
  # if not, use the station ID to get the filename:
  # filename <- paste0('http://www1.ncdc.noaa.gov/pub/data/ghcn/daily/all/', stn, '.dly')
  
  # file positions from http://www1.ncdc.noaa.gov/pub/data/ghcn/daily/readme.txt
  start_pos <- c(1, 12, 16, 18)
  end_pos <- c(11, 15, 17, 21)
  raw <- readr::read_fwf(filename,
                  readr::fwf_widths(widths = c(11, 4, 2, 4, rep(c(5, 1, 1, 1), 31))),
                  na = '-999',
                  col_types = paste(c('c', 'i', 'i', 'c', rep(c('i', 'c', 'c', 'c'), 31)), collapse = '')
  ) %>% data.table()

  # rename
  NOAA_names <- paste0('X', 1:128)
  new_names <- c('ID', 'YEAR', 'MONTH', 'ELEMENT', paste0(rep(c('VALUE', 'MFLAG', 'QFLAG', 'SFLAG'), 31), floor(1:(31 * 4 - 1) / 4) + 1))
  setnames(raw, NOAA_names, new_names)

  # only get rain
  raw <- raw[ELEMENT == 'PRCP']

  # melt to get tidy data
  rain <- data.table::melt(raw, id=c('ID', 'YEAR', 'MONTH', 'ELEMENT'), measure.vars = paste0('VALUE', 1:31))
  rain[value == -9999, value := NA]

  # get the date -- currently hidden
  # will get an error message because some months don't have 31 days -- OK just delete them!
  rain[, DAY := stringr::str_sub(variable, 6)]
  rain[, date := lubridate::ymd(paste(YEAR, MONTH, DAY, sep = '-'), quiet = T)]
  rain <- rain[!is.na(date)]

  # set final names to my liking
  out <- rain[, .(stnid = ID, date, prcp_mm = value / 10)]
  setkey(out, date)
  return(out)
}