import os
import datetime

print(f"Last updated on {datetime.date.today()}. __falwa__.version: {__import__('falwa').__version__}")

Last updated on 2026-03-15. __falwa__.version: 2.3.3

Using barotropic_qglat_lwa to compute equivalent latitude and local wave activity (total, cyclonic and anticyclonic)

Instructions

This sample code demonstrate how the wrapper function “barotropic_eqlat_lwa” in thepython package “hn2016_falwa” computes the finite-amplitude local wave activity (LWA) from absolute vorticity fields in a barotropic model with spherical geometry according to the definition in Huang & Nakamura (2016,JAS) equation (13). This sample code reproduces the LWA plots (Fig.4 in HN15) computed based on an absolute vorticity map.

Contact

Please make inquiries and report issues via Github: https://github.com/csyhuang/hn2016_falwa/issues

from falwa.wrapper import barotropic_eqlat_lwa # Module for plotting local wave activity (LWA) plots and
                        # the corresponding equivalent-latitude profile
from math import pi
import xarray as xr
import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline

# --- Parameters --- #
Earth_radius = 6.378e+6 # Earth's radius

# --- Load the absolute vorticity field [256x512] --- #
data_folder_path = os.getcwd() + "../../../tests/data/"
readFile = xr.open_dataset(f'{data_folder_path}barotropic_vorticity.nc', engine='netcdf4')


# --- Read in longitude and latitude arrays --- #
xlon = readFile.longitude.values
ylat = readFile.latitude.values
clat = np.abs(np.cos(ylat*pi/180.)) # cosine latitude
nlon = xlon.size
nlat = ylat.size

# --- Parameters needed to use the module HN2015_LWA --- #
dphi = (ylat[2]-ylat[1])*pi/180. # Equal spacing between latitude grid points, in radian
area = 2.*pi*Earth_radius**2 *(np.cos(ylat[:,np.newaxis]*pi/180.)*dphi)/float(nlon) * np.ones((nlat,nlon))
area = np.abs(area) # To make sure area element is always positive (given floating point errors).

# --- Read in the absolute vorticity field from the netCDF file --- #
absVorticity = readFile.absolute_vorticity.values

Obtain equivalent-latitude relationship and also the LWA from an absolute vorticity snapshot

# --- Obtain equivalent-latitude relationship and also the LWA from the absolute vorticity snapshot ---
qref, lwa = barotropic_eqlat_lwa(ylat,absVorticity,area,Earth_radius*clat*dphi,nlat) # Full domain included

Plotting the results

# --- Color axis for plotting LWA --- #
lwa_caxis = np.linspace(0,lwa.max(),31,endpoint=True)

# --- Plot the abs. vorticity field, LWA and equivalent-latitude relationship and LWA --- #
fig, (ax1, ax2, ax3) = plt.subplots(1, 3, figsize=(14,4))

# Absolute vorticity map
c = ax1.contourf(xlon,ylat,absVorticity,31)
cb = plt.colorbar(c)
cb.formatter.set_powerlimits((0, 0))
cb.ax.yaxis.set_offset_position('right')
cb.update_ticks()
ax1.set_title('Absolute vorticity [1/s]')
ax1.set_xlabel('Longitude (degree)')
ax1.set_ylabel('Latitude (degree)')

# LWA (full domain)
c2 = ax2.contourf(xlon,ylat,lwa,lwa_caxis)
plt.colorbar(c2)
ax2.set_title('Local Wave Activity [m/s]')
ax2.set_xlabel('Longitude (degree)')
ax2.set_ylabel('Latitude (degree)')

# Equivalent-latitude relationship Q(y)
ax3.plot(qref, ylat, 'b', label='Equivalent-latitude relationship')
ax3.plot(np.mean(absVorticity,axis=1),ylat,'g',label='zonal mean abs. vorticity')
plt.ticklabel_format(style='sci', axis='x', scilimits=(0,0))
plt.ylim(-90,90)
plt.legend(loc=4,fontsize=10)
plt.title('Equivalent-latitude profile')
plt.ylabel('Latitude (degree)')
plt.xlabel('Q(y) [1/s] | y = latitude')
plt.tight_layout()
plt.show()

Compute local wave activity partitioned into anticyclonic and cyclonic components

This can be done by setting the input parameter return_partitioned_lwa=True

qref, lwa_partitioned = barotropic_eqlat_lwa(ylat,absVorticity,area,Earth_radius*clat*dphi,nlat, return_partitioned_lwa=True) # Full domain included

fig2, (ax4, ax5, ax6) = plt.subplots(1, 3, figsize=(14,4))
c4 = ax4.contourf(xlon,ylat,lwa_partitioned.sum(axis=0),lwa_caxis)
plt.colorbar(c4)
ax4.set_title('Local Wave Activity [m/s]')
ax4.set_xlabel('Longitude (degree)')
ax4.set_ylabel('Latitude (degree)')

# Anti-cyclonic LWA (full domain)
antycyclonic_lwa = np.concatenate((lwa_partitioned[1, :128, :], lwa_partitioned[0, -128:, :]), axis=0)
c5 = ax5.contourf(xlon,ylat,antycyclonic_lwa,lwa_caxis)
plt.colorbar(c5)
ax5.set_title('Local Wave Activity (anticyclonic) [m/s]')
ax5.set_xlabel('Longitude (degree)')
ax5.set_ylabel('Latitude (degree)')

# Cyclonic LWA (full domain)
cyclonic_lwa = np.concatenate((lwa_partitioned[0, :128, :], lwa_partitioned[1, -128:, :]), axis=0)
c6 = ax6.contourf(xlon,ylat,cyclonic_lwa,lwa_caxis)
plt.colorbar(c6)
ax6.set_title('Local Wave Activity (cyclonic) [m/s]')
ax6.set_xlabel('Longitude (degree)')
ax6.set_ylabel('Latitude (degree)')

plt.tight_layout()
plt.show()