Wave Activity Horizontal Flux

Rossby waves are like giant atmospheric ripples that guide the movement of weather systems across the midlatitudes. When these waves amplify or propagate downstream, they can help explain why blocking highs, troughs, and persistent cold or warm events appear in preferred regions.

Wave activity flux (WAF) offers a compact way to diagnose both the strength and the direction of that propagation. In practice, the flux vectors behave a bit like arrows showing where stationary wave energy is traveling.

This example demonstrates two commonly used horizontal WAF diagnostics:

In all examples below, the geopotential height anomaly must be provided in meters, not in \(\mathrm{m^2\ s^{-2}}\).

See also

First, import the plotting and analysis packages. We also configure a scientific-notation formatter so the colorbar labels remain clean when plotting large or small atmospheric quantities.

import easyclimate as ecl
import cartopy.crs as ccrs
import cmaps
import matplotlib.pyplot as plt
import matplotlib.ticker as ticker
import numpy as np

formatter = ticker.ScalarFormatter(useMathText=True, useOffset=True)
formatter.set_powerlimits((0, 0))

Next, load the anomaly field and the climatological basic-state winds.

  • z500_prime_data1 and z500_prime_data2 are 500 hPa geopotential height anomalies for two different Novembers.

  • u500_climatology_data and v500_climatology_data provide the monthly climatological background flow required by the TN01 formulation.

  • We divide geopotential height by 9.8 so the anomaly is expressed in meters.

Sorting by latitude keeps the map grids ordered from south to north before plotting.

z500_prime_data1 = ecl.open_tutorial_dataset("era5_daily_z500_prime_201411_N15.nc").z /9.8
z500_prime_data2 = ecl.open_tutorial_dataset("era5_daily_z500_prime_202411_N15.nc").z /9.8
u500_climatology_data = ecl.open_tutorial_dataset("era5_ymean_monthly_u500_199101_202012_N15.nc").u.sortby("lat")
v500_climatology_data = ecl.open_tutorial_dataset("era5_ymean_monthly_v500_199101_202012_N15.nc").v.sortby("lat")
z500_prime_data1 = z500_prime_data1.sortby("lat")
z500_prime_data2 = z500_prime_data2.sortby("lat")
era5_daily_z500_prime_201411_N15.nc ━━━━━━━ 100.0% • 229.2/… • 9.6     • 0:00:00
                                                     kB        MB/s
era5_daily_z500_prime_202411_N15.nc ━━━━━━━ 100.0% • 227.9/… • 14.6    • 0:00:00
                                                     kB        MB/s
era5_ymean_monthly_u500_199101_202012_N15.nc ━━━━ 100.0% • 125.… • 10.5 • 0:00:…
                                                           kB      MB/s
era5_ymean_monthly_v500_199101_202012_N15.nc ━━━━ 100.0% • 126.… • 8.2  • 0:00:…
                                                           kB      MB/s

TN01 Wave Activity Horizontal Flux

We start with the Takaya-Nakamura horizontal wave activity flux.

The function easyclimate.calc_TN_wave_activity_horizontal_flux returns an xarray.Dataset containing:

  • psi_p: perturbation streamfunction,

  • fx: zonal component of the flux,

  • fy: meridional component of the flux.

Here we diagnose the November 2014 case and then average over time to obtain the monthly-mean propagation pattern.

tn01_result1 = ecl.calc_TN_wave_activity_horizontal_flux(
    z_prime_data = z500_prime_data1,
    u_climatology_data = u500_climatology_data,
    v_climatology_data = v500_climatology_data,
    vertical_dim = "level",
    vertical_dim_units = "hPa",
)
tn01_result1
<xarray.Dataset> Size: 649kB
Dimensions:  (time: 30, lat: 30, lon: 60)
Coordinates:
  * time     (time) datetime64[ns] 240B 2014-11-01T09:00:00 ... 2014-11-30T09...
  * lat      (lat) float64 240B -85.48 -79.63 -73.74 ... 73.74 79.63 85.48
  * lon      (lon) float64 480B 0.0 6.0 12.0 18.0 ... 336.0 342.0 348.0 354.0
    level    int32 4B 500
Data variables:
    psi_p    (time, lat, lon) float32 216kB 1.763e+06 1.389e+06 ... 1.808e+07
    fx       (time, lat, lon) float32 216kB -0.1388 -0.1091 ... 20.95 23.59
    fy       (time, lat, lon) float32 216kB -0.06418 -0.09851 ... 14.58 20.23
Attributes:
    title:         Takaya and Nakamura Wave Activity Flux
    description:   Horizontal components of quasi-geostrophic wave activity f...
    reference:     Takaya, K. and H. Nakamura, 2001: A Formulation of a Phase...
    created_with:  easyclimate: calc_TN_wave_activity_horizontal_flux function


<xarray.Dataset> Size: 22kB
Dimensions:  (lat: 30, lon: 60)
Coordinates:
  * lat      (lat) float64 240B -85.48 -79.63 -73.74 ... 73.74 79.63 85.48
  * lon      (lon) float64 480B 0.0 6.0 12.0 18.0 ... 336.0 342.0 348.0 354.0
    level    int32 4B 500
Data variables:
    psi_p    (lat, lon) float32 7kB 1.838e+06 1.889e+06 ... -1.599e+06
    fx       (lat, lon) float32 7kB 0.276 0.4461 0.4728 ... 2.661 2.836 2.983
    fy       (lat, lon) float32 7kB -0.3099 -0.3664 -0.4537 ... 3.762 3.713
Attributes:
    title:         Takaya and Nakamura Wave Activity Flux
    description:   Horizontal components of quasi-geostrophic wave activity f...
    reference:     Takaya, K. and H. Nakamura, 2001: A Formulation of a Phase...
    created_with:  easyclimate: calc_TN_wave_activity_horizontal_flux function


To make a continuous global contour map, we add a cyclic longitude point to the scalar field. The vector field is kept on the original grid and will be thinned during plotting.

This first figure uses a polar stereographic projection, which is especially helpful for viewing Rossby-wave trains across the Northern Hemisphere.

Shading shows the perturbation streamfunction, while the arrows show the horizontal WAF vectors. Larger arrows indicate stronger downstream propagation of wave activity.

fig, ax = plt.subplots(
    figsize = (6, 6),
    subplot_kw={"projection": ccrs.NorthPolarStereo(central_longitude=140)}
)

ax.coastlines(edgecolor="black", linewidths=0.5)
gl, meta = ecl.plot.draw_polar_basemap(
    ax = ax,
    lon_step=30,
    lat_step=20,
    lat_range=[10, 90],
    draw_labels=True,
    gridlines_kwargs={"color": "grey", "alpha": 0.5, "linestyle": "--"},
    lat_label_lon=-60
)

fg1 = draw_shaded1.plot.contourf(
    ax = ax,
    levels = np.linspace(-1.5e7, 1.5e7, 21),
    cbar_kwargs = {'location': 'bottom', 'aspect': 40, 'format': formatter},
    transform=ccrs.PlateCarree(),
    cmap = "RdBu_r",
    zorder = 1,
)

q = draw_quiver1.sel(lat = slice(10, None)).plot.quiver(
    ax = ax,
    x = "lon", y = "lat", u = "fx", v = "fy",
    transform = ccrs.PlateCarree(),
    scale = 600,
    regrid_shape = 18,
    headwidth = 5,
    width = 0.0045,
    headlength = 8,
    headaxislength = 8,
    minlength = 3,
    add_guide = False,
    zorder = 2,
)
qk = ax.quiverkey(
    q, 0, -0.1,
    50, "50", labelpos = "N",
    color = "k",
    coordinates = "axes",
    fontproperties = {"size": 12},
    labelsep = 0.05,
    transform = ccrs.PlateCarree(),
)
qk.set_zorder(3)

ax.set_title("")
ecl.plot.set_polar_title("TN01 WAF Analysis (500hPa, Nov. 2014)", meta, size = 15)
plot waf
Text(0.5, 1.0864137921696746, 'TN01 WAF Analysis (500hPa, Nov. 2014)')

The same calculation can be repeated for a different year. Because the TN01 scheme automatically matches each time step to the corresponding climatological month, it is convenient for composite and event analyses.

tn01_result2 = ecl.calc_TN_wave_activity_horizontal_flux(
    z_prime_data = z500_prime_data2,
    u_climatology_data = u500_climatology_data,
    v_climatology_data = v500_climatology_data,
    vertical_dim = "level",
    vertical_dim_units = "hPa",
)
tn01_mean_result2 = tn01_result2.mean(dim = "time").sortby("lat")
tn01_mean_result2
<xarray.Dataset> Size: 22kB
Dimensions:  (lat: 30, lon: 60)
Coordinates:
  * lat      (lat) float64 240B -85.48 -79.63 -73.74 ... 73.74 79.63 85.48
  * lon      (lon) float64 480B 0.0 6.0 12.0 18.0 ... 336.0 342.0 348.0 354.0
    level    int32 4B 500
Data variables:
    psi_p    (lat, lon) float32 7kB 5.594e+06 5.615e+06 ... 7.892e+05 7.138e+05
    fx       (lat, lon) float32 7kB -0.01532 0.1764 0.399 ... 1.667 1.663 1.561
    fy       (lat, lon) float32 7kB -0.3833 -0.3339 -0.3641 ... 2.298 2.0 1.417
Attributes:
    title:         Takaya and Nakamura Wave Activity Flux
    description:   Horizontal components of quasi-geostrophic wave activity f...
    reference:     Takaya, K. and H. Nakamura, 2001: A Formulation of a Phase...
    created_with:  easyclimate: calc_TN_wave_activity_horizontal_flux function


This second TN01 map uses a Mercator projection to emphasize the downstream extension of the wave train over Eurasia and the North Pacific.

fig, ax = plt.subplots(
    figsize = (12, 6),
    subplot_kw={"projection": ccrs.Mercator(central_longitude=140)}
)
ax.set_extent([0, 360, 10, 85], crs = ccrs.PlateCarree())
ax.coastlines(edgecolor="black", linewidths=0.5)
ax.gridlines(draw_labels=["bottom", "left"], alpha=0)

fg1 = draw_shaded2.plot.contourf(
    ax = ax,
    levels = np.linspace(-1e7, 1e7, 21),
    cbar_kwargs = {'location': 'bottom', 'aspect': 80, 'shrink': 0.7, 'pad': 0.1, 'format': formatter},
    transform=ccrs.PlateCarree(),
    cmap = "RdBu_r",
    zorder = 1,
)

q = draw_quiver2.sel(lat = slice(5, None)).plot.quiver(
    ax = ax,
    x = "lon", y = "lat", u = "fx", v = "fy",
    transform = ccrs.PlateCarree(),
    scale = 800,
    regrid_shape = 15,
    add_guide = False,
)
qk = ax.quiverkey(
    q, 0.95, 1.03,
    50, "50", labelpos = "N",
    color = "k",
    coordinates = "axes",
    fontproperties = {"size": 12},
    labelsep = 0.05,
    transform = ccrs.PlateCarree(),
)

ax.set_title("TN01 WAF Analysis (500hPa, Nov. 2024)")
TN01 WAF Analysis (500hPa, Nov. 2024)
Text(0.5, 1.0, 'TN01 WAF Analysis (500hPa, Nov. 2024)')

Plumb Wave Activity Horizontal Flux

Now switch to the Plumb horizontal flux.

easyclimate.calc_Plumb_wave_activity_horizontal_flux uses only the geopotential height anomaly field and is especially suitable for diagnosing stationary wave propagation.

pwaf_result = ecl.calc_Plumb_wave_activity_horizontal_flux(
    z_prime_data = z500_prime_data2,
    vertical_dim = "level",
    vertical_dim_units = "hPa",
)
pwaf_mean_result = pwaf_result.mean(dim = "time").sortby("lat")
pwaf_mean_result
<xarray.Dataset> Size: 22kB
Dimensions:  (lat: 30, lon: 60)
Coordinates:
  * lat      (lat) float64 240B -85.48 -79.63 -73.74 ... 73.74 79.63 85.48
  * lon      (lon) float64 480B 0.0 6.0 12.0 18.0 ... 336.0 342.0 348.0 354.0
    level    int32 4B 500
Data variables:
    psi_p    (lat, lon) float32 7kB 5.594e+06 5.615e+06 ... 7.892e+05 7.138e+05
    fx       (lat, lon) float32 7kB 0.06288 0.1975 0.3501 ... 2.258 2.289 2.281
    fy       (lat, lon) float32 7kB -0.0001407 -0.205 -0.4062 ... -0.1934 -0.613
Attributes:
    title:         Plumb Wave Activity Flux
    description:   Horizontal components of quasi-geostrophic wave activity f...
    reference:     Plumb, R. A., 1985: On the Three-Dimensional Propagation o...
    created_with:  easyclimate: calc_Plumb_wave_activity_horizontal_flux func...


The plotting procedure is identical, which makes comparison between the TN01 and Plumb diagnostics straightforward.

If the two flux patterns are broadly similar, that often indicates a robust stationary Rossby-wave signal. Differences can reflect the role of the background climatological flow included in the TN01 formulation.

fig, ax = plt.subplots(
    figsize = (12, 6),
    subplot_kw={"projection": ccrs.Mercator(central_longitude=140)}
)
ax.set_extent([0, 360, 10, 85], crs = ccrs.PlateCarree())
ax.coastlines(edgecolor="black", linewidths=0.5)
ax.gridlines(draw_labels=["bottom", "left"], alpha=0)

fg1 = draw_shaded3.plot.contourf(
    ax = ax,
    levels = np.linspace(-1e7, 1e7, 21),
    cbar_kwargs = {'location': 'bottom', 'aspect': 80, 'shrink': 0.7, 'pad': 0.1, 'format': formatter},
    transform=ccrs.PlateCarree(),
    cmap = "RdBu_r",
    zorder = 1,
)

q = draw_quiver3.sel(lat = slice(5, None)).plot.quiver(
    ax = ax,
    x = "lon", y = "lat", u = "fx", v = "fy",
    transform = ccrs.PlateCarree(),
    scale = 800,
    regrid_shape = 15,
    add_guide = False,
)
qk = ax.quiverkey(
    q, 0.95, 1.03,
    50, "50", labelpos = "N",
    color = "k",
    coordinates = "axes",
    fontproperties = {"size": 12},
    labelsep = 0.05,
    transform = ccrs.PlateCarree(),
)

ax.set_title("Plumb WAF Analysis (500hPa, Nov. 2024)")
Plumb WAF Analysis (500hPa, Nov. 2024)
Text(0.5, 1.0, 'Plumb WAF Analysis (500hPa, Nov. 2024)')

These figures provide a compact visual summary of wave propagation:

  • Shading highlights the perturbation streamfunction pattern associated with the geopotential height anomaly.

  • Vectors indicate the direction and relative strength of downstream wave activity propagation.

  • TN01 flux is generally preferred when a spatially varying background flow is important.

  • Plumb flux is a classic diagnostic for large-scale stationary Rossby waves.

Total running time of the script: (0 minutes 5.959 seconds)