First, I would like to thank Kevin Goebbert and the MetPy team for providing this very helpful Python code for creating a DIFAX chart replica: https://unidata.github.io/python-gallery/examples/Upperair_Obs.html. I just have a quick question regarding this code: Can values of height rises/falls be added to these DIFAX replicas, similar to the legacy products themselves? If so, how would I go about doing that? Here is an example DIFAX chart to illustrate what I am referring to: http://archive.atmos.colostate.edu/data/misc/QHTA11/0701/07013113QHTA11.png. The height rises/falls are displayed for some radiosonde stations in the bottom right-hand corner of the station plots.
Thanks for your time and assistance!
Jordan McLeod
So the changes to make are:
timedelta(hours=12)
to figure out what the right time is for the previous dataDataFrame
s using set_index('station')
so that you can subtract the two height
columns to get the changeplot_parameter
using a custom formatter to plot in dm with +/- signThe whole product for the station plotting then works out as
from datetime import datetime, timedelta
import cartopy.crs as ccrs
import cartopy.feature as cfeature
import matplotlib.pyplot as plt
from metpy.io import add_station_lat_lon
from metpy.plots import StationPlot
from siphon.simplewebservice.iastate import IAStateUpperAir
level = 500
# Set date for desired UPA data
today = datetime.utcnow()
# Go back one day to ensure data availability
date = datetime(today.year, today.month, today.day, 0) - timedelta(days=1)
# Previous is 12 hours before that time
prev = date - timedelta(hours=12)
# Request data using Siphon request for data from Iowa State Archive
data = IAStateUpperAir.request_all_data(date)
prev_data = IAStateUpperAir.request_all_data(prev)
# Given our data, add station information and drop any stations that are missing lat/lon. Then set
# 'station' as the index column so that operations between datasets will align on that column
data = add_station_lat_lon(data).dropna(how='any', subset=('longitude', 'latitude')).set_index('station')
prev_data = add_station_lat_lon(prev_data).dropna(how='any', subset=('longitude', 'latitude')).set_index('station')
# Create subset of all data for a given level
df = data[data.pressure == level].copy()
prev_df = prev_data[prev_data.pressure == level]
# Calculate the change on the aligned data frames--needs to be
# added to original dataframe to maintain consistent ordering
# with later use
df['height_change'] = df.height - prev_df.height
# Set up map coordinate reference system
mapcrs = ccrs.LambertConformal(
central_latitude=45, central_longitude=-100, standard_parallels=(30, 60))
# Set up station locations for plotting observations
point_locs = mapcrs.transform_points(
ccrs.PlateCarree(), df['longitude'].values, df['latitude'].values)
# Start figure and set graphics extent
fig = plt.figure(1, figsize=(17, 15))
ax = plt.subplot(111, projection=mapcrs)
ax.set_extent([-125, -70, 20, 55])
# Add map features for geographic reference
ax.add_feature(cfeature.COASTLINE.with_scale('50m'), edgecolor='grey')
ax.add_feature(cfeature.LAND.with_scale('50m'), facecolor='white')
ax.add_feature(cfeature.STATES.with_scale('50m'), edgecolor='grey')
# Start the station plot by specifying the axes to draw on, as well as the
# lon/lat of the stations (with transform). We also the fontsize to 10 pt.
stationplot = StationPlot(ax, df['longitude'].values, df['latitude'].values, clip_on=True,
transform=ccrs.PlateCarree(), fontsize=10)
# Plot the temperature and dew point to the upper and lower left, respectively, of
# the center point.
stationplot.plot_parameter('NW', df['temperature'], color='black')
stationplot.plot_parameter('SW', df['dewpoint'], color='black')
# A more complex example uses a custom formatter to control how the geopotential height
# values are plotted. This is set in an earlier if-statement to work appropriate for
# different levels.
def hght_format(v):
return format(v, '.0f')[:3]
stationplot.plot_parameter('NE', df['height'], formatter=hght_format)
# Add wind barbs
stationplot.plot_barb(df['u_wind'], df['v_wind'], length=7, pivot='tip')
# Add height falls. This converts the change to decameters with /10, then
# formats with a +/- and 0 padding
def height_fall_formatter(v):
return f'{int(v / 10):+03d}'
# Plot the parameter with an italic font
stationplot.plot_parameter('SE', df['height_change'], formatter=height_fall_formatter,
fontstyle='italic')