Just a short update on ICOADS.

I started out an ICOADS project with grand plan to reprocess the observation files into a “raster” format. That ends up being a really big job. Along the way I had to learn a new package Rcurl and relearn some old Unix skills for compiling RCurl. Along the way I found more R tricks for doing things but I hit a brick wall on downloading the Icoads observation data. ( my problem, not the data supplier ). During the course of getting the observation data i hit some problems with the FTP and so I wrote to Scott Woodruff of NOAA, Scott was very helpful and pointed me at some great resources and documentation. Currently, some work is being completed around the various datasets, so it makes sense for me to delay my project, but I was able to use some of the assets Scott pointed me to to put together a quick animation of ICOADS data, all in R.

I thought it would be an easy task, maybe 1 hour, if that. Let’s go to the code and I’ll explain some of the issues. It ended up taking 3 days.

# first we get the processed  Monthly means of the SST data. Its in 2 degree bins. There is another

# product with 1 degree bins that starts in 1960. When you go through all the selection menus at
# NOAA you will find this url:
url_Icoads <-"ftp://ftp.cdc.noaa.gov/Datasets/icoads/2degree/enh/sst.mean.nc"

icoadsFile  <- basename(url_Icoads)

if(!file.exists(icoadsFile))download.file(url_Icoads,icoadsFile,mode="wb")

# now try to load it as a brick

Icoads <- brick(icoadsFile,varname="sst")

dirIcoads <- "Icoads"

if(!file.exists(dirIcoads))dir.create(dirIcoads)

#  now we set up an animation.

All that is pretty straightforward. What’s hidden is the details. ICOADS is so large that in R it has to stay on disk. That means to get certain values we require for plotting we have to perform functions on the disk file which are painfully slow in R on a MAC. Functions like “setMinMax” and the function “rotate.” Because ICOADS uses 0-360 in longitude that has to be transformed to -180 to 180 to work with world map plots. That process is painfully slow, and once the file is read into memory, everything in R on a MAC grinds to a halt. There are some ways around this which I’ll explore, but for now I’ll just  use one approach to get the data out and build an animation from it. Note, the data goes from 0-360 in longitude. That can be changed with a “rotate” command.

Ani.options(nmax=50,ani.width=1080,ani.height=540,withprompt="> ",

ani.type ="png",outdir=dirIcoads,imgdir="IMAGES",

filename = "world.htm",ani.dev=png,title="Icoads SST",

description ="world animation",autobrowse=TRUE,interval=.1)

Ani.start()

i <- 1000

while (i <= 1050) {

y <- raster(Icoads,layer=i)

plot(y, zlim=c(-2,40),col=rainbow(n=100,start=4/6,end=0),main="iii")

i <- i + 1

}

Ani.stop()

First off, I’ve had to modify the sourcecode for the animation package. Since this is R that’s easy, I just download the source. The animation package works like this. First, a call is made to ani.options() where options are set for the animation, including the names of output directories. Then a call to ani.start() is made. In this call the standard package will create those directories and do some other set up work… INCLUDING changing the working directory. Next the package expects you to call a plot routine. But when we are plotting a brick that is on disk, our access path is tied to the working directory, which the call to ani.start() has changed. So the plot will fail to find the file. There is a way to get around this problem in raster that is currently being addressed, but I decided to change the animation code as well because I’m not too keen about a design that changes directories and imposes API ordering. Anyway, ani.stop() changes directories back to finish the work. ani.stop() basically uses the ani.options() to output a html file. You click on that file and  you will get the animation to play. I modified the ani() functions to remove the changing of directories and made a few other minor changes and I renamed them Ani.option(),Ani.start() and Ani.stop(). the code needs some clean up, and credit given, but it functions. Autobrowse has been removed as an option. When you execute that code you get 50 frames. In the entire brick there are over 2500 frames of data so if you want to animate the whole thing, it takes a while. Then we will need to add “credits” for the ICOADS team and we have a final version..A couple screen caps below