pythonmetpy

Can height rises/falls be added to the station plots in a DIFAX chart replica?


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


Solution

  • So the changes to make are:

    1. Use a timedelta(hours=12) to figure out what the right time is for the previous data
    2. Align the DataFrames using set_index('station') so that you can subtract the two height columns to get the change
    3. Plot with plot_parameter using a custom formatter to plot in dm with +/- sign

    The 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')