gridded.pyugrid package

Subpackages

Submodules

gridded.pyugrid.read_netcdf module

code to read the netcdf unstructured grid standard:

https://github.com/ugrid-conventions/ugrid-conventions/

This code is called by the UGrid class to load into a UGRID object.

gridded.pyugrid.read_netcdf.find_mesh_names(nc)[source]

find all the meshes in an open netcCDF4.DataSet

Parameters

nc – the netCDF4 Dataset object to look for mesh names in

NOTE: checks for 2-d topology_dimension

gridded.pyugrid.read_netcdf.is_valid_mesh(nc, varname)[source]

determine if the given variable name is a valid mesh definition

Parameters
  • nc – a netCDF4 Dataset to check

  • varname – name of the candidate mesh variable

gridded.pyugrid.read_netcdf.load_grid_from_nc_dataset(nc, grid, mesh_name=None)[source]

loads UGrid object from a netCDF4.DataSet object, adding the data to the passed-in grid object.

It will load the mesh specified, or look for the first one it finds if none is specified

Parameters
  • nc (netCDF4 Dataset object) – netcdf Dataset to be loaded up

  • grid (UGrid object.) – the grid object to put the mesh and data into.

  • mesh_name=None – name of the mesh to load

# :param load_data=False: flag to indicate whether you want to load the # associated data or not. The mesh will be loaded # in any case. If False, only the mesh will be # loaded. If True, then all the data associated # with the mesh will be loaded. This could be huge! # :type load_data: boolean

NOTE: passing the UGrid object in to avoid circular references, while keeping the netcdf reading code in its own file.

gridded.pyugrid.read_netcdf.load_grid_from_ncfilename(filename, grid, mesh_name=None)[source]

loads UGrid object from a netcdf file, adding the data to the passed-in grid object.

It will load the mesh specified, or look for the first one it finds if none is specified

Parameters
  • filename – filename or OpenDAP url of dataset.

  • grid (UGrid object.) – the grid object to put the mesh and data into.

  • mesh_name=None – name of the mesh to load

# :param load_data=False: flag to indicate whether you want to load the # associated data or not. The mesh will be loaded # in any case. If False, only the mesh will be # loaded. If True, then all the data associated # with the mesh will be loaded. This could be huge! # :type load_data: boolean

gridded.pyugrid.ugrid module

ugrid classes

set of classes for working with unstructured model grids

The “ugrid” class is the base class: it stores everything in memory

It can read from and write to netcdf files in the UGRID format.

It may be able to reference a netcdf file at some point, rather than storing directly in memory.

NOTE: only tested for triangular and quad mesh grids at the moment.

class gridded.pyugrid.ugrid.UGrid(nodes=None, node_lon=None, node_lat=None, faces=None, edges=None, boundaries=None, face_face_connectivity=None, face_edge_connectivity=None, edge_coordinates=None, face_coordinates=None, boundary_coordinates=None, data=None, grid_topology=None, mesh_name='mesh')[source]

Bases: object

A basic class to hold an unstructured grid as defined in the UGrid convention.

The internal structure mirrors the netcdf data standard.

property boundaries
build_boundaries()[source]

Builds the boundary segments from the cell array.

It is assumed that -1 means no neighbor, which indicates a boundary

This will over-write the existing boundaries array if there is one.

This is a not-very-smart just loop through all the faces method.

build_boundary_coordinates()[source]

Builds the boundary_coordinates array, using the average of the nodes defining each boundary segment.

Note that you may want a different definition of the boundary coordinates than this computes, but this is here to have an easy default.

This will write-over an existing face_coordinates array

Useful if you want this in the output file

build_celltree()[source]

Tries to build the celltree for the current UGrid. Will fail if nodes or faces is not defined.

build_edge_coordinates()[source]

Builds the edge_coordinates array, using the average of the nodes defining each edge.

Note that you may want a different definition of the edge coordinates than this computes, but this is here to have an easy default.

This will write-over an existing edge_coordinates array

Useful if you want this in the output file

build_edges()[source]

Builds the edges array: all the edges defined by the faces

This will replace the existing edge array, if there is one.

NOTE: arbitrary order – should the order be preserved?

build_face_coordinates()[source]

Builds the face_coordinates array, using the average of the nodes defining each face.

Note that you may want a different definition of the face coordinates than this computes, but this is here to have an easy default.

This will write-over an existing face_coordinates array.

Useful if you want this in the output file.

build_face_edge_connectivity()[source]

Builds the face-edge connectivity array

Not implemented yet.

build_face_face_connectivity()[source]

Builds the face_face_connectivity array: giving the neighbors of each cell.

Note: arbitrary order and CW vs CCW may not be consistent.

check_consistent()[source]

Check if the various data is consistent: the edges and faces reference existing nodes, etc.

property edges
property face_edge_connectivity
property face_face_connectivity
property faces
classmethod from_nc_dataset(nc, mesh_name=None)[source]

create a UGrid object from a netcdf file (or opendap url)

Parameters
  • nc (netCDF4.DataSet) – An already open Dataset object

  • mesh_name=None – the name of the mesh you want. If None, then you’ll get the only mesh in the file. If there is more than one mesh in the file, a ValueError Will be raised

# :param load_data=False: flag to indicate whether you want to load the # associated data or not. The mesh will be # loaded in any case. If False, only the mesh # will be loaded. If True, then all the data # associated with the mesh will be loaded. # This could be huge!

# :type load_data: boolean

classmethod from_ncfile(nc_url, mesh_name=None)[source]

create a UGrid object from a netcdf file name (or opendap url)

Parameters
  • nc_url – the filename or OpenDap url you want to load

  • mesh_name=None – the name of the mesh you want. If None, then you’ll get the only mesh in the file. If there is more than one mesh in the file, a ValueError Will be raised

get_lines()[source]
infer_location(data)[source]
Parameters

data

Returns

‘nodes’ if data will fit to the nodes, ‘faces’ if the data will fit to the faces, ‘boundaries’ if the data will fit the boundaries. None otherwise.

If data is a netcdf variable, the “location” attribute is checked.

property info

summary of information about the grid

interpolate(points, variable, location=None, fill_value=0, indices=None, alphas=None, slices=None, _copy=False, _memo=True, _hash=None)

Interpolates a variable on one of the grids to an array of points.

Parameters
  • points – Nx2 Array of lon/lat coordinates to be interpolated to.

  • variable – Array-like of values to associate at location on grid (node, center, edge1, edge2). This may be more than a 2-dimensional array, but you must pass ‘slices’ kwarg with appropriate slice collection to reduce it to 2 dimensions.

  • location – One of (‘node’, ‘center’, ‘edge1’, ‘edge2’) ‘edge1’ is conventionally associated with the ‘vertical’ edges and likewise ‘edge2’ with the ‘horizontal’

  • fill_value – If masked values are encountered in interpolation, this value takes the place of the masked value

  • indices – If computed already, array of Nx2 cell indices can be passed in to increase speed.

  • alphas – If computed already, array of alphas can be passed in to increase speed.

With a numpy array:

sgrid.interpolate_var_to_points(points, sgrid.u[time_idx, depth_idx])

With a raw netCDF Variable:

sgrid.interpolate_var_to_points(points, nc.variables[‘u’], slices=[time_idx, depth_idx])

If you have pre-computed information, you can pass it in to avoid unnecessary computation and increase performance.

  • ind = # precomputed indices of points

  • alphas = # precomputed alphas (useful if interpolating to the same points frequently)

interpolate_var_to_points(points, variable, location=None, fill_value=0, indices=None, alphas=None, slices=None, _copy=False, _memo=True, _hash=None)[source]

Interpolates a variable on one of the grids to an array of points.

Parameters
  • points – Nx2 Array of lon/lat coordinates to be interpolated to.

  • variable – Array-like of values to associate at location on grid (node, center, edge1, edge2). This may be more than a 2-dimensional array, but you must pass ‘slices’ kwarg with appropriate slice collection to reduce it to 2 dimensions.

  • location – One of (‘node’, ‘center’, ‘edge1’, ‘edge2’) ‘edge1’ is conventionally associated with the ‘vertical’ edges and likewise ‘edge2’ with the ‘horizontal’

  • fill_value – If masked values are encountered in interpolation, this value takes the place of the masked value

  • indices – If computed already, array of Nx2 cell indices can be passed in to increase speed.

  • alphas – If computed already, array of alphas can be passed in to increase speed.

With a numpy array:

sgrid.interpolate_var_to_points(points, sgrid.u[time_idx, depth_idx])

With a raw netCDF Variable:

sgrid.interpolate_var_to_points(points, nc.variables[‘u’], slices=[time_idx, depth_idx])

If you have pre-computed information, you can pass it in to avoid unnecessary computation and increase performance.

  • ind = # precomputed indices of points

  • alphas = # precomputed alphas (useful if interpolating to the same points frequently)

interpolation_alphas(points, indices=None, _copy=False, _memo=True, _hash=None)[source]

Given an array of points, this function will return the bilinear interpolation alphas for each of the three nodes of the face that the point is located in. If the point is not located on the grid, the alphas are set to 0 :param points: Nx2 numpy array of lat/lon coordinates

Parameters

indices – If the face indices of the points is already known, it can be passed in to save repeating the effort.

Returns

Nx3 numpy array of interpolation factors

TODO: mask the indices that aren’t on the grid properly.

locate_faces(points, method='celltree', _copy=False, _memo=True, _hash=None)[source]

Returns the face indices, one per point.

Points that are not in the mesh will have an index of -1

If a single point is passed in, a single index will be returned If a sequence of points is passed in an array of indexes will be returned.

Parameters
  • points – The points that you want to locate – (lon, lat). If the shape of point is 1D, function will return a scalar index. If it is 2D, it will return a 1D array of indices

  • method='celltree' – method to use. Options are ‘celltree’, ‘simple’. for ‘celltree’ the celltree2d pacakge must be installed: https://github.com/NOAA-ORR-ERD/cell_tree2d/ ‘simple’ is very, very slow for large grids.

This version utilizes the CellTree data structure.

locate_nodes(points)[source]

Returns the index of the closest nodes to the input locations.

Parameters

points – the lons/lats of locations you want the nodes closest to.

Returns

the index of the closest node.

property node_lat
property node_lon
property nodes
property num_vertices

Maximum number of vertices in a face.

save(filepath, format='netcdf4', variables={})[source]

Save the ugrid object as a netcdf file.

Parameters
  • filepath – path to file you want o save to. An existing one will be clobbered if it already exists.

  • variables – dict of gridded.Variable objects to save to file

Follows the convention established by the netcdf UGRID working group:

http://ugrid-conventions.github.io/ugrid-conventions

NOTE: Variables are saved here, because different conventions do it

differently.

save_as_netcdf(filename, format='netcdf4')[source]

save the dataset to a file

Parameters
  • filename – full path to file to save to.

  • format – format to save – ‘netcdf3’ or ‘netcdf4’ are the only options at this point.

gridded.pyugrid.util module

Miscellaneous util functions.

gridded.pyugrid.util.asarraylike(obj)[source]

If it satisfies the requirements of pyugrid the object is returned as is. If not, then numpy’s array() will be called on it.

Parameters

obj – The object to check if it’s like an array

gridded.pyugrid.util.isarraylike(obj)[source]

tests if obj acts enough like an array to be used in pyugrid.

This should catch netCDF4 variables, etc.

Note: these won’t check if the attributes required actually work right.

gridded.pyugrid.util.point_in_tri(face_points, point, return_weights=False)[source]

Calculates whether point is internal/external to element by comparing summed area of sub triangles with area of triangle element.

NOTE: this seems like a pretty expensive way to do it

– compared to three side of line tests..

gridded.pyugrid.uvar module

UVar object, used to hold variables that are associated with a ugrid

FixMe: should we enable direct attribute acces via python’s attribute access?

i.e. like netcdf variables – would use overloading __setattr__ and __getattr__

class gridded.pyugrid.uvar.UMVar(name, location='none', data=None, attributes=None)[source]

Bases: object

A class to group multiple UVars (or other data sources) and retrieve common information. All the variables grouped in this class must have the same shape, location, and unique names.

TODO: Add attribues that all grouped variables have in common to the UMVar?

add_var(var)[source]
class gridded.pyugrid.uvar.UVar(name, location, data=None, attributes=None)[source]

Bases: object

A class to hold a variable associated with the UGrid. Data can be on the nodes, edges, etc. – “UGrid Variable”

It holds an array of the data, as well as the attributes associated with that data – this is mapped to a netcdf variable with attributes(attributes get stored in the netcdf file)

property data

The data associated with the UVar – a numpy array-like object with the actual data

property dtype

data type of the variable

property max

maximum value of the variable

property min

minimum value of the variable

property ndim

number of dimensions of the variable

property shape

The shape of the data array associated with the UVar.

Module contents