API

Filter Class

class gcm_filters.Filter(filter_scale: float, dx_min: float, filter_shape: gcm_filters.filter.FilterShape = FilterShape.GAUSSIAN, transition_width: float = 3.141592653589793, ndim: int = 2, n_steps: int = 0, grid_type: gcm_filters.kernels.GridType = GridType.REGULAR, grid_vars: dict = <factory>)[source]

A class for applying diffusion-based smoothing filters to gridded data.

filter_scalefloat

The filter scale, which has different meaning depending on filter shape

dx_minfloat

The smallest grid spacing. Should have same units as filter_scale

n_stepsint, optional

Number of total steps in the filter (A biharmonic step counts as two steps) n_steps == 0 means the number of steps is chosen automatically

filter_shapeFilterShape
  • FilterShape.GAUSSIAN: The target filter has shape \(e^{-(k filter_scale)^2/24}\)

  • FilterShape.TAPER: The target filter has target grid scale Lf. Smaller scales are zeroed out. Scales larger than pi * filter_scale / 2 are left as-is. In between is a smooth transition.

transition_widthfloat, optional

Width of the transition region in the “Taper” filter. Theoretical minimum is 1; not recommended.

ndimint, optional

Laplacian is applied on a grid of dimension ndim

grid_typeGridType

what sort of grid we are dealing with

grid_varsdict

dictionary of extra parameters used to initialize the grid Laplacian

filter_spec
Type

FilterSpec

Methods

apply(field_or_dataset, dims)

Filter a field or xarray dataset with scalar Laplacian across the dimensions specified by dims.

apply_to_field(field, dims)

Filter a field with scalar Laplacian across the dimensions specified by dims.

apply_to_vector(ufield, vfield, dims)

Filter a vector field with vector Laplacian across the dimensions specified by dims.

plot_shape([ax])

Plot the shape of the target filter and approximation.

apply(field_or_dataset, dims)[source]

Filter a field or xarray dataset with scalar Laplacian across the dimensions specified by dims.

apply_to_field(field, dims)[source]

Filter a field with scalar Laplacian across the dimensions specified by dims.

apply_to_vector(ufield, vfield, dims)[source]

Filter a vector field with vector Laplacian across the dimensions specified by dims.

plot_shape(ax=None)[source]

Plot the shape of the target filter and approximation.

Kernels

Core smoothing routines that operate on 2D arrays.

class gcm_filters.kernels.AreaWeightedMixin(area: numpy.ndarray)[source]

Bases: abc.ABC

Mixin to weight and deweight a field by the cell area.

area
Type

cell area

Methods

finalize

prepare

class gcm_filters.kernels.BaseScalarLaplacian[source]

Bases: abc.ABC

̵Base class for scalar Laplacians.

Methods

__call__(field)

Call self as a function.

required_grid_args()

finalize

prepare

class gcm_filters.kernels.BaseVectorLaplacian[source]

Bases: abc.ABC

Base class for vector Laplacians.

Methods

__call__(ufield, vfield)

Call self as a function.

required_grid_args()

finalize

prepare

class gcm_filters.kernels.CgridVectorLaplacian(wet_mask_t: numpy.ndarray, wet_mask_q: numpy.ndarray, dxT: numpy.ndarray, dyT: numpy.ndarray, dxCu: numpy.ndarray, dyCu: numpy.ndarray, dxCv: numpy.ndarray, dyCv: numpy.ndarray, dxBu: numpy.ndarray, dyBu: numpy.ndarray, area_u: numpy.ndarray, area_v: numpy.ndarray, kappa_iso: numpy.ndarray, kappa_aniso: numpy.ndarray)[source]

Bases: gcm_filters.kernels.BaseVectorLaplacian

̵Vector Laplacian on C-Grid. Follows The implementation for viscosity operators on C-grids suggested by Griffies and Hallberg, 2000.

wet_mask_t
Type

Mask array for t points, 1 for ocean, 0 for land

wet_mask_q
Type

Mask array for q (vorticity) points, 1 for ocean, 0 for land

dxT
Type

x-spacing centered at t points

dyT
Type

y-spacing centered at t points

dxCu
Type

x-spacing centered at u points

dyCu
Type

y-spacing centered at u points

dxCv
Type

x-spacing centered at v points

dyCv
Type

y-spacing centered at v points

dxBu
Type

x-spacing centered at q points

dyBu
Type

y-spacing centered at q points

area_u
Type

U-cell area

area_v
Type

V-cell area

kappa_iso
Type

isotropic viscosity

kappa_aniso
Type

additive anisotropic viscosity aligned with x-direction

Methods

__call__(ufield, vfield)

Call self as a function.

required_grid_args()

finalize

prepare

class gcm_filters.kernels.GridType(value)

Bases: enum.Enum

An enumeration.

class gcm_filters.kernels.IrregularLaplacianWithLandMask(wet_mask: numpy.ndarray, dxw: numpy.ndarray, dyw: numpy.ndarray, dxs: numpy.ndarray, dys: numpy.ndarray, area: numpy.ndarray, kappa_w: numpy.ndarray, kappa_s: numpy.ndarray)[source]

Bases: gcm_filters.kernels.BaseScalarLaplacian

Scalar Laplacian for locally orthogonal grids with land mask.

It is possible to vary the filter scale over the domain by introducing a nondimensional “diffusivity” (attributes kappa_w and kappa_s). For reasons given in Grooms et al. (2021) https://doi.org/10.1002/essoar.10506591.1, we require that both kappa_w and kappa_s values must be <= 1 and that at least one of them is set to 1 somewhere in the domain. Otherwise the scale of the filter will not be equal to filter_scale anywhere in the domain.

wet_mask
Type

Mask array, 1 for ocean, 0 for land

dxw
Type

x-spacing centered at western cell edge

dyw
Type

y-spacing centered at western cell edge

dxs
Type

x-spacing centered at southern cell edge

dys
Type

y-spacing centered at southern cell edge

area
Type

cell area

kappa_w

least one place in the domain must have kappa_w = 1 if kappa_s < 1.

Type

zonal diffusivity centered at western cell edge, values must be <= 1, and at

kappa_s

least one place in the domain must have kappa_s = 1 if kappa_w < 1.

Type

meridional diffusivity centered at southern cell edge, values must be <= 1, and at

Methods

__call__(field)

Call self as a function.

required_grid_args()

finalize

prepare

class gcm_filters.kernels.MOM5LaplacianT(wet_mask: numpy.ndarray, dxt: numpy.ndarray, dyt: numpy.ndarray, dxu: numpy.ndarray, dyu: numpy.ndarray, area_t: numpy.ndarray)[source]

Bases: gcm_filters.kernels.BaseScalarLaplacian

Laplacian for MOM5 (tracer points). MOM5 uses a Northeast convention B-grid, where velocity point U(i,j) is NE of tracer point T(i,j). .. attribute:: For information on MOM5 discretization see

wet_mask
Type

Mask array, 1 for ocean, 0 for land

dxt
Type

width in x of T-cell, model diagnostic dxt

dyt
Type

height in y of T-cell, model diagnostic dyt

dxu
Type

width in x of U-cell, model diagnostic dxu

dyu
Type

height in y of U-cell, model diagnostic dyu

area_t
Type

area of T-cell, dxt*dyt

Methods

__call__(field)

Call self as a function.

required_grid_args()

finalize

prepare

class gcm_filters.kernels.MOM5LaplacianU(wet_mask: numpy.ndarray, dxt: numpy.ndarray, dyt: numpy.ndarray, dxu: numpy.ndarray, dyu: numpy.ndarray, area_u: numpy.ndarray)[source]

Bases: gcm_filters.kernels.BaseScalarLaplacian

Laplacian for MOM5 (velocity points). MOM5 uses a Northeast convention B-grid, where velocity point U(i,j) is NE of tracer point T(i,j). For information on MOM5 discretization see: https://mom-ocean.github.io/assets/pdfs/MOM5_manual.pdf .. attribute:: wet_mask

type

Mask array, 1 for ocean, 0 for land

dxt
Type

width in x of T-cell, model diagnostic dxt

dyt
Type

height in y of T-cell, model diagnostic dyt

dxu
Type

width in x of U-cell, model diagnostic dxu

dyu
Type

height in y of U-cell, model diagnostic dyu

area_u
Type

area of U-cell, dxu*dyu

Methods

__call__(field)

Call self as a function.

required_grid_args()

finalize

prepare

class gcm_filters.kernels.POPTripolarLaplacianTpoint(wet_mask: numpy.ndarray, dxe: numpy.ndarray, dye: numpy.ndarray, dxn: numpy.ndarray, dyn: numpy.ndarray, tarea: numpy.ndarray)[source]

Bases: gcm_filters.kernels.BaseScalarLaplacian

̵Scalar Laplacian for locally orthogonal grid with land mask and tripole boundary condition, as for example used in the global POP configuration. This Laplacian works for scalar fields located at T-points. .. attribute:: wet_mask

type

Mask array, 1 for ocean, 0 for land; can be obtained via xr.where(KMT>0, 1, 0)

dxe
Type

x-spacing centered at eastern T-cell edge, provided by POP model diagnostic HUS(nlat, nlon)

dye
Type

y-spacing centered at eastern T-cell edge, provided by POP model diagnostic HTE(nlat, nlon)

dxn
Type

x-spacing centered at northern T-cell edge, provided by POP model diagnostic HTN(nlat, nlon)

dyn
Type

y-spacing centered at northern T-cell edge, provided by POP model diagnostic HUW(nlat, nlon)

tarea
Type

cell area, provided by POP model diagnostic TAREA(nlat, nlon)

Methods

__call__(field)

Call self as a function.

required_grid_args()

finalize

prepare

class gcm_filters.kernels.RegularLaplacian[source]

Bases: gcm_filters.kernels.BaseScalarLaplacian

̵Scalar Laplacian for regularly spaced Cartesian grids.

Methods

__call__(field)

Call self as a function.

required_grid_args()

finalize

prepare

class gcm_filters.kernels.RegularLaplacianWithArea(area: numpy.ndarray)[source]

Bases: gcm_filters.kernels.AreaWeightedMixin, gcm_filters.kernels.RegularLaplacian

̵Scalar Laplacian operating on a locally orthogonal grid in three steps:

  1. prepare(): Field is multiplied by the cell area. This corresponds to transforming the field from the original locally orthogonal grid to a regularly spaced Cartesian grid with dx = dy = 1.

  2. __call__(): Laplacian acts on regular Cartesian grid.

  3. finalize(): Diffused field is divided by the cell area of the original grid. This corresponds to transforming the field from the regular Cartesian grid back to the original grid.

area
Type

cell area

Methods

__call__(field)

Call self as a function.

required_grid_args()

finalize

prepare

class gcm_filters.kernels.RegularLaplacianWithLandMask(wet_mask: numpy.ndarray)[source]

Bases: gcm_filters.kernels.BaseScalarLaplacian

̵Scalar Laplacian for regularly spaced Cartesian grids with land mask.

wet_mask
Type

Mask array, 1 for ocean, 0 for land

Methods

__call__(field)

Call self as a function.

required_grid_args()

finalize

prepare

class gcm_filters.kernels.RegularLaplacianWithLandMaskAndArea(wet_mask: numpy.ndarray, area: numpy.ndarray)[source]

Bases: gcm_filters.kernels.AreaWeightedMixin, gcm_filters.kernels.RegularLaplacianWithLandMask

̵Scalar Laplacian operating on a locally orthogonal grid with land mask in three steps:

  1. prepare(): Field is multiplied by the cell area. This corresponds to transforming the field from the original locally orthogonal grid to a regularly spaced Cartesian grid with dx = dy = 1.

  2. __call__(): Laplacian acts on regular Cartesian grid.

  3. finalize(): Diffused field is divided by the cell area of the original grid. This corresponds to transforming the field from the regular Cartesian grid back to the original grid.

area
Type

cell area

wet_mask
Type

Mask array, 1 for ocean, 0 for land

Methods

__call__(field)

Call self as a function.

required_grid_args()

finalize

prepare

class gcm_filters.kernels.TripolarRegularLaplacianTpoint(area: numpy.ndarray, wet_mask: numpy.ndarray)[source]

Bases: gcm_filters.kernels.AreaWeightedMixin, gcm_filters.kernels.BaseScalarLaplacian

Scalar Laplacian operating on a locally orthogonal grid with land mask and a tripole boundary. There are three steps:

  1. prepare(): Field is multiplied by the cell area. This corresponds to transforming the field from the original locally orthogonal grid to a regularly spaced Cartesian grid with dx = dy = 1.

  2. __call__(): Laplacian acts on regular Cartesian grid.

  3. finalize(): Diffused field is divided by the cell area of the original grid. This corresponds to transforming the field from the regular Cartesian grid back to the original grid.

area
Type

cell area

wet_mask
Type

Mask array, 1 for ocean, 0 for land

Methods

__call__(field)

Call self as a function.

required_grid_args()

finalize

prepare

gcm_filters.kernels.required_grid_vars(grid_type: gcm_filters.kernels.GridType)[source]

Utility function for figuring out the required grid variables needed by each grid type.

Parameters

grid_type (GridType) – The grid type

Returns

grid_vars – A list of names of required grid variables.

Return type

list