Downloading weather, sea ice, and wave model data with the rNOMADS package in R

The NOAA Operational Model Archive and Distribution System is a treasure trove of near real time and archived model outputs describing global and regional weather, sea ice, and wave data.  I developed the rNOMADS package about a year ago to make this available to R users. In this post, I’ll present some source code and a couple of figures showing a few of the useful things you can do with rNOMADS.

For detailed examples showing rNOMADS with GRIB file support (Linux only) see the vignette here, for a cross platform version of the same, see here.

1.  Getting wind speed at a specific point

library(rNOMADS)

#A location near my house
lat <- 35.828304
lon <- -79.107467

#Find the latest Global Forecast System model run
model.urls <- GetDODSDates("gfs_0p50")
latest.model <- tail(model.urls$url, 1)
model.runs <- GetDODSModelRuns(latest.model)
latest.model.run <- tail(model.runs$model.run, 1)

#Get nearest model nodes
lons <- seq(0, 359.5, by = 0.5)
lats <- seq(-90, 90, by = 0.5)
lon.diff <- abs(lon + 360 - lons)
lat.diff <- abs(lat - lats)
model.lon.ind <- which(lon.diff == min(lon.diff)) - 1 #NOMADS indexes at 0
model.lat.ind <- which(lat.diff == min(lat.diff)) - 1

#Subset model
time <- c(0,0) #Model status at initialization
lon.inds <- c(model.lon.ind - 2, model.lon.ind + 2)
lat.inds <- c(model.lat.ind - 2, model.lat.ind + 2)
variables <- c("ugrd10m", "vgrd10m") #E-W and N-S wind

wind.data <- DODSGrab(latest.model, latest.model.run, variables,
time, lon.inds, lat.inds)

profile <- BuildProfile(wind.data, lon, lat, spatial.average = TRUE, points = 4)

#Present results!
print(paste("At", profile[[1]]$forecast.date, "the East-West winds at Briar Chapel were going", sprintf("%.2f", profile[[1]]$profile.data[1, which(profile[[1]]$variables == "ugrd10m"), 1]),
"meters per second, and the north-south winds were going", sprintf("%.2f", profile[[1]]$profile.data[1, which(profile[[1]]$variables == "vgrd10m"), 1]),
"meters per second."))

#How did I know all these strange parameter names?
info <- GetDODSModelRunInfo(latest.model, latest.model.run)
print(info)

2. Getting a temperature profile from 0 to 40 km above a specific point

library(rNOMADS)

#A location near my house
lat <- 35.828304
lon <- -79.107467

#Find the latest Global Forecast System model run
model.urls <- GetDODSDates("gfs_0p50")
latest.model <- tail(model.urls$url, 1)
model.runs <- GetDODSModelRuns(latest.model)
latest.model.run <- tail(model.runs$model.run, 1)

#Get nearest model nodes
lons <- seq(0, 359.5, by = 0.5)
lats <- seq(-90, 90, by = 0.5)
lon.diff <- abs(lon + 360 - lons)
lat.diff <- abs(lat - lats)
model.lon.ind <- which(lon.diff == min(lon.diff)) - 1 #NOMADS indexes at 0
model.lat.ind <- which(lat.diff == min(lat.diff)) - 1

#Subset model
time <- c(0,0) #Model status at initialization
lon.inds <- c(model.lon.ind - 2, model.lon.ind + 2)
lat.inds <- c(model.lat.ind - 2, model.lat.ind + 2)
levels <- c(0, 46) #All pressure levels
variables <- c("tmpprs", "hgtprs") #First get temperature

model.data <- DODSGrab(latest.model, latest.model.run, variables,
time, lon.inds, lat.inds, levels)

#Interpolate to the point of interest
profile <- BuildProfile(model.data, lon, lat,
spatial.average = TRUE, points = 4)

#Plot it!
tmp <- profile[[1]]$profile.data[,which(profile[[1]]$variables == "tmpprs"),] - 272.15
hgt <- profile[[1]]$profile.data[,which(profile[[1]]$variables == "hgtprs"),]

plot(tmp, hgt, xlab = "Temperature (C)",
ylab = "Geopotential Height", main = paste("Temperature above Chapel Hill, NC, at",
profile[[1]]$forecast.date))
Temperature profile of troposphere and stratosphere above Chapel Hill, NC, on May 28.

Temperature profile of troposphere and stratosphere above Chapel Hill, NC, on May 28.

2. A world map of surface temperature

library(GEOmap)
library(rNOMADS)

model.urls <- GetDODSDates("gfs_0p50")
latest.model <- tail(model.urls$url, 1)
model.runs <- GetDODSModelRuns(latest.model)
latest.model.run <- tail(model.runs$model.run, 1)

time <- c(0,0) #Analysis model
lon <- c(0, 719) #All 720 longitude points
lat <- c(0, 360) #All 361 latitude points

tmp2m.data <- DODSGrab(latest.model, latest.model.run,
"tmp2m", time, lon, lat)
atmos <- ModelGrid(tmp2m.data, c(0.5, 0.5))
colormap <- rev(rainbow(500, start = 0 , end = 5/6))
image(atmos$x, sort(atmos$y), atmos$z[1,1,,], col = colormap,
xlab = "Longitude", ylab = "Latitude",
main = paste("World Temperature at Ground Level:", atmos$model.run.date))
plotGEOmap(coastmap, border = "black", add = TRUE,
MAPcol = NA)

World temperature on May 28.

World temperature on May 28.

2. Wave heights in the northwest Atlantic ocean

library(rNOMADS)
library(GEOmap)

model.urls <- GetDODSDates("wave")
latest.model <- tail(model.urls$url, 1)
model.runs <- GetDODSModelRuns(latest.model) #West Atlantic waves
latest.model.run <- tail(model.runs$model.run, 1)
time <- c(0,0)
lon <- c(0, 274)
lat <- c(0, 202)

wave.data <- DODSGrab(latest.model, latest.model.run,
"htsgwsfc", time, lon, lat)
wave.grid <- ModelGrid(wave.data, c(0.25, 0.25))
#Remove "no data" values
wave.grid$z[which(wave.grid$z>1e10, arr.ind=TRUE)] <- NA
colormap <- rainbow(500, start=0, end=5/6)
image(wave.grid$x, sort(wave.grid$y), wave.grid$z[1,1,,], col = colormap,
xlab = "Longitude", ylab = "Latitude",
main = paste("Wave Height:", wave.grid$model.run.date))
plotGEOmap(coastmap, border = "black", add = TRUE,
MAPcol = "black")

Wave heights in the north Atlantic on July 4th, 2014. High seas from Hurricane Arthur are visible just off the US East Coast.

Advertisements

One thought on “Downloading weather, sea ice, and wave model data with the rNOMADS package in R

  1. Pingback: Downloading weather, sea ice, and wave model data with the rNOMADS package in R ← Patient 2 Earn

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s