# 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
from pylab import *
import netCDF4
import pandas as pd
import datetime as dt
import bearcart
bearcart.initialize_notebook()
# 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
# NetCDF4-Python can open OPeNDAP dataset just like a local NetCDF file
nc = netCDF4.Dataset(url)
nc.variables.keys()
[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']
lat = nc.variables['lat'][:]
lon = nc.variables['lon'][:]
time_var = nc.variables['time']
dtime = netCDF4.num2date(time_var[:],time_var.units)
# determine what longitude convention is being used
print lon.min(),lon.max()
-99.0 -60.0
# 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
# Function to find index to nearest point
def near(array,value):
idx=(abs(array-value)).argmin()
return idx
# Find nearest point to desired location (no interpolation)
ix = near(lon, loni)
iy = near(lat, lati)
# 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
# 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')
# 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]
# Create Pandas time series object
ts = pd.Series(hs,index=tim)
# 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)
<matplotlib.text.Text at 0x38811d0>
vis = bearcart.Chart(ts)
vis
/home/local/python27_epd/lib/python2.7/site-packages/IPython/core/formatters.py:239: FormatterWarning: Exception in text/html formatter: 0.89999998 is not JSON serializable FormatterWarning,
<bearcart.bearcart.Chart at 0x38aa750>