Sunday, October 7, 2007

public access large-scale model data

This may not be new news, but it was a good find. Here is a very extensive archive of publicly available model output. This could be used, for example, for regional model boundary conditions:

[http://apdrc.soest.hawaii.edu/dods/public_data/](http://apdrc.soest.hawaii.edu/dods/public_data/)

An example python script for making a movie of ssh in the Indian ocean from the Navy Layer Ocean Model is found after the fold


plot\_nlom\_ssh.py:

from numpy import *
from dap.client import open as opendap
from matplotlib.toolkits.basemap import Basemap
import pylab as pl

dapserver = 'http://apdrc.soest.hawaii.edu/dods/public_data/SODA/soda_pop2.0.3'

nc = opendap(dapserver)

t = nc.time[:]
dates = asarray(pl.num2date(t))

lon_idx = slice(100, 1070)
lat_idx = slice(2350, 3300)

lon = nc.lon[lon_idx]
lat = nc.lat[lat_idx]

proj = Basemap(projection='lcc',
resolution='i',
llcrnrlon=44.5,
llcrnrlat= 1.5,
urcrnrlon=83.0,
urcrnrlat=33.0,
lat_0=20.0,
lon_0=60.0)

x, y = proj(*meshgrid(lon, lat))

fig = pl.figure(figsize=(5.5, 5.0))
for n in range(len(t)):
ssh = nc.ssh[n, 0, lat_idx, lon_idx]
ssh = ma.masked_where(ssh==ssh.min(), ssh)
if not all(ssh.mask):
fig.clear()
ax = fig.add_axes([0, 0, 1, 1])
proj.drawcoastlines(ax=ax)
#proj.drawrivers(color=(0.5, 0.5, 1.0))
pc = ax.pcolormesh(x, y, ssh-35.0, vmin=-50.0, vmax=50.0, cmap=pl.cm.RdBu_r)
proj.fillcontinents(color=3*(0.9,), ax=ax)
proj.drawcountries(color=(0.5, 0.5, 0.5), ax=ax)
txt = pl.text(2.5e5, 2.1e6, dates[n].strftime('%b %d, %Y'), family='monospace')
axcb = fig.add_axes([0.87, 0.48, 0.02, 0.35])
pl.colorbar(pc, axcb).set_label('ssh [cm]')
pl.savefig('frame_%d.png' % n)
print dates[n]

No comments:

Post a Comment