In [31]:
import netCDF4
%matplotlib inline
from scipy.stats import binned_statistic_2d
import datetime as dt
import pandas as pd
import oceans
import matplotlib.pyplot as plt
import numpy as np
import numpy.ma as ma

In [32]:
url='http://dods.ndbc.noaa.gov/thredds/dodsC/data/stdmet/brbn4/brbn4.ncml'
nc = netCDF4.Dataset(url)
ncv = nc.variables

In [33]:
ncv.keys()


Out[33]:
[u'latitude',
 u'longitude',
 u'time',
 u'wind_dir',
 u'wind_spd',
 u'gust',
 u'wave_height',
 u'dominant_wpd',
 u'average_wpd',
 u'mean_wave_dir',
 u'air_pressure',
 u'air_temperature',
 u'sea_surface_temperature',
 u'dewpt_temperature',
 u'visibility',
 u'water_level']

In [34]:
time_var = ncv['time']
dtime = netCDF4.num2date(time_var[:],time_var.units)

In [35]:
# Extract desired times.  Here we select a specific time of interest
start = dt.datetime(2010,8,1,0,0,0)
istart = netCDF4.date2index(start,time_var,select='nearest')
stop = dt.datetime(2010,9,1,0,0,0)
istop = netCDF4.date2index(stop,time_var,select='nearest')

In [36]:
# Get all time records of variable [vname] at indices [iy,ix]
vname = 'wind_spd'
var = ncv[vname]
var.shape
wspd = var[istart:istop,:,:].ravel()

In [37]:
vname = 'wind_dir'
var = ncv[vname]
var.shape
wdir = var[istart:istop,:,:].ravel()

In [38]:
tim = dtime[istart:istop]

In [39]:
# Create Pandas time series object
ts = pd.Series(wspd,index=tim)

In [40]:
ts.plot()


Out[40]:
<matplotlib.axes._subplots.AxesSubplot at 0x7fd7bb219e90>

In [41]:
u,v = oceans.spdir2uv(wspd,wdir,deg=True)

In [42]:
fig=plt.figure(figsize=(12,4))
plt.plot(tim,u,tim,v)


Out[42]:
[<matplotlib.lines.Line2D at 0x7fd7bad56550>,
 <matplotlib.lines.Line2D at 0x7fd7bad14650>]

In [43]:
f1,x,y,bin1 = binned_statistic_2d(u, v, u, statistic='mean', bins=10, range=None)
f2,x,y,bin2 = binned_statistic_2d(u, v, v, statistic='mean', bins=10, range=None)

In [44]:
u


Out[44]:
masked_array(data = [-- 2.4492935982947064e-16 -- ..., -- -- --],
             mask = [ True False  True ...,  True  True  True],
       fill_value = 99.0)

In [65]:
url='http://tds.marine.rutgers.edu/thredds/dodsC/met/ncdc-nam-3hour/Uwind_nam_3hourly_MAB_and_GoM_2009.nc'
nc = netCDF4.Dataset(url)
ncv = nc.variables
print ncv.keys()


[u'lat', u'lon', u'time', u'Uwind']

In [84]:
url='http://tds.marine.rutgers.edu/thredds/dodsC/met/ncdc-nam-3hour/Uwind_nam_3hourly_MAB_and_GoM_2009.nc'
vname = 'Uwind'
nc = netCDF4.Dataset(url)
ncv = nc.variables
time_var = ncv['time']
dtime = netCDF4.num2date(time_var[:],time_var.units)
istart = netCDF4.date2index(start,time_var,select='nearest')
istop = netCDF4.date2index(stop,time_var,select='nearest')
var = ncv[vname]
j=78
i=69
v = var[istart:istop,j,i]
tim = dtime[istart:istop]

In [83]:
tim


Out[83]:
array([ 1308.   ,  1308.125,  1308.25 ,  1308.375,  1308.5  ,  1308.625,
        1308.75 ,  1308.875,  1309.   ,  1309.125,  1309.25 ,  1309.375,
        1309.5  ,  1309.625,  1309.75 ,  1309.875,  1310.   ,  1310.125,
        1310.25 ,  1310.375,  1310.5  ,  1310.625,  1310.75 ,  1310.875,
        1311.   ,  1311.125,  1311.25 ,  1311.375,  1311.5  ,  1311.625,
        1311.75 ,  1311.875,  1312.   ,  1312.125,  1312.25 ,  1312.375,
        1312.5  ,  1312.625,  1312.75 ,  1312.875,  1313.   ,  1313.125,
        1313.25 ,  1313.375,  1313.5  ,  1313.625,  1313.75 ,  1313.875,
        1314.   ,  1314.125,  1314.25 ,  1314.375,  1314.5  ,  1314.625,
        1314.75 ,  1314.875,  1315.   ,  1315.125,  1315.25 ,  1315.375,
        1315.5  ,  1315.625,  1315.75 ,  1315.875,  1316.   ,  1316.125,
        1316.25 ,  1316.375,  1316.5  ,  1316.625,  1316.75 ,  1316.875,
        1317.   ,  1317.125,  1317.25 ,  1317.375,  1317.5  ,  1317.625,
        1317.75 ,  1317.875,  1318.   ,  1318.125,  1318.25 ,  1318.375,
        1318.5  ,  1318.625,  1318.75 ,  1318.875,  1319.   ,  1319.125,
        1319.25 ,  1319.375,  1319.5  ,  1319.625,  1319.75 ,  1319.875,
        1320.   ,  1320.125,  1320.25 ,  1320.375,  1320.5  ,  1320.625,
        1320.75 ,  1320.875,  1321.   ,  1321.125,  1321.25 ,  1321.375,
        1321.5  ,  1321.625,  1321.75 ,  1321.875,  1322.   ,  1322.125,
        1322.25 ,  1322.375,  1322.5  ,  1322.625,  1322.75 ,  1322.875,
        1323.   ,  1323.125,  1323.25 ,  1323.375,  1323.5  ,  1323.625,
        1323.75 ,  1323.875,  1324.   ,  1324.125,  1324.25 ,  1324.375,
        1324.5  ,  1324.625,  1324.75 ,  1324.875,  1325.   ,  1325.125,
        1325.25 ,  1325.375,  1325.5  ,  1325.625,  1325.75 ,  1325.875,
        1326.   ,  1326.125,  1326.25 ,  1326.375,  1326.5  ,  1326.625,
        1326.75 ,  1326.875,  1327.   ,  1327.125,  1327.25 ,  1327.375,
        1327.5  ,  1327.625,  1327.75 ,  1327.875,  1328.   ,  1328.125,
        1328.25 ,  1328.375,  1328.5  ,  1328.625,  1328.75 ,  1328.875,
        1329.   ,  1329.125,  1329.25 ,  1329.375,  1329.5  ,  1329.625,
        1329.75 ,  1329.875,  1330.   ,  1330.125,  1330.25 ,  1330.375,
        1330.5  ,  1330.625,  1330.75 ,  1330.875,  1331.   ,  1331.125,
        1331.25 ,  1331.375,  1331.5  ,  1331.625,  1331.75 ,  1331.875,
        1332.   ,  1332.125,  1332.25 ,  1332.375,  1332.5  ,  1332.625,
        1332.75 ,  1332.875,  1333.   ,  1333.125,  1333.25 ,  1333.375,
        1333.5  ,  1333.625,  1333.75 ,  1333.875,  1334.   ,  1334.125,
        1334.25 ,  1334.375,  1334.5  ,  1334.625,  1334.75 ,  1334.875,
        1335.   ,  1335.125,  1335.25 ,  1335.375,  1335.5  ,  1335.625,
        1335.75 ,  1335.875,  1336.   ,  1336.125,  1336.25 ,  1336.375,
        1336.5  ,  1336.625,  1336.75 ,  1336.875,  1337.   ,  1337.125,
        1337.25 ,  1337.375,  1337.5  ,  1337.625,  1337.75 ,  1337.875,
        1338.   ,  1338.125,  1338.25 ,  1338.375,  1338.5  ,  1338.625,
        1338.75 ,  1338.875])

In [73]:
istart


Out[73]:
2948

In [85]:
def get_nam_ts(url,vname,start=None,stop=None,j=None,i=None):
    nc = netCDF4.Dataset(url)
    ncv = nc.variables
    time_var = ncv['time']
    dtime = netCDF4.num2date(time_var[:],time_var.units)
    istart = netCDF4.date2index(start,time_var,select='nearest')
    istop = netCDF4.date2index(stop,time_var,select='nearest')
    var = ncv[vname]
    v = var[istart:istop,j,i]
    tim = dtime[istart:istop]
    return v,tim

In [59]:
ncv.keys()


Out[59]:
[u'lat', u'lon', u'time', u'swrad']

In [53]:
time_var = ncv['time']
dtime = netCDF4.num2date(time_var[:],time_var.units)

In [156]:
# Extract desired times.  Here we select a specific time of interest
start = dt.datetime(2009,1,1,0,0,0)
stop = dt.datetime(2009,12,1,0,0,0)


url='http://tds.marine.rutgers.edu/thredds/dodsC/met/ncdc-nam-3hour/Uwind_nam_3hourly_MAB_and_GoM_2009.nc'
vname = 'Uwind'
u,tim = get_nam_ts(url,vname=vname,start=start,stop=stop,j=78,i=69)

url='http://tds.marine.rutgers.edu/thredds/dodsC/met/ncdc-nam-3hour/Vwind_nam_3hourly_MAB_and_GoM_2009.nc'
vname = 'Vwind'
v,tim = get_nam_ts(url,vname=vname,start=start,stop=stop,j=78,i=69)

fig = plt.figure(figsize=(16,4))
plt.plot(tim,u,tim,v)
plt.grid()



In [157]:
# since we are binning both u and v the bins will be the same for u and v
ubin,x,y,bin1 = binned_statistic_2d(u, v, u, statistic='mean', bins=10, range=None)
vbin,x,y,bin2 = binned_statistic_2d(u, v, v, statistic='mean', bins=10, range=None)

In [158]:
ubin = ma.masked_invalid(ubin)
vbin = ma.masked_invalid(vbin)

In [159]:
plt.pcolormesh(x,y,vbin)


Out[159]:
<matplotlib.collections.QuadMesh at 0x7fd7b522c810>

In [160]:
dirbin,spdbin = oceans.uv2spdir(ubin,vbin)
spdbin = ma.masked_invalid(spdbin)
plt.pcolormesh(x,y,spdbin)


Out[160]:
<matplotlib.collections.QuadMesh at 0x7fd7b58c3cd0>

In [184]:
dir, spd = oceans.uv2spdir(u,v)

In [185]:
ind = np.where((spd>=3.6) & (spd<=10.8))[0]

In [186]:
len(ind)*1.0/len(spd)


Out[186]:
0.7039670658682635

In [189]:
hist, bin_edges = np.histogram(dir[ind],18)
bottom = bin_edges[:-1] 
heights = np.diff(bin_edges)
plt.barh(bottom,hist,height=heights)
#plt.barh(bin_edges,hist)


Out[189]:
<Container object of 18 artists>

In [190]:
d, s = oceans.uv2spdir(0,1)
print d,s


0.0 1.0

In [181]:
d,s = oceans.uv2spdir(0,0)
print d,s


90.0 0.0

In [133]:
heights


Out[133]:
array([ 35.64590912,  35.64590912,  35.64590912,  35.64590912,
        35.64590912,  35.64590912,  35.64590912,  35.64590912,
        35.64590912,  35.64590912])

In [134]:
centers


Out[134]:
array([ 179.70557404,   19.298983  ,   54.94489212,   90.59080124,
        126.23671036,  161.88261948,  197.52852859,  233.17443771,
        268.82034683,  304.46625595])

In [135]:
hist


Out[135]:
array([26, 17, 11,  9,  7, 14,  6,  5,  3, 22])

In [ ]:
7-21