In [1]:

    

# Extract a forecast time-series of wave height from Global NCEP WaveWatch 3
# and plot the time-series using Pandas

# This is just one of the many datasets in Unidata's THREDDS Catalog Collection:
# http://thredds.ucar.edu/thredds/catalog.html
# and this approach would work for any CF-Compliant NetCDF or OPeNDAP dataset


    

In [2]:

    

from pylab import *
import netCDF4
import pandas as pd
import datetime as dt


    

In [3]:

    

# NCEP WaveWatch 3 wave model: New England Coastal Grid:
url='http://thredds.ucar.edu/thredds/dodsC/grib/NCEP/WW3/Coastal_US_East_Coast/best'  # forecast
url='http://geoport.whoi.edu/thredds/dodsC/fmrc/NCEP/ww3/cfsr/4m/best'   # 30 year hindcast


    

In [4]:

    

# NetCDF4-Python can open OPeNDAP dataset just like a local NetCDF file
nc = netCDF4.Dataset(url)
nc.variables.keys()


    

    
Out[4]:





[u'lat',
 u'lon',
 u'time',
 u'time1',
 u'u-component_of_wind_surface',
 u'v-component_of_wind_surface',
 u'Significant_height_of_combined_wind_waves_and_swell_surface',
 u'Primary_wave_direction_degree_true_surface',
 u'Primary_wave_mean_period_surface']

In [5]:

    

lat = nc.variables['lat'][:]
lon = nc.variables['lon'][:]
time_var = nc.variables['time']
dtime = netCDF4.num2date(time_var[:],time_var.units)


    

In [6]:

    

# determine what longitude convention is being used
print lon.min(),lon.max()


    

    




-99.0 -60.0

In [7]:

    

# Specify desired station time series location
# note we add 360 because of the lon convention in this dataset
lati = 41.01; loni = -70.2  # Provincetown, MA


    

In [8]:

    

# Function to find index to nearest point
def near(array,value):
    idx=(abs(array-value)).argmin()
    return idx


    

In [9]:

    

# Find nearest point to desired location (no interpolation)
ix = near(lon, loni)
iy = near(lat, lati)


    

In [10]:

    

# Extract desired times.  Here we select 3 days before the present time through the end of the forecast
start = dt.datetime.utcnow()- dt.timedelta(days=3)
istart = netCDF4.date2index(start,time_var,select='nearest')
istop = len(time_var)-1


    

In [11]:

    

# Extract desired times.  Here we select a specific time of interest
start = dt.datetime(2000,1,1,0,0,0)
istart = netCDF4.date2index(start,time_var,select='nearest')
stop = dt.datetime(2000,1,7,0,0,0)
istop = netCDF4.date2index(stop,time_var,select='nearest')


    

In [12]:

    

# Get all time records of variable [vname] at indices [iy,ix]
vname = 'Significant_height_of_combined_wind_waves_and_swell_surface'
var = nc.variables[vname]
hs = var[istart:istop,iy,ix]
tim = dtime[istart:istop]


    

In [13]:

    

# Create Pandas time series object
ts = pd.Series(hs,index=tim)


    

In [14]:

    

# Use Pandas time series plot method
figure(figsize=(16,4))
ts.plot(title='%s (Lon=%.2f, Lat=%.2f)' % (vname, lon[ix], lat[iy]))
ylabel(var.units)


    

    
Out[14]:





<matplotlib.text.Text at 0x48ae850>






    

In [14]:

    




    

In [14]: