Skip to content

Commit

Permalink
version 0.7.5
Browse files Browse the repository at this point in the history
  • Loading branch information
davidcarslaw authored and cran-robot committed Jan 26, 2017
0 parents commit 672c3d0
Show file tree
Hide file tree
Showing 19 changed files with 1,326 additions and 0 deletions.
26 changes: 26 additions & 0 deletions DESCRIPTION
@@ -0,0 +1,26 @@
Package: worldmet
Type: Package
Title: Import Surface Meteorological Data from NOAA Integrated Surface
Database (ISD)
Version: 0.7.5
Date: 2017-01-25
Authors@R: c(person("David", "Carslaw", role = c("aut", "cre"), email =
"david.carslaw@york.ac.uk"))
ByteCompile: true
Depends: R (>= 3.2.0)
Imports: openair, plyr, dplyr, leaflet, RCurl, readr, zoo
Maintainer: David Carslaw <david.carslaw@york.ac.uk>
Description: Functions to import data from more than 30,000 surface
meteorological sites around the world managed by the National Oceanic and Atmospheric Administration (NOAA) Integrated Surface
Database (ISD, see <https://www.ncdc.noaa.gov/isd>).
License: GPL (>= 2)
URL: http://github.com/davidcarslaw/worldmet
BugReports: http://github.com/davidcarslaw/worldmet/issues
LazyLoad: true
LazyData: true
RoxygenNote: 5.0.1
NeedsCompilation: no
Packaged: 2017-01-25 10:36:19 UTC; DC23
Author: David Carslaw [aut, cre]
Repository: CRAN
Date/Publication: 2017-01-26 12:08:20
18 changes: 18 additions & 0 deletions MD5
@@ -0,0 +1,18 @@
ee56e89a45cdd75fe6663dcd66a83a20 *DESCRIPTION
41b51e8724b02957fc9006e05727029c *NAMESPACE
c907682c4ac72c8f0ae11146430a592a *R/exportADMS.R
91c5e2e7955c451fbacfdfa2a7439320 *R/getMeta.R
24d05868485997988bf40b7dee359ae2 *R/metNOAA.R
14fcd05358778bdc3207bed2db1ddc83 *R/meta.R
821f31580226b09a684f3c8349b19b9d *R/sysdata.rda
cf4532f8d10db6d94136193ccf0f1cdb *R/weatherCodes.R
1ff037620039bef68eeeab3e4a5e68d1 *R/worldmet-package.R
5444a97c8b920565ebe0544f3eb43170 *README.md
821f31580226b09a684f3c8349b19b9d *data/meta.rda
8699a0e8e38705634dd0ae8f861ea0c4 *data/weatherCodes.rda
265b7bf2f6e2263c74b7cec6d621e43a *man/exportADMS.Rd
21c14d6ea101c06126bf0a9c2a96fcc1 *man/getMeta.Rd
c99887a99c9c765940df377b4c136a18 *man/importNOAA.Rd
79517c932c8b1d331c648f7db528a02a *man/meta.Rd
a741b9fa5a114c5f5e6845a507639014 *man/weatherCodes.Rd
303a9d0e1a2cf5672f113a09b37b9d81 *man/worldmet.Rd
18 changes: 18 additions & 0 deletions NAMESPACE
@@ -0,0 +1,18 @@
# Generated by roxygen2: do not edit by hand

export(exportADMS)
export(getMeta)
export(importNOAA)
import(RCurl)
import(openair)
import(plyr)
import(readr)
importFrom(dplyr,"%>%")
importFrom(leaflet,addCircles)
importFrom(leaflet,addMarkers)
importFrom(leaflet,addTiles)
importFrom(leaflet,leaflet)
importFrom(utils,head)
importFrom(utils,read.csv)
importFrom(utils,write.table)
importFrom(zoo,na.approx)
109 changes: 109 additions & 0 deletions R/exportADMS.R
@@ -0,0 +1,109 @@
#' Export a meteorological data frame in ADMS format
#'
#' @param dat A data frame imported by \code{\link{importNOAA}}.
#' @param out A file name for the ADMS file. The file is written to the working
#' directory by default.
#' @param interp Should interpolation of missing values be undertaken? If
#' \code{TRUE} linear interpolation is carried out for gaps of up to and
#' including \code{maxgap}.
#' @param maxgap The maximum gap in hours that should be interpolated where
#' there are missing data when \code{interp = TRUE.} Data with gaps more than
#' \code{maxgap} are left as missing.
#'
#' @return Writes a text file to a location of the user's choosing.
#' @export
#' @importFrom zoo na.approx
#' @examples
#'
#' \dontrun{
#' ## import some data then export it
#' dat <- importNOAA(year = 2012)
#' exportADMS(dat, file = "~/temp/adms_met.MET")
#' }
exportADMS <- function(dat, out = "./ADMS_met.MET", interp = FALSE, maxgap = 2) {

# keep R check quiet
wd = u = v = NULL

## make sure the data do not have gaps
all.dates <- data.frame(
date = seq(ISOdatetime(year = as.numeric(format(dat$date[1], "%Y")),
month = 1, day = 1, hour = 0, min = 0,
sec = 0, tz = "GMT"),
ISOdatetime(year = as.numeric(format(dat$date[1], "%Y")),
month = 12, day = 31, hour = 23, min = 0,
sec = 0, tz = "GMT"), by = "hour")
)

dat <- merge(dat, all.dates, all = TRUE)

## make sure precipitation is available
if (!"precip" %in% names(dat))
dat$precip <- NA

if (interp) {

## variables to interpolate
## note need to deal with wd properly
dat <- transform(dat, u = sin(pi * wd / 180), v = cos(pi * wd / 180))

varInterp <- c("ws", "u", "v", "air_temp", "RH", "cl", "precip")

# don't want to try and interpret fields that are all missing
ids <- sapply(varInterp, function (x) !all(is.na(dat[[x]])))
varInterp <- varInterp[ids]

dat[varInterp] <- zoo::na.approx(dat[varInterp], maxgap = maxgap, na.rm = FALSE)

## now put wd back
dat <- transform(dat, wd = as.vector(atan2(u, v) * 360 / 2 / pi))

## correct for negative wind directions
ids <- which(dat$wd < 0) ## ids where wd < 0
dat$wd[ids] <- dat$wd[ids] + 360
}

## exports met data to ADMS format file
year <- as.numeric(format(dat$date, "%Y"))
day <- as.numeric(format(dat$date, "%j"))
hour <- as.numeric(format(dat$date, "%H"))
station <- "0000"

## data frame of met data needed
adms <- data.frame(station, year, day, hour,
round(dat$air_temp, 1), round(dat$ws, 1),
round(dat$wd, 1), round(dat$RH, 1),
round(dat$cl), round(dat$precip, 1),
stringsAsFactors = FALSE)

## print key data capture rates to the screen
dc <- round(100 - 100 * (length(which(is.na(dat$ws))) / length(dat$ws)), 1)
print(paste("Data capture for wind speed:", dc, "%"))

dc <- round(100 - 100 * (length(which(is.na(dat$wd))) / length(dat$wd)), 1)
print(paste("Data capture for wind direction:", dc, "%"))

dc <- round(100 - 100 * (length(which(is.na(dat$air_temp))) / length(dat$air_temp)), 1)
print(paste("Data capture for temperature:", dc, "%"))

dc <- round(100 - 100 * (length(which(is.na(dat$cl))) / length(dat$cl)), 1)
print(paste("Data capture for cloud cover:", dc, "%"))

## replace NA with -999
adms[] <- lapply(adms, function(x) replace(x, is.na(x), -999))

## write the data file
write.table(adms, file = out, col.names = FALSE, row.names = FALSE,
sep = ",", quote = FALSE)

## add the header lines
fConn <- file(out, 'r+')
Lines <- readLines(fConn)
writeLines(c("VARIABLES:\n10\nSTATION DCNN\nYEAR\nTDAY\nTHOUR\nT0C\nU\nPHI\nRHUM\nCL\nP\nDATA:",
Lines), con = fConn)
close(fConn)

}



192 changes: 192 additions & 0 deletions R/getMeta.R
@@ -0,0 +1,192 @@
##' Get information on meteorological sites
##'
##' This function is primarily used to find a site code that can be
##' used to access data using \code{\link{importNOAA}}. Sites searches
##' of approximately 30,000 sites can be carried out based on the site
##' name and based on the nearest locations based on user-supplied
##' latitude and logitude.
##' @title Find a ISD site code and other meta data
##' @param site A site name search string e.g. \code{site =
##' "heathrow"}. The search strings and be partial and can be upper
##' or lower case e.g. \code{site = "HEATHR"}.
##' @param lat A latitude in decimal degrees to search. Takes the
##' values -90 to 90.
##' @param lon A longitude in decimal degrees to search. Takes values
##' -180 to 180. Negative numbers are west of the Greenwich
##' meridian.
##' @param country The country code. This is a two letter code. For a
##' full listing see
##' \url{ftp://ftp.ncdc.noaa.gov/pub/data/noaa/country-list.txt}.
##' @param state The state code. This is a two letter code.
##' @param n The number of nearest sites to search based on
##' \code{latitude} and \code{longitude}.
##' @param end.year To help filter sites based on how recent the
##' available data are. \code{end.year} can be "current", "any" or a
##' numeric year such as 2016, or a range of years e.g. 1990:2016
##' (which would select any site that had an end date in that range.
##' \strong{By default only sites that have some data for the
##' current year are returned}.
##' @param plot If \code{TRUE} will plot sites on an interactive
##' leaflet map.
##' @param fresh Should the meta data be read from the NOAA server or
##' the \code{worldmet} package?. If \code{FALSE} it is read from
##' the package version, which is fast but might be out of date. If
##' \code{TRUE} the data are read from the NOAA server.
##' @param returnMap Should the leaflet map be returned instead of the
##' meta data? Default is \code{FALSE}.
##' @return A data frame is returned with all available meta data,
##' mostly importantly including a \code{code} that can be supplied
##' to \code{\link{importNOAA}}. If latitude and longitude searches
##' are made an approximate distance, \code{dist} in km is also
##' returned.
##' @export
##' @author David Carslaw
##' @examples
##'
##' \dontrun{
##' ## search for sites with name beijing
##' getMeta(site = "beijing")
##' }
##'
##' \dontrun{
##' ## search for near a specified lat/lon - near Beijing airport
##' ## returns 'n' nearest by default
##' getMeta(lat = 40, lon = 116.9)
##' }
getMeta <- function(site = "heathrow", lat = NA, lon = NA,
country = NA, state = NA, n = 10, end.year = "current",
plot = TRUE, fresh = TRUE, returnMap = FALSE) {
## read the meta data

# check year
if (!any(end.year %in% c("current", "all"))) {
if (!is.numeric(end.year)) {
stop("end.year should be one of 'current', 'all' or a numeric 4-digit year such as 2016.")
}
}

# we base the current year as the max available in the meta data
if ("current" %in% end.year) end.year <- max(as.numeric(format(meta$END, "%Y")), na.rm = TRUE)
if ("all" %in% end.year) end.year <- 1900:2100


## download the file, else use the package version
if (fresh) meta <- getMetaLive()

## search based on name of site
if (!missing(site)) {
## search for station
meta <- meta[grep(site, meta$STATION, ignore.case = TRUE), ]

}

## search based on country codes
if (!missing(country) && !is.na(country)) {
## search for country
id <- which(meta$CTRY %in% toupper(country))
meta <- meta[id, ]

}

## search based on state codes
if (!missing(state)) {
## search for state
id <- which(meta$ST %in% toupper(state))
meta <- meta[id, ]

}

# make sure no missing lat / lon
id <- which(is.na(meta$LON))
if (length(id) > 0) meta <- meta[-id, ]

id <- which(is.na(meta$LAT))
if (length(id) > 0) meta <- meta[-id, ]

# filter by end year
id <- which(format(meta$END, "%Y") %in% end.year)
meta <- meta[id, ]

## approximate distance to site
if (!missing(lat) && !missing(lon)) {
r <- 6371 # radius of the Earth

## Coordinates need to be in radians
meta$longR <- meta$LON * pi / 180
meta$latR <- meta$LAT * pi / 180
LON <- lon * pi / 180
LAT <- lat * pi / 180
meta$dist <- acos(sin(LAT) * sin(meta$latR) + cos(LAT) *
cos(meta$latR) * cos(meta$longR - LON)) * r

## sort and retrun top n nearest
meta <- head(openair:::sortDataFrame(meta, key = "dist"), n)

}



if (plot) {

dat <- dplyr::rename(meta, latitude = LAT, longitude = LON)

if (!"dist" %in% names(dat)) dat$dist <- NA

content <- paste(paste(dat$STATION,
paste("Code:", dat$code),
paste("Start:", dat$BEGIN),
paste("End:", dat$END),
paste("Distance (km)", round(dat$dist, 1)),
sep = "<br/>"))


m <- leaflet(dat) %>% addTiles() %>%
addMarkers(~longitude, ~latitude, popup = content)

if (!is.na(lat) && !is.na(lon))
m <- m %>% addCircles(lng = lon, lat = lat, weight = 20, radius = 200,
stroke = TRUE, color = "red",
popup = paste("Search location",
paste("Lat =", lat),
paste("Lon =", lon), sep = "<br/>"))
print(m)

}


if (returnMap) return(m) else return(dat)

}

getMetaLive <- function(...) {

## downloads the whole thing fresh

## use RCurl
bin <- getBinaryURL("ftp://ftp.ncdc.noaa.gov/pub/data/noaa/isd-history.csv",
ssl.verifypeer = FALSE)
tmp <- tempfile()
writeBin(bin, tmp)
meta <- read.csv(tmp, header = FALSE, skip = 21)

## names in the meta file
names(meta) <- c("USAF", "WBAN","STATION", "CTRY", "ST", "CALL", "LAT",
"LON", "ELEV(M)", "BEGIN", "END")

## full character string of site id
meta$USAF <- formatC(meta$USAF, width = 6, format = "d", flag = "0")

## start/end date of measurements
meta$BEGIN <- as.Date(as.character(meta$BEGIN), format = "%Y%m%d")
meta$END <- as.Date(as.character(meta$END), format = "%Y%m%d")

## code used to query data
meta$code <- paste(meta$USAF, meta$WBAN, sep = "-")

return(meta)
}

# how to update meta data
# meta <- getMeta(fresh = TRUE)
# devtools::use_data(meta, overwrite = TRUE, internal = TRUE)
# devtools::use_data(meta, overwrite = TRUE)

0 comments on commit 672c3d0

Please sign in to comment.