metpy

name	metpy
description	Calculate derived atmospheric quantities using metpy.calc — thermodynamics, soundings, kinematics, boundary layer, smoothing, and more.
user-invocable	true

Help the user calculate derived atmospheric quantities with MetPy. Reference: https://unidata.github.io/MetPy/latest/api/generated/metpy.calc.html

Installation

conda install -c conda-forge metpy

Units — Critical Concept

MetPy requires quantities to carry units via pint. Always attach units before calling metpy.calc functions.

import metpy.calc as mpcalc
from metpy.units import units
import numpy as np

T = np.array([20., 15., 10.]) * units.degC
p = np.array([1000., 850., 700.]) * units.hPa
rh = np.array([80., 70., 60.]) * units.percent

With xarray, use .metpy.quantify() to attach units from CF attributes:

import xarray as xr
ds = xr.open_dataset("file.nc").metpy.quantify()
T = ds["temperature"]

---

Dry Thermodynamics

# Potential temperature
theta = mpcalc.potential_temperature(p, T)

# Temperature from potential temperature
T = mpcalc.temperature_from_potential_temperature(p, theta)

# Dry lapse rate (parcel)
T_parcel = mpcalc.dry_lapse(p, T[0])

# Air density
rho = mpcalc.density(p, T, mixing_ratio)

# Dry/moist static energy
dse = mpcalc.dry_static_energy(height, T)
mse = mpcalc.moist_static_energy(height, T, specific_humidity)

# Layer thickness (hypsometric)
dz = mpcalc.thickness_hydrostatic(p, T)

# Geopotential <-> height
z = mpcalc.geopotential_to_height(geopotential)
phi = mpcalc.height_to_geopotential(z)

# Pressure weighted mean
T_mean = mpcalc.mean_pressure_weighted(p, T, height=z, bottom=z[0], depth=100*units.hPa)

---

Moist Thermodynamics

# Saturation vapor pressure
es = mpcalc.saturation_vapor_pressure(T)

# Dewpoint
Td = mpcalc.dewpoint(vapor_pressure)
Td = mpcalc.dewpoint_from_relative_humidity(T, rh)
Td = mpcalc.dewpoint_from_specific_humidity(specific_humidity, p)

# Relative humidity
rh = mpcalc.relative_humidity_from_dewpoint(T, Td)
rh = mpcalc.relative_humidity_from_mixing_ratio(mixing_ratio, T, p)
rh = mpcalc.relative_humidity_from_specific_humidity(specific_humidity, T, p)

# Mixing ratio
w  = mpcalc.mixing_ratio(vapor_pressure, p)
w  = mpcalc.mixing_ratio_from_relative_humidity(T, p, rh)
w  = mpcalc.saturation_mixing_ratio(p, T)

# Specific humidity
q  = mpcalc.specific_humidity_from_mixing_ratio(w)
q  = mpcalc.specific_humidity_from_dewpoint(p, Td)

# Equivalent / virtual potential temperature
theta_e  = mpcalc.equivalent_potential_temperature(p, T, Td)
theta_es = mpcalc.saturation_equivalent_potential_temperature(p, T, Td)
theta_v  = mpcalc.virtual_potential_temperature(p, T, w)

# Virtual temperature
Tv = mpcalc.virtual_temperature(T, w)
Tv = mpcalc.virtual_temperature_from_dewpoint(p, T, Td)

# Wet bulb temperature
Tw = mpcalc.wet_bulb_temperature(p, T, Td)

# Precipitable water
pw = mpcalc.precipitable_water(p, Td)

# Vertical velocity conversion
w  = mpcalc.vertical_velocity(omega, p, T)           # omega -> w (m/s)
om = mpcalc.vertical_velocity_pressure(w, p, T)       # w -> omega (Pa/s)

---

Soundings & Stability Indices

# Parcel profile
prof = mpcalc.parcel_profile(p, T[0], Td[0])
prof, p_lcl, T_lcl = mpcalc.parcel_profile_with_lcl(p, T, Td)

# Key levels
lcl_p, lcl_T = mpcalc.lcl(p[0], T[0], Td[0])
lfc_p, lfc_T = mpcalc.lfc(p, T, Td)
el_p,  el_T  = mpcalc.el(p, T, Td)

# CAPE / CIN
cape, cin = mpcalc.cape_cin(p, T, Td, prof)
cape, cin = mpcalc.surface_based_cape_cin(p, T, Td)
cape, cin = mpcalc.mixed_layer_cape_cin(p, T, Td)
cape, cin = mpcalc.most_unstable_cape_cin(p, T, Td)
dcape     = mpcalc.downdraft_cape(p, T, Td)

# Mixed layer parcel
p_ml, T_ml, Td_ml = mpcalc.mixed_parcel(p, T, Td)
T_ml = mpcalc.mixed_layer(p, T)

# Stability indices
k   = mpcalc.k_index(p, T, Td)
tt  = mpcalc.total_totals_index(p, T, Td)
ct  = mpcalc.cross_totals(p, T, Td)
vt  = mpcalc.vertical_totals(p, T)
sw  = mpcalc.sweat_index(p, T, Td, u, v)
ssi = mpcalc.showalter_index(p, T, Td)
li  = mpcalc.lifted_index(p, T, prof)

# Wind / shear
u, v   = mpcalc.wind_components(speed, direction)
shear  = mpcalc.bulk_shear(p, u, v, height=z, depth=6*units.km)
srh, _, _ = mpcalc.storm_relative_helicity(z, u, v, depth=3*units.km)
rm, lm, _ = mpcalc.bunkers_storm_motion(p, u, v, z)
mcs_up, mcs_dn = mpcalc.corfidi_storm_motion(p, u, v)

# Composite parameters
stp = mpcalc.significant_tornado(sbcape, lcl_height, srh, shear)
scp = mpcalc.supercell_composite(mucape, srh, shear)

---

Kinematics & Dynamics

Requires grid spacing (use mpcalc.lat_lon_grid_deltas for lat/lon grids):

dx, dy = mpcalc.lat_lon_grid_deltas(lon, lat)

# Wind
spd  = mpcalc.wind_speed(u, v)
dirn = mpcalc.wind_direction(u, v)

# Vorticity / divergence
zeta = mpcalc.vorticity(u, v, dx=dx, dy=dy)
abvort = mpcalc.absolute_vorticity(u, v, dx=dx, dy=dy, latitude=lat)
div  = mpcalc.divergence(u, v, dx=dx, dy=dy)

# Deformation
shear_def   = mpcalc.shearing_deformation(u, v, dx=dx, dy=dy)
stretch_def = mpcalc.stretching_deformation(u, v, dx=dx, dy=dy)
total_def   = mpcalc.total_deformation(u, v, dx=dx, dy=dy)

# Geostrophic / ageostrophic wind
ug, vg = mpcalc.geostrophic_wind(height, dx=dx, dy=dy, latitude=lat)
uag, vag = mpcalc.ageostrophic_wind(height, u, v, dx=dx, dy=dy)

# Advection
T_adv = mpcalc.advection(T, u=u, v=v, dx=dx, dy=dy)

# Frontogenesis
fronto = mpcalc.frontogenesis(theta, u, v, dx=dx, dy=dy)

# Potential vorticity
pv = mpcalc.potential_vorticity_baroclinic(theta, u, v, p)
pv = mpcalc.potential_vorticity_barotropic(height, u, v, dx=dx, dy=dy)

# Q-vectors
qx, qy = mpcalc.q_vector(u, v, T, p, dx=dx, dy=dy)

# Coriolis
f = mpcalc.coriolis_parameter(lat)

---

Boundary Layer / Turbulence

# Brunt-Vaisala frequency
N  = mpcalc.brunt_vaisala_frequency(z, T, p)
N2 = mpcalc.brunt_vaisala_frequency_squared(z, T, p)

# Richardson number
Ri = mpcalc.gradient_richardson_number(z, T, u, v)

# TKE / friction velocity (from time series)
tke  = mpcalc.tke(u, v, w, perturbation=True)
ustar = mpcalc.friction_velocity(u, w)

---

Apparent Temperature

hi   = mpcalc.heat_index(T, rh)
wc   = mpcalc.windchill(T, speed)
at   = mpcalc.apparent_temperature(T, rh, speed)

---

Standard Atmosphere

p_std = mpcalc.height_to_pressure_std(z)
z_std = mpcalc.pressure_to_height_std(p)
p_sl  = mpcalc.altimeter_to_sea_level_pressure(altimeter, z)

---

Smoothing

T_smooth = mpcalc.smooth_gaussian(T_2d, n=3)
T_smooth = mpcalc.smooth_n_point(T_2d, n=9, passes=1)
T_smooth = mpcalc.smooth_rectangular(T_2d, size=(3,3), passes=2)
T_smooth = mpcalc.smooth_circular(T_2d, radius=2, passes=1)

---

Isentropic Interpolation

theta_levels = np.array([300, 310, 320]) * units.K
ds_isentropic = mpcalc.isentropic_interpolation_as_dataset(
    theta_levels, T, p, u, v
)

---

Utility

# Grid spacing from lat/lon
dx, dy = mpcalc.lat_lon_grid_deltas(lon, lat)

# Layer extraction
p_layer, T_layer = mpcalc.get_layer(p, T, height=z, bottom=z[0], depth=3*units.km)

# Find intersections (e.g. parcel/environment profile)
x_int, y_int = mpcalc.find_intersections(p, T, prof)

# Perturbation from mean
T_prime = mpcalc.get_perturbation(T_ts)