Search This Blog

Friday, October 28, 2011

Global RTOFS

Sea surface temperature, 2011-10-27
Thanks to John Kelley for pointing this out to me.  There's a Global RTOFS ocean model running now with pretty easy access via python.  I spent last night trying to get it up and running on my MacBook Pro, here are a few notes.  Before diving into it, check out the website for Global RTOFS here.

Step 1: get some python modules

 # unpack the tarball
cd ~/installs_Mac/python_modules/
tar -zxvf netCDF4-0.9.7.tar.gz
cd netCDF4-0.9.7


# install some hdf5 libs
fink install hdf5-18

# install some netcdf4 libs, this installs hdf5 stuff also so I'm not sure the above step was needed?
fink install netcdf7

# I had to figure out where fink installed the shared libraries for netcdf4 (they're in /sw/opt/netcdf7/)

otool -L /sw/bin/ncdump
/sw/bin/ncdump:
/sw/opt/netcdf7/lib/libnetcdf.7.dylib (compatibility version 9.0.0, current version 9.1.0)
/sw/lib/libsz.2.dylib (compatibility version 3.0.0, current version 3.0.0)
/sw/lib/libhdf5_hl.7.dylib (compatibility version 8.0.0, current version 8.1.0)
/sw/lib/libhdf5.7.dylib (compatibility version 8.0.0, current version 8.1.0)
/usr/lib/libSystem.B.dylib (compatibility version 1.0.0, current version 125.2.11)
/usr/lib/libz.1.dylib (compatibility version 1.0.0, current version 1.2.3)
/sw/lib/libcurl.4.dylib (compatibility version 7.0.0, current version 7.0.0)

# set the NETCDF4_DIR based on output from otool above
export NETCDF4_DIR=/sw/opt/netcdf7/

# okay, now ready tot build the netCDF4 for python bits
# Set some enviroenment variablesexport CFLAGS=-m32
python2.6 setup.py install --home=~

# Now do the same for the basemap module
cd ~/installs_Mac/python_modules/
tar -zxvf basemap-1.0.1.tar.gz
cd basemap-1.0.1

# Need to use the default 'python' instead of the version installed through fink...wish I knew why
PATH="/Library/Frameworks/Python.framework/Versions/2.6/bin:${PATH}"
export PATH

# can't do the usual CFLAGS=-m32 (which is usually required for building against python as installed
# by fink?
unset CFLAGS
export GEOS_DIR=/sw/opt/libgeos3.2.2/
python2.6 setup.py install --home=~



Step 2.  Plot some data

Now that all goodies are installed, I tried out the script from here (note that I had to change the date to 20111027, the date provided in the example didn't work for me).

Here's my version of the script, mine does a salinity plot as well.


#!/usr/bin/env python2.6

from mpl_toolkits.basemap import Basemap
import numpy as np
import matplotlib.pyplot as plt
from pylab import *
import netCDF4


mydate='20111027'

url_temp='http://nomads.ncep.noaa.gov:9090/dods/rtofs/rtofs_global'+mydate+'/rtofs_glo_3dz_forecast_daily_temp'
file = netCDF4.Dataset(url_temp)
lat  = file.variables['lat'][:]
lon  = file.variables['lon'][:]
data_temp      = file.variables['temperature'][1,1,:,:]
file.close()


m=Basemap(projection='mill',lat_ts=10, llcrnrlon=lon.min(),urcrnrlon=lon.max(), llcrnrlat=lat.min(),urcrnrlat=lat.max(), resolution='c')
Lon, Lat = meshgrid(lon,lat)
x, y = m(Lon,Lat)

plt.figure()
cs = m.pcolormesh(x,y,data_temp,shading='flat', cmap=plt.cm.jet)

m.drawcoastlines()
m.fillcontinents()
m.drawmapboundary()
m.drawparallels(np.arange(-90.,120.,30.), labels=[1,0,0,0])
m.drawmeridians(np.arange(-180.,180.,60.), labels=[0,0,0,1])

colorbar(cs)
plt.title('Example 1: RTOFS Global Temperature')
plt.show()

# Now do the salinity
url_salt='http://nomads.ncep.noaa.gov:9090/dods/rtofs/rtofs_global'+mydate+'/rtofs_glo_3dz_forecast_daily_salt'
file = netCDF4.Dataset(url_salt)
lat  = file.variables['lat'][:]
lon  = file.variables['lon'][:]
data_salt      = file.variables['salinity'][1,1,:,:]
file.close()

m=Basemap(projection='mill',lat_ts=10, llcrnrlon=lon.min(),urcrnrlon=lon.max(), llcrnrlat=lat.min(),urcrnrlat=lat.max(), resolution='c')
Lon, Lat = meshgrid(lon,lat)
x, y = m(Lon,Lat)

plt.figure()
cs = m.pcolormesh(x,y,data_salt,shading='flat', cmap=plt.cm.jet)

m.drawcoastlines()
m.fillcontinents()
m.drawmapboundary()
m.drawparallels(np.arange(-90.,120.,30.), labels=[1,0,0,0])
m.drawmeridians(np.arange(-180.,180.,60.), labels=[0,0,0,1])

colorbar(cs)
plt.title('Example 1: RTOFS Global Salinity')
plt.show()



Results?

Here's a few screen shots of the plots I managed to generate using my script.  Next up?  Extracting vertical profiles of temperature and salinity.


No comments:

Post a Comment