NCL_dataonmap_10.py

This script illustrates the following concepts:

• Plotting WRF data on native grid

• Plotting data using wrf python functions

• Overlaying continent outlines on a map

• Following best practices when choosing a colormap. More information on colormap best practices can be found here.

See following URLs to see the reproduced NCL plot & script:

• Original NCL script: https://www.ncl.ucar.edu/Applications/Scripts/dataonmap_10.ncl

• Original NCL plot: https://www.ncl.ucar.edu/Applications/Images/dataonmap_10_1_lg.png

Import packages

from netCDF4 import Dataset
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.ticker as mticker
import cartopy.crs as ccrs
from wrf import (getvar, to_np, latlon_coords)

import geocat.datafiles as gdf
import geocat.viz as gv

Read in the data

wrfin = Dataset(gdf.get("netcdf_files/wrfout_d01_2003-07-15_00_00_00"),
                decode_times=True)
q2 = getvar(wrfin, "Q2")

Plot the data

# Get the latitude and longitude coordinate. This is usually needed for plotting.
lats, lons = latlon_coords(q2)

# Generate figure (set its size (width, height) in inches)
fig = plt.figure(figsize=(10, 10))

# Generate axes using Cartopy
ax = plt.axes(projection=ccrs.PlateCarree())

# Add filled contours
plt.contourf(to_np(lons),
             to_np(lats),
             q2,
             levels=np.linspace(0.01125, 0.05, 32),
             cmap="magma",
             vmin=0,
             vmax=0.05,
             zorder=4)

# Add a colorbar
cbar = plt.colorbar(ax=ax,
                    orientation="vertical",
                    ticks=np.arange(0.0125, 0.0476, 0.0025),
                    drawedges=True,
                    extendrect=True,
                    shrink=0.65)

# Format colorbar ticks and labels
cbar.ax.tick_params(size=0, labelsize=10)

# Draw gridlines
gl = ax.gridlines(crs=ccrs.PlateCarree(),
                  draw_labels=True,
                  dms=False,
                  x_inline=False,
                  y_inline=False,
                  linewidth=1,
                  color="k",
                  alpha=0.25,
                  zorder=4)

# Manipulate latitude and longitude gridline numbers and spacing
gl.top_labels = False
gl.right_labels = False
gl.xlocator = mticker.FixedLocator(np.arange(-105, -80, 5))
gl.ylocator = mticker.FixedLocator(np.arange(18, 35, 2))
gl.xlabel_style = {"rotation": 0, "size": 10}
gl.ylabel_style = {"rotation": 0, "size": 10}
gl.xlines = True
gl.ylines = True

# Add titles and labels to projection
gv.set_titles_and_labels(ax,
                         maintitle="WRF data on native grid",
                         lefttitle="QV at 2 M",
                         maintitlefontsize=16,
                         lefttitlefontsize=14)

plt.show()