Extract snapshot from 3D [time,lat,lon] dataset


In [178]:
# Extract subsetted snapshots of wind from 4km NAM model

# This approach would work for any CF-Compliant NetCDF or OPeNDAP dataset

In [261]:
from pylab import *
import netCDF4
import datetime as dt
from mpl_toolkits.basemap import Basemap

In [262]:
# wind now
dnow = dt.datetime.utcnow() - dt.timedelta(hours=17)
# or wind at specific time (UTC)
#dnow = dt.datetime(2013,6,13,21)
daystr = dnow.strftime('%Y%m%d')
print daystr


20130703

In [263]:
cycle=max(0,(dnow.hour/6-1)*6)  # rounding by 6 hours to last cycle
print cycle


0

In [264]:
#NAM 4km
url='http://nomads.ncep.noaa.gov:9090/dods/nam/nam%s/nam_conusnest_%2.2dz' % (daystr,cycle)
url='http://nomads.ncep.noaa.gov:9090/dods/nam/nam%s/nam_conusnest_%2.2dz' % ('20130613',6)

lon_range=[-74.5, -67.];lat_range=[39., 44.];  # 4km lon is [-180 180]

# NAM 12km
#url=sprintf('http://nomads.ncep.noaa.gov:9090/dods/nam/nam%s/nam_%2.2dz',str,d,h);
#lon_range=[-71.5 -63];lat_range=[41 46];  % 12km lon is [-180 180]

# NAM 32km
#url=sprintf('http://nomads.ncep.noaa.gov:9090/dods/nam/nam%s/nam_na_%2.2dz',str,d,h);
#lon_range=[-71.5 -63]+360;lat_range=[41 46]; % 32km lon is [0 360]

print url
bbox=[lon_range[0], lat_range[0], lon_range[1], lat_range[1]]
print bbox


http://nomads.ncep.noaa.gov:9090/dods/nam/nam20130613/nam_conusnest_06z
[-74.5, 39.0, -67.0, 44.0]

In [265]:
# NetCDF4-Python can open OPeNDAP dataset just like a local NetCDF file
nc = netCDF4.Dataset(url)

In [266]:
nc.variables.keys()


Out[266]:
[u'time',
 u'lev',
 u'lat',
 u'lon',
 u'absvprs',
 u'no4lftx180_0mb',
 u'acpcpsfc',
 u'albdosfc',
 u'apcpsfc',
 u'bgrunsfc',
 u'bmixlhlev1',
 u'brtmpsfc',
 u'brtmptoa',
 u'capesfc',
 u'cape180_0mb',
 u'cape90_0mb',
 u'cape255_0mb',
 u'ccondsfc',
 u'cdsfc',
 u'cdconclm',
 u'cfrzrsfc',
 u'ciceprs',
 u'cicehlev1',
 u'cicepsfc',
 u'cinsfc',
 u'cin180_0mb',
 u'cin90_0mb',
 u'cin255_0mb',
 u'clwmrprs',
 u'clwmrhlev1',
 u'cnvhrclm',
 u'cnwatsfc',
 u'cpofpsfc',
 u'cpratsfc',
 u'crainsfc',
 u'csdsfsfc',
 u'csnowsfc',
 u'cueficlm',
 u'dlwrfavesfc',
 u'dlwrfsfc',
 u'dptprs',
 u'dpt2m',
 u'dpt30_0mb',
 u'dswrfavesfc',
 u'dswrfsfc',
 u'dzdtprs',
 u'evpsfc',
 u'fricvsfc',
 u'gfluxavesfc',
 u'gfluxsfc',
 u'gustsfc',
 u'hcdchcll',
 u'hgtsfc',
 u'hgtprs',
 u'hgthlev1',
 u'hgtclb',
 u'hgttop0c',
 u'hgtceil',
 u'hgtpbl',
 u'hgtlwb0',
 u'hgtlblsw',
 u'hgthtlsw',
 u'hgtclt',
 u'hgt0c',
 u'hgtl5',
 u'hgtmwl',
 u'hlcy3000_0m',
 u'hlcy1000_0m',
 u'hpblsfc',
 u'icecsfc',
 u'landsfc',
 u'lcdclcll',
 u'lftxl100_100',
 u'lhtflavesfc',
 u'lhtflsfc',
 u'lrghrclm',
 u'ltngsfc',
 u'lwhrclm',
 u'maxrh2m',
 u'mcdcmcll',
 u'mconv850mb',
 u'mconv950mb',
 u'mconv30_0mb',
 u'mconvclm',
 u'minrh2m',
 u'msletmsl',
 u'mstav0_100cm',
 u'mxsalbsfc',
 u'ncpcpsfc',
 u'pevapsfc',
 u'pli30_0mb',
 u'plpl255_0mb',
 u'porossfc',
 u'pot10m',
 u'pot30_0mb',
 u'pratesfc',
 u'pressfc',
 u'pres80m',
 u'preshlev1',
 u'pres30_0mb',
 u'pres60_30mb',
 u'pres90_60mb',
 u'pres120_90mb',
 u'pres150_120mb',
 u'pres180_150mb',
 u'presclb',
 u'presgclb',
 u'presgclt',
 u'prescclb',
 u'prescclt',
 u'presscclb',
 u'presscclt',
 u'presdcclb',
 u'presdcclt',
 u'presclt',
 u'presl5',
 u'presmwl',
 u'prestrop',
 u'prmslmsl',
 u'pwat30_0mb',
 u'pwatclm',
 u'rcqsfc',
 u'rcssfc',
 u'rcsolsfc',
 u'rctsfc',
 u'refcclm',
 u'refd1000m',
 u'refd4000m',
 u'refdhlev1',
 u'retopclm',
 u'rhprs',
 u'rh2m',
 u'rhhlev1',
 u'rh30_0mb',
 u'rh60_30mb',
 u'rh90_60mb',
 u'rh120_90mb',
 u'rh150_120mb',
 u'rh180_150mb',
 u'rh0c',
 u'rimeprs',
 u'rimehlev1',
 u'rlyrssfc',
 u'rsminsfc',
 u'rwmrprs',
 u'rwmrhlev1',
 u'sfcrsfc',
 u'sfexcsfc',
 u'shtflavesfc',
 u'shtflsfc',
 u'smdrysfc',
 u'smrefsfc',
 u'snfalbsfc',
 u'snmrprs',
 u'snmrhlev1',
 u'snodsfc',
 u'snohfsfc',
 u'snomsfc',
 u'snowcsfc',
 u'soill0_10cm',
 u'soill10_40cm',
 u'soill40_100cm',
 u'soill100_200cm',
 u'soilm0_200cm',
 u'soilw0_10cm',
 u'soilw10_40cm',
 u'soilw40_100cm',
 u'soilw100_200cm',
 u'sotypsfc',
 u'spfhsfc',
 u'spfhprs',
 u'spfh2m',
 u'spfh10m',
 u'spfh80m',
 u'spfhhlev1',
 u'spfh30_0mb',
 u'spfh60_30mb',
 u'spfh90_60mb',
 u'spfh120_90mb',
 u'spfh150_120mb',
 u'spfh180_150mb',
 u'ssrunsfc',
 u'strm250mb',
 u'swhrclm',
 u'tcdcclm',
 u'tclswclm',
 u'tcolcclm',
 u'tcoliclm',
 u'tcolmclm',
 u'tcolrclm',
 u'tcolsclm',
 u'tcolwclm',
 u'tcondhlev1',
 u'tkeprs',
 u'tkehlev1',
 u'tkehlev2',
 u'tmax2m',
 u'tmin2m',
 u'tmpsfc',
 u'tmpprs',
 u'tmp_914m',
 u'tmp_1524m',
 u'tmp_1829m',
 u'tmp_2134m',
 u'tmp_2743m',
 u'tmp_3658m',
 u'tmp_305m',
 u'tmp_457m',
 u'tmp_610m',
 u'tmp_4572m',
 u'tmp2m',
 u'tmp80m',
 u'tmphlev1',
 u'tmp30_0mb',
 u'tmp60_30mb',
 u'tmp90_60mb',
 u'tmp120_90mb',
 u'tmp150_120mb',
 u'tmp180_150mb',
 u'tmpclt',
 u'tmptrop',
 u'tsoil0_10cm',
 u'tsoil10_40cm',
 u'tsoil40_100cm',
 u'tsoil100_200cm',
 u'tsoill106',
 u'uflxsfc',
 u'ugrdprs',
 u'ugrd_914m',
 u'ugrd_1524m',
 u'ugrd_1829m',
 u'ugrd_2134m',
 u'ugrd_2743m',
 u'ugrd_3658m',
 u'ugrd_305m',
 u'ugrd_457m',
 u'ugrd_610m',
 u'ugrd_4572m',
 u'ugrd10m',
 u'ugrd80m',
 u'ugrdhlev1',
 u'ugrdhlev2',
 u'ugrd30_0mb',
 u'ugrd60_30mb',
 u'ugrd90_60mb',
 u'ugrd120_90mb',
 u'ugrd150_120mb',
 u'ugrd180_150mb',
 u'ugrdpbl',
 u'ugrdmwl',
 u'ugrdtrop',
 u'ulwrfavesfc',
 u'ulwrfavetoa',
 u'ulwrfsfc',
 u'ulwrftoa',
 u'uphl5000_2000m',
 u'ustm6000_0m',
 u'uswrfavesfc',
 u'uswrftoa',
 u'uswrfsfc',
 u'vegsfc',
 u'vflxsfc',
 u'vgrdprs',
 u'vgrd_914m',
 u'vgrd_1524m',
 u'vgrd_1829m',
 u'vgrd_2134m',
 u'vgrd_2743m',
 u'vgrd_3658m',
 u'vgrd_305m',
 u'vgrd_457m',
 u'vgrd_610m',
 u'vgrd_4572m',
 u'vgrd10m',
 u'vgrd80m',
 u'vgrdhlev1',
 u'vgrdhlev2',
 u'vgrd30_0mb',
 u'vgrd60_30mb',
 u'vgrd90_60mb',
 u'vgrd120_90mb',
 u'vgrd150_120mb',
 u'vgrd180_150mb',
 u'vgrdpbl',
 u'vgrdmwl',
 u'vgrdtrop',
 u'vgtypsfc',
 u'vissfc',
 u'vstm6000_0m',
 u'vvelprs',
 u'vvelhlev1',
 u'vvel30_0mb',
 u'vvel90_60mb',
 u'vvel180_150mb',
 u'vwshtrop',
 u'weasdaccsfc',
 u'weasdsfc',
 u'wiltsfc',
 u'wtmpsfc',
 u'var0161981000m',
 u'var02220l100_100',
 u'var02221l100_100',
 u'var0222210m',
 u'var0222310m',
 u'var02224pbl',
 u'var071995000_2000m']

In [267]:
uvar=nc.variables['ugrd10m']
vvar=nc.variables['vgrd10m']

In [268]:
def nam_subset(nc,bbox):
    lat = nc.variables['lat'][:]
    lon = nc.variables['lon'][:]
#    lon = lon - 360.0 
    ibox =(lon>=bbox[0])&(lon<=bbox[1])
    jbox =(lat>=bbox[2])&(lat<=bbox[3])
    return lon[ibox],lat[jbox],ibox,jbox

In [269]:
lon,lat,ibox,jbox = nam_subset(nc,bbox)

time_var = nc.variables['time']

In [270]:
# Extract desired time.  
start = dnow + dt.timedelta(hours=0)
istart = netCDF4.date2index(start,time_var,select='nearest')
dtime = netCDF4.num2date(time_var[istart],time_var.units)
print start.strftime('%Y/%m/%d %H:%M UTC')
print dtime.strftime('%Y/%m/%d %H:%M UTC')


2013/07/03 03:06 UTC
2013/06/15 18:00 UTC

In [271]:
# Get snapshot of variable [vname] at time = istart
u = uvar[istart,jbox,ibox]
v = vvar[istart,jbox,ibox]
spd = sqrt(u*u + v*v)
#hs1 = ma.masked_where(isnan(hs1),hs1)
# tim = dtime[istart:]

In [272]:
nsub=4 # subsample vectors
scale=0.002 # length of vectors

In [273]:
fig=figure(figsize=(12,12)) 
m = Basemap(resolution='i', llcrnrlon=bbox[0],llcrnrlat=bbox[1],
            urcrnrlon=bbox[2],urcrnrlat=bbox[3],lat_0=bbox[1],lon_0=bbox[0],lat_ts=bbox[1])
m.drawcoastlines()
#ax=fig.add_subplot(111,aspect=1.0/cos(min(lat) * pi / 180.0))
pc3=pcolormesh(lon,lat,spd, vmin=0, vmax=25)
quiver(lon[::nsub],lat[::nsub],u[::nsub,::nsub],v[::nsub,::nsub],scale=1.0/scale, zorder=1e35, width=0.002)
title('NAM Wind Speed in the Gulf of Maine: ' + dtime.strftime('%Y/%m/%d %H:%M UTC'))
#parallels = np.arange(-75.,-63,5.)
#    meridians = np.arange(-100.,-50.,5.)
#    m.drawparallels(parallels,labels=[1,0,0,0],fontsize=10);
#    m.drawmeridians(meridians,labels=[0,0,0,1],fontsize=10);
cbax3 = fig.add_axes([0.95, 0.3, 0.03, 0.4])
cb3 = plt.colorbar(pc3, cax=cbax3,  orientation='vertical')
cb3.set_label('m/s')
plt.savefig(dtime.strftime('%Y%m%dT%H%M00Z')+'.png')



In [245]: