How to install rNOMADS with GRIB file support on Windows

Two years ago, I wrote a software package for R called “rNOMADS” that interfaces with online weather and sea ice model repositories to gather data in real time, for free.  The data are delivered in two ways: a simple, pure R, cross platform interface using GrADS-DODS, and binary files in GRIB format.  The one issue with GRIB is that this format can’t be read directly into R; it requires the external program “wgrib2.”  Installing rNOMADS with GRIB support for Linux is covered in this post (Mac OS is probably similar).  I thought GRIB support for Windows was impossible until the guy who runs this blog casually told me he’d figured it out.  So, I finally got around to trying it myself, and I’m happy to say I got it to work!  Here’s how:

Step 1:

Download the most recent wgrib2.exe version and all DLL files from Wesley Ebisuzaki’s web site:  Click here.

Step 2:

Make a directory somewhere on your computer.   I chose:  C:\Program Files\wgrib2

Copy the .exe and the .dll files to this directory.

Step 3:

Append the directory path to the Windows PATH variable.  Find out how to do this here.

Edit: a user points out that sometimes this doesn’t work.  If not, change the path variable in R using Sys.setenv():

Sys.setenv("PATH" = "C:/Program Files/wgrib2")

 

Step 4:  If R is open, close and reopen it.  Then try the following command in the R interpreter:
system("wgrib2")

If you get a bunch of text that looks like the image below, you’ve succeeded.  You can start reading GRIB files with rNOMADS!

Output of wgrib2 call from R.

Output of wgrib2 call from R.

New updates to the rNOMADS package and big changes in the GFS model

I rolled out a big update to the rNOMADS package in R about two weeks ago.  Now, the list of real time weather, ocean, and sea ice models available through rNOMADS updates automatically by scraping the NOMADS web site.  This way, changes in model inventories will be instantly reflected in rNOMADS without the need for a new version release.

Keep abreast of future updates to rNOMADS by subscribing to the mailing list here.  Feel free to ask for help or make comments on this list as well.

In other news, NOAA just updated the Global Forecast System to provide 0.25 x 0.25 degree output – doubling the resolution of the model!  Check out this crystal clear views of surface temperatures across the planet (source code below the image):

World temperature at 2 m above ground using the 0.25 x 0.25 degree output of the Global Forecast System model.

World temperature at 2 m above ground using the 0.25 x 0.25 degree output of the Global Forecast System model.

 


library(rNOMADS)
library(GEOmap)

#Get dates of model output
model.urls <- GetDODSDates(“gfs_0p25”)

#Find day of most recent model run
latest.model <- tail(model.urls$url, 1)

#Find most recent model run on that day
model.runs <- GetDODSModelRuns(latest.model)

#Get the most recent model (excluding analysis only)
latest.model.run <- tail(model.runs$model.run[which(grepl(“z$”, model.runs$model.run))], 1)

#Define model domain
time <- c(0,0) #Analysis model
lon <- c(0, 1439) #All longitude points
lat <- c(0, 720) #All latitude points
variables <- c(“tmp2m”) #Temperature 2 m above ground

#Get data from NOMADS real time server
tmp.data <- DODSGrab(latest.model, latest.model.run,
variables, time, lon, lat, display.url = FALSE)

#Reformat it
tmp.grid <- ModelGrid(tmp.data, c(0.25, 0.25))

#Define color scale
colormap <- rev(rainbow(500, start = 0 , end = 5/6))

#Plot it
image(x = tmp.grid$x, y = sort(tmp.grid$y), z = tmp.grid$z[1,1,,], col = colormap,
xlab = “Longitude”, ylab = “Latitude”,
main = paste(“World Temperature at Ground Level:”,
tmp.grid$fcst.date))

plotGEOmap(coastmap, border = “black”, add = TRUE,
MAPcol = NA)

rNOMADS 2.0.2 released

I uploaded the newest version of rNOMADS to CRAN yesterday. This one has a new plotting function for wind altitude azimuth, and magnitude (see below for plot and source code). I also added a function for reading GRIB inventories, fixed a few typos, and tweaked a few functions. In other news, the rNOMADS subversion repository has migrated to R-Forge: http://rnomads.r-forge.r-project.org/.

I’ve also set up an rNOMADS mailing list; subscribe here if you need to ask for help or want to hear about new releases.

Wind profiles at each infrasound station in the Transportable Array (currently on the US East Coast). Ground surface is at the center of the plot, the top of the stratosphere is at the outer radius. Wind magnitudes are denoted by color, and azimuth by position around the circle.

Wind profiles at each infrasound station in the Transportable Array (currently on the US East Coast). Ground surface is at the center of the plot, the top of the stratosphere is at the outer radius. Wind magnitudes are denoted by color, and azimuth by position around the circle.


library(rNOMADS)

download.file("http://www.unc.edu/~haksaeng/rNOMADS/myTA.RDATA",
   destfile = "myTA.RDATA")
load("myTA.RDATA")
#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 model nodes

lons <- seq(0, 359.5, by = 0.5)
lats <- seq(-90, 90, by = 0.5)
lon.ind <- which(lons <= (max(myTA$lon + 360) + 1) & lons >= (min(myTA$lon + 360) - 1))
lat.ind <- which(lats <= (max(myTA$lat) + 1) & lats >= (min(myTA$lat) - 1))
levels <- c(0, 46)
time <- c(0, 0)

#Get data 
variables <- c("hgtprs", "ugrdprs", "vgrdprs")
model.data <- DODSGrab(latest.model, latest.model.run,
   variables, time, c(min(lon.ind), max(lon.ind)),
   c(min(lat.ind), max(lat.ind)), levels)

#Build profiles
profile <- BuildProfile(model.data, myTA$lon, myTA$lat,
    spatial.average = FALSE)

#Build profiles
zonal.wind      <- NULL
meridional.wind <- NULL
height          <- NULL

for(k in 1:length(profile)) {
   hgt  <- profile[[k]]$profile.data[, which(profile[[k]]$variables == "hgtprs"),]
   ugrd <- profile[[k]]$profile.data[, which(profile[[k]]$variables == "ugrdprs"),]
   vgrd <- profile[[k]]$profile.data[, which(profile[[k]]$variables == "vgrdprs"),]

   synth.hgt <- seq(min(hgt),
       max(hgt), length.out = 1000)
   ugrd.spline <- splinefun(hgt, ugrd, method = "natural")
   vgrd.spline <- splinefun(hgt, vgrd, method = "natural")
   zonal.wind[[k]] <- ugrd.spline(synth.hgt)
   meridional.wind[[k]] <- vgrd.spline(synth.hgt)
   height[[k]] <- synth.hgt
}

#Plot them all
PlotWindProfile(zonal.wind, meridional.wind, height, lines = TRUE,
    points = FALSE, elev.circles = c(0, 25000, 50000), elev.labels = c(0, 25, 50),
    radial.lines = seq(45, 360, by = 45), colorbar = TRUE, invert = FALSE,
    point.cex = 2, pch = 19, lty = 1, lwd = 1,
    height.range = c(0, 50000), colorbar.label = "Wind Speed (m/s)")


rNOMADS 2.0.1 released

The latest version of rNOMADS is now available on CRAN.  This update resolves several minor bugs and one major one involving multiple variable/level selections when using the ModelGrid function.  I have also added support for two more models on NOMADS:  Climate Forecast System Flux Products and Climate Forecast System 3D Pressure Products.  This brings the total number of real time model products supported by rNOMADS to 57.  A quick test of the Climate Forecast System Flux Products produced the following image of surface albedo for 1200 GMT on July 4th, 2014:

Albedo of the Earth's surface at 1200 GMT, July 4th, 2014.

Albedo of the Earth’s surface at 1200 GMT, July 4th, 2014.

In other news, data curators at the University Corporation for Atmospheric Research (UCAR) are currently investigating rNOMADS as an open source software package to read their GRIB archive.  UCAR hosts a series of model archive data, some of which is freely available to the public.  Their specific interest is in using rNOMADS to read their archive of the FNL model.  They have gotten wgrib2 to work on MAC OS X, allowing rNOMADS GRIB functionality to be extended to Macintosh machines.

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.

Accessing global weather model data using the rNOMADS package in R

*** UPDATE ***

rNOMADS 2.0.0 was released yesterday (5-15-2014), and it now has cross-platform support (meaning you can use its data retrieval capabilities in Windows and Mac as well as Linux).

Stay tuned for more tutorials and usage examples in the coming days.

***

The rNOMADS package interfaces with the NOAA Operational Model Archive and Distribution System to provide access to 60 operational (i.e. real time and prediction) models describing the state of the ocean and the atmosphere. rNOMADS has been used to get wind and wave data for a real time sailing game, to quantify solar energy available for power plants in Europe, and to predict helium balloon flights. We look forward to continuing to hear about novel and unexpected uses for this spectacular data set.

In this post, we show how to use rNOMADS to do the following:
1. Plot global temperature and wind at different atmospheric pressure levels
2. Produce a precise atmospheric temperature and wind profile for a specific point at a specific time

Links to source code are provided below each set of figures. A link to rNOMADS installation instructions is provided at the end of the post.

Global Temperature and Wind Maps

The Global Forecast System (GFS) model provides weather data on a 0.5 x 0.5 degree grid for the entire planet.  It is run 4 times daily, and produces a prediction every three hours out to 180 hours.  The zero hour “analysis” forecast is the state of the atmosphere at the model run time, and it uses some observational data to increase its accuracy.  Here, we plot the temperature at 2 m above the surface, the wind at 10 m above the surface, and the winds at 300 mb (around 9 kilometers elevation).  The 300 mb plot often shows the northern and southern jet streams quite well.

Temperature at the Earth's surface determined using the Global Forecast System model.

Temperature at the Earth’s surface determined using the Global Forecast System model.

Winds at the surface of the Earth from the GFS model.  Note the little spot of high winds - that's Tropical Cyclone Gillian, a Category 3 storm when this image was generated.

Winds at the surface of the Earth from the GFS model. Note the little spot of high winds south of Indonesia – that’s Tropical Cyclone Gillian, a Category 3 storm when this image was generated.

Jet streams and Rossby waves are clearly visible in this image of the wind speeds of the upper troposphere/lower stratosphere.

Jet streams and Rossby waves are clearly visible in this image of the wind speeds of the upper troposphere/lower stratosphere.

Download the source code for these images here.

 

Instantaneous Atmospheric Profile over Sakura-Jima volcano, Japan

It’s important to know which direction the winds are going directly above active volcanoes, because eruptions can carry ash into air space and over inhabited areas.  One impetus for the continued development of rNOMADS was to provide a one-stop solution for generating high precision wind profiles over specific points, allowing ash distribution prediction as soon as an eruption begins.  Here, we have generated a spatially and temporally interpolated wind and temperature profile over Sakura-Jima volcano, Japan.  The profile is calculated for the exact time when the source code is run.

Instantaneous temperature profile above Sakura-jima volcano, Japan.

Instantaneous temperature profile above Sakura-jima volcano, Japan.

An eruption at this instant would produce ashfalls east of the volcano for a plume height of 15 km.  However, if a truly massive 30 km high plume was produced, ashfalls would occur both east and west of the volcano.

An eruption at this instant would produce ashfalls east of the volcano for a plume height of 15 km. However, if a truly massive 30 km high plume was produced, ashfalls would occur both east and west of the volcano.

Download the source code for these images here.

 

Some of these scripts require the aqfig package in R to generate the colorbar legends.

Instructions on installing R, rNOMADS, and their dependencies are here.

How to install rNOMADS with GRIB support on Linux

About 6 months ago, I wrote a package for the R programming language called rNOMADS.  It interfaces with the National Oceanic and Atmospheric Association’s Operational Archive and Distribution System (NOMADS), comprising everything from global weather forecast models, to high resolution regional models, to wave and sea ice models.  rNOMADS uses R to download requested data from the NOMADS servers using two different methods: direct ascii data retrieval and GRIB file retrieval.  GRIB is a binary file format that can’t be read into R directly – it requires external programs that must be installed separately. This post describes how to install rNOMADS with GRIB support, which is available for Linux systems only.  Windows and Mac users should use the ascii download routines; I will write a tutorial and post it shortly.

rNOMADS allows R users to make some spectacular pictures (see below) but I also hope that it will lead to greater availability of these valuable resources to data analysts all over the world.  Already, rNOMADS is being used in an online sailboat racing game and to get data for predicting solar and wind energy output.

World temperature at 2 m above ground, generated using rNOMADS.

World temperature at 2 m above ground, generated using rNOMADS.

Step 1

The first step is to make sure R is configured correctly. You can skip the reinstallation (step 2) if need be, but I recommend getting r-base-dev anyway for all the libraries it has.

su
apt-get install r-base-dev
apt-get build-dep r-base
exit

Step 2

Then, I build R from the source. You can skip this method if you’ve already done this.

You can download the latest version  here.Get the tar.gz file next to the bullet point that says “Latest Release.”
Navigate to the file and type:

tar -xvf R-3.0.3.tar.gz
cd R-3.0.3
./configure
make
su
make install
exit

Step 3

Now we get dependencies for specific packages that rNOMADS requires.

su
apt-get install libxml2-dev
apt-get install libcurl4-openssl-dev
apt-get install libboost-dev
exit

Step 4

Now, we install wgrib2 and wgrib and add them to our path.
Download wgrib2 here.
Download wgrib here.
Then:

tar -xvf wgrib2.tgz
cd grib2/
make
su
mv wgrib2/ /usr/local
exit

mkdir wgrib
mv wgrib.tar wgrib
cd wgrib
tar -xvf wgrib.tar
make
su
mkdir /usr/local/wgrib
mv * /usr/local/wgrib
exit


su
gedit /etc/environment

add “/usr/local/wgrib2” and “/usr/local/wgrib” to the end of the PATH variable, then close /etc/environment and close the superuser terminal.
Next:

source /etc/environment

Step 5

Now we install rNOMADS and another package called GEOmap so that we can build and test the example vignette.

su
R

Now we are in the R interpreter, as root. Type:

install.packages('rNOMADS')
install.packages('GEOmap')

Step 6

If everything went well, you now have rNOMADS! Test rNOMADS (and learn how it works) by building the example vignette. First, download the source here and navigate to the saved file. Then, type:

R CMD Sweave rNOMADS_grib_examples.Rnw
pdflatex rNOMADS_grib_examples.tex

The latest R release has been causing issues with Sweave; if you don’t see any figures in your PDF or if the document compilation fails, try running Sweave using the commands described here.

This should create a PDF file called “rNOMADS_grib_examples.pdf” with code and images. It should look similar to this one.