Windspharm

Contents

Windspharm#

Spherical harmonic vector wind analysis for Skyborn.

This module is based on the ajdawson/windspharm project (ajdawson/windspharm), originally authored by Andrew Dawson. The current version is maintained by Qianye Su and is licensed under the BSD-3-Clause.

Main Classes:

VectorWind: Enhanced interface for wind analysis on regular grids ReducedVectorWind: Packed reduced-Gaussian wind analysis

Example

>>> from skyborn.windspharm import VectorWind
>>> import numpy as np
>>>
>>> # Create sample wind data
>>> nlat, nlon = 73, 144
>>> u = np.random.randn(nlat, nlon)
>>> v = np.random.randn(nlat, nlon)
>>>
>>> # Initialize VectorWind with type hints and modern interface
>>> vw = VectorWind(u, v, gridtype='gaussian')
>>>
>>> # Calculate various fields with improved documentation
>>> vorticity = vw.vorticity()
>>> divergence = vw.divergence()
>>> psi, chi = vw.sfvp()

Notes

The NumPy-facing VectorWind interface expects latitude to run north-to-south. If your array is south-to-north, reorder it first with skyborn.windspharm.tools.order_latdim or use the xarray interface, which handles latitude reordering automatically.

class VectorWind(u, v, gridtype='regular', rsphere=6.3712e6, legfunc='stored', precision='auto')[source]#

Bases: object

Vector wind analysis using spherical harmonic transforms.

This class provides methods for analyzing atmospheric wind fields through spherical harmonic decomposition. It can compute vorticity, divergence, streamfunction, velocity potential, and perform Helmholtz decomposition.

Parameters:
  • u (array_like) – Zonal wind component. Shape should be (nlat, nlon) or (nlat, nlon, *extra_dims), where nlat is latitude points and nlon is longitude points. Latitude dimension should be north-to-south. This follows the underlying SPHEREPACK convention that the first latitude row is the northernmost point and latitude indices increase southward.

  • v (array_like) – Meridional wind component. Must have same shape as u.

  • gridtype ({'regular', 'gaussian'}, default 'regular') – Type of input grid: - ‘regular’: evenly-spaced lat/lon grid - ‘gaussian’: Gaussian latitude grid

  • rsphere (float, default 6.3712e6) – Earth radius in meters for spherical harmonic computations.

  • legfunc ({'stored', 'computed'}, default 'stored') – Legendre function computation method: - ‘stored’: precompute and store (faster, more memory) - ‘computed’: compute on-the-fly (slower, less memory)

  • precision ({'auto', 'single', 'double'}, default 'auto') – Public output precision mode. ‘auto’ preserves the promoted input floating precision, ‘single’ returns float32/complex64 outputs, and ‘double’ returns float64/complex128 outputs.

u#

Zonal wind component

Type:

ndarray

v#

Meridional wind component

Type:

ndarray

gridtype#

Grid type used

Type:

str

rsphere#

Sphere radius in meters used by the underlying transform

Type:

float

legfunc#

Legendre function handling mode used by the underlying transform

Type:

str

precision#

Public output precision mode

Type:

str

s#

Spherical harmonic transform object

Type:

Spharmt

Examples

>>> import numpy as np
>>> from skyborn.windspharm import VectorWind
>>>
>>> # Create sample wind field
>>> nlat, nlon = 73, 144
>>> u = np.random.randn(nlat, nlon)
>>> v = np.random.randn(nlat, nlon)
>>>
>>> # Initialize VectorWind
>>> vw = VectorWind(u, v, gridtype='gaussian')
>>>
>>> # Calculate dynamical quantities
>>> vorticity = vw.vorticity()
>>> divergence = vw.divergence()
>>> streamfunction = vw.streamfunction()
>>> velocity_potential = vw.velocitypotential()

Notes

The NumPy standard.VectorWind interface does not automatically reorder latitude. If your data is south-to-north, reverse it first or use skyborn.windspharm.tools.order_latdim(). The xarray interface handles this reordering automatically.

__init__(u, v, gridtype='regular', rsphere=6.3712e6, legfunc='stored', precision='auto')[source]#

Initialize VectorWind instance with comprehensive data validation.

This method performs thorough validation of input wind components including checks for missing values, infinite values, shape compatibility, and dimensional requirements.

Parameters:
  • u (ndarray | MaskedArray) – Zonal wind component with latitude on axis 0 and longitude on axis 1.

  • v (ndarray | MaskedArray) – Meridional wind component with the same shape as u.

  • gridtype (Literal['regular', 'gaussian']) – Grid type - "regular" or "gaussian".

  • rsphere (float) – Sphere radius in meters.

  • legfunc (Literal['stored', 'computed']) – Legendre function handling mode - "stored" or "computed".

  • precision (Literal['auto', 'single', 'double']) – Public output precision mode. "auto" preserves the promoted input floating precision, "single" returns float32 outputs, and "double" returns float64 outputs.

Raises:

ValueError – If the input arrays fail shape, dimensional, or finite value validation.

Return type:

None

property grid_info: dict#

Get information about the grid configuration.

Returns:

Dictionary containing grid information

Return type:

dict

magnitude()[source]#

Calculate wind speed (magnitude of vector wind).

Returns:

Wind speed field computed as sqrt(u² + v²)

Return type:

ndarray

Examples

>>> wind_speed = vw.magnitude()
vrtdiv(truncation=None)[source]#

Calculate relative vorticity and horizontal divergence.

Parameters:

truncation (int, optional) – Triangular truncation limit for spherical harmonic computation. If None, uses the default truncation based on grid resolution.

Returns:

  • vorticity (ndarray) – Relative vorticity field

  • divergence (ndarray) – Horizontal divergence field

Return type:

Tuple[ndarray, ndarray]

See also

vorticity

Calculate only vorticity

divergence

Calculate only divergence

Examples

>>> vrt, div = vw.vrtdiv()
>>> vrt_t13, div_t13 = vw.vrtdiv(truncation=13)
vorticity(truncation=None)[source]#

Calculate relative vorticity.

Parameters:

truncation (int, optional) – Triangular truncation limit for spherical harmonic computation.

Returns:

Relative vorticity field (∇ × v)

Return type:

ndarray

See also

vrtdiv

Calculate both vorticity and divergence

absolutevorticity

Calculate absolute vorticity

Examples

>>> vrt = vw.vorticity()
>>> vrt_t13 = vw.vorticity(truncation=13)
divergence(truncation=None)[source]#

Calculate horizontal divergence.

Parameters:

truncation (int, optional) – Triangular truncation limit for spherical harmonic computation.

Returns:

Horizontal divergence field (∇ · v)

Return type:

ndarray

See also

vrtdiv

Calculate both vorticity and divergence

Examples

>>> div = vw.divergence()
>>> div_t13 = vw.divergence(truncation=13)
planetaryvorticity(omega=None)[source]#

Calculate planetary vorticity (Coriolis parameter).

Parameters:

omega (float, optional) – Earth’s angular velocity in rad/s. Default is 7.292e-5 s⁻¹.

Returns:

Planetary vorticity field (2Ω sin φ)

Return type:

ndarray

See also

absolutevorticity

Calculate absolute vorticity

Examples

>>> f = vw.planetaryvorticity()
>>> f_custom = vw.planetaryvorticity(omega=7.2921150e-5)
absolutevorticity(omega=None, truncation=None)[source]#

Calculate absolute vorticity (relative + planetary vorticity).

Parameters:
  • omega (float, optional) – Earth’s angular velocity in rad/s. Default is 7.292e-5 s⁻¹.

  • truncation (int, optional) – Triangular truncation limit for spherical harmonic computation.

Returns:

Absolute vorticity field (ζ + f)

Return type:

ndarray

See also

vorticity

Calculate relative vorticity

planetaryvorticity

Calculate planetary vorticity

Examples

>>> abs_vrt = vw.absolutevorticity()
>>> abs_vrt_t13 = vw.absolutevorticity(omega=7.2921150e-5, truncation=13)
sfvp(truncation=None)[source]#

Calculate streamfunction and velocity potential.

Parameters:

truncation (int, optional) – Triangular truncation limit for spherical harmonic computation.

Returns:

  • streamfunction (ndarray) – Streamfunction field (ψ)

  • velocity_potential (ndarray) – Velocity potential field (χ)

Return type:

Tuple[ndarray, ndarray]

See also

streamfunction

Calculate only streamfunction

velocitypotential

Calculate only velocity potential

Examples

>>> psi, chi = vw.sfvp()
>>> psi_t13, chi_t13 = vw.sfvp(truncation=13)
streamfunction(truncation=None)[source]#

Calculate streamfunction.

Parameters:

truncation (int, optional) – Triangular truncation limit for spherical harmonic computation.

Returns:

Streamfunction field (ψ)

Return type:

ndarray

See also

sfvp

Calculate both streamfunction and velocity potential

Examples

>>> psi = vw.streamfunction()
>>> psi_t13 = vw.streamfunction(truncation=13)
velocitypotential(truncation=None)[source]#

Calculate velocity potential.

Parameters:

truncation (int, optional) – Triangular truncation limit for spherical harmonic computation.

Returns:

Velocity potential field (χ)

Return type:

ndarray

See also

sfvp

Calculate both streamfunction and velocity potential

Examples

>>> chi = vw.velocitypotential()
>>> chi_t13 = vw.velocitypotential(truncation=13)
helmholtz(truncation=None)[source]#

Perform Helmholtz decomposition of vector wind.

Decomposes the wind field into irrotational (divergent) and non-divergent (rotational) components.

Parameters:

truncation (int, optional) – Triangular truncation limit for spherical harmonic computation.

Returns:

  • u_chi (ndarray) – Zonal component of irrotational wind

  • v_chi (ndarray) – Meridional component of irrotational wind

  • u_psi (ndarray) – Zonal component of non-divergent wind

  • v_psi (ndarray) – Meridional component of non-divergent wind

Return type:

Tuple[ndarray, ndarray, ndarray, ndarray]

See also

irrotationalcomponent

Get only irrotational component

nondivergentcomponent

Get only non-divergent component

Examples

>>> u_chi, v_chi, u_psi, v_psi = vw.helmholtz()
>>> u_chi_t13, v_chi_t13, u_psi_t13, v_psi_t13 = vw.helmholtz(truncation=13)
irrotationalcomponent(truncation=None)[source]#

Calculate irrotational (divergent) component of vector wind.

Note

If both irrotational and non-divergent components are needed, use helmholtz() method for efficiency.

Parameters:

truncation (int, optional) – Triangular truncation limit for spherical harmonic computation.

Returns:

  • u_chi (ndarray) – Zonal component of irrotational wind

  • v_chi (ndarray) – Meridional component of irrotational wind

Return type:

Tuple[ndarray, ndarray]

See also

helmholtz

Complete Helmholtz decomposition

nondivergentcomponent

Non-divergent component

Examples

>>> u_chi, v_chi = vw.irrotationalcomponent()
>>> u_chi_t13, v_chi_t13 = vw.irrotationalcomponent(truncation=13)
nondivergentcomponent(truncation=None)[source]#

Calculate non-divergent (rotational) component of vector wind.

Note

If both non-divergent and irrotational components are needed, use helmholtz() method for efficiency.

Parameters:

truncation (int, optional) – Triangular truncation limit for spherical harmonic computation.

Returns:

  • u_psi (ndarray) – Zonal component of non-divergent wind

  • v_psi (ndarray) – Meridional component of non-divergent wind

Return type:

Tuple[ndarray, ndarray]

See also

helmholtz

Complete Helmholtz decomposition

irrotationalcomponent

Irrotational component

Examples

>>> u_psi, v_psi = vw.nondivergentcomponent()
>>> u_psi_t13, v_psi_t13 = vw.nondivergentcomponent(truncation=13)
gradient(chi, truncation=None)[source]#

Calculate vector gradient of a scalar field on the sphere.

Parameters:
  • chi (array_like) – Scalar field with shape (nlat, nlon) or (nlat, nlon, nfields) matching the dimensions of the VectorWind instance.

  • truncation (int, optional) – Triangular truncation limit for spherical harmonic computation.

Returns:

  • u_chi (ndarray) – Zonal component of vector gradient (∂χ/∂λ)

  • v_chi (ndarray) – Meridional component of vector gradient (∂χ/∂φ)

Return type:

Tuple[ndarray, ndarray]

Examples

>>> abs_vrt = vw.absolutevorticity()
>>> avrt_u, avrt_v = vw.gradient(abs_vrt)
>>> avrt_u_t13, avrt_v_t13 = vw.gradient(abs_vrt, truncation=13)
truncate(field, truncation=None)[source]#

Apply spectral truncation to a scalar field.

This is useful to represent fields consistently with the output of other VectorWind methods.

Parameters:
  • field (array_like) – Scalar field with shape (nlat, nlon) or (nlat, nlon, nfields) matching the dimensions of the VectorWind instance.

  • truncation (int, optional) – Triangular truncation limit. If None, defaults to nlat-1.

Returns:

Field with spectral truncation applied

Return type:

ndarray

Examples

>>> field_trunc = vw.truncate(scalar_field)
>>> field_t21 = vw.truncate(scalar_field, truncation=21)
rossbywavesource(truncation=None, omega=None)[source]#

Calculate Rossby wave source.

The Rossby wave source is defined as: S = -ζₐ∇·v - v_χ·∇ζₐ

where: - ζₐ is absolute vorticity (relative + planetary) - ∇·v is horizontal divergence - v_χ is the irrotational (divergent) wind component - ∇ζₐ is the gradient of absolute vorticity

Parameters:
  • truncation (int, optional) – Triangular truncation limit for spherical harmonic computation. If None, uses the default truncation based on grid resolution.

  • omega (float, optional) – Earth’s angular velocity in rad/s. Default is 7.292e-5 s⁻¹.

Returns:

Rossby wave source field (s⁻²)

Return type:

ndarray

See also

absolutevorticity

Calculate absolute vorticity

divergence

Calculate horizontal divergence

irrotationalcomponent

Calculate irrotational wind component

gradient

Calculate vector gradient

Notes

The Rossby wave source quantifies the generation of Rossby wave activity in the atmosphere. Positive values indicate wave generation, while negative values indicate wave absorption or dissipation.

The calculation involves several steps: 1. Compute absolute vorticity (relative + planetary) 2. Compute horizontal divergence 3. Compute irrotational wind components 4. Compute gradients of absolute vorticity 5. Combine terms according to the RWS formula

Examples

>>> rws = vw.rossbywavesource()
>>> rws_t21 = vw.rossbywavesource(truncation=21)
>>> rws_custom_omega = vw.rossbywavesource(omega=7.2921150e-5)

References

Sardeshmukh, P. D., & Hoskins, B. J. (1988). The generation of global rotational flow by steady idealized tropical heating. Journal of the Atmospheric Sciences, 45(7), 1228-1251.

class ReducedVectorWind(u, v, pl, rsphere=6.3712e6, legfunc='stored', precision='auto')[source]#

Bases: object

Vector wind analysis on packed reduced Gaussian grids.

Parameters:
  • u (ArrayLike) – Packed zonal and meridional wind components with shape (sum(pl), ...).

  • v (ArrayLike) – Packed zonal and meridional wind components with shape (sum(pl), ...).

  • pl (ArrayLike) – Reduced Gaussian row lengths, ordered north-to-south.

  • rsphere (float) – Sphere radius in meters.

  • legfunc (str) – Legendre function mode passed to grouped Gaussian Spharmt backends.

  • precision (str) – Public output precision mode. "auto" preserves the promoted input floating precision, "single" returns float32 outputs, and "double" returns float64 outputs.

u, v

Packed zonal and meridional wind components.

pl#

Reduced Gaussian row lengths ordered north-to-south.

gridtype#

Fixed grid type label "reduced_gaussian".

rsphere#

Sphere radius in meters.

legfunc#

Legendre function handling mode used by the underlying transform.

precision#

Public output precision mode.

s#

Underlying skyborn.spharm.ReducedGaussianSpharmt instance.

__init__(u, v, pl, rsphere=6.3712e6, legfunc='stored', precision='auto')[source]#

Create a reduced-Gaussian vector wind analysis wrapper.

Parameters:
  • u (ndarray | MaskedArray) – Packed zonal wind component with shape (sum(pl), ...).

  • v (ndarray | MaskedArray) – Packed meridional wind component with the same shape as u.

  • pl (ndarray | MaskedArray) – Reduced Gaussian row lengths ordered north-to-south.

  • rsphere (float) – Sphere radius in meters.

  • legfunc (str) – Legendre function handling mode - "stored" or "computed".

  • precision (str) – Public output precision mode. "auto" preserves the promoted input floating precision, "single" returns float32 outputs, and "double" returns float64 outputs.

Raises:

ValueError – If u, v, or pl fails validation.

Return type:

None

property grid_info: dict#
magnitude()[source]#

Return wind speed.

Return type:

ndarray

vrtdiv(truncation=None)[source]#

Return relative vorticity and horizontal divergence.

Parameters:

truncation (int | None)

Return type:

Tuple[ndarray, ndarray]

vorticity(truncation=None)[source]#

Return relative vorticity.

Parameters:

truncation (int | None)

Return type:

ndarray

divergence(truncation=None)[source]#

Return horizontal divergence.

Parameters:

truncation (int | None)

Return type:

ndarray

planetaryvorticity(omega=None)[source]#

Return planetary vorticity on the packed reduced grid.

Parameters:

omega (float | None)

Return type:

ndarray

absolutevorticity(omega=None, truncation=None)[source]#

Return absolute vorticity.

Parameters:
  • omega (float | None)

  • truncation (int | None)

Return type:

ndarray

sfvp(truncation=None)[source]#

Return streamfunction and velocity potential.

Parameters:

truncation (int | None)

Return type:

Tuple[ndarray, ndarray]

streamfunction(truncation=None)[source]#

Return streamfunction.

Parameters:

truncation (int | None)

Return type:

ndarray

velocitypotential(truncation=None)[source]#

Return velocity potential.

Parameters:

truncation (int | None)

Return type:

ndarray

helmholtz(truncation=None)[source]#

Return irrotational and non-divergent wind components.

Parameters:

truncation (int | None)

Return type:

Tuple[ndarray, ndarray, ndarray, ndarray]

irrotationalcomponent(truncation=None)[source]#

Return irrotational wind component.

Parameters:

truncation (int | None)

Return type:

Tuple[ndarray, ndarray]

nondivergentcomponent(truncation=None)[source]#

Return non-divergent wind component.

Parameters:

truncation (int | None)

Return type:

Tuple[ndarray, ndarray]

gradient(chi, truncation=None)[source]#

Return vector gradient of a packed scalar field.

Parameters:
Return type:

Tuple[ndarray, ndarray]

truncate(field, truncation=None)[source]#

Apply spectral truncation to a packed scalar field.

Parameters:
Return type:

ndarray

rossbywavesource(truncation=None, omega=None)[source]#

Return Rossby wave source on the packed reduced grid.

Parameters:
  • truncation (int | None)

  • omega (float | None)

Return type:

ndarray

Overview#

The skyborn.windspharm module provides spherical harmonic vector wind analysis capabilities for atmospheric and oceanic sciences. This module is based on the ajdawson/windspharm project and has been enhanced with modern Python features including type hints, comprehensive error handling, and improved documentation.

Key Use Case: Divergent Wind Component Analysis The module is particularly useful for calculating divergent wind components from 4D atmospheric data (time, level, lat, lon), which is essential for understanding atmospheric circulation patterns and energy transport mechanisms.

Interface Selection#

Skyborn provides two public VectorWind interfaces:

  • from skyborn.windspharm.standard import VectorWind Standard NumPy interface for array inputs.

  • from skyborn.windspharm.xarray import VectorWind xarray-aware interface for xarray.DataArray inputs that preserves dimensions, coordinates, and metadata.

At present, from skyborn.windspharm import VectorWind also resolves to skyborn.windspharm.standard.VectorWind. For clarity in user code and examples, the documentation uses the explicit standard import for NumPy inputs.

Use the standard interface for plain arrays:

import numpy as np
from skyborn.windspharm.standard import VectorWind

u = np.random.randn(73, 144)
v = np.random.randn(73, 144)
vw = VectorWind(u, v, gridtype="gaussian")

Use the xarray interface for DataArray inputs:

import xarray as xr
from skyborn.windspharm.xarray import VectorWind

u = xr.open_dataarray("u_wind.nc")
v = xr.open_dataarray("v_wind.nc")
vw = VectorWind(u, v)

Latitude Ordering#

The NumPy standard.VectorWind interface expects latitude to be ordered north-to-south. This is not an arbitrary Skyborn preference: it follows the underlying SPHEREPACK convention used by spharm, where the first latitude row is interpreted as the northernmost point and latitude indices increase southward.

If a south-to-north array is passed directly into the standard interface, the backend will still interpret the first row as the north side of the grid, which leads to incorrect spherical harmonic analysis and large diagnostic errors.

Recommended usage:

  • NumPy arrays: reorder first with skyborn.windspharm.tools.order_latdim.

  • xarray interface: latitude is reordered automatically before calling the standard backend.

Main Features#

  • Vector Wind Analysis: Compute vorticity, divergence, streamfunction, and velocity potential

  • Helmholtz Decomposition: Separate wind fields into rotational and divergent components

  • Data Preparation Tools: Robust prep_data and recover_data utilities for handling multi-dimensional arrays

  • Parallel Processing Support: Efficient computation for large time-series datasets

  • Multiple Grid Support: Works with regular and Gaussian latitude grids

  • Efficient Computation: Uses optimized spherical harmonic transforms

  • Comprehensive Validation: Thorough input validation and error handling

Real-World Application Examples#

Computing Divergent Wind Components from 4D Data

A common use case involves processing atmospheric wind data with shape (time, level, lat, lon) to extract divergent wind components for circulation analysis. This requires careful data preparation and can benefit from parallel processing for large datasets.

Core Classes#

VectorWind#

class VectorWind(u, v, gridtype='regular', rsphere=6.3712e6, legfunc='stored', precision='auto')[source]#

Bases: object

Vector wind analysis using spherical harmonic transforms.

This class provides methods for analyzing atmospheric wind fields through spherical harmonic decomposition. It can compute vorticity, divergence, streamfunction, velocity potential, and perform Helmholtz decomposition.

Parameters:
  • u (array_like) – Zonal wind component. Shape should be (nlat, nlon) or (nlat, nlon, *extra_dims), where nlat is latitude points and nlon is longitude points. Latitude dimension should be north-to-south. This follows the underlying SPHEREPACK convention that the first latitude row is the northernmost point and latitude indices increase southward.

  • v (array_like) – Meridional wind component. Must have same shape as u.

  • gridtype ({'regular', 'gaussian'}, default 'regular') – Type of input grid: - ‘regular’: evenly-spaced lat/lon grid - ‘gaussian’: Gaussian latitude grid

  • rsphere (float, default 6.3712e6) – Earth radius in meters for spherical harmonic computations.

  • legfunc ({'stored', 'computed'}, default 'stored') – Legendre function computation method: - ‘stored’: precompute and store (faster, more memory) - ‘computed’: compute on-the-fly (slower, less memory)

  • precision ({'auto', 'single', 'double'}, default 'auto') – Public output precision mode. ‘auto’ preserves the promoted input floating precision, ‘single’ returns float32/complex64 outputs, and ‘double’ returns float64/complex128 outputs.

u#

Zonal wind component

Type:

ndarray

v#

Meridional wind component

Type:

ndarray

gridtype#

Grid type used

Type:

str

rsphere#

Sphere radius in meters used by the underlying transform

Type:

float

legfunc#

Legendre function handling mode used by the underlying transform

Type:

str

precision#

Public output precision mode

Type:

str

s#

Spherical harmonic transform object

Type:

Spharmt

Examples

>>> import numpy as np
>>> from skyborn.windspharm import VectorWind
>>>
>>> # Create sample wind field
>>> nlat, nlon = 73, 144
>>> u = np.random.randn(nlat, nlon)
>>> v = np.random.randn(nlat, nlon)
>>>
>>> # Initialize VectorWind
>>> vw = VectorWind(u, v, gridtype='gaussian')
>>>
>>> # Calculate dynamical quantities
>>> vorticity = vw.vorticity()
>>> divergence = vw.divergence()
>>> streamfunction = vw.streamfunction()
>>> velocity_potential = vw.velocitypotential()

Notes

The NumPy standard.VectorWind interface does not automatically reorder latitude. If your data is south-to-north, reverse it first or use skyborn.windspharm.tools.order_latdim(). The xarray interface handles this reordering automatically.

Methods

Field Computation

VectorWind.vorticity([truncation])

Calculate relative vorticity.

VectorWind.divergence([truncation])

Calculate horizontal divergence.

VectorWind.streamfunction([truncation])

Calculate streamfunction.

VectorWind.velocitypotential([truncation])

Calculate velocity potential.

VectorWind.sfvp([truncation])

Calculate streamfunction and velocity potential.

Vector Field Operations

VectorWind.helmholtz([truncation])

Perform Helmholtz decomposition of vector wind.

VectorWind.nondivergentcomponent([truncation])

Calculate non-divergent (rotational) component of vector wind.

VectorWind.irrotationalcomponent([truncation])

Calculate irrotational (divergent) component of vector wind.

VectorWind.gradient(chi[, truncation])

Calculate vector gradient of a scalar field on the sphere.

Utility Methods

VectorWind.truncate(field[, truncation])

Apply spectral truncation to a scalar field.

__init__(u, v, gridtype='regular', rsphere=6.3712e6, legfunc='stored', precision='auto')[source]#

Initialize VectorWind instance with comprehensive data validation.

This method performs thorough validation of input wind components including checks for missing values, infinite values, shape compatibility, and dimensional requirements.

Parameters:
  • u (ndarray | MaskedArray) – Zonal wind component with latitude on axis 0 and longitude on axis 1.

  • v (ndarray | MaskedArray) – Meridional wind component with the same shape as u.

  • gridtype (Literal['regular', 'gaussian']) – Grid type - "regular" or "gaussian".

  • rsphere (float) – Sphere radius in meters.

  • legfunc (Literal['stored', 'computed']) – Legendre function handling mode - "stored" or "computed".

  • precision (Literal['auto', 'single', 'double']) – Public output precision mode. "auto" preserves the promoted input floating precision, "single" returns float32 outputs, and "double" returns float64 outputs.

Raises:

ValueError – If the input arrays fail shape, dimensional, or finite value validation.

Return type:

None

property grid_info: dict#

Get information about the grid configuration.

Returns:

Dictionary containing grid information

Return type:

dict

magnitude()[source]#

Calculate wind speed (magnitude of vector wind).

Returns:

Wind speed field computed as sqrt(u² + v²)

Return type:

ndarray

Examples

>>> wind_speed = vw.magnitude()
vrtdiv(truncation=None)[source]#

Calculate relative vorticity and horizontal divergence.

Parameters:

truncation (int, optional) – Triangular truncation limit for spherical harmonic computation. If None, uses the default truncation based on grid resolution.

Returns:

  • vorticity (ndarray) – Relative vorticity field

  • divergence (ndarray) – Horizontal divergence field

Return type:

Tuple[ndarray, ndarray]

See also

vorticity

Calculate only vorticity

divergence

Calculate only divergence

Examples

>>> vrt, div = vw.vrtdiv()
>>> vrt_t13, div_t13 = vw.vrtdiv(truncation=13)
vorticity(truncation=None)[source]#

Calculate relative vorticity.

Parameters:

truncation (int, optional) – Triangular truncation limit for spherical harmonic computation.

Returns:

Relative vorticity field (∇ × v)

Return type:

ndarray

See also

vrtdiv

Calculate both vorticity and divergence

absolutevorticity

Calculate absolute vorticity

Examples

>>> vrt = vw.vorticity()
>>> vrt_t13 = vw.vorticity(truncation=13)
divergence(truncation=None)[source]#

Calculate horizontal divergence.

Parameters:

truncation (int, optional) – Triangular truncation limit for spherical harmonic computation.

Returns:

Horizontal divergence field (∇ · v)

Return type:

ndarray

See also

vrtdiv

Calculate both vorticity and divergence

Examples

>>> div = vw.divergence()
>>> div_t13 = vw.divergence(truncation=13)
planetaryvorticity(omega=None)[source]#

Calculate planetary vorticity (Coriolis parameter).

Parameters:

omega (float, optional) – Earth’s angular velocity in rad/s. Default is 7.292e-5 s⁻¹.

Returns:

Planetary vorticity field (2Ω sin φ)

Return type:

ndarray

See also

absolutevorticity

Calculate absolute vorticity

Examples

>>> f = vw.planetaryvorticity()
>>> f_custom = vw.planetaryvorticity(omega=7.2921150e-5)
absolutevorticity(omega=None, truncation=None)[source]#

Calculate absolute vorticity (relative + planetary vorticity).

Parameters:
  • omega (float, optional) – Earth’s angular velocity in rad/s. Default is 7.292e-5 s⁻¹.

  • truncation (int, optional) – Triangular truncation limit for spherical harmonic computation.

Returns:

Absolute vorticity field (ζ + f)

Return type:

ndarray

See also

vorticity

Calculate relative vorticity

planetaryvorticity

Calculate planetary vorticity

Examples

>>> abs_vrt = vw.absolutevorticity()
>>> abs_vrt_t13 = vw.absolutevorticity(omega=7.2921150e-5, truncation=13)
sfvp(truncation=None)[source]#

Calculate streamfunction and velocity potential.

Parameters:

truncation (int, optional) – Triangular truncation limit for spherical harmonic computation.

Returns:

  • streamfunction (ndarray) – Streamfunction field (ψ)

  • velocity_potential (ndarray) – Velocity potential field (χ)

Return type:

Tuple[ndarray, ndarray]

See also

streamfunction

Calculate only streamfunction

velocitypotential

Calculate only velocity potential

Examples

>>> psi, chi = vw.sfvp()
>>> psi_t13, chi_t13 = vw.sfvp(truncation=13)
streamfunction(truncation=None)[source]#

Calculate streamfunction.

Parameters:

truncation (int, optional) – Triangular truncation limit for spherical harmonic computation.

Returns:

Streamfunction field (ψ)

Return type:

ndarray

See also

sfvp

Calculate both streamfunction and velocity potential

Examples

>>> psi = vw.streamfunction()
>>> psi_t13 = vw.streamfunction(truncation=13)
velocitypotential(truncation=None)[source]#

Calculate velocity potential.

Parameters:

truncation (int, optional) – Triangular truncation limit for spherical harmonic computation.

Returns:

Velocity potential field (χ)

Return type:

ndarray

See also

sfvp

Calculate both streamfunction and velocity potential

Examples

>>> chi = vw.velocitypotential()
>>> chi_t13 = vw.velocitypotential(truncation=13)
helmholtz(truncation=None)[source]#

Perform Helmholtz decomposition of vector wind.

Decomposes the wind field into irrotational (divergent) and non-divergent (rotational) components.

Parameters:

truncation (int, optional) – Triangular truncation limit for spherical harmonic computation.

Returns:

  • u_chi (ndarray) – Zonal component of irrotational wind

  • v_chi (ndarray) – Meridional component of irrotational wind

  • u_psi (ndarray) – Zonal component of non-divergent wind

  • v_psi (ndarray) – Meridional component of non-divergent wind

Return type:

Tuple[ndarray, ndarray, ndarray, ndarray]

See also

irrotationalcomponent

Get only irrotational component

nondivergentcomponent

Get only non-divergent component

Examples

>>> u_chi, v_chi, u_psi, v_psi = vw.helmholtz()
>>> u_chi_t13, v_chi_t13, u_psi_t13, v_psi_t13 = vw.helmholtz(truncation=13)
irrotationalcomponent(truncation=None)[source]#

Calculate irrotational (divergent) component of vector wind.

Note

If both irrotational and non-divergent components are needed, use helmholtz() method for efficiency.

Parameters:

truncation (int, optional) – Triangular truncation limit for spherical harmonic computation.

Returns:

  • u_chi (ndarray) – Zonal component of irrotational wind

  • v_chi (ndarray) – Meridional component of irrotational wind

Return type:

Tuple[ndarray, ndarray]

See also

helmholtz

Complete Helmholtz decomposition

nondivergentcomponent

Non-divergent component

Examples

>>> u_chi, v_chi = vw.irrotationalcomponent()
>>> u_chi_t13, v_chi_t13 = vw.irrotationalcomponent(truncation=13)
nondivergentcomponent(truncation=None)[source]#

Calculate non-divergent (rotational) component of vector wind.

Note

If both non-divergent and irrotational components are needed, use helmholtz() method for efficiency.

Parameters:

truncation (int, optional) – Triangular truncation limit for spherical harmonic computation.

Returns:

  • u_psi (ndarray) – Zonal component of non-divergent wind

  • v_psi (ndarray) – Meridional component of non-divergent wind

Return type:

Tuple[ndarray, ndarray]

See also

helmholtz

Complete Helmholtz decomposition

irrotationalcomponent

Irrotational component

Examples

>>> u_psi, v_psi = vw.nondivergentcomponent()
>>> u_psi_t13, v_psi_t13 = vw.nondivergentcomponent(truncation=13)
gradient(chi, truncation=None)[source]#

Calculate vector gradient of a scalar field on the sphere.

Parameters:
  • chi (array_like) – Scalar field with shape (nlat, nlon) or (nlat, nlon, nfields) matching the dimensions of the VectorWind instance.

  • truncation (int, optional) – Triangular truncation limit for spherical harmonic computation.

Returns:

  • u_chi (ndarray) – Zonal component of vector gradient (∂χ/∂λ)

  • v_chi (ndarray) – Meridional component of vector gradient (∂χ/∂φ)

Return type:

Tuple[ndarray, ndarray]

Examples

>>> abs_vrt = vw.absolutevorticity()
>>> avrt_u, avrt_v = vw.gradient(abs_vrt)
>>> avrt_u_t13, avrt_v_t13 = vw.gradient(abs_vrt, truncation=13)
truncate(field, truncation=None)[source]#

Apply spectral truncation to a scalar field.

This is useful to represent fields consistently with the output of other VectorWind methods.

Parameters:
  • field (array_like) – Scalar field with shape (nlat, nlon) or (nlat, nlon, nfields) matching the dimensions of the VectorWind instance.

  • truncation (int, optional) – Triangular truncation limit. If None, defaults to nlat-1.

Returns:

Field with spectral truncation applied

Return type:

ndarray

Examples

>>> field_trunc = vw.truncate(scalar_field)
>>> field_t21 = vw.truncate(scalar_field, truncation=21)
rossbywavesource(truncation=None, omega=None)[source]#

Calculate Rossby wave source.

The Rossby wave source is defined as: S = -ζₐ∇·v - v_χ·∇ζₐ

where: - ζₐ is absolute vorticity (relative + planetary) - ∇·v is horizontal divergence - v_χ is the irrotational (divergent) wind component - ∇ζₐ is the gradient of absolute vorticity

Parameters:
  • truncation (int, optional) – Triangular truncation limit for spherical harmonic computation. If None, uses the default truncation based on grid resolution.

  • omega (float, optional) – Earth’s angular velocity in rad/s. Default is 7.292e-5 s⁻¹.

Returns:

Rossby wave source field (s⁻²)

Return type:

ndarray

See also

absolutevorticity

Calculate absolute vorticity

divergence

Calculate horizontal divergence

irrotationalcomponent

Calculate irrotational wind component

gradient

Calculate vector gradient

Notes

The Rossby wave source quantifies the generation of Rossby wave activity in the atmosphere. Positive values indicate wave generation, while negative values indicate wave absorption or dissipation.

The calculation involves several steps: 1. Compute absolute vorticity (relative + planetary) 2. Compute horizontal divergence 3. Compute irrotational wind components 4. Compute gradients of absolute vorticity 5. Combine terms according to the RWS formula

Examples

>>> rws = vw.rossbywavesource()
>>> rws_t21 = vw.rossbywavesource(truncation=21)
>>> rws_custom_omega = vw.rossbywavesource(omega=7.2921150e-5)

References

Sardeshmukh, P. D., & Hoskins, B. J. (1988). The generation of global rotational flow by steady idealized tropical heating. Journal of the Atmospheric Sciences, 45(7), 1228-1251.

Standard Interface#

Spherical harmonic vector wind computations.

This module provides the main VectorWind class for analyzing atmospheric wind fields using spherical harmonic transforms. It enables calculation of vorticity, divergence, streamfunction, velocity potential, and other dynamical quantities.

Xarray Interface#

Spherical harmonic vector wind computations with xarray interface.

This module provides a VectorWind class that works with xarray DataArrays, preserving coordinate information and metadata throughout the computation process. It serves as a high-level interface to the standard VectorWind implementation.

Main Class:

VectorWind: xarray-aware interface for wind field analysis

Example

>>> import xarray as xr
>>> from skyborn.windspharm.xarray import VectorWind
>>>
>>> # Load wind data as xarray DataArrays
>>> u = xr.open_dataarray('u_wind.nc')
>>> v = xr.open_dataarray('v_wind.nc')
>>>
>>> # Create VectorWind instance
>>> vw = VectorWind(u, v)
>>>
>>> # Compute with preserved metadata
>>> vorticity = vw.vorticity()
>>> streamfunction = vw.streamfunction()
class ReducedVectorWind(u, v, pl, rsphere=6.3712e6, legfunc='stored', precision='auto')[source]#

Bases: object

Xarray interface for packed reduced-Gaussian vector wind analysis.

The first dimension of u and v is interpreted as the packed reduced-Gaussian point dimension and must have length sum(pl). Any remaining dimensions are carried through unchanged.

Parameters:
  • u (DataArray)

  • v (DataArray)

  • pl (Any)

  • rsphere (float)

  • legfunc (LegFunc)

  • precision (Literal['auto', 'single', 'double'])

__init__(u, v, pl, rsphere=6.3712e6, legfunc='stored', precision='auto')[source]#
Parameters:
Return type:

None

u()[source]#

Return the packed zonal wind component.

Return type:

DataArray

v()[source]#

Return the packed meridional wind component.

Return type:

DataArray

magnitude()[source]#

Return wind speed on the packed reduced grid.

Return type:

DataArray

vrtdiv(truncation=None)[source]#

Return relative vorticity and horizontal divergence.

Parameters:

truncation (int | None)

Return type:

Tuple[DataArray, DataArray]

vorticity(truncation=None)[source]#

Return relative vorticity.

Parameters:

truncation (int | None)

Return type:

DataArray

divergence(truncation=None)[source]#

Return horizontal divergence.

Parameters:

truncation (int | None)

Return type:

DataArray

planetaryvorticity(omega=None)[source]#

Return planetary vorticity.

Parameters:

omega (float | None)

Return type:

DataArray

absolutevorticity(omega=None, truncation=None)[source]#

Return absolute vorticity.

Parameters:
  • omega (float | None)

  • truncation (int | None)

Return type:

DataArray

sfvp(truncation=None)[source]#

Return streamfunction and velocity potential.

Parameters:

truncation (int | None)

Return type:

Tuple[DataArray, DataArray]

streamfunction(truncation=None)[source]#

Return streamfunction.

Parameters:

truncation (int | None)

Return type:

DataArray

velocitypotential(truncation=None)[source]#

Return velocity potential.

Parameters:

truncation (int | None)

Return type:

DataArray

helmholtz(truncation=None)[source]#

Return irrotational and non-divergent wind components.

Parameters:

truncation (int | None)

Return type:

Tuple[DataArray, DataArray, DataArray, DataArray]

irrotationalcomponent(truncation=None)[source]#

Return irrotational wind component.

Parameters:

truncation (int | None)

Return type:

Tuple[DataArray, DataArray]

nondivergentcomponent(truncation=None)[source]#

Return non-divergent wind component.

Parameters:

truncation (int | None)

Return type:

Tuple[DataArray, DataArray]

gradient(chi, truncation=None)[source]#

Return vector gradient of a packed scalar field.

Parameters:
Return type:

Tuple[DataArray, DataArray]

truncate(field, truncation=None)[source]#

Apply spectral truncation to a packed scalar field.

Parameters:
Return type:

DataArray

rossbywavesource(truncation=None, omega=None)[source]#

Return Rossby wave source on the packed reduced grid.

Parameters:
  • truncation (int | None)

  • omega (float | None)

Return type:

DataArray

class VectorWind(u, v, rsphere=6.3712e6, legfunc='stored', gridtype=None, precision='auto')[source]#

Bases: object

Vector wind analysis using xarray DataArrays.

This class provides a high-level interface for spherical harmonic wind analysis that preserves xarray coordinate information and metadata. It wraps the standard VectorWind implementation while maintaining CF-compliant attributes.

Parameters:
  • u (xarray.DataArray) – Zonal and meridional wind components. Must have the same dimensions, coordinates, and contain no missing values. Should include latitude and longitude dimensions with appropriate coordinate information.

  • v (xarray.DataArray) – Zonal and meridional wind components. Must have the same dimensions, coordinates, and contain no missing values. Should include latitude and longitude dimensions with appropriate coordinate information.

  • rsphere (float, default 6.3712e6) – Earth radius in meters for spherical harmonic computations.

  • legfunc ({'stored', 'computed'}, default 'stored') – Legendre function computation method: - ‘stored’: precompute and store (faster, more memory) - ‘computed’: compute on-the-fly (slower, less memory)

  • gridtype ({'regular', 'gaussian'}, optional) – Explicit grid type override. If omitted, the grid type is inferred from latitude coordinates.

  • precision (Literal['auto', 'single', 'double'])

_api#

Underlying standard VectorWind instance

Type:

standard.VectorWind

_reorder#

Original dimension ordering for output reconstruction

Type:

tuple

_ishape#

Original data shape

Type:

tuple

_coords#

Original coordinate information

Type:

list

Examples

>>> import xarray as xr
>>> from skyborn.windspharm.xarray import VectorWind
>>>
>>> # Load wind components
>>> u = xr.open_dataarray('u850.nc')
>>> v = xr.open_dataarray('v850.nc')
>>>
>>> # Create VectorWind instance
>>> vw = VectorWind(u, v)
>>>
>>> # Compute vorticity with preserved metadata
>>> vorticity = vw.vorticity()
>>> print(vorticity.attrs)  # CF-compliant attributes
>>>
>>> # Helmholtz decomposition
>>> u_chi, v_chi, u_psi, v_psi = vw.helmholtz()
__init__(u, v, rsphere=6.3712e6, legfunc='stored', gridtype=None, precision='auto')[source]#

Initialize VectorWind instance with comprehensive validation.

Parameters:
Return type:

None

u()[source]#

Get zonal component of vector wind.

Returns:

Zonal (eastward) wind component with CF-compliant attributes

Return type:

DataArray

Examples

>>> u_wind = vw.u()
>>> print(u_wind.attrs['standard_name'])  # 'eastward_wind'
v()[source]#

Get meridional component of vector wind.

Returns:

Meridional (northward) wind component with CF-compliant attributes

Return type:

DataArray

Examples

>>> v_wind = vw.v()
>>> print(v_wind.attrs['standard_name'])  # 'northward_wind'
magnitude()[source]#

Calculate wind speed (magnitude of vector wind).

Returns:

Wind speed with CF-compliant attributes

Return type:

DataArray

Examples

>>> wind_speed = vw.magnitude()
>>> print(wind_speed.attrs['standard_name'])  # 'wind_speed'
vrtdiv(truncation=None)[source]#

Calculate relative vorticity and horizontal divergence.

Parameters:

truncation (int, optional) – Triangular truncation limit for spherical harmonic computation

Returns:

  • vorticity (DataArray) – Relative vorticity with CF-compliant attributes

  • divergence (DataArray) – Horizontal divergence with CF-compliant attributes

Return type:

Tuple[DataArray, DataArray]

See also

vorticity

Calculate only vorticity

divergence

Calculate only divergence

Examples

>>> vrt, div = vw.vrtdiv()
>>> vrt_t13, div_t13 = vw.vrtdiv(truncation=13)
vorticity(truncation=None)[source]#

Calculate relative vorticity.

Parameters:

truncation (int, optional) – Triangular truncation limit for spherical harmonic computation

Returns:

Relative vorticity field with CF-compliant attributes

Return type:

DataArray

See also

vrtdiv

Calculate both vorticity and divergence

absolutevorticity

Calculate absolute vorticity

Examples

>>> vrt = vw.vorticity()
>>> vrt_t13 = vw.vorticity(truncation=13)
divergence(truncation=None)[source]#

Calculate horizontal divergence.

Parameters:

truncation (int, optional) – Triangular truncation limit for spherical harmonic computation

Returns:

Horizontal divergence field with CF-compliant attributes

Return type:

DataArray

See also

vrtdiv

Calculate both vorticity and divergence

Examples

>>> div = vw.divergence()
>>> div_t13 = vw.divergence(truncation=13)
planetaryvorticity(omega=None)[source]#

Calculate planetary vorticity (Coriolis parameter).

Parameters:

omega (float, optional) – Earth’s angular velocity in rad/s. Default is 7.292e-5 s⁻¹

Returns:

Planetary vorticity (Coriolis parameter) with CF-compliant attributes

Return type:

DataArray

See also

absolutevorticity

Calculate absolute vorticity

Examples

>>> f = vw.planetaryvorticity()
>>> f_custom = vw.planetaryvorticity(omega=7.2921150e-5)
absolutevorticity(omega=None, truncation=None)[source]#

Absolute vorticity (sum of relative and planetary vorticity).

Optional arguments:

omega

Earth’s angular velocity. The default value if not specified is 7.292x10**-5 s**-1.

truncation

Truncation limit (triangular truncation) for the spherical harmonic computation.

Returns:

avorticity

The absolute (relative + planetary) vorticity.

See also:

~VectorWind.vorticity, ~VectorWind.planetaryvorticity.

Examples:

Compute absolute vorticity:

avrt = w.absolutevorticity()

Compute absolute vorticity and apply spectral truncation at triangular T13, also override the default value for Earth’s angular velocity:

avrt = w.absolutevorticity(omega=7.2921150, truncation=13)
sfvp(truncation=None)[source]#

Calculate streamfunction and velocity potential.

Parameters:

truncation (int, optional) – Triangular truncation limit for spherical harmonic computation

Returns:

  • streamfunction (DataArray) – Streamfunction field with CF-compliant attributes

  • velocity_potential (DataArray) – Velocity potential field with CF-compliant attributes

Return type:

Tuple[DataArray, DataArray]

See also

streamfunction

Calculate only streamfunction

velocitypotential

Calculate only velocity potential

Examples

>>> psi, chi = vw.sfvp()
>>> psi_t13, chi_t13 = vw.sfvp(truncation=13)
streamfunction(truncation=None)[source]#

Streamfunction.

Optional argument:

truncation

Truncation limit (triangular truncation) for the spherical harmonic computation.

Returns:

sf

The streamfunction.

See also:

~VectorWind.sfvp.

Examples:

Compute streamfunction:

sf = w.streamfunction()

Compute streamfunction and apply spectral truncation at triangular T13:

sfT13 = w.streamfunction(truncation=13)
velocitypotential(truncation=None)[source]#

Velocity potential.

Optional argument:

truncation

Truncation limit (triangular truncation) for the spherical harmonic computation.

Returns:

vp

The velocity potential.

See also:

~VectorWind.sfvp.

Examples:

Compute velocity potential:

vp = w.velocity potential()

Compute velocity potential and apply spectral truncation at triangular T13:

vpT13 = w.velocity potential(truncation=13)
helmholtz(truncation=None)[source]#

Irrotational and non-divergent components of the vector wind.

Optional argument:

truncation

Truncation limit (triangular truncation) for the spherical harmonic computation.

Returns:

uchi, vchi, upsi, vpsi

Zonal and meridional components of irrotational and non-divergent wind components respectively.

See also:

~VectorWind.irrotationalcomponent, ~VectorWind.nondivergentcomponent.

Examples:

Compute the irrotational and non-divergent components of the vector wind:

uchi, vchi, upsi, vpsi = w.helmholtz()

Compute the irrotational and non-divergent components of the vector wind and apply spectral truncation at triangular T13:

uchiT13, vchiT13, upsiT13, vpsiT13 = w.helmholtz(truncation=13)
irrotationalcomponent(truncation=None)[source]#

Irrotational (divergent) component of the vector wind.

Note

If both the irrotational and non-divergent components are required then ~VectorWind.helmholtz should be used instead.

Optional argument:

truncation

Truncation limit (triangular truncation) for the spherical harmonic computation.

Returns:

uchi, vchi

The zonal and meridional components of the irrotational wind respectively.

See also:

~VectorWind.helmholtz.

Examples:

Compute the irrotational component of the vector wind:

uchi, vchi = w.irrotationalcomponent()

Compute the irrotational component of the vector wind and apply spectral truncation at triangular T13:

uchiT13, vchiT13 = w.irrotationalcomponent(truncation=13)
nondivergentcomponent(truncation=None)[source]#

Non-divergent (rotational) component of the vector wind.

Note

If both the non-divergent and irrotational components are required then ~VectorWind.helmholtz should be used instead.

Optional argument:

truncation

Truncation limit (triangular truncation) for the spherical harmonic computation.

Returns:

upsi, vpsi

The zonal and meridional components of the non-divergent wind respectively.

See also:

~VectorWind.helmholtz.

Examples:

Compute the non-divergent component of the vector wind:

upsi, vpsi = w.nondivergentcomponent()

Compute the non-divergent component of the vector wind and apply spectral truncation at triangular T13:

upsiT13, vpsiT13 = w.nondivergentcomponent(truncation=13)
gradient(chi, truncation=None)[source]#

Calculate vector gradient of a scalar field on the sphere.

Parameters:
  • chi (DataArray) – Scalar field with same latitude/longitude dimensions as wind components

  • truncation (int, optional) – Triangular truncation limit for spherical harmonic computation

Returns:

  • u_gradient (DataArray) – Zonal component of vector gradient

  • v_gradient (DataArray) – Meridional component of vector gradient

Return type:

Tuple[DataArray, DataArray]

Examples

>>> abs_vrt = vw.absolutevorticity()
>>> avrt_u, avrt_v = vw.gradient(abs_vrt)
>>> avrt_u_t13, avrt_v_t13 = vw.gradient(abs_vrt, truncation=13)
rossbywavesource(truncation=None, omega=None)[source]#

Calculate Rossby wave source.

The Rossby wave source quantifies the generation of Rossby wave activity in the atmosphere through the interaction of divergence with absolute vorticity and the advection of absolute vorticity by the irrotational wind.

Parameters:
  • truncation (int, optional) – Triangular truncation limit for spherical harmonic computation. If None, uses the default truncation based on grid resolution.

  • omega (float, optional) – Earth’s angular velocity in rad/s. Default is 7.292e-5 s⁻¹.

Returns:

Rossby wave source field with CF-compliant attributes

Return type:

DataArray

See also

absolutevorticity

Calculate absolute vorticity

divergence

Calculate horizontal divergence

irrotationalcomponent

Calculate irrotational wind component

gradient

Calculate vector gradient

Notes

The Rossby wave source is defined as: S = -ζₐ∇·v - v_χ·∇ζₐ

where: - ζₐ is absolute vorticity (relative + planetary) - ∇·v is horizontal divergence - v_χ is the irrotational (divergent) wind component - ∇ζₐ is the gradient of absolute vorticity

Positive values indicate Rossby wave generation, while negative values indicate wave absorption or dissipation.

Examples

>>> rws = vw.rossbywavesource()
>>> rws_t21 = vw.rossbywavesource(truncation=21)
>>> rws_custom_omega = vw.rossbywavesource(omega=7.2921150e-5)

# Create a plot of Rossby wave source >>> import matplotlib.pyplot as plt >>> import cartopy.crs as ccrs >>> >>> ax = plt.axes(projection=ccrs.PlateCarree()) >>> rws.plot.contourf(ax=ax, transform=ccrs.PlateCarree(), … levels=20, cmap=’RdBu_r’) >>> ax.coastlines() >>> ax.gridlines() >>> plt.title(‘Rossby Wave Source’) >>> plt.show()

References

Sardeshmukh, P. D., & Hoskins, B. J. (1988). The generation of global rotational flow by steady idealized tropical heating. Journal of the Atmospheric Sciences, 45(7), 1228-1251.

truncate(field, truncation=None)[source]#

Apply spectral truncation to a scalar field.

Parameters:
  • field (DataArray) – Scalar field with same latitude/longitude dimensions as wind components

  • truncation (int, optional) – Triangular truncation limit. If None, defaults to nlat-1

Returns:

Field with spectral truncation applied

Return type:

DataArray

Examples

>>> field_trunc = vw.truncate(scalar_field)
>>> field_t21 = vw.truncate(scalar_field, truncation=21)

Tools and Utilities#

Data management tools for windspharm.VectorWind and spharm.Spharmt.

This module provides utilities for preparing data for use with spherical harmonic wind analysis, including dimension reordering, shape management, and coordinate system handling.

Main Functions:

prep_data: Prepare data for VectorWind input recover_data: Restore original data shape/order get_recovery: Create recovery function for multiple arrays reverse_latdim: Reverse latitude dimension order order_latdim: Ensure north-to-south latitude ordering

Key Functions for Multi-dimensional Data Processing

prep_data(data, dim_order)[source]#

Prepare data for input to VectorWind or spherical harmonic transforms.

This function reorders and reshapes input data to the format required by VectorWind and spherical harmonic analysis routines. It returns the prepared data along with information needed to recover the original format.

Parameters:
  • data (array_like) – Input data array with at least 2 dimensions. Must contain latitude and longitude dimensions.

  • dim_order (str) – String specifying dimension order using ‘x’ for longitude and ‘y’ for latitude. Other characters represent additional dimensions.

Returns:

  • prepared_data (ndarray) – Data reshaped to (latitude, longitude, other) format

  • recovery_info (dict) – Dictionary containing information needed to recover original format. Contains keys: ‘intermediate_shape’, ‘intermediate_order’, ‘original_order’

Raises:

ValueError – If dim_order doesn’t contain both ‘x’ and ‘y’

Return type:

Tuple[ndarray, Dict[str, Any]]

See also

recover_data

Recover original data format

get_recovery

Create recovery function for multiple arrays

Examples

>>> import numpy as np
>>> # Prepare 4D data: (time=12, level=17, lat=73, lon=144)
>>> data = np.random.randn(12, 17, 73, 144)
>>> pdata, info = prep_data(data, 'tzyx')
>>> print(pdata.shape)  # (73, 144, 204)
>>>
>>> # Prepare data with arbitrary dimension labels
>>> data = np.random.randn(144, 16, 73, 21)
>>> pdata, info = prep_data(data, 'xayb')
recover_data(prepared_data, recovery_info)[source]#

Recover original shape and dimension order of processed data.

This function reverses the operations performed by prep_data, restoring the original shape and dimension ordering of data that has been processed by VectorWind or spherical harmonic methods.

Parameters:
  • prepared_data (ndarray) – Data array with 2 or 3 dimensions where first two are (lat, lon)

  • recovery_info (dict) – Recovery information dictionary from prep_data

Returns:

Data restored to original shape and dimension order

Return type:

ndarray

See also

prep_data

Prepare data for processing

get_recovery

Create recovery function for multiple arrays

Examples

>>> # Recover data processed with prep_data
>>> original_data = recover_data(processed_data, info)
>>>
>>> # Use with VectorWind output
>>> pdata, info = prep_data(wind_data, 'tzyx')
>>> vw = VectorWind(u_prepared, v_prepared)
>>> vorticity = vw.vorticity()
>>> vort_original = recover_data(vorticity, info)
get_recovery(recovery_info)[source]#

Create a recovery function for multiple arrays.

Returns a function that can recover the original shape and dimension order of multiple arrays using a single recovery information dictionary. This is useful when processing multiple related arrays.

Parameters:

recovery_info (dict) – Recovery information dictionary from prep_data

Returns:

recover_func – Function that takes variable number of arrays and returns them with original shape and dimension order restored

Return type:

callable

See also

prep_data

Prepare data for processing

recover_data

Recover single array

Examples

>>> # Prepare multiple wind components
>>> u_prep, info = prep_data(u, 'tzyx')
>>> v_prep, _ = prep_data(v, 'tzyx')
>>>
>>> # Process with VectorWind
>>> vw = VectorWind(u_prep, v_prep)
>>> sf, vp = vw.sfvp()
>>>
>>> # Create recovery function and use it
>>> recover = get_recovery(info)
>>> u_orig, v_orig, sf_orig, vp_orig = recover(u_prep, v_prep, sf, vp)
reverse_latdim(u, v, axis=0)[source]#

Reverse the latitude dimension of wind components.

Creates copies of the input wind components with the latitude dimension reversed. This is useful for converting between different latitude ordering conventions (e.g., north-to-south vs south-to-north).

Parameters:
  • u (array_like) – Zonal and meridional wind components

  • v (array_like) – Zonal and meridional wind components

  • axis (int, default 0) – Index of the latitude dimension to reverse

Returns:

  • u_reversed (ndarray) – Zonal wind component with latitude dimension reversed (copy)

  • v_reversed (ndarray) – Meridional wind component with latitude dimension reversed (copy)

Return type:

Tuple[ndarray, ndarray]

See also

order_latdim

Ensure north-to-south latitude ordering

Examples

>>> # Reverse first dimension (default)
>>> u_rev, v_rev = reverse_latdim(u, v)
>>>
>>> # Reverse third dimension
>>> u_rev, v_rev = reverse_latdim(u, v, axis=2)
order_latdim(lat_coords, u, v, axis=0)[source]#

Ensure latitude dimension is ordered north-to-south.

Checks the latitude coordinate array and ensures it runs from north to south. If necessary, reverses the latitude coordinates and the corresponding dimension in the wind components. Always returns copies of the input arrays.

This helper is mainly for the NumPy standard.VectorWind interface, whose underlying spherical harmonic backend follows the SPHEREPACK convention that the first latitude row is the northernmost point and latitude indices increase southward.

Parameters:
  • lat_coords (array_like) – Array of latitude coordinate values

  • u (array_like) – Zonal and meridional wind components

  • v (array_like) – Zonal and meridional wind components

  • axis (int, default 0) – Index of the latitude dimension in the wind components

Returns:

  • lat_ordered (ndarray) – Latitude coordinates ordered north-to-south (copy)

  • u_ordered (ndarray) – Zonal wind component with latitude ordered north-to-south (copy)

  • v_ordered (ndarray) – Meridional wind component with latitude ordered north-to-south (copy)

Return type:

Tuple[ndarray, ndarray, ndarray]

See also

reverse_latdim

Reverse latitude dimension

Examples

>>> # Order latitude when it's the first dimension
>>> lat_ord, u_ord, v_ord = order_latdim(lat, u, v)
>>>
>>> # Order latitude when it's the third dimension
>>> lat_ord, u_ord, v_ord = order_latdim(lat, u, v, axis=2)
  • prep_data: Prepares wind data for spherical harmonic analysis by handling dimension reordering

  • recover_data: Recovers original data structure after spherical harmonic processing

These tools are essential when working with 4D atmospheric data (time, level, lat, lon) and need to process individual time slices while preserving the original data structure.

Real-World Usage Pattern#

For processing large atmospheric datasets with parallel computation:

from skyborn.windspharm.tools import prep_data, recover_data
from skyborn.windspharm.standard import VectorWind
from concurrent.futures import ThreadPoolExecutor, as_completed
from tqdm import tqdm
import numpy as np

def calculate_divergent_wind_component(time_idx: int,
                                     zonal_wind: np.ndarray,
                                     meridional_wind: np.ndarray) -> tuple:
    """
    Calculate divergent wind components for a single time step.

    Args:
        time_idx: Time index for processing
        zonal_wind: Zonal wind component array (time, level, lat, lon)
        meridional_wind: Meridional wind component array (time, level, lat, lon)

    Returns:
        tuple: (time_idx, divergent_zonal_wind, divergent_meridional_wind)
    """
    try:
        # Prepare data for spherical harmonic analysis
        # Convert from (level, lat, lon) to analysis format
        zonal_prepared, wind_info = prep_data(zonal_wind[time_idx], 'zyx')
        meridional_prepared, _ = prep_data(meridional_wind[time_idx], 'zyx')

        # Create vector wind object and compute irrotational component
        vector_wind = VectorWind(zonal_prepared, meridional_prepared)
        divergent_u, divergent_v = vector_wind.irrotationalcomponent()

        # Recover original data structure
        divergent_u_recovered = recover_data(divergent_u, wind_info)
        divergent_v_recovered = recover_data(divergent_v, wind_info)

        return time_idx, divergent_u_recovered, divergent_v_recovered

    except Exception as e:
        print(f"Error processing time step {time_idx}: {e}")
        raise

def compute_divergent_wind_components(zonal_wind: np.ndarray,
                                    meridional_wind: np.ndarray,
                                    max_workers: int = 4) -> tuple:
    """
    Compute divergent wind components using parallel processing.

    Args:
        zonal_wind: 4D array of zonal wind (time, level, lat, lon)
        meridional_wind: 4D array of meridional wind (time, level, lat, lon)
        max_workers: Maximum number of worker threads

    Returns:
        tuple: (divergent_zonal_wind, divergent_meridional_wind)
    """
    input_shape = zonal_wind.shape
    n_timesteps = input_shape[0]

    # Initialize output arrays with same shape as input
    divergent_zonal = np.zeros_like(zonal_wind)
    divergent_meridional = np.zeros_like(meridional_wind)

    # Process all time steps in parallel
    with ThreadPoolExecutor(max_workers=max_workers) as executor:
        # Submit all tasks
        future_to_time = {
            executor.submit(calculate_divergent_wind_component, t,
                          zonal_wind, meridional_wind): t
            for t in range(n_timesteps)
        }

        # Collect results with progress bar
        for future in tqdm(as_completed(future_to_time),
                         total=n_timesteps,
                         desc='Computing divergent wind components'):
            try:
                time_idx, div_u, div_v = future.result()
                divergent_zonal[time_idx] = div_u
                divergent_meridional[time_idx] = div_v

            except Exception as e:
                print(f"Error in future result: {e}")
                raise

    return divergent_zonal, divergent_meridional

Usage Example:

# Load your 4D wind data (time, level, lat, lon)
# zonal_wind.shape = (365, 37, 181, 360)  # Daily data, 37 levels
# meridional_wind.shape = (365, 37, 181, 360)

# Import skyborn windspharm components
from skyborn.windspharm.tools import prep_data, recover_data
from skyborn.windspharm.standard import VectorWind

# Compute divergent components efficiently
div_u, div_v = compute_divergent_wind_components(
    zonal_wind, meridional_wind, max_workers=8
)

print(f"Processed {zonal_wind.shape[0]} time steps")
print(f"Output shape: {div_u.shape}")

This approach is particularly effective for:

  • Long time series atmospheric data analysis

  • Climate model output processing

  • Reanalysis data divergent circulation studies

  • Large-scale atmospheric dynamics research

Common Functions#

Common functionality shared across windspharm interfaces.

This module provides utility functions that are used by multiple windspharm interfaces (standard, xarray, etc.) for grid type detection, dimension ordering, and data reshaping operations.

Main Functions:

get_apiorder: Calculate dimension reordering for API compatibility inspect_gridtype: Determine grid type from latitude coordinates to3d: Reshape array to 3D format for spherical harmonic analysis

get_apiorder(ndim, latitude_dim, longitude_dim)[source]#

Calculate dimension ordering for API-compatible array transposition.

This function determines the proper dimension ordering to move latitude and longitude dimensions to the first two positions, which is required by the spherical harmonic analysis routines.

Parameters:
  • ndim (int) – Total number of dimensions in the array

  • latitude_dim (int) – Index of the latitude dimension

  • longitude_dim (int) – Index of the longitude dimension

Returns:

  • apiorder (list of int) – List of dimension indices in the order required for API compatibility. Latitude and longitude dimensions are moved to positions 0 and 1.

  • reorder (list of int) – Inverse mapping to restore original dimension order after processing

Return type:

Tuple[List[int], List[int]]

Examples

>>> # For 4D array (time, level, lat, lon)
>>> apiorder, reorder = get_apiorder(4, 2, 3)
>>> print(apiorder)  # [2, 3, 0, 1] - lat, lon, time, level
>>> print(reorder)   # [2, 3, 0, 1] - inverse mapping
inspect_gridtype(latitudes)[source]#

Determine grid type by examining latitude coordinate values.

Analyzes the spacing and distribution of latitude points to determine whether they correspond to a regular (equally-spaced) or Gaussian grid. This is essential for selecting the appropriate spherical harmonic analysis method.

Parameters:

latitudes (array_like) – Array of latitude coordinate values in degrees. Should be ordered from north to south (90° to -90°).

Returns:

gridtype – Detected grid type: - ‘regular’: Equally-spaced latitude grid - ‘gaussian’: Gaussian latitude grid (optimal for spectral methods)

Return type:

{‘regular’, ‘gaussian’}

Raises:

ValueError – If the grid type cannot be determined or if latitudes don’t match either regular or Gaussian patterns within tolerance

Examples

>>> # Regular grid with 73 points including poles
>>> lats_regular = np.linspace(90, -90, 73)
>>> gridtype = inspect_gridtype(lats_regular)
>>> print(gridtype)  # 'regular'
>>>
>>> # Gaussian grid
>>> from skyborn.spharm import gaussian_lats_wts
>>> lats_gauss, _ = gaussian_lats_wts(64)
>>> gridtype = inspect_gridtype(lats_gauss)
>>> print(gridtype)  # 'gaussian'
to3d(array)[source]#

Reshape array to 3D format for spherical harmonic analysis.

Converts input arrays to the (nlat, nlon, nfields) format required by spherical harmonic routines, where any additional dimensions beyond latitude and longitude are flattened into a single ‘fields’ dimension.

Parameters:

array (array_like) – Input array with shape (nlat, nlon, …) where latitude and longitude are the first two dimensions

Returns:

Reshaped array with shape (nlat, nlon, nfields) where nfields is the product of all dimensions beyond the first two

Return type:

ndarray

Examples

>>> # 4D array: (nlat=73, nlon=144, ntime=12, nlevel=17)
>>> data_4d = np.random.randn(73, 144, 12, 17)
>>> data_3d = to3d(data_4d)
>>> print(data_3d.shape)  # (73, 144, 204)
>>>
>>> # Already 3D array remains unchanged
>>> data_3d_input = np.random.randn(73, 144, 10)
>>> result = to3d(data_3d_input)
>>> print(result.shape)  # (73, 144, 10)

xarray Interface#

Spherical harmonic vector wind computations with xarray interface.

This module provides a VectorWind class that works with xarray DataArrays, preserving coordinate information and metadata throughout the computation process. It serves as a high-level interface to the standard VectorWind implementation.

Main Class:

VectorWind: xarray-aware interface for wind field analysis

Example

>>> import xarray as xr
>>> from skyborn.windspharm.xarray import VectorWind
>>>
>>> # Load wind data as xarray DataArrays
>>> u = xr.open_dataarray('u_wind.nc')
>>> v = xr.open_dataarray('v_wind.nc')
>>>
>>> # Create VectorWind instance
>>> vw = VectorWind(u, v)
>>>
>>> # Compute with preserved metadata
>>> vorticity = vw.vorticity()
>>> streamfunction = vw.streamfunction()
class ReducedVectorWind(u, v, pl, rsphere=6.3712e6, legfunc='stored', precision='auto')[source]#

Bases: object

Xarray interface for packed reduced-Gaussian vector wind analysis.

The first dimension of u and v is interpreted as the packed reduced-Gaussian point dimension and must have length sum(pl). Any remaining dimensions are carried through unchanged.

Parameters:
  • u (DataArray)

  • v (DataArray)

  • pl (Any)

  • rsphere (float)

  • legfunc (LegFunc)

  • precision (Literal['auto', 'single', 'double'])

__init__(u, v, pl, rsphere=6.3712e6, legfunc='stored', precision='auto')[source]#
Parameters:
Return type:

None

u()[source]#

Return the packed zonal wind component.

Return type:

DataArray

v()[source]#

Return the packed meridional wind component.

Return type:

DataArray

magnitude()[source]#

Return wind speed on the packed reduced grid.

Return type:

DataArray

vrtdiv(truncation=None)[source]#

Return relative vorticity and horizontal divergence.

Parameters:

truncation (int | None)

Return type:

Tuple[DataArray, DataArray]

vorticity(truncation=None)[source]#

Return relative vorticity.

Parameters:

truncation (int | None)

Return type:

DataArray

divergence(truncation=None)[source]#

Return horizontal divergence.

Parameters:

truncation (int | None)

Return type:

DataArray

planetaryvorticity(omega=None)[source]#

Return planetary vorticity.

Parameters:

omega (float | None)

Return type:

DataArray

absolutevorticity(omega=None, truncation=None)[source]#

Return absolute vorticity.

Parameters:
  • omega (float | None)

  • truncation (int | None)

Return type:

DataArray

sfvp(truncation=None)[source]#

Return streamfunction and velocity potential.

Parameters:

truncation (int | None)

Return type:

Tuple[DataArray, DataArray]

streamfunction(truncation=None)[source]#

Return streamfunction.

Parameters:

truncation (int | None)

Return type:

DataArray

velocitypotential(truncation=None)[source]#

Return velocity potential.

Parameters:

truncation (int | None)

Return type:

DataArray

helmholtz(truncation=None)[source]#

Return irrotational and non-divergent wind components.

Parameters:

truncation (int | None)

Return type:

Tuple[DataArray, DataArray, DataArray, DataArray]

irrotationalcomponent(truncation=None)[source]#

Return irrotational wind component.

Parameters:

truncation (int | None)

Return type:

Tuple[DataArray, DataArray]

nondivergentcomponent(truncation=None)[source]#

Return non-divergent wind component.

Parameters:

truncation (int | None)

Return type:

Tuple[DataArray, DataArray]

gradient(chi, truncation=None)[source]#

Return vector gradient of a packed scalar field.

Parameters:
Return type:

Tuple[DataArray, DataArray]

truncate(field, truncation=None)[source]#

Apply spectral truncation to a packed scalar field.

Parameters:
Return type:

DataArray

rossbywavesource(truncation=None, omega=None)[source]#

Return Rossby wave source on the packed reduced grid.

Parameters:
  • truncation (int | None)

  • omega (float | None)

Return type:

DataArray

xarray VectorWind Class#

class VectorWind(u, v, rsphere=6.3712e6, legfunc='stored', gridtype=None, precision='auto')[source]#

Bases: object

Vector wind analysis using xarray DataArrays.

This class provides a high-level interface for spherical harmonic wind analysis that preserves xarray coordinate information and metadata. It wraps the standard VectorWind implementation while maintaining CF-compliant attributes.

Parameters:
  • u (xarray.DataArray) – Zonal and meridional wind components. Must have the same dimensions, coordinates, and contain no missing values. Should include latitude and longitude dimensions with appropriate coordinate information.

  • v (xarray.DataArray) – Zonal and meridional wind components. Must have the same dimensions, coordinates, and contain no missing values. Should include latitude and longitude dimensions with appropriate coordinate information.

  • rsphere (float, default 6.3712e6) – Earth radius in meters for spherical harmonic computations.

  • legfunc ({'stored', 'computed'}, default 'stored') – Legendre function computation method: - ‘stored’: precompute and store (faster, more memory) - ‘computed’: compute on-the-fly (slower, less memory)

  • gridtype ({'regular', 'gaussian'}, optional) – Explicit grid type override. If omitted, the grid type is inferred from latitude coordinates.

  • precision (Literal['auto', 'single', 'double'])

_api#

Underlying standard VectorWind instance

Type:

standard.VectorWind

_reorder#

Original dimension ordering for output reconstruction

Type:

tuple

_ishape#

Original data shape

Type:

tuple

_coords#

Original coordinate information

Type:

list

Examples

>>> import xarray as xr
>>> from skyborn.windspharm.xarray import VectorWind
>>>
>>> # Load wind components
>>> u = xr.open_dataarray('u850.nc')
>>> v = xr.open_dataarray('v850.nc')
>>>
>>> # Create VectorWind instance
>>> vw = VectorWind(u, v)
>>>
>>> # Compute vorticity with preserved metadata
>>> vorticity = vw.vorticity()
>>> print(vorticity.attrs)  # CF-compliant attributes
>>>
>>> # Helmholtz decomposition
>>> u_chi, v_chi, u_psi, v_psi = vw.helmholtz()
__init__(u, v, rsphere=6.3712e6, legfunc='stored', gridtype=None, precision='auto')[source]#

Initialize VectorWind instance with comprehensive validation.

Parameters:
Return type:

None

u()[source]#

Get zonal component of vector wind.

Returns:

Zonal (eastward) wind component with CF-compliant attributes

Return type:

DataArray

Examples

>>> u_wind = vw.u()
>>> print(u_wind.attrs['standard_name'])  # 'eastward_wind'
v()[source]#

Get meridional component of vector wind.

Returns:

Meridional (northward) wind component with CF-compliant attributes

Return type:

DataArray

Examples

>>> v_wind = vw.v()
>>> print(v_wind.attrs['standard_name'])  # 'northward_wind'
magnitude()[source]#

Calculate wind speed (magnitude of vector wind).

Returns:

Wind speed with CF-compliant attributes

Return type:

DataArray

Examples

>>> wind_speed = vw.magnitude()
>>> print(wind_speed.attrs['standard_name'])  # 'wind_speed'
vrtdiv(truncation=None)[source]#

Calculate relative vorticity and horizontal divergence.

Parameters:

truncation (int, optional) – Triangular truncation limit for spherical harmonic computation

Returns:

  • vorticity (DataArray) – Relative vorticity with CF-compliant attributes

  • divergence (DataArray) – Horizontal divergence with CF-compliant attributes

Return type:

Tuple[DataArray, DataArray]

See also

vorticity

Calculate only vorticity

divergence

Calculate only divergence

Examples

>>> vrt, div = vw.vrtdiv()
>>> vrt_t13, div_t13 = vw.vrtdiv(truncation=13)
vorticity(truncation=None)[source]#

Calculate relative vorticity.

Parameters:

truncation (int, optional) – Triangular truncation limit for spherical harmonic computation

Returns:

Relative vorticity field with CF-compliant attributes

Return type:

DataArray

See also

vrtdiv

Calculate both vorticity and divergence

absolutevorticity

Calculate absolute vorticity

Examples

>>> vrt = vw.vorticity()
>>> vrt_t13 = vw.vorticity(truncation=13)
divergence(truncation=None)[source]#

Calculate horizontal divergence.

Parameters:

truncation (int, optional) – Triangular truncation limit for spherical harmonic computation

Returns:

Horizontal divergence field with CF-compliant attributes

Return type:

DataArray

See also

vrtdiv

Calculate both vorticity and divergence

Examples

>>> div = vw.divergence()
>>> div_t13 = vw.divergence(truncation=13)
planetaryvorticity(omega=None)[source]#

Calculate planetary vorticity (Coriolis parameter).

Parameters:

omega (float, optional) – Earth’s angular velocity in rad/s. Default is 7.292e-5 s⁻¹

Returns:

Planetary vorticity (Coriolis parameter) with CF-compliant attributes

Return type:

DataArray

See also

absolutevorticity

Calculate absolute vorticity

Examples

>>> f = vw.planetaryvorticity()
>>> f_custom = vw.planetaryvorticity(omega=7.2921150e-5)
absolutevorticity(omega=None, truncation=None)[source]#

Absolute vorticity (sum of relative and planetary vorticity).

Optional arguments:

omega

Earth’s angular velocity. The default value if not specified is 7.292x10**-5 s**-1.

truncation

Truncation limit (triangular truncation) for the spherical harmonic computation.

Returns:

avorticity

The absolute (relative + planetary) vorticity.

See also:

~VectorWind.vorticity, ~VectorWind.planetaryvorticity.

Examples:

Compute absolute vorticity:

avrt = w.absolutevorticity()

Compute absolute vorticity and apply spectral truncation at triangular T13, also override the default value for Earth’s angular velocity:

avrt = w.absolutevorticity(omega=7.2921150, truncation=13)
sfvp(truncation=None)[source]#

Calculate streamfunction and velocity potential.

Parameters:

truncation (int, optional) – Triangular truncation limit for spherical harmonic computation

Returns:

  • streamfunction (DataArray) – Streamfunction field with CF-compliant attributes

  • velocity_potential (DataArray) – Velocity potential field with CF-compliant attributes

Return type:

Tuple[DataArray, DataArray]

See also

streamfunction

Calculate only streamfunction

velocitypotential

Calculate only velocity potential

Examples

>>> psi, chi = vw.sfvp()
>>> psi_t13, chi_t13 = vw.sfvp(truncation=13)
streamfunction(truncation=None)[source]#

Streamfunction.

Optional argument:

truncation

Truncation limit (triangular truncation) for the spherical harmonic computation.

Returns:

sf

The streamfunction.

See also:

~VectorWind.sfvp.

Examples:

Compute streamfunction:

sf = w.streamfunction()

Compute streamfunction and apply spectral truncation at triangular T13:

sfT13 = w.streamfunction(truncation=13)
velocitypotential(truncation=None)[source]#

Velocity potential.

Optional argument:

truncation

Truncation limit (triangular truncation) for the spherical harmonic computation.

Returns:

vp

The velocity potential.

See also:

~VectorWind.sfvp.

Examples:

Compute velocity potential:

vp = w.velocity potential()

Compute velocity potential and apply spectral truncation at triangular T13:

vpT13 = w.velocity potential(truncation=13)
helmholtz(truncation=None)[source]#

Irrotational and non-divergent components of the vector wind.

Optional argument:

truncation

Truncation limit (triangular truncation) for the spherical harmonic computation.

Returns:

uchi, vchi, upsi, vpsi

Zonal and meridional components of irrotational and non-divergent wind components respectively.

See also:

~VectorWind.irrotationalcomponent, ~VectorWind.nondivergentcomponent.

Examples:

Compute the irrotational and non-divergent components of the vector wind:

uchi, vchi, upsi, vpsi = w.helmholtz()

Compute the irrotational and non-divergent components of the vector wind and apply spectral truncation at triangular T13:

uchiT13, vchiT13, upsiT13, vpsiT13 = w.helmholtz(truncation=13)
irrotationalcomponent(truncation=None)[source]#

Irrotational (divergent) component of the vector wind.

Note

If both the irrotational and non-divergent components are required then ~VectorWind.helmholtz should be used instead.

Optional argument:

truncation

Truncation limit (triangular truncation) for the spherical harmonic computation.

Returns:

uchi, vchi

The zonal and meridional components of the irrotational wind respectively.

See also:

~VectorWind.helmholtz.

Examples:

Compute the irrotational component of the vector wind:

uchi, vchi = w.irrotationalcomponent()

Compute the irrotational component of the vector wind and apply spectral truncation at triangular T13:

uchiT13, vchiT13 = w.irrotationalcomponent(truncation=13)
nondivergentcomponent(truncation=None)[source]#

Non-divergent (rotational) component of the vector wind.

Note

If both the non-divergent and irrotational components are required then ~VectorWind.helmholtz should be used instead.

Optional argument:

truncation

Truncation limit (triangular truncation) for the spherical harmonic computation.

Returns:

upsi, vpsi

The zonal and meridional components of the non-divergent wind respectively.

See also:

~VectorWind.helmholtz.

Examples:

Compute the non-divergent component of the vector wind:

upsi, vpsi = w.nondivergentcomponent()

Compute the non-divergent component of the vector wind and apply spectral truncation at triangular T13:

upsiT13, vpsiT13 = w.nondivergentcomponent(truncation=13)
gradient(chi, truncation=None)[source]#

Calculate vector gradient of a scalar field on the sphere.

Parameters:
  • chi (DataArray) – Scalar field with same latitude/longitude dimensions as wind components

  • truncation (int, optional) – Triangular truncation limit for spherical harmonic computation

Returns:

  • u_gradient (DataArray) – Zonal component of vector gradient

  • v_gradient (DataArray) – Meridional component of vector gradient

Return type:

Tuple[DataArray, DataArray]

Examples

>>> abs_vrt = vw.absolutevorticity()
>>> avrt_u, avrt_v = vw.gradient(abs_vrt)
>>> avrt_u_t13, avrt_v_t13 = vw.gradient(abs_vrt, truncation=13)
rossbywavesource(truncation=None, omega=None)[source]#

Calculate Rossby wave source.

The Rossby wave source quantifies the generation of Rossby wave activity in the atmosphere through the interaction of divergence with absolute vorticity and the advection of absolute vorticity by the irrotational wind.

Parameters:
  • truncation (int, optional) – Triangular truncation limit for spherical harmonic computation. If None, uses the default truncation based on grid resolution.

  • omega (float, optional) – Earth’s angular velocity in rad/s. Default is 7.292e-5 s⁻¹.

Returns:

Rossby wave source field with CF-compliant attributes

Return type:

DataArray

See also

absolutevorticity

Calculate absolute vorticity

divergence

Calculate horizontal divergence

irrotationalcomponent

Calculate irrotational wind component

gradient

Calculate vector gradient

Notes

The Rossby wave source is defined as: S = -ζₐ∇·v - v_χ·∇ζₐ

where: - ζₐ is absolute vorticity (relative + planetary) - ∇·v is horizontal divergence - v_χ is the irrotational (divergent) wind component - ∇ζₐ is the gradient of absolute vorticity

Positive values indicate Rossby wave generation, while negative values indicate wave absorption or dissipation.

Examples

>>> rws = vw.rossbywavesource()
>>> rws_t21 = vw.rossbywavesource(truncation=21)
>>> rws_custom_omega = vw.rossbywavesource(omega=7.2921150e-5)

# Create a plot of Rossby wave source >>> import matplotlib.pyplot as plt >>> import cartopy.crs as ccrs >>> >>> ax = plt.axes(projection=ccrs.PlateCarree()) >>> rws.plot.contourf(ax=ax, transform=ccrs.PlateCarree(), … levels=20, cmap=’RdBu_r’) >>> ax.coastlines() >>> ax.gridlines() >>> plt.title(‘Rossby Wave Source’) >>> plt.show()

References

Sardeshmukh, P. D., & Hoskins, B. J. (1988). The generation of global rotational flow by steady idealized tropical heating. Journal of the Atmospheric Sciences, 45(7), 1228-1251.

truncate(field, truncation=None)[source]#

Apply spectral truncation to a scalar field.

Parameters:
  • field (DataArray) – Scalar field with same latitude/longitude dimensions as wind components

  • truncation (int, optional) – Triangular truncation limit. If None, defaults to nlat-1

Returns:

Field with spectral truncation applied

Return type:

DataArray

Examples

>>> field_trunc = vw.truncate(scalar_field)
>>> field_t21 = vw.truncate(scalar_field, truncation=21)

Examples#

Basic Wind Analysis#

import numpy as np
from skyborn.windspharm.standard import VectorWind

# Create sample wind data on a Gaussian grid
nlat, nlon = 73, 144
u = np.random.randn(nlat, nlon)  # Zonal wind
v = np.random.randn(nlat, nlon)  # Meridional wind

# Initialize VectorWind
vw = VectorWind(u, v, gridtype='gaussian')

# Calculate dynamical quantities
vorticity = vw.vorticity()
divergence = vw.divergence()
streamfunction = vw.streamfunction()
velocity_potential = vw.velocitypotential()

Helmholtz Decomposition#

# Decompose wind into rotational and divergent components
u_rot, v_rot, u_div, v_div = vw.helmholtz()

# Or get components separately
u_nondiv, v_nondiv = vw.nondivergentcomponent()
u_irrot, v_irrot = vw.irrotationalcomponent()

Time Series Analysis#

# For time series data
nt = 100  # number of time steps
u_time = np.random.randn(nlat, nlon, nt)
v_time = np.random.randn(nlat, nlon, nt)

vw_time = VectorWind(u_time, v_time, gridtype='gaussian')
vorticity_time = vw_time.vorticity()  # Shape: (nlat, nlon, nt)

Error Handling#

The module provides comprehensive error handling for common issues:

try:
    # This will raise an error if shapes don't match
    vw = VectorWind(u[:, :50], v, gridtype='gaussian')
except ValueError as e:
    print(f"Shape mismatch error: {e}")

try:
    # This will raise an error for invalid grid type
    vw = VectorWind(u, v, gridtype='invalid')
except ValueError as e:
    print(f"Invalid grid type: {e}")

Performance Considerations#

  • Grid Type: Gaussian grids are generally more efficient for spherical harmonic transforms

  • Legendre Functions: Use ‘stored’ mode for repeated computations, ‘computed’ for memory-constrained environments

  • Array Layout: Ensure latitude dimension decreases from north to south

  • Data Type: Use float64 for best numerical accuracy

Mathematical Background#

The module implements spherical harmonic analysis for vector fields on the sphere. Key mathematical concepts include:

Vorticity (ζ):
\[\zeta = \frac{1}{a\cos\phi}\left(\frac{\partial v}{\partial \lambda} - \frac{\partial (u\cos\phi)}{\partial \phi}\right)\]
Divergence (δ):
\[\delta = \frac{1}{a\cos\phi}\left(\frac{\partial u}{\partial \lambda} + \frac{\partial (v\cos\phi)}{\partial \phi}\right)\]
Streamfunction (ψ) and Velocity Potential (χ):
\[u = -\frac{1}{a}\frac{\partial \psi}{\partial \phi} + \frac{1}{a\cos\phi}\frac{\partial \chi}{\partial \lambda}\]
\[v = \frac{1}{a\cos\phi}\frac{\partial \psi}{\partial \lambda} + \frac{1}{a}\frac{\partial \chi}{\partial \phi}\]

where: - a is Earth’s radius - φ is latitude - λ is longitude - u, v are zonal and meridional wind components

See Also#