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.

Anyone seen a red lunch box with cameras attached?

Catastrophe!  Yesterday, while launching a hydrogen balloon in Chapel Hill, North Carolina, our satellite tracker fell out of the payload.  The camera box, however, cheerfully continued its ascent into the wild blue yonder.  The result?  There is a red lunch box with three cameras in it somewhere in North Carolina.

Detailed description of payload as it will appear on the ground:

Shreds of mylar and a balloon nozzle, followed by a 40′ string, then a plastic parachute (may look like a plastic bag from far away).  Another length of string approximately 10′ long, followed by a red lunch box with “Solar Balloon Payload” and my phone number written on it.  The lunch box has three cameras, one pointing down, one pointing out the side, and a smaller one pointing up.  A milk jug with a sticker saying “Soleil Multimedia” may also still be attached.

The projected landing zone was in north/northwest Harnett County.  I believe that is reasonably accurate.  So our payload, with spectacular footage, is somewhere out there.  I want it back.

Here’s a couple of potential flight trajectories:

Trajectory 1

Trajectory 2

Trajectory 3

If you have any information or know anyone in the region who might be able to help, please contact us.