wrf-python/test/ipynb/SLS_Example.ipynb

In [ ]:

%matplotlib inline

from __future__ import (absolute_import, division, print_function, unicode_literals)

import numpy as np
import matplotlib.pyplot as plt
from matplotlib.cm import get_cmap
import cartopy.crs as crs
import cartopy.feature as cfeature
from netCDF4 import Dataset as nc

from wrf import to_np, getvar, smooth2d

# Open the output NetCDF with netcdf-python
filename = "/Users/ladwig/Documents/wrf_files/wrfout_d01_2016-10-07_00_00_00"
ncfile = nc(filename)

# Get sea level pressure and cloud top temperature
slp = getvar(ncfile, "slp")
ctt = getvar(ncfile, "ctt")

print (ctt)

# Smooth the SLP
smooth_slp = smooth2d(slp, 3)

# Extract the latitude and longitude coordinate arrays as regular numpy array instead of xarray.DataArray
lons = ctt.coords["XLONG"]
lats = ctt.coords["XLAT"]

# Get the cartopy projection class
wrf_proj = slp.attrs["projection"]
cart_proj = wrf_proj.cartopy()

# Create the figure
fig = plt.figure(figsize=(8,8))
ax = plt.axes([0.1,0.1,0.8,0.8], projection=cart_proj)

# Download and create the states, land, and oceans using cartopy features
states = cfeature.NaturalEarthFeature(category='cultural', scale='50m', facecolor='none',
                             name='admin_1_states_provinces_shp')
land = cfeature.NaturalEarthFeature(category='physical', name='land', scale='50m', 
                                    facecolor=cfeature.COLORS['land'])
ocean = cfeature.NaturalEarthFeature(category='physical', name='ocean', scale='50m', 
                                     facecolor=cfeature.COLORS['water'])

# Make the pressure contours.
#contour_levels = [960, 962, 965, 967, 970, 975, 980, 985, 990, 995]
#c1 = plt.contour(lons, lats, to_np(smooth_slp), levels=contour_levels, colors="white", 
#            transform=crs.PlateCarree(), zorder=3, linewidth=1.3)

# Add pressure contour labels
#plt.clabel(c1, contour_levels, inline=True, fmt='%1.1f', fontsize=8)

# Create the filled cloud top temperature contours
#contour_levels = [-80, -70, -60, -50, -40, -30, -20, -10, 0, 10, 20, 30]
contour_levels = np.arange(-80.0, 10.0, 10.0)
plt.contourf(to_np(lons), to_np(lats), to_np(ctt), contour_levels, cmap=get_cmap("Greys"),
             transform=crs.PlateCarree(), zorder=2)

# Create the label bar for cloud top temperature
cb2 = plt.colorbar(ax=ax, shrink=.80)

# Draw the oceans, land, and states
ax.add_feature(ocean)
ax.add_feature(land)
ax.add_feature(states, linewidth=.5, edgecolor="black")

# Crop the domain to the region around the hurricane
ax.set_extent([-85., -75.0, np.amin(lats), 30.0], crs=crs.PlateCarree())
# Set the map limits
#ax.set_xlim(wrf_proj.cartopy_xlim())
#ax.set_ylim(wrf_proj.cartopy_ylim())

ax.gridlines(crs=crs.PlateCarree(), draw_labels=False)

# Add the title and show the image
#plt.title("Hurricane Matthew Cloud Top Temperature (degC) and Sea Level Pressure (hPa)")
#plt.title("Hurricane Matthew Cloud Top Temperature (degC)")
plt.show()

In [ ]: