opendrift.readers.reader_FVCOM_xarray
Attributes
Classes
A reader for unstructured (irregularily gridded) `CF compliant |
Functions
|
Works around xarray limitation to load FVCOM grid data |
Convert from iterable of days since epoch to a list of pandas Timestamps |
|
|
name: siglay or siglev |
Module Contents
- opendrift.readers.reader_FVCOM_xarray.logger
- class opendrift.readers.reader_FVCOM_xarray.Reader(filename, proj4, filetype='zarr', cache=None, name=None)[source]
Bases:
opendrift.readers.basereader.BaseReader,opendrift.readers.basereader.UnstructuredReaderA reader for unstructured (irregularily gridded) CF compliant netCDF files.
- Args:
- param filename:
A pattern of netCDF files or a zarr archive. The netCDF file can also be an URL to an OPeNDAP server.
- type filename:
string, required.
- param filetype:
‘netcdf’ or ‘zarr’.
- type filetype:
string, required.
- param name:
Name of reader
- type name:
string, optional
- param proj4:
PROJ.4 string describing projection of data.
- type proj4:
string, optional
See also
py:mod:opendrift.readers.basereader.unstructured.
filename: zarr store or netcdf file proj4: PROJ4 string filetype: “zarr” or “netcdf” cache: either None (no caching) or int (size of cache in bytes).
Only respected for zarr stores.
name: name of Reader.
- variable_aliases
- node_variables = ['sea_floor_depth_below_sea_level', 'sea_water_salinity', 'sea_water_temperature',...
- face_variables = ['x_wind', 'y_wind', 'x_sea_water_velocity', 'y_sea_water_velocity', 'upward_sea_water_velocity']
- dataset = None
- siglay = None
- siglev = None
- siglay_center = None
- siglev_center = None
- ocean_depth_nele = None
- ocean_depth_node = None
- proj4
- times
Setting this to True overrides temporal and spatial bounds checks. Also useful for readers that are constant and do not have a temporal dimension.
- start_time
- end_time
- xmin
- xmax
- ymin
- ymax
- variable_mapping
- variables = []
- boundary
- nodes_idx
- faces_idx
- get_variables(requested_variables, time=None, x=None, y=None, z=None)[source]
FVCOM Grid:
FVCOM uses ‘triangular prisms’ for gridding. Some variables are defined on the faces of the triangles, while others at the node.
x and y holds the positions of the node, while xc and yc holds the positions on the centroids/faces. The centroids/faces are also known as ‘zonal’, or elements (presumably as in finite element).
Note
Currently this reader does not really interpolate. It looks up the closest point in time and space.
Each element has a lookup-table of its surrounding elements, this list can be used when looking up elements for the interpolator of an arbitrary point on the grid. The same goes for the nodes.
Let E be number of elements and N be number of nodes.
Relevant lookup-tables:
nbe: [3 x E] elements surround each element nbve: [9 x N] elements surrounding each node, minimum 3 nbsn: [11 x N] nodes surrounding each node
Variables:
Face: * u * v
Node: * temperature * salinity
- static _vector_nearest_(X, xp)[source]
Find nearest element in vector of vectors X for each xp.
Args:
X NxM matrix of levels xp M vector of positions
- Returns:
- i M vector of indices [0..N] of closest element in
X[0..N, i] to xp[i]
- __nearest_node_sigma__(var, nodes, z)[source]
Find nearest depth at node (sigma layer or sigma level depending on where the variable is defined).
- __nearest_face_sigma__(var, el, z)[source]
TBD: This is extremely hacky…
Find nearest depth at element (sigma layer or sigma level depending on where the variable is defined).
- static z_from_sigma(sigma, depth, elevation=0)[source]
Calculate z-depth from sigma constant.
https://rdrr.io/github/edwardlavender/fvcom.tbx/src/R/depth_from_known.R
Args:
sigma Sigma coefficient(s) depth Depth below mean sea-surface elevation Elevation of sea-surface (e.g. tidal)
Returns: z, depth below sea-surface in meters.
- opendrift.readers.reader_FVCOM_xarray.open_fvcom_files_as_xarray(fname, load=False, **kwarg)[source]
Works around xarray limitation to load FVCOM grid data