Skip to content

Oceanography

Physical oceanography diagnostics on gridded fields. The Coriolis parameter \(f=2\Omega\sin\phi\) sets the rotation; geostrophic balance gives the surface velocity from sea-surface height \(\eta\):

\[ u_g=-\frac{g}{f}\frac{\partial\eta}{\partial y},\qquad v_g=\;\frac{g}{f}\frac{\partial\eta}{\partial x}. \]

Vorticity, strain, kinetic energy, the Okubo–Weiss parameter, and stratification (Brunt–Väisälä frequency, mixed-layer depth) follow from the velocity and density fields.

Cartesian / β-plane grids

By default every horizontal diagnostic works on a lon/lat grid in degrees (dims=("lon", "lat"), geometry="spherical") and derives \(f\) from latitude. Idealised model output on an \(f\)- or β-plane — x/y in metres, no lon/lat — passes the dims, the geometry and an explicit Coriolis parameter instead:

from xrops.ocn import GeostrophicVelocities, OkuboWeiss, relative_vorticity

cart = {"dims": ("x", "y"), "geometry": "cartesian"}  # or "rectilinear"

zeta = relative_vorticity(ds, "u", "v", **cart)
ug = GeostrophicVelocities("eta", **cart, f=(1e-4, 1.6e-11))(ds)
ow = OkuboWeiss("u", "v", **cart)(ds)

f= accepts a scalar \(f_0\) (f-plane), an (f0, beta) tuple for the β-plane \(f = f_0 + \beta\,(y - y_0)\) with \(y_0\) the centre of the y coordinate, (f0, beta, y0) to pin \(y_0\), or a DataArray used as-is. It is required on non-spherical grids, where there is no latitude to derive it from. Fields must be collocated on one grid; destagger C-grid output before calling these diagnostics.

Rotation & geostrophy

CoriolisNormalized

CoriolisNormalized(variable: str, f0: float | None = None)

Bases: Operator

Normalise a field by the Coriolis parameter f₀.

Divides variable by f (e.g. to form the Rossby number ζ/f).

Parameters:

Name Type Description Default
variable str

Name of the variable to normalise.

required
f0 float | None

Reference Coriolis parameter (s⁻¹). None derives f from latitude.

None

Returns:

Type Description

The input dataset with variable replaced by variable / f.

Source code in packages/xrtoolz-ops/src/xrops/ocn/operators.py
def __init__(self, variable: str, f0: float | None = None):
    self.variable = variable
    self.f0 = f0

Streamfunction

Streamfunction(variable: str = 'ssh', g: float | None = None, f0: float | None = None)

Bases: Operator

Geostrophic stream function ψ = (g / f₀) · η from SSH.

Parameters:

Name Type Description Default
variable str

Name of the sea-surface-height (η) variable.

'ssh'
g float | None

Gravitational acceleration (m s⁻²). None uses the package default.

None
f0 float | None

Reference Coriolis parameter (s⁻¹). None derives f from latitude on the f-plane reference.

None

Returns:

Type Description

The input dataset with a psi stream-function variable (m² s⁻¹).

Source code in packages/xrtoolz-ops/src/xrops/ocn/operators.py
def __init__(
    self,
    variable: str = "ssh",
    g: float | None = None,
    f0: float | None = None,
):
    self.variable = variable
    self.g = g
    self.f0 = f0

GeostrophicVelocities

GeostrophicVelocities(variable: str = 'ssh', *, dims: tuple[str, str] = ('lon', 'lat'), geometry: Geometry = 'spherical', f: Coriolis = None)

Bases: Operator

Geostrophic velocities (u_g, v_g) from a height field.

Computes u_g = −(g/f) ∂η/∂y and v_g = (g/f) ∂η/∂x, on the lon/lat sphere by default.

Parameters:

Name Type Description Default
variable str

Name of the sea-surface-height (η) variable.

'ssh'
dims tuple[str, str]

(x_dim, y_dim) horizontal dims. Default ("lon", "lat").

('lon', 'lat')
geometry Geometry

"spherical" (default), "cartesian" or "rectilinear".

'spherical'
f Coriolis

Coriolis parameter. None derives it from latitude (spherical only); otherwise a scalar f₀, an (f₀, β) / (f₀, β, y₀) β-plane tuple, or a DataArray.

None

Returns:

Type Description

The input dataset with the geostrophic u and v velocity

variables (m s⁻¹).

Source code in packages/xrtoolz-ops/src/xrops/ocn/operators.py
def __init__(
    self,
    variable: str = "ssh",
    *,
    dims: tuple[str, str] = ("lon", "lat"),
    geometry: Geometry = "spherical",
    f: Coriolis = None,
):
    self.variable = variable
    self.dims = dims
    self.geometry = geometry
    self.f = f

AgeostrophicVelocities

AgeostrophicVelocities(variable: str = 'ssh', u: str = 'u', v: str = 'v', *, dims: tuple[str, str] = ('lon', 'lat'), geometry: Geometry = 'spherical', f: Coriolis = None)

Bases: Operator

Ageostrophic velocities (u_a, v_a) = (u, v) − u_g(η).

Subtracts the geostrophic velocities derived from variable (SSH) from the total velocities.

Parameters:

Name Type Description Default
variable str

Name of the sea-surface-height (η) variable.

'ssh'
u str

Name of the total eastward velocity variable.

'u'
v str

Name of the total northward velocity variable.

'v'
dims tuple[str, str]

(x_dim, y_dim) horizontal dims. Default ("lon", "lat").

('lon', 'lat')
geometry Geometry

"spherical" (default), "cartesian" or "rectilinear".

'spherical'
f Coriolis

Coriolis parameter. None derives it from latitude (spherical only); otherwise a scalar f₀, an (f₀, β) / (f₀, β, y₀) β-plane tuple, or a DataArray.

None

Returns:

Type Description

The input dataset with u_a and v_a ageostrophic-velocity

variables (m s⁻¹).

Source code in packages/xrtoolz-ops/src/xrops/ocn/operators.py
def __init__(
    self,
    variable: str = "ssh",
    u: str = "u",
    v: str = "v",
    *,
    dims: tuple[str, str] = ("lon", "lat"),
    geometry: Geometry = "spherical",
    f: Coriolis = None,
):
    self.variable = variable
    self.u = u
    self.v = v
    self.dims = dims
    self.geometry = geometry
    self.f = f

Vorticity & strain

RelativeVorticity

RelativeVorticity(u: str = 'u', v: str = 'v', *, dims: tuple[str, str] = ('lon', 'lat'), geometry: Geometry = 'spherical')

Bases: _UVGradOperator

Relative (vertical) vorticity ζ = ∂v/∂x − ∂u/∂y.

Parameters:

Name Type Description Default
u str

Name of the eastward velocity variable.

'u'
v str

Name of the northward velocity variable.

'v'
dims tuple[str, str]

(x_dim, y_dim) horizontal dims. Default ("lon", "lat").

('lon', 'lat')
geometry Geometry

"spherical" (default), "cartesian" or "rectilinear".

'spherical'

Returns:

Type Description

The input dataset with a vort_r variable (s⁻¹).

Source code in packages/xrtoolz-ops/src/xrops/ocn/operators.py
def __init__(
    self,
    u: str = "u",
    v: str = "v",
    *,
    dims: tuple[str, str] = ("lon", "lat"),
    geometry: Geometry = "spherical",
):
    super().__init__(u=u, v=v)
    self.dims = dims
    self.geometry = geometry

AbsoluteVorticity

AbsoluteVorticity(u: str = 'u', v: str = 'v', *, dims: tuple[str, str] = ('lon', 'lat'), geometry: Geometry = 'spherical', f: Coriolis = None)

Bases: _UVGradOperator

Absolute vorticity ζ_a = ζ + f (relative plus planetary).

Parameters:

Name Type Description Default
u str

Name of the eastward velocity variable.

'u'
v str

Name of the northward velocity variable.

'v'
dims tuple[str, str]

(x_dim, y_dim) horizontal dims. Default ("lon", "lat").

('lon', 'lat')
geometry Geometry

"spherical" (default), "cartesian" or "rectilinear".

'spherical'
f Coriolis

Coriolis parameter. None derives it from latitude (spherical only); otherwise a scalar f₀, an (f₀, β) / (f₀, β, y₀) β-plane tuple, or a DataArray.

None

Returns:

Type Description

The input dataset with a vort_a variable (s⁻¹).

Source code in packages/xrtoolz-ops/src/xrops/ocn/operators.py
def __init__(
    self,
    u: str = "u",
    v: str = "v",
    *,
    dims: tuple[str, str] = ("lon", "lat"),
    geometry: Geometry = "spherical",
    f: Coriolis = None,
):
    super().__init__(u=u, v=v, dims=dims, geometry=geometry)
    self.f = f

ShearVorticity

ShearVorticity(u: str = 'u', v: str = 'v', *, dims: tuple[str, str] = ('lon', 'lat'), geometry: Geometry = 'spherical')

Bases: _UVGradOperator

Along-flow shear component of the relative vorticity.

The part of ζ due to cross-stream changes in flow speed (as opposed to streamline curvature).

Parameters:

Name Type Description Default
u str

Name of the eastward velocity variable.

'u'
v str

Name of the northward velocity variable.

'v'
dims tuple[str, str]

(x_dim, y_dim) horizontal dims. Default ("lon", "lat").

('lon', 'lat')
geometry Geometry

"spherical" (default), "cartesian" or "rectilinear".

'spherical'

Returns:

Type Description

The input dataset with a vort_shear variable (s⁻¹).

Source code in packages/xrtoolz-ops/src/xrops/ocn/operators.py
def __init__(
    self,
    u: str = "u",
    v: str = "v",
    *,
    dims: tuple[str, str] = ("lon", "lat"),
    geometry: Geometry = "spherical",
):
    super().__init__(u=u, v=v)
    self.dims = dims
    self.geometry = geometry

CurvatureVorticity

CurvatureVorticity(u: str = 'u', v: str = 'v', *, dims: tuple[str, str] = ('lon', 'lat'), geometry: Geometry = 'spherical')

Bases: _UVGradOperator

Cross-flow curvature component of the relative vorticity.

The part of ζ due to streamline curvature; vort_shear and vort_curv sum to the relative vorticity vort_r.

Parameters:

Name Type Description Default
u str

Name of the eastward velocity variable.

'u'
v str

Name of the northward velocity variable.

'v'
dims tuple[str, str]

(x_dim, y_dim) horizontal dims. Default ("lon", "lat").

('lon', 'lat')
geometry Geometry

"spherical" (default), "cartesian" or "rectilinear".

'spherical'

Returns:

Type Description

The input dataset with a vort_curv variable (s⁻¹).

Source code in packages/xrtoolz-ops/src/xrops/ocn/operators.py
def __init__(
    self,
    u: str = "u",
    v: str = "v",
    *,
    dims: tuple[str, str] = ("lon", "lat"),
    geometry: Geometry = "spherical",
):
    super().__init__(u=u, v=v)
    self.dims = dims
    self.geometry = geometry

Divergence

Divergence(u: str = 'u', v: str = 'v', *, dims: tuple[str, str] = ('lon', 'lat'), geometry: Geometry = 'spherical')

Bases: _UVGradOperator

Horizontal divergence ∇·u = ∂u/∂x + ∂v/∂y.

Parameters:

Name Type Description Default
u str

Name of the eastward velocity variable.

'u'
v str

Name of the northward velocity variable.

'v'
dims tuple[str, str]

(x_dim, y_dim) horizontal dims. Default ("lon", "lat").

('lon', 'lat')
geometry Geometry

"spherical" (default), "cartesian" or "rectilinear".

'spherical'

Returns:

Type Description

The input dataset with a div variable (s⁻¹).

Source code in packages/xrtoolz-ops/src/xrops/ocn/operators.py
def __init__(
    self,
    u: str = "u",
    v: str = "v",
    *,
    dims: tuple[str, str] = ("lon", "lat"),
    geometry: Geometry = "spherical",
):
    super().__init__(u=u, v=v)
    self.dims = dims
    self.geometry = geometry

ShearStrain

ShearStrain(u: str = 'u', v: str = 'v', *, dims: tuple[str, str] = ('lon', 'lat'), geometry: Geometry = 'spherical')

Bases: _UVGradOperator

Shear strain rate S_s = ∂v/∂x + ∂u/∂y.

Parameters:

Name Type Description Default
u str

Name of the eastward velocity variable.

'u'
v str

Name of the northward velocity variable.

'v'
dims tuple[str, str]

(x_dim, y_dim) horizontal dims. Default ("lon", "lat").

('lon', 'lat')
geometry Geometry

"spherical" (default), "cartesian" or "rectilinear".

'spherical'

Returns:

Type Description

The input dataset with a shear_strain variable (s⁻¹).

Source code in packages/xrtoolz-ops/src/xrops/ocn/operators.py
def __init__(
    self,
    u: str = "u",
    v: str = "v",
    *,
    dims: tuple[str, str] = ("lon", "lat"),
    geometry: Geometry = "spherical",
):
    super().__init__(u=u, v=v)
    self.dims = dims
    self.geometry = geometry

TensorStrain

TensorStrain(u: str = 'u', v: str = 'v', *, dims: tuple[str, str] = ('lon', 'lat'), geometry: Geometry = 'spherical')

Bases: _UVGradOperator

Normal (tensor) strain rate S_n = ∂u/∂x − ∂v/∂y.

Parameters:

Name Type Description Default
u str

Name of the eastward velocity variable.

'u'
v str

Name of the northward velocity variable.

'v'
dims tuple[str, str]

(x_dim, y_dim) horizontal dims. Default ("lon", "lat").

('lon', 'lat')
geometry Geometry

"spherical" (default), "cartesian" or "rectilinear".

'spherical'

Returns:

Type Description

The input dataset with a tensor_strain variable (s⁻¹).

Source code in packages/xrtoolz-ops/src/xrops/ocn/operators.py
def __init__(
    self,
    u: str = "u",
    v: str = "v",
    *,
    dims: tuple[str, str] = ("lon", "lat"),
    geometry: Geometry = "spherical",
):
    super().__init__(u=u, v=v)
    self.dims = dims
    self.geometry = geometry

StrainMagnitude

StrainMagnitude(u: str = 'u', v: str = 'v', *, dims: tuple[str, str] = ('lon', 'lat'), geometry: Geometry = 'spherical')

Bases: _UVGradOperator

Total strain-rate magnitude √(S_n² + S_s²).

Parameters:

Name Type Description Default
u str

Name of the eastward velocity variable.

'u'
v str

Name of the northward velocity variable.

'v'
dims tuple[str, str]

(x_dim, y_dim) horizontal dims. Default ("lon", "lat").

('lon', 'lat')
geometry Geometry

"spherical" (default), "cartesian" or "rectilinear".

'spherical'

Returns:

Type Description

The input dataset with a strain variable (s⁻¹).

Source code in packages/xrtoolz-ops/src/xrops/ocn/operators.py
def __init__(
    self,
    u: str = "u",
    v: str = "v",
    *,
    dims: tuple[str, str] = ("lon", "lat"),
    geometry: Geometry = "spherical",
):
    super().__init__(u=u, v=v)
    self.dims = dims
    self.geometry = geometry

OkuboWeiss

OkuboWeiss(u: str = 'u', v: str = 'v', *, dims: tuple[str, str] = ('lon', 'lat'), geometry: Geometry = 'spherical')

Bases: _UVGradOperator

Okubo–Weiss parameter W = S_n² + S_s² − ζ².

Negative W marks vorticity-dominated (eddy-core) regions; positive W marks strain-dominated regions.

Parameters:

Name Type Description Default
u str

Name of the eastward velocity variable.

'u'
v str

Name of the northward velocity variable.

'v'
dims tuple[str, str]

(x_dim, y_dim) horizontal dims. Default ("lon", "lat").

('lon', 'lat')
geometry Geometry

"spherical" (default), "cartesian" or "rectilinear".

'spherical'

Returns:

Type Description

The input dataset with an ow variable (s⁻²).

Source code in packages/xrtoolz-ops/src/xrops/ocn/operators.py
def __init__(
    self,
    u: str = "u",
    v: str = "v",
    *,
    dims: tuple[str, str] = ("lon", "lat"),
    geometry: Geometry = "spherical",
):
    super().__init__(u=u, v=v)
    self.dims = dims
    self.geometry = geometry

Enstrophy

Enstrophy(variable: str = 'vort_r')

Bases: Operator

Enstrophy ½ ζ² from a relative-vorticity field.

Parameters:

Name Type Description Default
variable str

Name of the (relative) vorticity variable to square.

'vort_r'

Returns:

Type Description

The input dataset with an ens variable (s⁻²).

Source code in packages/xrtoolz-ops/src/xrops/ocn/operators.py
def __init__(self, variable: str = "vort_r"):
    self.variable = variable

PotentialVorticityBarotropic

PotentialVorticityBarotropic(height: str = 'h', u: str = 'u', v: str = 'v', *, dims: tuple[str, str] = ('lon', 'lat'), geometry: Geometry = 'spherical', f: Coriolis = None)

Bases: Operator

Single-layer barotropic potential vorticity (ζ + f) / h.

Parameters:

Name Type Description Default
height str

Name of the layer-thickness (h) variable.

'h'
u str

Name of the eastward velocity variable.

'u'
v str

Name of the northward velocity variable.

'v'
dims tuple[str, str]

(x_dim, y_dim) horizontal dims. Default ("lon", "lat").

('lon', 'lat')
geometry Geometry

"spherical" (default), "cartesian" or "rectilinear".

'spherical'
f Coriolis

Coriolis parameter. None derives it from latitude (spherical only); otherwise a scalar f₀, an (f₀, β) / (f₀, β, y₀) β-plane tuple, or a DataArray.

None

Returns:

Type Description

The input dataset with a pv_barotropic variable (m⁻¹ s⁻¹).

Source code in packages/xrtoolz-ops/src/xrops/ocn/operators.py
def __init__(
    self,
    height: str = "h",
    u: str = "u",
    v: str = "v",
    *,
    dims: tuple[str, str] = ("lon", "lat"),
    geometry: Geometry = "spherical",
    f: Coriolis = None,
):
    self.height = height
    self.u = u
    self.v = v
    self.dims = dims
    self.geometry = geometry
    self.f = f

Energy & transport

KineticEnergy

KineticEnergy(u: str = 'u', v: str = 'v')

Bases: _UVOperator

Kinetic energy KE = ½ (u² + v²).

Parameters:

Name Type Description Default
u str

Name of the eastward velocity variable.

'u'
v str

Name of the northward velocity variable.

'v'

Returns:

Type Description

The input dataset with a ke variable (m² s⁻²).

Source code in packages/xrtoolz-ops/src/xrops/ocn/operators.py
def __init__(self, u: str = "u", v: str = "v"):
    self.u = u
    self.v = v

EddyKineticEnergy

EddyKineticEnergy(u_anom: str = 'u_anom', v_anom: str = 'v_anom')

Bases: Operator

Eddy kinetic energy EKE = ½ (u'² + v'²) from velocity anomalies.

Parameters:

Name Type Description Default
u_anom str

Name of the eastward velocity-anomaly variable u'.

'u_anom'
v_anom str

Name of the northward velocity-anomaly variable v'.

'v_anom'

Returns:

Type Description

The input dataset with an eke variable (m² s⁻²).

Source code in packages/xrtoolz-ops/src/xrops/ocn/operators.py
def __init__(self, u_anom: str = "u_anom", v_anom: str = "v_anom"):
    self.u_anom = u_anom
    self.v_anom = v_anom

VelocityMagnitude

VelocityMagnitude(u: str = 'u', v: str = 'v', w: str | None = None)

Bases: Operator

Velocity magnitude |U| = √(u² + v² [+ w²]).

Parameters:

Name Type Description Default
u str

Name of the eastward velocity variable.

'u'
v str

Name of the northward velocity variable.

'v'
w str | None

Name of the vertical velocity variable, or None for the horizontal-only magnitude.

None

Returns:

Type Description

The input dataset with a speed variable (m s⁻¹).

Source code in packages/xrtoolz-ops/src/xrops/ocn/operators.py
def __init__(self, u: str = "u", v: str = "v", w: str | None = None):
    self.u = u
    self.v = v
    self.w = w

HorizontalVelocityMagnitude

HorizontalVelocityMagnitude(u: str = 'u', v: str = 'v')

Bases: _UVOperator

Horizontal current speed √(u² + v²).

Parameters:

Name Type Description Default
u str

Name of the eastward velocity variable.

'u'
v str

Name of the northward velocity variable.

'v'

Returns:

Type Description

The input dataset with a speed variable (m s⁻¹).

Source code in packages/xrtoolz-ops/src/xrops/ocn/operators.py
def __init__(self, u: str = "u", v: str = "v"):
    self.u = u
    self.v = v

Advection

Advection(scalar: str, components: tuple[str, ...] = ('u', 'v'), dim: tuple[str, ...] = ('lon', 'lat'), *, geometry: Geometry = 'spherical')

Bases: Operator

Tracer advection −u·∇c, on the lon/lat sphere by default.

Parameters:

Name Type Description Default
scalar str

Name of the advected tracer variable c.

required
components tuple[str, ...]

Names of the velocity components — (u, v) for the horizontal form, or (u, v, w) to include the vertical term −w ∂c/∂z.

('u', 'v')
dim tuple[str, ...]

Spatial dimension names the gradient is taken over: (lon, lat), optionally plus a vertical dim such as "depth", which is differentiated as a rectilinear axis.

('lon', 'lat')
geometry Geometry

Geometry of the horizontal pair — "spherical" (default, requires lon/lat), "cartesian" or "rectilinear".

'spherical'

Returns:

Type Description

The input dataset with a <scalar>_advection tendency variable

added.

Source code in packages/xrtoolz-ops/src/xrops/ocn/operators.py
def __init__(
    self,
    scalar: str,
    components: tuple[str, ...] = ("u", "v"),
    dim: tuple[str, ...] = ("lon", "lat"),
    *,
    geometry: Geometry = "spherical",
):
    self.scalar = scalar
    self.components = components
    self.dim = dim
    self.geometry = geometry

Frontogenesis

Frontogenesis(scalar: str, u: str = 'u', v: str = 'v', *, dims: tuple[str, str] = ('lon', 'lat'), geometry: Geometry = 'spherical')

Bases: Operator

Petterssen 2-D kinematic frontogenesis of a scalar.

The rate of change of the horizontal gradient magnitude of scalar following the flow; positive values indicate front sharpening.

Parameters:

Name Type Description Default
scalar str

Name of the tracer variable (e.g. temperature).

required
u str

Name of the eastward velocity variable.

'u'
v str

Name of the northward velocity variable.

'v'
dims tuple[str, str]

(x_dim, y_dim) horizontal dims. Default ("lon", "lat").

('lon', 'lat')
geometry Geometry

"spherical" (default), "cartesian" or "rectilinear".

'spherical'

Returns:

Type Description

The input dataset with a <scalar>_frontogenesis variable.

Source code in packages/xrtoolz-ops/src/xrops/ocn/operators.py
def __init__(
    self,
    scalar: str,
    u: str = "u",
    v: str = "v",
    *,
    dims: tuple[str, str] = ("lon", "lat"),
    geometry: Geometry = "spherical",
):
    self.scalar = scalar
    self.u = u
    self.v = v
    self.dims = dims
    self.geometry = geometry

Stratification

BruntVaisalaFrequency

BruntVaisalaFrequency(density: str = 'rho', depth: str = 'depth', rho0: float = 1025.0, g: float | None = None, positive: str = 'down')

Bases: Operator

Squared Brunt–Väisälä (buoyancy) frequency N² = −(g/ρ₀) ∂ρ/∂z.

Parameters:

Name Type Description Default
density str

Name of the (potential) density variable ρ.

'rho'
depth str

Name of the vertical coordinate variable.

'depth'
rho0 float

Reference density ρ₀ (kg m⁻³).

1025.0
g float | None

Gravitational acceleration (m s⁻²). None uses the package default.

None
positive str

Orientation of depth — "down" (depth increases downward) or "up".

'down'

Returns:

Type Description

The input dataset with an n_squared variable (s⁻²).

Source code in packages/xrtoolz-ops/src/xrops/ocn/operators.py
def __init__(
    self,
    density: str = "rho",
    depth: str = "depth",
    rho0: float = 1025.0,
    g: float | None = None,
    positive: str = "down",
):
    self.density = density
    self.depth = depth
    self.rho0 = rho0
    self.g = g
    self.positive = positive

MixedLayerDepth

MixedLayerDepth(density: str = 'rho', depth: str = 'depth', reference_depth: float = 10.0, threshold: float = 0.03)

Bases: Operator

Mixed-layer depth via the de Boyer Montégut density threshold.

Finds the depth at which density first exceeds the reference_depth value by threshold (kg m⁻³).

Parameters:

Name Type Description Default
density str

Name of the (potential) density variable.

'rho'
depth str

Name of the vertical coordinate variable.

'depth'
reference_depth float

Near-surface reference depth (m) the threshold is measured from.

10.0
threshold float

Density increase defining the mixed-layer base (kg m⁻³).

0.03

Returns:

Type Description

The input dataset with an mld variable (m).

Source code in packages/xrtoolz-ops/src/xrops/ocn/operators.py
def __init__(
    self,
    density: str = "rho",
    depth: str = "depth",
    reference_depth: float = 10.0,
    threshold: float = 0.03,
):
    self.density = density
    self.depth = depth
    self.reference_depth = reference_depth
    self.threshold = threshold

LapseRate

LapseRate(temperature: str = 'T', depth: str = 'depth', positive: str = 'down')

Bases: Operator

Vertical lapse rate Γ = −∂T/∂z (z positive upward).

Parameters:

Name Type Description Default
temperature str

Name of the temperature variable.

'T'
depth str

Name of the vertical coordinate variable.

'depth'
positive str

Orientation of depth — "down" (depth increases downward) or "up".

'down'

Returns:

Type Description

The input dataset with a lapse_rate variable (K m⁻¹).

Source code in packages/xrtoolz-ops/src/xrops/ocn/operators.py
def __init__(
    self,
    temperature: str = "T",
    depth: str = "depth",
    positive: str = "down",
):
    self.temperature = temperature
    self.depth = depth
    self.positive = positive

Sea-surface height

CalculateSSHAlongtrack

CalculateSSHAlongtrack(variable: str = 'ssh', sla: str = 'sla_filtered', mdt: str = 'mdt', lwe: str | None = 'lwe')

Bases: Operator

Compose along-track SSH from its altimetry components.

Reconstructs sea-surface height as ssh = sla + mdt − lwe: filtered sea-level anomaly plus the mean dynamic topography, minus an optional long-wavelength-error correction.

Parameters:

Name Type Description Default
variable str

Name of the output SSH variable to create.

'ssh'
sla str

Name of the (filtered) sea-level-anomaly variable.

'sla_filtered'
mdt str

Name of the mean-dynamic-topography variable.

'mdt'
lwe str | None

Name of the long-wavelength-error variable to subtract, or None to skip the correction.

'lwe'

Returns:

Type Description

The input dataset with the composed variable (SSH) added.

Source code in packages/xrtoolz-ops/src/xrops/ocn/operators.py
def __init__(
    self,
    variable: str = "ssh",
    sla: str = "sla_filtered",
    mdt: str = "mdt",
    lwe: str | None = "lwe",
):
    self.variable = variable
    self.sla = sla
    self.mdt = mdt
    self.lwe = lwe

Validation

ValidateSSH

ValidateSSH(variable: str = 'ssh')

Bases: Operator

Attach canonical CF metadata to a sea-surface-height field.

Assigns the canonical SSH attributes (units, standard_name, long_name) to variable via assign_attrs. No value checks (finiteness, range, shape) are performed.

Parameters:

Name Type Description Default
variable str

Name of the sea-surface-height variable to tag.

'ssh'

Returns:

Type Description

The input dataset with variable CF-tagged.

Source code in packages/xrtoolz-ops/src/xrops/ocn/operators.py
def __init__(self, variable: str = "ssh"):
    self.variable = variable

ValidateVelocity

ValidateVelocity(u: str = 'u', v: str = 'v')

Bases: Operator

Attach canonical CF metadata to a horizontal velocity pair (u, v).

Assigns the canonical zonal / meridional velocity attributes (units, standard_name, long_name) to both components via assign_attrs. No value checks (finiteness, shape consistency) are performed.

Parameters:

Name Type Description Default
u str

Name of the eastward (zonal) velocity variable.

'u'
v str

Name of the northward (meridional) velocity variable.

'v'

Returns:

Type Description

The input dataset with u and v CF-tagged.

Source code in packages/xrtoolz-ops/src/xrops/ocn/operators.py
def __init__(self, u: str = "u", v: str = "v"):
    self.u = u
    self.v = v

Functional primitives (Layer 0)

These pure functions back the operators above; each takes xr.DataArray/xr.Dataset and a dim: argument.

coriolis_parameter

coriolis_parameter(lat: DataArray | float | ndarray) -> DataArray | float | ndarray

Coriolis parameter f = 2 Ω sin(φ) in s⁻¹.

Parameters:

Name Type Description Default
lat DataArray | float | ndarray

Latitude(s) in degrees.

required

Returns:

Type Description
DataArray | float | ndarray

Coriolis parameter with the same shape as lat.

Source code in packages/xrtoolz-ops/src/xrops/ocn/_src/kinematics.py
def coriolis_parameter(
    lat: xr.DataArray | float | np.ndarray,
) -> xr.DataArray | float | np.ndarray:
    """Coriolis parameter ``f = 2 Ω sin(φ)`` in s⁻¹.

    Args:
        lat: Latitude(s) in **degrees**.

    Returns:
        Coriolis parameter with the same shape as ``lat``.
    """
    return 2.0 * grad.OMEGA * np.sin(np.deg2rad(lat))

coriolis_normalized

coriolis_normalized(ds: Dataset, variable: str, f0: float | None = None) -> Dataset

Normalize a variable by the Coriolis parameter f₀.

Common for plotting Rossby numbers (ζ / f, σ / f).

Parameters:

Name Type Description Default
ds Dataset

Dataset.

required
variable str

Variable to normalize.

required
f0 float | None

Coriolis parameter. Defaults to the mean over ds.lat.

None

Returns:

Type Description
Dataset

Dataset with variable replaced by variable / f0.

Source code in packages/xrtoolz-ops/src/xrops/ocn/_src/kinematics.py
def coriolis_normalized(
    ds: xr.Dataset,
    variable: str,
    f0: float | None = None,
) -> xr.Dataset:
    """Normalize a variable by the Coriolis parameter ``f₀``.

    Common for plotting Rossby numbers (``ζ / f``, ``σ / f``).

    Args:
        ds: Dataset.
        variable: Variable to normalize.
        f0: Coriolis parameter. Defaults to the mean over ``ds.lat``.

    Returns:
        Dataset with ``variable`` replaced by ``variable / f0``.
    """
    f0_value = float(coriolis_parameter(ds["lat"]).mean()) if f0 is None else float(f0)
    out = ds.copy()
    out[variable] = out[variable] / f0_value
    return out

streamfunction

streamfunction(ds: Dataset, variable: str = 'ssh', g: float | None = None, f0: float | None = None) -> Dataset

Stream function from sea-surface height.

Uses the linear geostrophic approximation ψ = (g / f₀) η, i.e. η = (f₀ / g) ψ.

Parameters:

Name Type Description Default
ds Dataset

Dataset containing variable (SSH) and a lat coord.

required
variable str

Name of the SSH variable.

'ssh'
g float | None

Gravity in m/s². Defaults to :data:xrops.grad.GRAVITY.

None
f0 float | None

Coriolis parameter in s⁻¹. Defaults to the mean over the latitude coordinate of ds.

None

Returns:

Type Description
Dataset

Dataset with a single variable "psi" (stream function).

Source code in packages/xrtoolz-ops/src/xrops/ocn/_src/kinematics.py
def streamfunction(
    ds: xr.Dataset,
    variable: str = "ssh",
    g: float | None = None,
    f0: float | None = None,
) -> xr.Dataset:
    """Stream function from sea-surface height.

    Uses the linear geostrophic approximation ``ψ = (g / f₀) η``, i.e.
    ``η = (f₀ / g) ψ``.

    Args:
        ds: Dataset containing ``variable`` (SSH) and a ``lat`` coord.
        variable: Name of the SSH variable.
        g: Gravity in m/s². Defaults to :data:`xrops.grad.GRAVITY`.
        f0: Coriolis parameter in s⁻¹. Defaults to the mean over the
            latitude coordinate of ``ds``.

    Returns:
        Dataset with a single variable ``"psi"`` (stream function).
    """
    ssh = ds[variable]
    g_value = grad.GRAVITY if g is None else g
    f0_value = float(coriolis_parameter(ssh["lat"]).mean()) if f0 is None else float(f0)
    psi = (g_value / f0_value) * ssh
    psi.attrs.update(
        long_name="Stream Function",
        standard_name="stream_function",
    )
    return psi.to_dataset(name="psi")

geostrophic_velocities

geostrophic_velocities(ds: Dataset, variable: str = 'ssh', *, dims: tuple[str, str] = ('lon', 'lat'), geometry: Geometry = 'spherical', f: Coriolis = None) -> Dataset

Geostrophic u and v from a height field (SSH).

Applies the standard linear geostrophic balance::

u = -(g / f) ∂η/∂y
v =  (g / f) ∂η/∂x

with ∂/∂x and ∂/∂y taken by :func:xrops.grad.partial under geometry (the lon/lat sphere by default).

Do not pass a stream-function field here — applying the same formula to ψ = (g/f) η would double the scaling and give unit-inconsistent velocities. :func:streamfunction is provided as a separate diagnostic.

Parameters:

Name Type Description Default
ds Dataset

Dataset containing the height variable.

required
variable str

Name of the SSH (height) variable, in metres.

'ssh'
dims tuple[str, str]

(x_dim, y_dim) horizontal dims. Default ("lon", "lat").

('lon', 'lat')
geometry Geometry

"spherical" (default; lon/lat in degrees), "cartesian" (uniform metres) or "rectilinear".

'spherical'
f Coriolis

Coriolis parameter. None (default) derives it from the latitude dims[1] (spherical only); otherwise a scalar f₀, an (f₀, β) / (f₀, β, y₀) β-plane tuple, or a DataArray.

None

Returns:

Type Description
Dataset

Dataset with u and v (m/s).

Example

A Gaussian SSH bump on a β-plane in metres::

x = y = np.linspace(-5e5, 5e5, 101)
eta = 0.5 * np.exp(-(x[None, :] ** 2 + y[:, None] ** 2) / 1e10)
ds = xr.Dataset({"eta": (("y", "x"), eta)}, coords={"x": x, "y": y})
geostrophic_velocities(
    ds, "eta", dims=("x", "y"), geometry="cartesian", f=(1e-4, 1.6e-11)
)
Source code in packages/xrtoolz-ops/src/xrops/ocn/_src/kinematics.py
def geostrophic_velocities(
    ds: xr.Dataset,
    variable: str = "ssh",
    *,
    dims: tuple[str, str] = ("lon", "lat"),
    geometry: Geometry = "spherical",
    f: Coriolis = None,
) -> xr.Dataset:
    """Geostrophic ``u`` and ``v`` from a height field (SSH).

    Applies the standard linear geostrophic balance::

        u = -(g / f) ∂η/∂y
        v =  (g / f) ∂η/∂x

    with ``∂/∂x`` and ``∂/∂y`` taken by :func:`xrops.grad.partial` under
    ``geometry`` (the lon/lat sphere by default).

    Do **not** pass a stream-function field here — applying the same
    formula to ``ψ = (g/f) η`` would double the scaling and give
    unit-inconsistent velocities. :func:`streamfunction` is provided as
    a separate diagnostic.

    Args:
        ds: Dataset containing the height ``variable``.
        variable: Name of the SSH (height) variable, in metres.
        dims: ``(x_dim, y_dim)`` horizontal dims. Default ``("lon", "lat")``.
        geometry: ``"spherical"`` (default; lon/lat in degrees),
            ``"cartesian"`` (uniform metres) or ``"rectilinear"``.
        f: Coriolis parameter. ``None`` (default) derives it from the
            latitude ``dims[1]`` (spherical only); otherwise a scalar
            ``f₀``, an ``(f₀, β)`` / ``(f₀, β, y₀)`` β-plane tuple, or a
            DataArray.

    Returns:
        Dataset with ``u`` and ``v`` (m/s).

    Example:
        A Gaussian SSH bump on a β-plane in metres::

            x = y = np.linspace(-5e5, 5e5, 101)
            eta = 0.5 * np.exp(-(x[None, :] ** 2 + y[:, None] ** 2) / 1e10)
            ds = xr.Dataset({"eta": (("y", "x"), eta)}, coords={"x": x, "y": y})
            geostrophic_velocities(
                ds, "eta", dims=("x", "y"), geometry="cartesian", f=(1e-4, 1.6e-11)
            )
    """
    ssh = ds[variable]
    f_field = _coriolis_field(ssh, dims, geometry, f)
    deta_dx = _partial(ssh, dims[0], dims, geometry)
    deta_dy = _partial(ssh, dims[1], dims, geometry)
    u = -(grad.GRAVITY / f_field) * deta_dy
    v = (grad.GRAVITY / f_field) * deta_dx
    u.attrs.update(long_name="Zonal Velocity", standard_name="zonal_velocity")
    v.attrs.update(long_name="Meridional Velocity", standard_name="meridional_velocity")
    return xr.Dataset({"u": u, "v": v})

ageostrophic_velocities

ageostrophic_velocities(ds: Dataset, variable: str = 'ssh', u: str = 'u', v: str = 'v', *, dims: tuple[str, str] = ('lon', 'lat'), geometry: Geometry = 'spherical', f: Coriolis = None) -> Dataset

Ageostrophic velocity (u_a, v_a) = (u, v) − u_g(η).

Parameters:

Name Type Description Default
ds Dataset

Dataset containing the SSH variable and the total velocities u/v.

required
variable str

Name of the SSH (height) variable from which the geostrophic velocities are derived.

'ssh'
u, v

Names of the total horizontal velocity components.

required
dims tuple[str, str]

(x_dim, y_dim) horizontal dims. Default ("lon", "lat").

('lon', 'lat')
geometry Geometry

"spherical" (default; lon/lat in degrees), "cartesian" (uniform metres) or "rectilinear".

'spherical'
f Coriolis

Coriolis parameter. None (default) derives it from the latitude dims[1] (spherical only); otherwise a scalar f₀, an (f₀, β) / (f₀, β, y₀) β-plane tuple, or a DataArray.

None

Returns:

Type Description
Dataset

Dataset with u_a and v_a (ageostrophic components).

Source code in packages/xrtoolz-ops/src/xrops/ocn/_src/kinematics.py
def ageostrophic_velocities(
    ds: xr.Dataset,
    variable: str = "ssh",
    u: str = "u",
    v: str = "v",
    *,
    dims: tuple[str, str] = ("lon", "lat"),
    geometry: Geometry = "spherical",
    f: Coriolis = None,
) -> xr.Dataset:
    """Ageostrophic velocity ``(u_a, v_a) = (u, v) − u_g(η)``.

    Args:
        ds: Dataset containing the SSH ``variable`` and the total
            velocities ``u``/``v``.
        variable: Name of the SSH (height) variable from which the
            geostrophic velocities are derived.
        u, v: Names of the total horizontal velocity components.
        dims: ``(x_dim, y_dim)`` horizontal dims. Default ``("lon", "lat")``.
        geometry: ``"spherical"`` (default; lon/lat in degrees),
            ``"cartesian"`` (uniform metres) or ``"rectilinear"``.
        f: Coriolis parameter. ``None`` (default) derives it from the
            latitude ``dims[1]`` (spherical only); otherwise a scalar
            ``f₀``, an ``(f₀, β)`` / ``(f₀, β, y₀)`` β-plane tuple, or a
            DataArray.

    Returns:
        Dataset with ``u_a`` and ``v_a`` (ageostrophic components).
    """
    geo = geostrophic_velocities(
        ds, variable=variable, dims=dims, geometry=geometry, f=f
    )
    u_a = ds[u] - geo["u"]
    v_a = ds[v] - geo["v"]
    u_a.attrs.update(
        long_name="Ageostrophic Zonal Velocity",
        standard_name="ageostrophic_zonal_velocity",
    )
    v_a.attrs.update(
        long_name="Ageostrophic Meridional Velocity",
        standard_name="ageostrophic_meridional_velocity",
    )
    return xr.Dataset({"u_a": u_a, "v_a": v_a})

relative_vorticity

relative_vorticity(ds: Dataset, u: str = 'u', v: str = 'v', *, dims: tuple[str, str] = ('lon', 'lat'), geometry: Geometry = 'spherical') -> Dataset

Relative vorticity ζ = ∂v/∂x - ∂u/∂y.

On the default spherical geometry the u tan φ / R curvature term is included.

Parameters:

Name Type Description Default
ds Dataset

Dataset with the velocity components.

required
u, v

Velocity component names.

required
dims tuple[str, str]

(x_dim, y_dim) horizontal dims. Default ("lon", "lat").

('lon', 'lat')
geometry Geometry

"spherical" (default; lon/lat in degrees), "cartesian" (uniform metres) or "rectilinear".

'spherical'

Returns:

Type Description
Dataset

Dataset with a single variable "vort_r" (s⁻¹).

Example

Solid-body rotation u = −y, v = x on a Cartesian grid has ζ = 2 everywhere::

x = y = np.linspace(0.0, 1e4, 11)
X, Y = np.meshgrid(x, y)
ds = xr.Dataset(
    {"u": (("y", "x"), -Y), "v": (("y", "x"), X)},
    coords={"x": x, "y": y},
)
relative_vorticity(ds, dims=("x", "y"), geometry="cartesian")
Source code in packages/xrtoolz-ops/src/xrops/ocn/_src/kinematics.py
def relative_vorticity(
    ds: xr.Dataset,
    u: str = "u",
    v: str = "v",
    *,
    dims: tuple[str, str] = ("lon", "lat"),
    geometry: Geometry = "spherical",
) -> xr.Dataset:
    """Relative vorticity ``ζ = ∂v/∂x - ∂u/∂y``.

    On the default spherical geometry the ``u tan φ / R`` curvature term
    is included.

    Args:
        ds: Dataset with the velocity components.
        u, v: Velocity component names.
        dims: ``(x_dim, y_dim)`` horizontal dims. Default ``("lon", "lat")``.
        geometry: ``"spherical"`` (default; lon/lat in degrees),
            ``"cartesian"`` (uniform metres) or ``"rectilinear"``.

    Returns:
        Dataset with a single variable ``"vort_r"`` (s⁻¹).

    Example:
        Solid-body rotation ``u = −y, v = x`` on a Cartesian grid has
        ``ζ = 2`` everywhere::

            x = y = np.linspace(0.0, 1e4, 11)
            X, Y = np.meshgrid(x, y)
            ds = xr.Dataset(
                {"u": (("y", "x"), -Y), "v": (("y", "x"), X)},
                coords={"x": x, "y": y},
            )
            relative_vorticity(ds, dims=("x", "y"), geometry="cartesian")
    """
    zeta = grad.curl(
        ds, (u, v), dims=dims, geometry=geometry, **_geom_kw(dims, geometry)
    )
    zeta.attrs.update(
        long_name="Relative Vorticity", standard_name="relative_vorticity"
    )
    return zeta.to_dataset(name="vort_r")

absolute_vorticity

absolute_vorticity(ds: Dataset, u: str = 'u', v: str = 'v', *, dims: tuple[str, str] = ('lon', 'lat'), geometry: Geometry = 'spherical', f: Coriolis = None) -> Dataset

Absolute vorticity η = ζ + f.

Parameters:

Name Type Description Default
ds Dataset

Dataset with the velocity components.

required
u, v

Velocity component names.

required
dims tuple[str, str]

(x_dim, y_dim) horizontal dims. Default ("lon", "lat").

('lon', 'lat')
geometry Geometry

"spherical" (default; lon/lat in degrees), "cartesian" (uniform metres) or "rectilinear".

'spherical'
f Coriolis

Coriolis parameter. None (default) derives it from the latitude dims[1] (spherical only); otherwise a scalar f₀, an (f₀, β) / (f₀, β, y₀) β-plane tuple, or a DataArray.

None

Returns:

Type Description
Dataset

Dataset with a single variable "vort_a" (s⁻¹).

Source code in packages/xrtoolz-ops/src/xrops/ocn/_src/kinematics.py
def absolute_vorticity(
    ds: xr.Dataset,
    u: str = "u",
    v: str = "v",
    *,
    dims: tuple[str, str] = ("lon", "lat"),
    geometry: Geometry = "spherical",
    f: Coriolis = None,
) -> xr.Dataset:
    """Absolute vorticity ``η = ζ + f``.

    Args:
        ds: Dataset with the velocity components.
        u, v: Velocity component names.
        dims: ``(x_dim, y_dim)`` horizontal dims. Default ``("lon", "lat")``.
        geometry: ``"spherical"`` (default; lon/lat in degrees),
            ``"cartesian"`` (uniform metres) or ``"rectilinear"``.
        f: Coriolis parameter. ``None`` (default) derives it from the
            latitude ``dims[1]`` (spherical only); otherwise a scalar
            ``f₀``, an ``(f₀, β)`` / ``(f₀, β, y₀)`` β-plane tuple, or a
            DataArray.

    Returns:
        Dataset with a single variable ``"vort_a"`` (s⁻¹).
    """
    zeta = relative_vorticity(ds, u=u, v=v, dims=dims, geometry=geometry)["vort_r"]
    eta = zeta + _coriolis_field(zeta, dims, geometry, f)
    eta.attrs.update(long_name="Absolute Vorticity", standard_name="absolute_vorticity")
    return eta.to_dataset(name="vort_a")

shear_vorticity

shear_vorticity(ds: Dataset, u: str = 'u', v: str = 'v', *, dims: tuple[str, str] = ('lon', 'lat'), geometry: Geometry = 'spherical') -> Dataset

Vertical shear vorticity component (along-flow).

From [Majumdar2024]_::

ζ_s = (v u ∂u/∂x + v² ∂v/∂x − u² ∂u/∂y − u v ∂v/∂y) / (u² + v²)

The total relative vorticity decomposes as ζ = shear_vorticity + curvature_vorticity. The decomposition is undefined where the wind vanishes — calm points (u = v = 0) are returned as 0 rather than NaN/Inf.

Parameters:

Name Type Description Default
ds Dataset

Dataset with the velocity components.

required
u, v

Velocity component names.

required
dims tuple[str, str]

(x_dim, y_dim) horizontal dims. Default ("lon", "lat").

('lon', 'lat')
geometry Geometry

"spherical" (default; lon/lat in degrees), "cartesian" (uniform metres) or "rectilinear".

'spherical'

Returns:

Type Description
Dataset

Dataset with a single variable vort_shear.

Source code in packages/xrtoolz-ops/src/xrops/ocn/_src/kinematics.py
def shear_vorticity(
    ds: xr.Dataset,
    u: str = "u",
    v: str = "v",
    *,
    dims: tuple[str, str] = ("lon", "lat"),
    geometry: Geometry = "spherical",
) -> xr.Dataset:
    """Vertical shear vorticity component (along-flow).

    From [Majumdar2024]_::

        ζ_s = (v u ∂u/∂x + v² ∂v/∂x − u² ∂u/∂y − u v ∂v/∂y) / (u² + v²)

    The total relative vorticity decomposes as
    ``ζ = shear_vorticity + curvature_vorticity``. The decomposition
    is undefined where the wind vanishes — calm points
    (``u = v = 0``) are returned as ``0`` rather than NaN/Inf.

    Args:
        ds: Dataset with the velocity components.
        u, v: Velocity component names.
        dims: ``(x_dim, y_dim)`` horizontal dims. Default ``("lon", "lat")``.
        geometry: ``"spherical"`` (default; lon/lat in degrees),
            ``"cartesian"`` (uniform metres) or ``"rectilinear"``.

    Returns:
        Dataset with a single variable ``vort_shear``.
    """
    uu = ds[u]
    vv = ds[v]
    dudx, dudy, dvdx, dvdy = _vector_derivatives(ds, u, v, dims, geometry)
    speed_sq = uu**2 + vv**2
    numerator = vv * uu * dudx + vv * vv * dvdx - uu * uu * dudy - uu * vv * dvdy
    safe_speed_sq = speed_sq.where(speed_sq != 0.0, 1.0)
    zeta_s = xr.where(speed_sq != 0.0, numerator / safe_speed_sq, 0.0)
    zeta_s.attrs.update(long_name="Shear Vorticity", standard_name="shear_vorticity")
    return zeta_s.rename(None).to_dataset(name="vort_shear")

curvature_vorticity

curvature_vorticity(ds: Dataset, u: str = 'u', v: str = 'v', *, dims: tuple[str, str] = ('lon', 'lat'), geometry: Geometry = 'spherical') -> Dataset

Vertical curvature vorticity component (cross-flow).

From [Majumdar2024]_::

ζ_c = (u² ∂v/∂x − v² ∂u/∂y − v u ∂u/∂x + u v ∂v/∂y) / (u² + v²)

Calm points (u = v = 0) are returned as 0, matching the convention in :func:shear_vorticity.

Parameters:

Name Type Description Default
ds Dataset

Dataset with the velocity components.

required
u, v

Velocity component names.

required
dims tuple[str, str]

(x_dim, y_dim) horizontal dims. Default ("lon", "lat").

('lon', 'lat')
geometry Geometry

"spherical" (default; lon/lat in degrees), "cartesian" (uniform metres) or "rectilinear".

'spherical'

Returns:

Type Description
Dataset

Dataset with a single variable vort_curv.

Source code in packages/xrtoolz-ops/src/xrops/ocn/_src/kinematics.py
def curvature_vorticity(
    ds: xr.Dataset,
    u: str = "u",
    v: str = "v",
    *,
    dims: tuple[str, str] = ("lon", "lat"),
    geometry: Geometry = "spherical",
) -> xr.Dataset:
    """Vertical curvature vorticity component (cross-flow).

    From [Majumdar2024]_::

        ζ_c = (u² ∂v/∂x − v² ∂u/∂y − v u ∂u/∂x + u v ∂v/∂y) / (u² + v²)

    Calm points (``u = v = 0``) are returned as ``0``, matching the
    convention in :func:`shear_vorticity`.

    Args:
        ds: Dataset with the velocity components.
        u, v: Velocity component names.
        dims: ``(x_dim, y_dim)`` horizontal dims. Default ``("lon", "lat")``.
        geometry: ``"spherical"`` (default; lon/lat in degrees),
            ``"cartesian"`` (uniform metres) or ``"rectilinear"``.

    Returns:
        Dataset with a single variable ``vort_curv``.
    """
    uu = ds[u]
    vv = ds[v]
    dudx, dudy, dvdx, dvdy = _vector_derivatives(ds, u, v, dims, geometry)
    speed_sq = uu**2 + vv**2
    numerator = uu * uu * dvdx - vv * vv * dudy - vv * uu * dudx + uu * vv * dvdy
    safe_speed_sq = speed_sq.where(speed_sq != 0.0, 1.0)
    zeta_c = xr.where(speed_sq != 0.0, numerator / safe_speed_sq, 0.0)
    zeta_c.attrs.update(
        long_name="Curvature Vorticity", standard_name="curvature_vorticity"
    )
    return zeta_c.rename(None).to_dataset(name="vort_curv")

divergence

divergence(ds: Dataset, u: str = 'u', v: str = 'v', *, dims: tuple[str, str] = ('lon', 'lat'), geometry: Geometry = 'spherical') -> Dataset

Horizontal divergence ∂u/∂x + ∂v/∂y.

On the default spherical geometry the −v tan φ / R curvature term is included.

Parameters:

Name Type Description Default
ds Dataset

Dataset with the velocity components.

required
u, v

Velocity component names.

required
dims tuple[str, str]

(x_dim, y_dim) horizontal dims. Default ("lon", "lat").

('lon', 'lat')
geometry Geometry

"spherical" (default; lon/lat in degrees), "cartesian" (uniform metres) or "rectilinear".

'spherical'

Returns:

Type Description
Dataset

Dataset with a single variable "div" (s⁻¹).

Source code in packages/xrtoolz-ops/src/xrops/ocn/_src/kinematics.py
def divergence(
    ds: xr.Dataset,
    u: str = "u",
    v: str = "v",
    *,
    dims: tuple[str, str] = ("lon", "lat"),
    geometry: Geometry = "spherical",
) -> xr.Dataset:
    """Horizontal divergence ``∂u/∂x + ∂v/∂y``.

    On the default spherical geometry the ``−v tan φ / R`` curvature term
    is included.

    Args:
        ds: Dataset with the velocity components.
        u, v: Velocity component names.
        dims: ``(x_dim, y_dim)`` horizontal dims. Default ``("lon", "lat")``.
        geometry: ``"spherical"`` (default; lon/lat in degrees),
            ``"cartesian"`` (uniform metres) or ``"rectilinear"``.

    Returns:
        Dataset with a single variable ``"div"`` (s⁻¹).
    """
    div = grad.divergence(
        ds, (u, v), dims=dims, geometry=geometry, **_geom_kw(dims, geometry)
    )
    div.attrs.update(long_name="Divergence", standard_name="divergence")
    return div.to_dataset(name="div")

shear_strain

shear_strain(ds: Dataset, u: str = 'u', v: str = 'v', *, dims: tuple[str, str] = ('lon', 'lat'), geometry: Geometry = 'spherical') -> Dataset

Shear strain Sₛ = ∂v/∂x + ∂u/∂y.

Parameters:

Name Type Description Default
ds Dataset

Dataset with the velocity components.

required
u, v

Velocity component names.

required
dims tuple[str, str]

(x_dim, y_dim) horizontal dims. Default ("lon", "lat").

('lon', 'lat')
geometry Geometry

"spherical" (default; lon/lat in degrees), "cartesian" (uniform metres) or "rectilinear".

'spherical'

Returns:

Type Description
Dataset

Dataset with a single variable "shear_strain" (s⁻¹).

Source code in packages/xrtoolz-ops/src/xrops/ocn/_src/kinematics.py
def shear_strain(
    ds: xr.Dataset,
    u: str = "u",
    v: str = "v",
    *,
    dims: tuple[str, str] = ("lon", "lat"),
    geometry: Geometry = "spherical",
) -> xr.Dataset:
    """Shear strain ``Sₛ = ∂v/∂x + ∂u/∂y``.

    Args:
        ds: Dataset with the velocity components.
        u, v: Velocity component names.
        dims: ``(x_dim, y_dim)`` horizontal dims. Default ``("lon", "lat")``.
        geometry: ``"spherical"`` (default; lon/lat in degrees),
            ``"cartesian"`` (uniform metres) or ``"rectilinear"``.

    Returns:
        Dataset with a single variable ``"shear_strain"`` (s⁻¹).
    """
    dvdx = _partial(ds[v], dims[0], dims, geometry)
    dudy = _partial(ds[u], dims[1], dims, geometry)
    sh = dvdx + dudy
    sh.attrs.update(long_name="Shear Strain", standard_name="shear_strain")
    return sh.to_dataset(name="shear_strain")

tensor_strain

tensor_strain(ds: Dataset, u: str = 'u', v: str = 'v', *, dims: tuple[str, str] = ('lon', 'lat'), geometry: Geometry = 'spherical') -> Dataset

Normal / tensor strain Sₙ = ∂u/∂x - ∂v/∂y.

Parameters:

Name Type Description Default
ds Dataset

Dataset with the velocity components.

required
u, v

Velocity component names.

required
dims tuple[str, str]

(x_dim, y_dim) horizontal dims. Default ("lon", "lat").

('lon', 'lat')
geometry Geometry

"spherical" (default; lon/lat in degrees), "cartesian" (uniform metres) or "rectilinear".

'spherical'

Returns:

Type Description
Dataset

Dataset with a single variable "tensor_strain" (s⁻¹).

Source code in packages/xrtoolz-ops/src/xrops/ocn/_src/kinematics.py
def tensor_strain(
    ds: xr.Dataset,
    u: str = "u",
    v: str = "v",
    *,
    dims: tuple[str, str] = ("lon", "lat"),
    geometry: Geometry = "spherical",
) -> xr.Dataset:
    """Normal / tensor strain ``Sₙ = ∂u/∂x - ∂v/∂y``.

    Args:
        ds: Dataset with the velocity components.
        u, v: Velocity component names.
        dims: ``(x_dim, y_dim)`` horizontal dims. Default ``("lon", "lat")``.
        geometry: ``"spherical"`` (default; lon/lat in degrees),
            ``"cartesian"`` (uniform metres) or ``"rectilinear"``.

    Returns:
        Dataset with a single variable ``"tensor_strain"`` (s⁻¹).
    """
    dudx = _partial(ds[u], dims[0], dims, geometry)
    dvdy = _partial(ds[v], dims[1], dims, geometry)
    st = dudx - dvdy
    st.attrs.update(long_name="Tensor Strain", standard_name="tensor_strain")
    return st.to_dataset(name="tensor_strain")

strain_magnitude

strain_magnitude(ds: Dataset, u: str = 'u', v: str = 'v', *, dims: tuple[str, str] = ('lon', 'lat'), geometry: Geometry = 'spherical') -> Dataset

Total strain magnitude sqrt(Sₙ² + Sₛ²).

Parameters:

Name Type Description Default
ds Dataset

Dataset with the velocity components.

required
u, v

Velocity component names.

required
dims tuple[str, str]

(x_dim, y_dim) horizontal dims. Default ("lon", "lat").

('lon', 'lat')
geometry Geometry

"spherical" (default; lon/lat in degrees), "cartesian" (uniform metres) or "rectilinear".

'spherical'

Returns:

Type Description
Dataset

Dataset with a single variable "strain" (s⁻¹).

Source code in packages/xrtoolz-ops/src/xrops/ocn/_src/kinematics.py
def strain_magnitude(
    ds: xr.Dataset,
    u: str = "u",
    v: str = "v",
    *,
    dims: tuple[str, str] = ("lon", "lat"),
    geometry: Geometry = "spherical",
) -> xr.Dataset:
    """Total strain magnitude ``sqrt(Sₙ² + Sₛ²)``.

    Args:
        ds: Dataset with the velocity components.
        u, v: Velocity component names.
        dims: ``(x_dim, y_dim)`` horizontal dims. Default ``("lon", "lat")``.
        geometry: ``"spherical"`` (default; lon/lat in degrees),
            ``"cartesian"`` (uniform metres) or ``"rectilinear"``.

    Returns:
        Dataset with a single variable ``"strain"`` (s⁻¹).
    """
    sn = tensor_strain(ds, u=u, v=v, dims=dims, geometry=geometry)["tensor_strain"]
    ss = shear_strain(ds, u=u, v=v, dims=dims, geometry=geometry)["shear_strain"]
    total = np.sqrt(sn**2 + ss**2)
    total.attrs.update(long_name="Strain Magnitude", standard_name="strain")
    return total.to_dataset(name="strain")

okubo_weiss

okubo_weiss(ds: Dataset, u: str = 'u', v: str = 'v', *, dims: tuple[str, str] = ('lon', 'lat'), geometry: Geometry = 'spherical') -> Dataset

Okubo–Weiss parameter Sₙ² + Sₛ² − ζ².

Positive in strain-dominated regions, negative in vortical regions.

Parameters:

Name Type Description Default
ds Dataset

Dataset with the velocity components.

required
u, v

Velocity component names.

required
dims tuple[str, str]

(x_dim, y_dim) horizontal dims. Default ("lon", "lat").

('lon', 'lat')
geometry Geometry

"spherical" (default; lon/lat in degrees), "cartesian" (uniform metres) or "rectilinear".

'spherical'

Returns:

Type Description
Dataset

Dataset with a single variable "ow" (s⁻²).

Source code in packages/xrtoolz-ops/src/xrops/ocn/_src/kinematics.py
def okubo_weiss(
    ds: xr.Dataset,
    u: str = "u",
    v: str = "v",
    *,
    dims: tuple[str, str] = ("lon", "lat"),
    geometry: Geometry = "spherical",
) -> xr.Dataset:
    """Okubo–Weiss parameter ``Sₙ² + Sₛ² − ζ²``.

    Positive in strain-dominated regions, negative in vortical regions.

    Args:
        ds: Dataset with the velocity components.
        u, v: Velocity component names.
        dims: ``(x_dim, y_dim)`` horizontal dims. Default ``("lon", "lat")``.
        geometry: ``"spherical"`` (default; lon/lat in degrees),
            ``"cartesian"`` (uniform metres) or ``"rectilinear"``.

    Returns:
        Dataset with a single variable ``"ow"`` (s⁻²).
    """
    kw: dict[str, Any] = {"u": u, "v": v, "dims": dims, "geometry": geometry}
    sn = tensor_strain(ds, **kw)["tensor_strain"]
    ss = shear_strain(ds, **kw)["shear_strain"]
    zeta = relative_vorticity(ds, **kw)["vort_r"]
    ow = sn**2 + ss**2 - zeta**2
    ow.attrs.update(long_name="Okubo-Weiss Parameter", standard_name="okubo_weiss")
    return ow.to_dataset(name="ow")

enstrophy

enstrophy(ds: Dataset, variable: str = 'vort_r') -> Dataset

Enstrophy 0.5 ζ².

Source code in packages/xrtoolz-ops/src/xrops/ocn/_src/kinematics.py
def enstrophy(
    ds: xr.Dataset,
    variable: str = "vort_r",
) -> xr.Dataset:
    """Enstrophy ``0.5 ζ²``."""
    ens = 0.5 * (ds[variable] ** 2)
    ens.attrs.update(long_name="Enstrophy", standard_name="enstrophy")
    return ens.to_dataset(name="ens")

potential_vorticity_barotropic

potential_vorticity_barotropic(ds: Dataset, height: str = 'h', u: str = 'u', v: str = 'v', *, dims: tuple[str, str] = ('lon', 'lat'), geometry: Geometry = 'spherical', f: Coriolis = None) -> Dataset

Single-layer barotropic potential vorticity (ζ + f) / h.

For a shallow-water / single-layer model h is the layer thickness; for atmospheric height-coordinate use it is the geopotential height of the layer.

Parameters:

Name Type Description Default
ds Dataset

Dataset with u, v, and the layer thickness/height.

required
height str

Variable name of the layer thickness / height.

'h'
u, v

Velocity component names.

required
dims tuple[str, str]

(x_dim, y_dim) horizontal dims. Default ("lon", "lat").

('lon', 'lat')
geometry Geometry

"spherical" (default; lon/lat in degrees), "cartesian" (uniform metres) or "rectilinear".

'spherical'
f Coriolis

Coriolis parameter. None (default) derives it from the latitude dims[1] (spherical only); otherwise a scalar f₀, an (f₀, β) / (f₀, β, y₀) β-plane tuple, or a DataArray.

None

Returns:

Type Description
Dataset

Dataset with a single variable "pv_barotropic".

Source code in packages/xrtoolz-ops/src/xrops/ocn/_src/kinematics.py
def potential_vorticity_barotropic(
    ds: xr.Dataset,
    height: str = "h",
    u: str = "u",
    v: str = "v",
    *,
    dims: tuple[str, str] = ("lon", "lat"),
    geometry: Geometry = "spherical",
    f: Coriolis = None,
) -> xr.Dataset:
    """Single-layer barotropic potential vorticity ``(ζ + f) / h``.

    For a shallow-water / single-layer model ``h`` is the layer
    thickness; for atmospheric height-coordinate use it is the
    geopotential height of the layer.

    Args:
        ds: Dataset with ``u``, ``v``, and the layer thickness/height.
        height: Variable name of the layer thickness / height.
        u, v: Velocity component names.
        dims: ``(x_dim, y_dim)`` horizontal dims. Default ``("lon", "lat")``.
        geometry: ``"spherical"`` (default; lon/lat in degrees),
            ``"cartesian"`` (uniform metres) or ``"rectilinear"``.
        f: Coriolis parameter. ``None`` (default) derives it from the
            latitude ``dims[1]`` (spherical only); otherwise a scalar
            ``f₀``, an ``(f₀, β)`` / ``(f₀, β, y₀)`` β-plane tuple, or a
            DataArray.

    Returns:
        Dataset with a single variable ``"pv_barotropic"``.
    """
    eta = absolute_vorticity(ds, u=u, v=v, dims=dims, geometry=geometry, f=f)["vort_a"]
    pv = eta / ds[height]
    pv.attrs.update(
        long_name="Barotropic Potential Vorticity",
        standard_name="barotropic_potential_vorticity",
    )
    return pv.rename(None).to_dataset(name="pv_barotropic")

kinetic_energy

kinetic_energy(ds: Dataset, u: str = 'u', v: str = 'v') -> Dataset

Kinetic energy 0.5 (u² + v²).

Source code in packages/xrtoolz-ops/src/xrops/ocn/_src/kinematics.py
def kinetic_energy(
    ds: xr.Dataset,
    u: str = "u",
    v: str = "v",
) -> xr.Dataset:
    """Kinetic energy ``0.5 (u² + v²)``."""
    ke = 0.5 * (ds[u] ** 2 + ds[v] ** 2)
    ke.attrs.update(long_name="Kinetic Energy", standard_name="kinetic_energy")
    return ke.to_dataset(name="ke")

eddy_kinetic_energy

eddy_kinetic_energy(ds: Dataset, u_anom: str = 'u_anom', v_anom: str = 'v_anom') -> Dataset

Eddy kinetic energy EKE = ½ (u'² + v'²).

The caller supplies anomaly fields (typically a Reynolds decomposition: u' = u − ⟨u⟩). Choosing the temporal/spatial mean is a science decision left to the caller; we just compute the quadratic combination.

Parameters:

Name Type Description Default
ds Dataset

Dataset with the velocity anomaly fields.

required
u_anom, v_anom

Variable names of the anomalies.

required

Returns:

Type Description
Dataset

Dataset with a single variable "eke".

Source code in packages/xrtoolz-ops/src/xrops/ocn/_src/kinematics.py
def eddy_kinetic_energy(
    ds: xr.Dataset,
    u_anom: str = "u_anom",
    v_anom: str = "v_anom",
) -> xr.Dataset:
    """Eddy kinetic energy ``EKE = ½ (u'² + v'²)``.

    The caller supplies anomaly fields (typically a Reynolds
    decomposition: ``u' = u − ⟨u⟩``). Choosing the temporal/spatial mean
    is a science decision left to the caller; we just compute the
    quadratic combination.

    Args:
        ds: Dataset with the velocity anomaly fields.
        u_anom, v_anom: Variable names of the anomalies.

    Returns:
        Dataset with a single variable ``"eke"``.
    """
    eke = 0.5 * (ds[u_anom] ** 2 + ds[v_anom] ** 2)
    eke.attrs.update(
        long_name="Eddy Kinetic Energy", standard_name="eddy_kinetic_energy"
    )
    return eke.rename(None).to_dataset(name="eke")

velocity_magnitude

velocity_magnitude(ds: Dataset, u: str = 'u', v: str = 'v', w: str | None = None) -> Dataset

Velocity magnitude |U|.

2-D sqrt(u² + v²) if w is None; otherwise the full 3-D sqrt(u² + v² + w²).

Returns:

Type Description
Dataset

Dataset with a single variable "speed".

Source code in packages/xrtoolz-ops/src/xrops/ocn/_src/kinematics.py
def velocity_magnitude(
    ds: xr.Dataset,
    u: str = "u",
    v: str = "v",
    w: str | None = None,
) -> xr.Dataset:
    """Velocity magnitude ``|U|``.

    2-D ``sqrt(u² + v²)`` if ``w`` is ``None``; otherwise the full
    3-D ``sqrt(u² + v² + w²)``.

    Returns:
        Dataset with a single variable ``"speed"``.
    """
    sq = ds[u] ** 2 + ds[v] ** 2
    if w is not None:
        sq = sq + ds[w] ** 2
    speed = np.sqrt(sq)
    speed.attrs.update(
        long_name="Velocity Magnitude", standard_name="velocity_magnitude"
    )
    return speed.rename(None).to_dataset(name="speed")

horizontal_velocity_magnitude

horizontal_velocity_magnitude(ds: Dataset, u: str = 'u', v: str = 'v') -> Dataset

Horizontal velocity magnitude sqrt(u² + v²).

Equivalent to :func:velocity_magnitude with w=None; provided under the oceanspy-style name for callers that prefer the explicit horizontal-only form.

Source code in packages/xrtoolz-ops/src/xrops/ocn/_src/kinematics.py
def horizontal_velocity_magnitude(
    ds: xr.Dataset,
    u: str = "u",
    v: str = "v",
) -> xr.Dataset:
    """Horizontal velocity magnitude ``sqrt(u² + v²)``.

    Equivalent to :func:`velocity_magnitude` with ``w=None``; provided
    under the oceanspy-style name for callers that prefer the explicit
    horizontal-only form.
    """
    return velocity_magnitude(ds, u=u, v=v)

advection

advection(ds: Dataset, scalar: str, components: tuple[str, ...] = ('u', 'v'), dims: tuple[str, ...] = ('lon', 'lat'), *, geometry: Geometry = 'spherical') -> Dataset

Tracer advection −u·∇c, on the lon/lat sphere by default.

Sign convention follows :func:metpy.calc.advection: a positive value means the tracer is being increased at that point by the flow.

Parameters:

Name Type Description Default
ds Dataset

Dataset containing scalar and the wind/current components.

required
scalar str

Name of the advected scalar field (e.g. "sst", "ssh").

required
components tuple[str, ...]

Wind/current component variable names paired with dims — (u, v) for the horizontal form, or (u, v, w) to include vertical advection.

('u', 'v')
dims tuple[str, ...]

Two horizontal coordinate names, optionally followed by a vertical coordinate name. Under the default spherical geometry the horizontal pair must be lon/lat. The vertical axis is differentiated with geometry="rectilinear", so a stretched depth axis needs no special handling. w must live on the same grid as scalar — staggered vertical grids are not supported.

('lon', 'lat')
geometry Geometry

Geometry of the horizontal pair — "spherical" (default), "cartesian" or "rectilinear".

'spherical'

Returns:

Type Description
Dataset

Dataset with a single variable f"{scalar}_advection".

Raises:

Type Description
ValueError

if components and dims lengths differ, if there are not 2 or 3 of each, if the first two dims are not the lon/lat pair under spherical geometry, or if the vertical dim repeats one of them.

Source code in packages/xrtoolz-ops/src/xrops/ocn/_src/kinematics.py
def advection(
    ds: xr.Dataset,
    scalar: str,
    components: tuple[str, ...] = ("u", "v"),
    dims: tuple[str, ...] = ("lon", "lat"),
    *,
    geometry: Geometry = "spherical",
) -> xr.Dataset:
    """Tracer advection ``−u·∇c``, on the lon/lat sphere by default.

    Sign convention follows :func:`metpy.calc.advection`: a positive
    value means the tracer is being increased at that point by the flow.

    Args:
        ds: Dataset containing ``scalar`` and the wind/current
            components.
        scalar: Name of the advected scalar field (e.g. ``"sst"``,
            ``"ssh"``).
        components: Wind/current component variable names paired with
            ``dims`` — ``(u, v)`` for the horizontal form, or
            ``(u, v, w)`` to include vertical advection.
        dims: Two horizontal coordinate names, optionally followed by a
            vertical coordinate name. Under the default spherical
            geometry the horizontal pair must be ``lon``/``lat``. The
            vertical axis is differentiated with
            ``geometry="rectilinear"``, so a stretched depth axis needs no
            special handling. ``w`` must live on the same grid as
            ``scalar`` — staggered vertical grids are not supported.
        geometry: Geometry of the horizontal pair — ``"spherical"``
            (default), ``"cartesian"`` or ``"rectilinear"``.

    Returns:
        Dataset with a single variable ``f"{scalar}_advection"``.

    Raises:
        ValueError: if ``components`` and ``dims`` lengths differ, if
            there are not 2 or 3 of each, if the first two dims are not
            the lon/lat pair under spherical geometry, or if the vertical
            dim repeats one of them.
    """
    if len(components) != len(dims):
        raise ValueError(
            f"components ({len(components)}) and dims ({len(dims)}) must "
            f"have the same length; got components={components!r} "
            f"dims={dims!r}."
        )
    if len(dims) not in (2, 3):
        raise ValueError(
            "advection takes two horizontal dims, optionally plus one "
            f"vertical dim; got dims={dims!r}."
        )
    if geometry == "spherical" and sorted(dims[:2]) != ["lat", "lon"]:
        raise ValueError(
            f"the first two dims must be the horizontal lon/lat pair so the "
            f"spherical metric applies; got {tuple(dims[:2])!r}."
        )
    if len(dims) == 3 and dims[2] in dims[:2]:
        raise ValueError(
            f"the vertical dim must differ from the horizontal pair; got dims={dims!r}."
        )
    missing = [d for d in dims if d not in ds[scalar].dims]
    if missing:
        raise ValueError(
            f"scalar {scalar!r} has dims {tuple(ds[scalar].dims)!r} and does "
            f"not span {missing!r}; every advected dim must be on the tracer."
        )
    # Each component must live on the axis it is paired with. Without this a
    # staggered velocity (w on 'depth_w' while the tracer is on 'depth')
    # broadcasts instead of aligning, silently returning a field carrying
    # both vertical dims rather than the tracer's grid.
    for comp_name, dim in zip(components, dims, strict=True):
        if dim not in ds[comp_name].dims:
            raise ValueError(
                f"component {comp_name!r} has dims {tuple(ds[comp_name].dims)!r} "
                f"and does not span its paired dim {dim!r}. Staggered grids are "
                f"not supported — interpolate {comp_name!r} onto the {scalar!r} "
                "grid first."
            )
    # The spherical pair may arrive in either order; xrops.grad needs it by name.
    horizontal = ("lon", "lat") if geometry == "spherical" else (dims[0], dims[1])
    flux: xr.DataArray | None = None
    for comp_name, dim in zip(components, dims, strict=True):
        # Horizontal axes carry ``geometry``'s metric; the vertical axis is
        # a plain 1-D coordinate, typically stretched, hence rectilinear.
        if dim in horizontal:
            partial_c = _partial(ds[scalar], dim, horizontal, geometry)
        else:
            partial_c = grad.partial(ds[scalar], dim, geometry="rectilinear")
        term = ds[comp_name] * partial_c
        flux = term if flux is None else flux + term
    assert flux is not None
    out = -flux
    out.attrs.update(
        long_name=f"Advection of {scalar}",
        standard_name=f"{scalar}_advection",
    )
    return out.to_dataset(name=f"{scalar}_advection")

frontogenesis

frontogenesis(ds: Dataset, scalar: str, u: str = 'u', v: str = 'v', *, dims: tuple[str, str] = ('lon', 'lat'), geometry: Geometry = 'spherical') -> Dataset

Petterssen 2-D kinematic frontogenesis, on the lon/lat sphere by default.

Implements the form from [Bluestein1993]_ pp. 248–253::

F = ½ |∇c| [D cos(2β) − δ]

where c is the scalar (potential temperature in atmospheric use, or any other density-correlated tracer such as SST in oceanography), D is the total deformation, β is the angle between the axis of dilatation and the isentropes, and δ is the divergence.

Parameters:

Name Type Description Default
ds Dataset

Dataset with the scalar and the velocity components.

required
scalar str

Name of the scalar field (e.g. "sst", "theta").

required
u, v

Velocity component names.

required
dims tuple[str, str]

(x_dim, y_dim) horizontal dims. Default ("lon", "lat").

('lon', 'lat')
geometry Geometry

"spherical" (default; lon/lat in degrees), "cartesian" (uniform metres) or "rectilinear".

'spherical'

Returns:

Type Description
Dataset

Dataset with a single variable f"{scalar}_frontogenesis".

Source code in packages/xrtoolz-ops/src/xrops/ocn/_src/kinematics.py
def frontogenesis(
    ds: xr.Dataset,
    scalar: str,
    u: str = "u",
    v: str = "v",
    *,
    dims: tuple[str, str] = ("lon", "lat"),
    geometry: Geometry = "spherical",
) -> xr.Dataset:
    """Petterssen 2-D kinematic frontogenesis, on the lon/lat sphere by default.

    Implements the form from [Bluestein1993]_ pp. 248–253::

        F = ½ |∇c| [D cos(2β) − δ]

    where ``c`` is the scalar (potential temperature in atmospheric use,
    or any other density-correlated tracer such as SST in oceanography),
    ``D`` is the total deformation, ``β`` is the angle between the axis
    of dilatation and the isentropes, and ``δ`` is the divergence.

    Args:
        ds: Dataset with the scalar and the velocity components.
        scalar: Name of the scalar field (e.g. ``"sst"``, ``"theta"``).
        u, v: Velocity component names.
        dims: ``(x_dim, y_dim)`` horizontal dims. Default ``("lon", "lat")``.
        geometry: ``"spherical"`` (default; lon/lat in degrees),
            ``"cartesian"`` (uniform metres) or ``"rectilinear"``.

    Returns:
        Dataset with a single variable ``f"{scalar}_frontogenesis"``.
    """
    dcdx = _partial(ds[scalar], dims[0], dims, geometry)
    dcdy = _partial(ds[scalar], dims[1], dims, geometry)
    mag = np.sqrt(dcdx**2 + dcdy**2)
    kw: dict[str, Any] = {"u": u, "v": v, "dims": dims, "geometry": geometry}
    sh = shear_strain(ds, **kw)["shear_strain"]
    st = tensor_strain(ds, **kw)["tensor_strain"]
    total = strain_magnitude(ds, **kw)["strain"]
    div = divergence(ds, **kw)["div"]
    # ψ is the angle of the axis of dilatation; sin β projects the
    # scalar gradient onto that axis.
    psi = 0.5 * np.arctan2(sh, st)
    # ``cos(2β) = 1 − 2 sin²β``; keep the form metpy uses.
    sin_beta = xr.where(
        mag != 0.0,
        (dcdx * np.cos(psi) + dcdy * np.sin(psi)) / mag.where(mag != 0.0, 1.0),
        0.0,
    )
    front = 0.5 * mag * (total * (1.0 - 2.0 * sin_beta**2) - div)
    front.attrs.update(
        long_name=f"Frontogenesis of {scalar}",
        standard_name=f"{scalar}_frontogenesis",
    )
    return front.rename(None).to_dataset(name=f"{scalar}_frontogenesis")

brunt_vaisala_frequency

brunt_vaisala_frequency(ds: Dataset, density: str = 'rho', depth: str = 'depth', rho0: float = 1025.0, g: float | None = None, positive: str = 'down') -> Dataset

Squared Brunt–Väisälä (buoyancy) frequency.

Defined as N² = −(g/ρ₀) ∂ρ/∂z with z positive upward. In a stable column N² > 0. We return the same sign convention independent of how the input vertical coordinate is oriented; set positive to match the dataset.

Vertical FD uses the rectilinear geometry so non-uniform pressure / depth grids are supported.

Parameters:

Name Type Description Default
ds Dataset

Dataset with a density variable on a vertical coordinate.

required
density str

In-situ or potential density variable (kg m⁻³).

'rho'
depth str

Name of the vertical coordinate (metres).

'depth'
rho0 float

Reference density (kg m⁻³). Defaults to 1025.

1025.0
g float | None

Gravity (m s⁻²). Defaults to :data:xrops.grad.GRAVITY.

None
positive str

Either "down" (default — oceanographic depth coordinate, increasing into the column) or "up" (height above some reference, decreasing into the column).

'down'

Returns:

Type Description
Dataset

Dataset with a single variable "n_squared".

Source code in packages/xrtoolz-ops/src/xrops/ocn/_src/kinematics.py
def brunt_vaisala_frequency(
    ds: xr.Dataset,
    density: str = "rho",
    depth: str = "depth",
    rho0: float = 1025.0,
    g: float | None = None,
    positive: str = "down",
) -> xr.Dataset:
    """Squared Brunt–Väisälä (buoyancy) frequency.

    Defined as ``N² = −(g/ρ₀) ∂ρ/∂z`` with ``z`` positive **upward**.
    In a stable column ``N² > 0``. We return the same sign convention
    independent of how the input vertical coordinate is oriented;
    set ``positive`` to match the dataset.

    Vertical FD uses the ``rectilinear`` geometry so non-uniform
    pressure / depth grids are supported.

    Args:
        ds: Dataset with a ``density`` variable on a vertical
            coordinate.
        density: In-situ or potential density variable (kg m⁻³).
        depth: Name of the vertical coordinate (metres).
        rho0: Reference density (kg m⁻³). Defaults to 1025.
        g: Gravity (m s⁻²). Defaults to :data:`xrops.grad.GRAVITY`.
        positive: Either ``"down"`` (default — oceanographic depth
            coordinate, increasing into the column) or ``"up"`` (height
            above some reference, decreasing into the column).

    Returns:
        Dataset with a single variable ``"n_squared"``.
    """
    if positive not in ("down", "up"):
        raise ValueError(f"positive={positive!r} must be 'down' or 'up'.")
    g_value = grad.GRAVITY if g is None else g
    drho_dz = grad.partial(ds[density], depth, geometry="rectilinear")
    sign = 1.0 if positive == "down" else -1.0
    n_sq = sign * (g_value / rho0) * drho_dz
    n_sq.attrs.update(
        long_name="Squared Brunt-Väisälä Frequency",
        standard_name="square_of_brunt_vaisala_frequency_in_sea_water",
        units="s-2",
    )
    return n_sq.rename(None).to_dataset(name="n_squared")

mixed_layer_depth

mixed_layer_depth(ds: Dataset, density: str = 'rho', depth: str = 'depth', reference_depth: float = 10.0, threshold: float = 0.03) -> Dataset

Mixed-layer depth via the de Boyer Montégut density threshold.

Implements the criterion of [deBoyerMontegut2004]_::

MLD = first depth z below z_ref where ρ(z) − ρ(z_ref) > Δρ

with z_ref = 10 m and Δρ = 0.03 kg/m³ as the standard defaults. The reference density is taken as the value of density at the model level closest to reference_depth.

For columns where the criterion is never met (fully mixed water column), the returned value is the deepest depth in the profile.

Parameters:

Name Type Description Default
ds Dataset

Dataset with a density variable on a vertical coordinate named depth (positive downward, metres).

required
density str

Density variable name.

'rho'
depth str

Vertical coordinate name.

'depth'
reference_depth float

Depth at which ρ_ref is read (m).

10.0
threshold float

Density jump used to mark the base of the mixed layer (kg m⁻³).

0.03

Returns:

Type Description
Dataset

Dataset with a single variable "mld" collapsed over the

Dataset

vertical coordinate.

Source code in packages/xrtoolz-ops/src/xrops/ocn/_src/kinematics.py
def mixed_layer_depth(
    ds: xr.Dataset,
    density: str = "rho",
    depth: str = "depth",
    reference_depth: float = 10.0,
    threshold: float = 0.03,
) -> xr.Dataset:
    """Mixed-layer depth via the de Boyer Montégut density threshold.

    Implements the criterion of [deBoyerMontegut2004]_::

        MLD = first depth z below z_ref where ρ(z) − ρ(z_ref) > Δρ

    with ``z_ref = 10 m`` and ``Δρ = 0.03 kg/m³`` as the standard
    defaults. The reference density is taken as the value of ``density``
    at the model level closest to ``reference_depth``.

    For columns where the criterion is never met (fully mixed water
    column), the returned value is the deepest depth in the profile.

    Args:
        ds: Dataset with a ``density`` variable on a vertical coordinate
            named ``depth`` (positive downward, metres).
        density: Density variable name.
        depth: Vertical coordinate name.
        reference_depth: Depth at which ``ρ_ref`` is read (m).
        threshold: Density jump used to mark the base of the mixed layer
            (kg m⁻³).

    Returns:
        Dataset with a single variable ``"mld"`` collapsed over the
        vertical coordinate.
    """
    rho = ds[density]
    if depth not in rho.dims:
        raise ValueError(
            f"Density {density!r} is not defined on the {depth!r} dimension."
        )
    z = rho[depth]
    rho_ref = rho.sel({depth: reference_depth}, method="nearest")
    excess = rho - rho_ref - threshold
    # First index along ``depth`` (and below ``reference_depth``) where
    # ``excess > 0``. ``argmax`` over a boolean field returns the first
    # ``True`` index; if no level satisfies the criterion we fall back
    # to the deepest level.
    below_ref = (z >= reference_depth).astype(bool)
    crossing = (excess > 0.0) & below_ref
    found_any = crossing.any(dim=depth)
    first_idx = crossing.argmax(dim=depth)
    last_idx = xr.full_like(first_idx, fill_value=z.size - 1)
    pick_idx = xr.where(found_any, first_idx, last_idx)
    mld = z.isel({depth: pick_idx})
    mld.attrs.update(
        long_name="Mixed Layer Depth",
        standard_name="ocean_mixed_layer_thickness_defined_by_sigma_t",
        units="m",
    )
    return mld.rename(None).to_dataset(name="mld")

lapse_rate

lapse_rate(ds: Dataset, temperature: str = 'T', depth: str = 'depth', positive: str = 'down') -> Dataset

Vertical lapse rate Γ = −∂T/∂z (z positive upward).

The sign convention matches atmospheric usage: a column with temperature decreasing upward has a positive lapse rate. We flip the sign automatically when positive="down" (oceanographic depth convention) so the returned values still follow the upward-z definition.

Parameters:

Name Type Description Default
ds Dataset

Dataset with a temperature (or potential-temperature) variable on a vertical coordinate.

required
temperature str

Variable name (e.g. "T", "theta").

'T'
depth str

Vertical coordinate name.

'depth'
positive str

"down" (oceanographic depth, default) or "up" (height-above-surface).

'down'

Returns:

Type Description
Dataset

Dataset with a single variable "lapse_rate".

Source code in packages/xrtoolz-ops/src/xrops/ocn/_src/kinematics.py
def lapse_rate(
    ds: xr.Dataset,
    temperature: str = "T",
    depth: str = "depth",
    positive: str = "down",
) -> xr.Dataset:
    """Vertical lapse rate ``Γ = −∂T/∂z`` (z positive upward).

    The sign convention matches atmospheric usage: a column with
    temperature decreasing upward has a positive lapse rate. We flip
    the sign automatically when ``positive="down"`` (oceanographic
    depth convention) so the returned values still follow the
    upward-z definition.

    Args:
        ds: Dataset with a temperature (or potential-temperature)
            variable on a vertical coordinate.
        temperature: Variable name (e.g. ``"T"``, ``"theta"``).
        depth: Vertical coordinate name.
        positive: ``"down"`` (oceanographic depth, default) or ``"up"``
            (height-above-surface).

    Returns:
        Dataset with a single variable ``"lapse_rate"``.
    """
    if positive not in ("down", "up"):
        raise ValueError(f"positive={positive!r} must be 'down' or 'up'.")
    dT_dz = grad.partial(ds[temperature], depth, geometry="rectilinear")
    sign = 1.0 if positive == "down" else -1.0
    gamma = sign * dT_dz
    gamma.attrs.update(
        long_name="Lapse Rate",
        standard_name="air_temperature_lapse_rate",
    )
    return gamma.rename(None).to_dataset(name="lapse_rate")

density_from_ts

density_from_ts(ds: Dataset, *, salinity: str = 'so', temperature: str = 'thetao', pressure: str | float | DataArray = 0.0, lon: str = 'lon', lat: str = 'lat', eos: str = 'teos10') -> Dataset

In-situ seawater density via TEOS-10 (lazy gsw import).

Parameters:

Name Type Description Default
ds Dataset

Dataset with practical salinity salinity (PSU) and potential temperature temperature (°C).

required
salinity, temperature

Variable names.

required
pressure str | float | DataArray

Sea pressure in dbar — either a variable name in ds, a scalar (default 0.0 for surface density), or an :class:xr.DataArray that broadcasts against the salinity field.

0.0
lon, lat

Coordinate names used by gsw.SA_from_SP. The coords may be 1-D; they are broadcast to the salinity shape before calling :mod:gsw.

required
eos str

"teos10" (default, via :mod:gsw). Reserved for future "linear" / "jmd" alternatives.

'teos10'

Returns:

Type Description
Dataset

Dataset with a single variable "rho" (kg/m³).

Raises:

Type Description
ImportError

If :mod:gsw is not installed (it is an optional extra; install via pip install xrtoolz-ops[oceanography] or pip install gsw).

Source code in packages/xrtoolz-ops/src/xrops/ocn/_src/kinematics.py
def density_from_ts(
    ds: xr.Dataset,
    *,
    salinity: str = "so",
    temperature: str = "thetao",
    pressure: str | float | xr.DataArray = 0.0,
    lon: str = "lon",
    lat: str = "lat",
    eos: str = "teos10",
) -> xr.Dataset:
    """In-situ seawater density via TEOS-10 (lazy ``gsw`` import).

    Args:
        ds: Dataset with practical salinity ``salinity`` (PSU) and
            potential temperature ``temperature`` (°C).
        salinity, temperature: Variable names.
        pressure: Sea pressure in dbar — either a variable name in
            ``ds``, a scalar (default ``0.0`` for surface density), or
            an :class:`xr.DataArray` that broadcasts against the
            salinity field.
        lon, lat: Coordinate names used by ``gsw.SA_from_SP``. The
            coords may be 1-D; they are broadcast to the salinity
            shape before calling :mod:`gsw`.
        eos: ``"teos10"`` (default, via :mod:`gsw`). Reserved for future
            ``"linear"`` / ``"jmd"`` alternatives.

    Returns:
        Dataset with a single variable ``"rho"`` (kg/m³).

    Raises:
        ImportError: If :mod:`gsw` is not installed (it is an optional
            extra; install via ``pip install xrtoolz-ops[oceanography]``
            or ``pip install gsw``).
    """
    if eos != "teos10":
        raise NotImplementedError(
            f"eos={eos!r} is not implemented; only 'teos10' is supported. "
            f"Open an issue if you need a linear / JMD EOS."
        )
    import importlib

    try:
        gsw = importlib.import_module("gsw")
    except ImportError as exc:
        raise ImportError(
            "density_from_ts requires the optional 'gsw' dependency. "
            "Install with: pip install xrtoolz-ops[oceanography] (or pip install gsw)."
        ) from exc

    sp = ds[salinity]
    pt = ds[temperature]
    if isinstance(pressure, str):
        p_input: float | xr.DataArray = ds[pressure]
    else:
        p_input = pressure
    lon_da = ds[lon]
    lat_da = ds[lat]
    sp_b, lon_b, lat_b = xr.broadcast(sp, lon_da, lat_da)
    if isinstance(p_input, xr.DataArray):
        sp_b, p_b = xr.broadcast(sp_b, p_input)
    else:
        p_b = p_input
    sa = gsw.SA_from_SP(sp_b, p_b, lon_b, lat_b)
    ct = gsw.CT_from_pt(sa, pt)
    rho = gsw.rho(sa, ct, p_b)
    rho_da = xr.DataArray(rho, coords=sp_b.coords, dims=sp_b.dims, name="rho")
    rho_da.attrs.update(
        long_name="In-situ Seawater Density (TEOS-10)",
        standard_name="sea_water_density",
        units="kg m-3",
    )
    return rho_da.to_dataset()

calculate_ssh_alongtrack

calculate_ssh_alongtrack(ds: Dataset, variable: str = 'ssh', sla: str = 'sla_filtered', mdt: str = 'mdt', lwe: str | None = 'lwe') -> Dataset

Compose along-track SSH from SLA + MDT (- LWE, optional).

Equivalent altimetry-convention formula:

ssh = sla + mdt - lwe (when lwe is given)

ssh = sla + mdt (when lwe is None)

Upstream sla_to_ssh convention (no LWE correction) is reproduced by passing lwe=None::

calculate_ssh_alongtrack(ds, sla="sla", mdt="mdt", lwe=None)

When MDT lives in a separate file it must first be regridded onto the SLA grid before the datasets can be merged. The canonical pattern is::

from xrops.interpolate import regrid_like
from xrops.ocn import calculate_ssh_alongtrack

ds_mdt_on_sla = regrid_like(ds_mdt, ds_sla, dims=("lat", "lon"))
ds = ds_sla.assign(mdt=ds_mdt_on_sla["mdt"])
ds = calculate_ssh_alongtrack(ds, sla="sla", mdt="mdt", lwe=None)

Parameters:

Name Type Description Default
ds Dataset

Dataset containing the sla, mdt, and (optionally) lwe variables.

required
variable str

Name under which to store the SSH output.

'ssh'
sla str

Sea-level anomaly variable name.

'sla_filtered'
mdt str

Mean dynamic topography variable name.

'mdt'
lwe str | None

Land-water equivalent correction variable name. Pass None to skip the LWE correction (simple sla + mdt convention).

'lwe'

Returns:

Type Description
Dataset

ds with variable added.

Source code in packages/xrtoolz-ops/src/xrops/ocn/_src/ssh.py
def calculate_ssh_alongtrack(
    ds: xr.Dataset,
    variable: str = "ssh",
    sla: str = "sla_filtered",
    mdt: str = "mdt",
    lwe: str | None = "lwe",
) -> xr.Dataset:
    """Compose along-track SSH from SLA + MDT (- LWE, optional).

    Equivalent altimetry-convention formula:

    ``ssh = sla + mdt - lwe``  (when *lwe* is given)

    ``ssh = sla + mdt``        (when *lwe* is ``None``)

    Upstream ``sla_to_ssh`` convention (no LWE correction) is reproduced by
    passing ``lwe=None``::

        calculate_ssh_alongtrack(ds, sla="sla", mdt="mdt", lwe=None)

    When MDT lives in a separate file it must first be regridded onto the
    SLA grid before the datasets can be merged.  The canonical pattern is::

        from xrops.interpolate import regrid_like
        from xrops.ocn import calculate_ssh_alongtrack

        ds_mdt_on_sla = regrid_like(ds_mdt, ds_sla, dims=("lat", "lon"))
        ds = ds_sla.assign(mdt=ds_mdt_on_sla["mdt"])
        ds = calculate_ssh_alongtrack(ds, sla="sla", mdt="mdt", lwe=None)

    Args:
        ds: Dataset containing the ``sla``, ``mdt``, and (optionally) ``lwe``
            variables.
        variable: Name under which to store the SSH output.
        sla: Sea-level anomaly variable name.
        mdt: Mean dynamic topography variable name.
        lwe: Land-water equivalent correction variable name.  Pass ``None``
            to skip the LWE correction (simple ``sla + mdt`` convention).

    Returns:
        ``ds`` with ``variable`` added.
    """
    ds = ds.copy()
    ssh = ds[sla] + ds[mdt]
    if lwe is not None:
        ssh = ssh - ds[lwe]
    ds[variable] = ssh
    ds[variable].attrs.update(
        units="m",
        standard_name="sea_surface_height",
        long_name="Sea Surface Height",
    )
    return ds

calculate_ssh_unfiltered

calculate_ssh_unfiltered(ds: Dataset, variable: str = 'ssh', sla: str = 'sla_unfiltered', mdt: str = 'mdt', lwe: str | None = 'lwe') -> Dataset

Same as :func:calculate_ssh_alongtrack but from the unfiltered SLA.

Source code in packages/xrtoolz-ops/src/xrops/ocn/_src/ssh.py
def calculate_ssh_unfiltered(
    ds: xr.Dataset,
    variable: str = "ssh",
    sla: str = "sla_unfiltered",
    mdt: str = "mdt",
    lwe: str | None = "lwe",
) -> xr.Dataset:
    """Same as :func:`calculate_ssh_alongtrack` but from the unfiltered SLA."""
    return calculate_ssh_alongtrack(ds, variable=variable, sla=sla, mdt=mdt, lwe=lwe)

validate_ssh

validate_ssh(ds: Dataset, variable: str = 'ssh') -> Dataset

Attach CF-style attrs to an SSH variable.

Source code in packages/xrtoolz-ops/src/xrops/ocn/_src/validation.py
def validate_ssh(ds: xr.Dataset, variable: str = "ssh") -> xr.Dataset:
    """Attach CF-style attrs to an SSH variable."""
    ds = ds.copy()
    ds[variable] = ds[variable].assign_attrs(
        units="m",
        standard_name="sea_surface_height",
        long_name="Sea Surface Height",
    )
    return ds

validate_velocity

validate_velocity(ds: Dataset, u: str = 'u', v: str = 'v') -> Dataset

Attach CF-style attrs to zonal / meridional velocity variables.

Source code in packages/xrtoolz-ops/src/xrops/ocn/_src/validation.py
def validate_velocity(
    ds: xr.Dataset,
    u: str = "u",
    v: str = "v",
) -> xr.Dataset:
    """Attach CF-style attrs to zonal / meridional velocity variables."""
    ds = ds.copy()
    ds[u] = ds[u].assign_attrs(
        units="m s-1",
        standard_name="sea_water_x_velocity",
        long_name="Zonal Velocity",
    )
    ds[v] = ds[v].assign_attrs(
        units="m s-1",
        standard_name="sea_water_y_velocity",
        long_name="Meridional Velocity",
    )
    return ds