Interpolation on regular grids in 3D

Once you have constructed an ItpMetdata object p (e.g. using read_ocean_velocities, or by hand (see below), you can interpolate the underlying scalar/velocity field.

To do so, call fun(u,p,t) where fun is one of the functions below. Here, u are the spatial coordinates and t is the time.

Velocity fields

OceanTools.uv_tricubicFunction
uv_tricubic(u, p, t)

Tricubic interpolation (Lekien-Marsden + finite differences for values not specified in their paper) of velocity field at location u at time t. Velocity field components are assumed to be stored in p.data[1] and p.data[2].

source
OceanTools.uv_trilinearFunction
uv_trilinear(u, p, t)

Trilinear interpolation of velocity field at location u at time t. Velocity field stored in p.data[1] and p.data[2].

source

Scalar fields

OceanTools.scalar_trilinearFunction
scalar_trilinear(u, p, t)

Trilinear interpolation of scalar field at location u at time t. Scalar field is assumed to be stored in p.data[1].

source
OceanTools.scalar_tricubicFunction
scalar_tricubic(u, p, t)

Tricubic interpolation (Lekien-Marsden + finite differences for values not specified in their paper) of scalar field at location u at time t. Scalar field is assumed to be stored in p.data[1].

source

Gradients

There is support for calculating gradients of the interpolants in the spatial dimensions.

This function has not been extensively tested/used.

Velocity fields deriving from sea surface heights/Hamiltonians

OceanTools.ssh_rhsFunction
ssh_rhs(u, p, t)

Approximating geostrophic sea-surface velocities with the well-known formula

\[u = -A(y)\partial_y h(x,y,t)\ v = A(y)\partial_x h(x,y,t)\]

where:

  • u – longitudinal component of the velocity,
  • v – latitudinal component of the velocity,
  • x – longitude,
  • y – latitude,
  • h – sea-surface height,

and

\[A(y) = g/(2 R^2 \Omega \sin y).\]

Note

The sea surface height is assumed to be given in meters, the grid is assumed to be given in degrees (longitude/latitude). The output unit of the velocity is degree/day.

source

This function has not been extensively tested/used.

Solving the linearized flow map

This function has not been extensively tested/used.

Boundary Behaviour

The boundary behaviour in each direction is specified using p.boundaryX, p.boundaryY, and p.boundaryT. Conceptually speaking, the grid is extended to an infinite grid and interpolation is then performed on this grid.

The value flat means that grid points outside the original grid are assigned the value of the closest point of the original grid.

The value outofbounds means that values outside raise an error whenever they are accessed for interpolation. Note that this means that some values that lie within the spatial bounds of the original grid will raise an error because interpolation also requires nearby points.

The value periodic means that the grid is extended periodically in that direction.

The value semiperiodic means that the grid is extended periodically (with a given periodic), but that there can be "missing" values in between. For example, the case LL = [350,10,0] and UR = [10,20,10] with p.boundaryY=semiperiodic corresponds to the case where values are present for a small part of the globe. Evaluating at an x coordinate of 0 (or any other multiple of 360) works (provided that the values of y and t are inbounds. However, evaluating at an x-coordinate of 180 would raise an error.

Constructing an ItpMetadata object by hand

The easiest way to do this is to call the constructor

ItpMetadata(xspan, yspan, tspan, data, boundaryX, boundarY, boundaryT)

Here xspan, yspan and tspan are ranges that correspond to the values of the coordinates at the datapoints. For example, if you have a $2\pi$ periodic grid with $10$ points in each direction, these should all be range(0,stop=2\pi,length=11)[1:end-1].

The argument data is either (a) a 1-tuple containing a 3d array like object of the values at grid points or (b) a 2-tuple containing two 3d array like objects of values at grid points. Orderin should be in the form array_like_object[x_index,y_index,t_index] and with column-major layout. The case (a) is if you want to interpolate scalars, (b) is if you want to interpolate vectors (i.e. time-dependent spatial velocity fields).

Interpolation on regular grids in 4D

So far, only quad-linear interpolation is supported.

This function has not been extensively tested/used.