Basemap: contourf and quiver

Date: May 24th 2016
Last updated: May 24th 2016

Useful resources

• http://matplotlib.1069221.n5.nabble.com/plot-arrows-for-wind-direction-degrees-td13499.html
• http://earthpy.org/earthpy-basemap-amazon.html

plot using basemap
Assumes that a netCDF4 dataset has been imported and assigned to a variable named "wave_dataset".

# set plot size
fig=plt.figure(figsize=(5,11))

# collect data for East NSW
# drop the last 80 values 
lat  = wave_dataset.variables['lat'][:-80]
# drop the first 110 values and 
# dropping 120 from the end of the array
lon  = wave_dataset.variables['lon'][110:-120]
# use .shape function to check that arrays have 
# the correct size.
# e.g. lon.shape

# slice netCDF4 dataset
# Note both datasets come from wave_dataset
# subset swell direction dataset
wavedir = wave_dataset.variables['dirpwsfc'][1,:-80,110:-120]
# subset swell period dataset
swellperiod = wave_dataset.variables['perpwsfc'][1,:-80,110:-120]

# create basemap
# Note the projection - mercator
m = Basemap(llcrnrlon=145,llcrnrlat=-50,\
        urcrnrlon=160,urcrnrlat=-10,\
        projection='merc', resolution='i')

# add navigation features to basemap
m.drawcoastlines(linewidth=1.25)
m.fillcontinents(color='0.8')

# create numpy mesh grid
x, y = m(*np.meshgrid(lon,lat))

# flip the direction of arrows using minus 1
# adjust arrow size by converting the scale of 
# each arrow based on swellperiod
u = 1-(swellperiod**2/100)*np.sin(np.deg2rad(wavedir))
v = 1-(swellperiod**2/100)*np.cos(np.deg2rad(wavedir))

# transform variables to basemap projection
uproj,vproj,xx,yy = m.transform_vector(\
        u,v,lon,lat,5,10,returnxy=True,masked=True)

# add a filled gradient using countourf 
cs = m.contourf(x,y,swellperiod,
        # Note that levels changes the blending 
        # of the contour fill
        levels=np.arange(1,20,0.75),
        cmap=cm.coolwarm,
        vmin=1,
        vmax=20.5 )

# add a legend for the contour gradient
cbar = m.colorbar(cs,
        location='right',pad="5%")

# add arrows using quiver
Q = m.quiver(xx,yy,uproj,vproj,
        # changing the width changes the boldness of
        # each arrow. Pivot at tip so arrow head is
        # over water, not land.
        scale=10,width=0.005,pivot='tip')


# save with timestamp
date = datetime.now().strftime('%Y%m%d %H%M')
plt.title('Swell direction and\n wave period (%s)' % date)
plt.savefig('img/wavedirection_{}.png'.format(date))

Output