OTIS

  • Reads OTIS format tidal solutions provided by Ohio State University and ESR

    • multi-constituent binary

    • ATLAS-compact binary

    • single-constituent binary

  • Spatially interpolates tidal constituents to input coordinates

  • Writes OTIS-format tide files

Calling Sequence

import pyTMD.io
amp,ph,D,c = pyTMD.io.OTIS.extract_constants(ilon, ilat, grid_file, model_file, EPSG,
    type='z',
    method='spline',
    grid='OTIS')

Source code

pyTMD.io.OTIS.extract_constants(ilon: ndarray, ilat: ndarray, grid_file: str | pathlib.Path | None = None, model_file: str | pathlib.Path | list | None = None, EPSG: str | int | None = None, **kwargs)[source]

Reads files from tide models in OTIS and ATLAS-compact formats

Makes initial calculations to run the tide program

Spatially interpolates tidal constituents to input coordinates

Parameters
ilon: np.ndarray

longitude to interpolate

ilat: np.ndarray

latitude to interpolate

grid_file: str, pathlib.Path or NoneType, default None

grid file for model

model_file: str, pathlib.Path, list or NoneType, default None

model file containing each constituent

EPSG: str or NoneType, default None,

projection of tide model data

type: str, default ‘z’

Tidal variable to read

  • 'z': heights

  • 'u': horizontal transport velocities

  • 'U': horizontal depth-averaged transport

  • 'v': vertical transport velocities

  • 'V': vertical depth-averaged transport

method: str, default ‘spline’

Interpolation method

  • 'bilinear': quick bilinear interpolation

  • 'spline': scipy bivariate spline interpolation

  • 'linear', 'nearest': scipy regular grid interpolations

extrapolate: bool, default False

Extrapolate model using nearest-neighbors

cutoff: float, default 10.0

Extrapolation cutoff in kilometers

Set to np.inf to extrapolate for all points

grid: str, default ‘OTIS’

Tide model file type to read

  • 'ATLAS': reading a global solution with localized solutions

  • 'OTIS': combined global or local solution

  • 'TMD3': combined global or local netCDF4 solution

apply_flexure: bool, default False

Apply ice flexure scaling factor to height values

Returns
amplitude: np.ndarray

amplitudes of tidal constituents

phase: np.ndarray

phases of tidal constituents

D: np.ndarray

bathymetry of tide model

constituents: list

list of model constituents

pyTMD.io.OTIS.read_constants(grid_file: str | pathlib.Path | None = None, model_file: str | pathlib.Path | list | None = None, EPSG: str | int | None = None, **kwargs)[source]

Reads files from tide models in OTIS and ATLAS-compact formats

Parameters
grid_file: str, pathlib.Path or NoneType, default None

grid file for model

model_file: str, pathlib.Path, list or NoneType, default None

model file containing each constituent

EPSG: str or NoneType, default None,

projection of tide model data

type: str, default ‘z’

Tidal variable to read

  • 'z': heights

  • 'u': horizontal transport velocities

  • 'U': horizontal depth-averaged transport

  • 'v': vertical transport velocities

  • 'V': vertical depth-averaged transport

grid: str, default ‘OTIS’

Tide model file type to read

  • 'ATLAS': reading a global solution with localized solutions

  • 'OTIS': combined global or local solution

  • 'TMD3': combined global or local netCDF4 solution

apply_flexure: bool, default False

Apply ice flexure scaling factor to height values

Returns
constituents: obj

complex form of tide model constituents

pyTMD.io.OTIS.interpolate_constants(ilon: ndarray, ilat: ndarray, constituents, EPSG: str | int | None = None, **kwargs)[source]

Interpolate constants from OTIS/ATLAS-compact tidal models to input coordinates

Makes initial calculations to run the tide program

Parameters
ilon: np.ndarray

longitude to interpolate

ilat: np.ndarray

latitude to interpolate

constituents: obj

Tide model constituents (complex form)

EPSG: str or NoneType, default None,

projection of tide model data

type: str, default ‘z’

Tidal variable to read

  • 'z': heights

  • 'u': horizontal transport velocities

  • 'U': horizontal depth-averaged transport

  • 'v': vertical transport velocities

  • 'V': vertical depth-averaged transport

method: str, default ‘spline’

Interpolation method

  • 'bilinear': quick bilinear interpolation

  • 'spline': scipy bivariate spline interpolation

  • 'linear', 'nearest': scipy regular grid interpolations

extrapolate: bool, default False

Extrapolate model using nearest-neighbors

cutoff: float, default 10.0

Extrapolation cutoff in kilometers

Set to np.inf to extrapolate for all points

Returns
amplitude: np.ndarray

amplitudes of tidal constituents

phase: np.ndarray

phases of tidal constituents

D: np.ndarray

bathymetry of tide model

pyTMD.io.OTIS.read_otis_grid(input_file: str | pathlib.Path)[source]

Read grid file to extract model coordinates, bathymetry, masks and indices

Parameters
input_file: str or pathlib.Path

input grid file

Returns
x: np.ndarray

x-coordinates of input grid

y: np.ndarray

y-coordinates of input grid

hz: np.ndarray

model bathymetry

mz: np.ndarray

land/water mask

iob: np.ndarray

open boundary index

dt: np.ndarray

time step

pyTMD.io.OTIS.read_atlas_grid(input_file: str | pathlib.Path)[source]

Read ATLAS grid file to extract model coordinates, bathymetry, masks and indices for both global and local solutions

Parameters
input_file: str or pathlib.Path

input ATLAS grid file

Returns
x: np.ndarray

x-coordinates of input ATLAS grid

y: np.ndarray

y-coordinates of input ATLAS grid

hz: np.ndarray

model bathymetry

mz: np.ndarray

land/water mask

iob: np.ndarray

open boundary index

dt: float

time step

pmask: np.ndarray

global mask

local: dict

dictionary of local tidal solutions for grid variables

  • 'depth': model bathymetry

pyTMD.io.OTIS.read_netcdf_grid(input_file: str | pathlib.Path)[source]

Read netCDF4 grid file to extract model coordinates, bathymetry, masks and flexure scaling factors

Parameters
input_file: str or pathlib.Path

input grid file

Returns
x: np.ndarray

x-coordinates of input grid

y: np.ndarray

y-coordinates of input grid

hz: np.ndarray

model bathymetry

mz: np.ndarray

land/water mask

sf: np.ndarray

scaling factor for applying ice flexure

pyTMD.io.OTIS.read_constituents(input_file: str | pathlib.Path, grid: str = 'OTIS')[source]

Read the list of constituents from an elevation or transport file

Parameters
input_file: str or pathlib.Path

input tidal file

grid: str, default ‘OTIS’

Tide model file type to read

  • 'ATLAS': reading a global solution with localized solutions

  • 'OTIS': combined global or local solution

  • 'TMD3': combined global or local netCDF4 solution

Returns
constituents: list

list of tidal constituent IDs

nc: int

number of constituents

pyTMD.io.OTIS.read_otis_elevation(input_file: str | pathlib.Path, ic: int)[source]

Read elevation file to extract real and imaginary components for constituent

Parameters
input_file: str or pathlib.Path

input elevation file

ic: int

index of constituent

Returns
h: np.ndarray

tidal elevation

pyTMD.io.OTIS.read_atlas_elevation(input_file: str | pathlib.Path, ic: int, constituent: str)[source]

Read elevation file with localized solutions to extract real and imaginary components for constituent

Parameters
input_file: str or pathlib.Path

input ATLAS elevation file

ic: int

index of constituent

constituent: str

tidal constituent ID

Returns
h: float

global tidal elevation

local: dict

dictionary of local tidal solutions for elevation variables

  • 'z': tidal elevation

pyTMD.io.OTIS.read_otis_transport(input_file: str | pathlib.Path, ic: int)[source]

Read transport file to extract real and imaginary components for constituent

Parameters
input_file: str or pathlib.Path

input transport file

ic: int

index of constituent

Returns
u: float

zonal tidal transport

v: float

meridional zonal transport

pyTMD.io.OTIS.read_atlas_transport(input_file: str | pathlib.Path, ic: int, constituent: str)[source]

Read transport file with localized solutions to extract real and imaginary components for constituent

Parameters
input_file: str or pathlib.Path

input ATLAS transport file

ic: int

index of constituent

constituent: str

tidal constituent ID

Returns
u: np.ndarray

global zonal tidal transport

v: np.ndarray

global meridional zonal transport

local: dict

dictionary of local tidal solutions for transport variables

  • 'u': zonal tidal transport

  • 'v': meridional zonal transport

pyTMD.io.OTIS.create_atlas_mask(xi: ndarray, yi: ndarray, mz: ndarray, local: dict, variable: str | None = None)[source]

Creates a high-resolution grid mask from model variables

Parameters
xi: np.ndarray

input x-coordinates of global tide model

yi: np.ndarray

input y-coordinates of global tide model

mz: np.ndarray

global land/water mask

local: dict

dictionary of local tidal solutions

variable: str or NoneType, default None

key for variable within each local solution

  • 'depth': model bathymetry

  • 'z': tidal elevation

  • 'u': zonal tidal transport

  • 'v': meridional zonal transport

Returns
x30: np.ndarray

x-coordinates of high-resolution tide model

y30: np.ndarray

y-coordinates of high-resolution tide model

m30: np.ndarray

high-resolution land/water mask

pyTMD.io.OTIS.interpolate_atlas_model(xi: ndarray, yi: ndarray, zi: ndarray, spacing: float = 0.03333333333333333)[source]

Interpolates global ATLAS tidal solutions into a higher-resolution sampling

Parameters
xi: np.ndarray

input x-coordinates of global tide model

yi: np.ndarray

input y-coordinates of global tide model

zi: np.ndarray

global tide model data

spacing: float

output grid spacing

Returns
xs: np.ndarray

x-coordinates of high-resolution tide model

ys: np.ndarray

y-coordinates of high-resolution tide model

zs: np.ndarray

high-resolution tidal solution for variable

pyTMD.io.OTIS.combine_atlas_model(xi: ndarray, yi: ndarray, zi: ndarray, pmask: ndarray, local: dict, variable: str | None = None)[source]

Combines global and local ATLAS tidal solutions into a single high-resolution solution

Parameters
xi: np.ndarray

input x-coordinates of global tide model

yi: np.ndarray

input y-coordinates of global tide model

zi: np.ndarray

global tide model data

pmask: np.ndarray

global mask

local: dict

dictionary of local tidal solutions

variable: str or NoneType, default None

key for variable within each local solution

  • 'depth': model bathymetry

  • 'z': tidal elevation

  • 'u': zonal tidal transport

  • 'v': meridional zonal transport

Returns
x30: np.ndarray

x-coordinates of high-resolution tide model

y30: np.ndarray

y-coordinates of high-resolution tide model

z30: np.ndarray

combined high-resolution tidal solution for variable

pyTMD.io.OTIS.read_netcdf_file(input_file: str | pathlib.Path, ic: int, variable: str | None = None)[source]

Read netCDF4 file to extract real and imaginary components for constituent

Parameters
input_file: str or pathlib.Path

input transport file

ic: int

index of constituent

variable: str or NoneType, default None

Tidal variable to read

  • 'z': heights

  • 'u': horizontal transport velocities

  • 'U': horizontal depth-averaged transport

  • 'v': vertical transport velocities

  • 'V': vertical depth-averaged transport

Returns
hc: complex

complex form of tidal constituent oscillation

pyTMD.io.OTIS.output_otis_grid(FILE: str | pathlib.Path, xlim: numpy.ndarray | list, ylim: numpy.ndarray | list, hz: ndarray, mz: ndarray, iob: ndarray, dt: float)[source]

Writes OTIS-format grid files

Parameters
FILE: str or pathlib.Path

output OTIS grid file name

xlim: np.ndarray

x-coordinate grid-cell edges of output grid

ylim: np.ndarray

y-coordinate grid-cell edges of output grid

hz: np.ndarray

bathymetry

mz: np.ndarray

land/water mask

iob: np.ndarray

open boundary index

dt: float

time step

pyTMD.io.OTIS.output_otis_elevation(FILE: str | pathlib.Path, h: ndarray, xlim: numpy.ndarray | list, ylim: numpy.ndarray | list, constituents: list)[source]

Writes OTIS-format elevation files

Parameters
FILE: str or pathlib.Path

output OTIS elevation file name

h: np.ndarray

Eulerian form of tidal height oscillation

xlim: np.ndarray

x-coordinate grid-cell edges of output grid

ylim: np.ndarray

y-coordinate grid-cell edges of output grid

constituents: list

tidal constituent IDs

pyTMD.io.OTIS.output_otis_transport(FILE: str | pathlib.Path, u: ndarray, v: ndarray, xlim: numpy.ndarray | list, ylim: numpy.ndarray | list, constituents: list)[source]

Writes OTIS-format transport files

Parameters
FILE: str or pathlib.Path

output OTIS transport file name

u: complex

Eulerian form of tidal zonal transport oscillation

v: complex

Eulerian form of tidal meridional transport oscillation

xlim: float

x-coordinate grid-cell edges of output grid

ylim: float

y-coordinate grid-cell edges of output grid

constituents: list

tidal constituent IDs

pyTMD.io.OTIS._extend_array(input_array: ndarray, step_size: float)[source]

Extends a longitude array

Parameters
input_array: np.ndarray

array to extend

step_size: float

step size between elements of array

Returns
temp: np.ndarray

extended array

pyTMD.io.OTIS._extend_matrix(input_matrix: ndarray)[source]

Extends a global matrix

Parameters
input_matrix: np.ndarray

matrix to extend

Returns
temp: np.ndarray

extended matrix

pyTMD.io.OTIS._mask_nodes(hz: ndarray, is_global: bool = True)[source]

Construct masks for u and v nodes on a C-grid

Parameters
hz: np.ndarray

bathymetry of grid centers

is_global: bool, default True

input grid is global in terms of longitude

pyTMD.io.OTIS._interpolate_mask(mz: ndarray, is_global: bool = True)[source]

Interpolate mask from zeta nodes to u and v nodes on a C-grid

Parameters
mz: np.ndarray

mask at grid centers

is_global: bool, default True

input grid is global in terms of longitude

pyTMD.io.OTIS._interpolate_zeta(hz: ndarray, is_global: bool = True)[source]

Interpolate data from zeta nodes to u and v nodes on a C-grid

Parameters
hz: np.ndarray

data at grid centers

is_global: bool, default True

input grid is global in terms of longitude