diff --git a/.travis.yml b/.travis.yml index a9f6bf3bfe..9ff3371abb 100644 --- a/.travis.yml +++ b/.travis.yml @@ -28,7 +28,7 @@ git: install: - > - export IRIS_TEST_DATA_REF="1696ac3a823a06b95f430670f285ee97671d2cf2"; + export IRIS_TEST_DATA_REF="672dbb46c986038fa5d06a3d8aad691fd1951e07"; export IRIS_TEST_DATA_SUFFIX=$(echo "${IRIS_TEST_DATA_REF}" | sed "s/^v//"); # Install miniconda diff --git a/lib/iris/fileformats/cf.py b/lib/iris/fileformats/cf.py index 1db4e6c61e..f243fa845c 100644 --- a/lib/iris/fileformats/cf.py +++ b/lib/iris/fileformats/cf.py @@ -17,11 +17,9 @@ from abc import ABCMeta, abstractmethod from collections.abc import Iterable, MutableMapping -import os import re import warnings -import netCDF4 import numpy as np import numpy.ma as ma @@ -1008,8 +1006,12 @@ class CFReader: """ - def __init__(self, filename, warn=False, monotonic=False): - self._filename = os.path.expanduser(filename) + def __init__( + self, dataset, warn=False, monotonic=False, exclude_var_names=None + ): + self._dataset = dataset + self._filename = dataset.filepath() + # All CF variable types EXCEPT for the "special cases" of # CFDataVariable, CFCoordinateVariable and _CFFormulaTermsVariable. self._variable_types = ( @@ -1025,8 +1027,6 @@ def __init__(self, filename, warn=False, monotonic=False): #: Collection of CF-netCDF variables associated with this netCDF file self.cf_group = CFGroup() - self._dataset = netCDF4.Dataset(self._filename, mode="r") - # Issue load optimisation warning. if warn and self._dataset.file_format in [ "NETCDF3_CLASSIC", @@ -1039,6 +1039,7 @@ def __init__(self, filename, warn=False, monotonic=False): self._check_monotonic = monotonic + self.exclude_var_names = exclude_var_names or [] self._translate() self._build_cf_groups() self._reset() @@ -1049,14 +1050,20 @@ def __repr__(self): def _translate(self): """Classify the netCDF variables into CF-netCDF variables.""" - netcdf_variable_names = list(self._dataset.variables.keys()) + netcdf_variable_names = [ + var_name + for var_name in self._dataset.variables.keys() + if var_name not in self.exclude_var_names + ] # Identify all CF coordinate variables first. This must be done # first as, by CF convention, the definition of a CF auxiliary # coordinate variable may include a scalar CF coordinate variable, # whereas we want these two types of variables to be mutually exclusive. coords = CFCoordinateVariable.identify( - self._dataset.variables, monotonic=self._check_monotonic + self._dataset.variables, + ignore=self.exclude_var_names, + monotonic=self._check_monotonic, ) self.cf_group.update(coords) coordinate_names = list(self.cf_group.coordinates.keys()) @@ -1064,11 +1071,9 @@ def _translate(self): # Identify all CF variables EXCEPT for the "special cases". for variable_type in self._variable_types: # Prevent grid mapping variables being mis-identified as CF coordinate variables. - ignore = ( - None - if issubclass(variable_type, CFGridMappingVariable) - else coordinate_names - ) + ignore = self.exclude_var_names + if not issubclass(variable_type, CFGridMappingVariable): + ignore += coordinate_names self.cf_group.update( variable_type.identify(self._dataset.variables, ignore=ignore) ) @@ -1082,7 +1087,7 @@ def _translate(self): # Identify and register all CF formula terms. formula_terms = _CFFormulaTermsVariable.identify( - self._dataset.variables + self._dataset.variables, ignore=self.exclude_var_names ) for cf_var in formula_terms.values(): @@ -1125,10 +1130,9 @@ def _build(cf_variable): for variable_type in self._variable_types: # Prevent grid mapping variables being mis-identified as # CF coordinate variables. - if issubclass(variable_type, CFGridMappingVariable): - ignore = None - else: - ignore = coordinate_names + ignore = self.exclude_var_names + if not issubclass(variable_type, CFGridMappingVariable): + ignore += coordinate_names match = variable_type.identify( self._dataset.variables, ignore=ignore, @@ -1258,11 +1262,8 @@ def _build(cf_variable): def _reset(self): """Reset the attribute touch history of each variable.""" for nc_var_name in self._dataset.variables.keys(): - self.cf_group[nc_var_name].cf_attrs_reset() - - def __del__(self): - # Explicitly close dataset to prevent file remaining open. - self._dataset.close() + if nc_var_name not in self.exclude_var_names: + self.cf_group[nc_var_name].cf_attrs_reset() def _getncattr(dataset, attr, default=None): diff --git a/lib/iris/fileformats/netcdf.py b/lib/iris/fileformats/netcdf.py index 4d7ddedc61..24457caeea 100644 --- a/lib/iris/fileformats/netcdf.py +++ b/lib/iris/fileformats/netcdf.py @@ -44,6 +44,7 @@ import iris.exceptions import iris.fileformats.cf import iris.fileformats._pyke_rules +from iris.fileformats.ugrid_cf_reader import UGridCFReader import iris.io import iris.util from iris._lazy_data import as_lazy_data @@ -752,7 +753,7 @@ def coord_from_term(term): cube.add_aux_factory(factory) -def load_cubes(filenames, callback=None): +def load_cubes(filenames, callback=None, *args, **kwargs): """ Loads cubes from a list of NetCDF filenames/URLs. @@ -777,15 +778,20 @@ def load_cubes(filenames, callback=None): filenames = [filenames] for filename in filenames: - # Ingest the netCDF file. - cf = iris.fileformats.cf.CFReader(filename) + # Ingest the netCDF file, creating a reader which also checks for UGRID + # content. + reader = UGridCFReader(filename, *args, **kwargs) # Process each CF data variable. - data_variables = list(cf.cf_group.data_variables.values()) + list( - cf.cf_group.promoted.values() - ) + data_variables = list( + reader.cfreader.cf_group.data_variables.values() + ) + list(reader.cfreader.cf_group.promoted.values()) for cf_var in data_variables: - cube = _load_cube(engine, cf, cf_var, filename) + cube = _load_cube(engine, reader.cfreader, cf_var, filename) + + # Post-process each cube to attach information describing the + # unstructured mesh dimension, if any. + reader.complete_ugrid_cube(cube) # Process any associated formula terms and attach # the corresponding AuxCoordFactory. diff --git a/lib/iris/fileformats/ugrid_cf_reader.py b/lib/iris/fileformats/ugrid_cf_reader.py new file mode 100644 index 0000000000..f0d7ef62af --- /dev/null +++ b/lib/iris/fileformats/ugrid_cf_reader.py @@ -0,0 +1,201 @@ +# Copyright Iris contributors +# +# This file is part of Iris and is released under the LGPL license. +# See COPYING and COPYING.LESSER in the root of the repository for full +# licensing details. +""" +Adds a UGRID extension layer to netCDF file loading. + +""" +from collections import namedtuple +import os + +import netCDF4 + +from gridded.pyugrid.ugrid import UGrid +from gridded.pyugrid.read_netcdf import ( + find_mesh_names, + load_grid_from_nc_dataset, +) +from iris.fileformats.cf import CFReader + + +_UGRID_ELEMENT_TYPE_NAMES = ("node", "edge", "face", "volume") + +# Generate all possible UGRID structural property names. +# These are the UGRID mesh properties that contain variable names for linkage, +# which may appear as recognised properties of the main mesh variable. + +# Start with coordinate variables for each element type (aka "mesh_location"). +_UGRID_LINK_PROPERTIES = [ + "{}_coordinates".format(elem) for elem in _UGRID_ELEMENT_TYPE_NAMES +] + +# Add in all possible type-to-type_connectivity elements. +# NOTE: this actually generates extra unused names, such as +# "node_face_connectivity", because we are not bothering to distinguish +# between lower- and higher-order elements. +# For now just don't worry about that, as long as we get all the ones which +# *are* needed. +_UGRID_LINK_PROPERTIES += [ + "{}_{}_connectivity".format(e1, e2) + for e1 in _UGRID_ELEMENT_TYPE_NAMES + for e2 in _UGRID_ELEMENT_TYPE_NAMES +] + +# Also allow for boundary information. +_UGRID_LINK_PROPERTIES += ["boundary_node_connectivity"] + + +class CubeUgrid( + namedtuple("CubeUgrid", ["cube_dim", "grid", "mesh_location"]) +): + """ + Object recording the unstructured grid dimension of a cube. + + * cube_dim (int): + The cube dimension which maps the unstructured grid. + There can be only one. + + * grid (`gridded.pyugrid.UGrid`): + A 'gridded' description of a UGRID mesh. + + * mesh_location (str): + Which element of the mesh the cube is mapped to. + Can be 'face', 'edge' or 'node'. A 'volume' is not supported. + + """ + + def __str__(self): + result = "Cube unstructured-grid dimension:" + result += "\n cube dimension = {}".format(self.cube_dim) + result += '\n mesh_location = "{}"'.format(self.mesh_location) + result += '\n mesh "{}" :\n'.format(self.grid.mesh_name) + try: + mesh_str = str(self.grid.info) + except TypeError: + mesh_str = "" + result += "\n".join([" " + line for line in mesh_str.split("\n")]) + result += "\n" + return result + + +class UGridCFReader: + """ + A CFReader extension to add UGRID information to netcdf cube loading. + + Identifies UGRID-specific parts of a netcdf file, providing: + + * `self.cfreader` : a CFReader object to interpret the CF data from the + file for cube creation, while ignoring the UGRID mesh data. + + * `self.complete_ugrid_cube(cube)` a call to add the relevant UGRID + information to a cube created from the cfreader data. + + This allows us to decouple UGRID from CF support with minimal changes to + the existing `iris.fileformats.netcdf` code, which is intimately coupled to + both the CFReader class and the netCDF4 file interface. + + """ + + def __init__(self, filename, *args, **kwargs): + self.filename = os.path.expanduser(filename) + dataset = netCDF4.Dataset(self.filename, mode="r") + self.dataset = dataset + meshes = {} + for meshname in find_mesh_names(self.dataset): + mesh = UGrid() + load_grid_from_nc_dataset(dataset, mesh, mesh_name=meshname) + meshes[meshname] = mesh + self.meshes = meshes + + # Generate list of excluded variable names. + exclude_vars = list(meshes.keys()) + + temp_xios_fix = kwargs.pop("temp_xios_fix", False) + if not temp_xios_fix: + # This way *ought* to work, but maybe problems with the test file ? + for mesh in meshes.values(): + mesh_var = dataset.variables[mesh.mesh_name] + for attr in mesh_var.ncattrs(): + if attr in _UGRID_LINK_PROPERTIES: + exclude_vars.extend(mesh_var.getncattr(attr).split()) + else: + # A crude and XIOS-specific alternative .. + exclude_vars += [ + name + for name in dataset.variables.keys() + if any(name.startswith(meshname) for meshname in meshes.keys()) + ] + + # Identify possible mesh dimensions and make a map of them. + meshdims_map = {} # Maps {dimension-name: (mesh, mesh-location)} + for mesh in meshes.values(): + mesh_var = dataset.variables[mesh.mesh_name] + if mesh.faces is not None: + # Work out name of faces dimension and record it. + if "face_dimension" in mesh_var.ncattrs(): + faces_dim_name = mesh_var.getncattr("face_dimension") + else: + # Assume default dimension ordering, and get the dim name + # from dims of a non-optional connectivity variable. + faces_varname = mesh_var.face_node_connectivity + faces_var = dataset.variables[faces_varname] + faces_dim_name = faces_var.dimensions[0] + meshdims_map[faces_dim_name] = (mesh, "face") + if mesh.edges is not None: + # Work out name of edges dimension and record it. + if "edge_dimension" in mesh_var.ncattrs(): + edges_dim_name = mesh_var.getncattr("edge_dimension") + else: + # Assume default dimension ordering, and get the dim name + # from dims of a non-optional connectivity variable. + edges_varname = mesh_var.edge_node_connectivity + edges_var = dataset.variables[edges_varname] + edges_dim_name = edges_var.dimensions[0] + meshdims_map[edges_dim_name] = (mesh, "edge") + if mesh.nodes is not None: + # Work out name of nodes dimension and record it. + # Get it from a non-optional coordinate variable. + nodes_varname = mesh_var.node_coordinates.split()[0] + nodes_var = dataset.variables[nodes_varname] + nodes_dim_name = nodes_var.dimensions[0] + meshdims_map[nodes_dim_name] = (mesh, "node") + self.meshdims_map = meshdims_map + + # Create a CFReader object which skips the UGRID-related variables. + kwargs["exclude_var_names"] = exclude_vars + self.cfreader = CFReader(self.dataset, *args, **kwargs) + + def complete_ugrid_cube(self, cube): + """ + Add the ".ugrid" property to a cube loaded with the `self.cfreader`. + + We identify the unstructured-grid dimension of the cube (if any), and + attach a suitable CubeUgrid object, linking the cube mesh dimension to + an element-type (aka "mesh_location") of a mesh. + + """ + # Set a 'cube.ugrid' property. + data_var = self.dataset.variables[cube.var_name] + meshes_info = [ + (i_dim, self.meshdims_map.get(dim_name)) + for i_dim, dim_name in enumerate(data_var.dimensions) + if dim_name in self.meshdims_map + ] + if len(meshes_info) > 1: + msg = "Cube maps more than one mesh dimension: {}" + raise ValueError(msg.format(meshes_info)) + if meshes_info: + i_dim, (mesh, mesh_location) = meshes_info[0] + cube.ugrid = CubeUgrid( + cube_dim=i_dim, grid=mesh, mesh_location=mesh_location + ) + else: + # Add an empty 'cube.ugrid' to all cubes otherwise. + cube.ugrid = None + return + + def __del__(self): + # Explicitly close dataset to prevent file remaining open. + self.dataset.close() diff --git a/lib/iris/tests/integration/experimental/test_regrid_ProjectedUnstructured.py b/lib/iris/tests/integration/experimental/test_regrid_ProjectedUnstructured.py index 5b5f838f85..06b94f2384 100644 --- a/lib/iris/tests/integration/experimental/test_regrid_ProjectedUnstructured.py +++ b/lib/iris/tests/integration/experimental/test_regrid_ProjectedUnstructured.py @@ -27,7 +27,7 @@ class TestProjectedUnstructured(tests.IrisTest): def setUp(self): path = tests.get_data_path( - ("NetCDF", "unstructured_grid", "theta_nodal_xios.nc") + ("NetCDF", "unstructured_grid", "theta_nodal_not_ugrid.nc") ) self.src = iris.load_cube(path, "Potential Temperature") diff --git a/lib/iris/tests/integration/test_regridding.py b/lib/iris/tests/integration/test_regridding.py index f16b7f4ab5..2399364d4a 100644 --- a/lib/iris/tests/integration/test_regridding.py +++ b/lib/iris/tests/integration/test_regridding.py @@ -99,7 +99,7 @@ def test_nearest(self): class TestUnstructured(tests.IrisTest): def setUp(self): path = tests.get_data_path( - ("NetCDF", "unstructured_grid", "theta_nodal_xios.nc") + ("NetCDF", "unstructured_grid", "theta_nodal_not_ugrid.nc") ) self.src = iris.load_cube(path, "Potential Temperature") self.grid = simple_3d()[0, :, :] diff --git a/lib/iris/tests/integration/ugrid_cf_reader/__init__.py b/lib/iris/tests/integration/ugrid_cf_reader/__init__.py new file mode 100644 index 0000000000..5e2390f9a6 --- /dev/null +++ b/lib/iris/tests/integration/ugrid_cf_reader/__init__.py @@ -0,0 +1,10 @@ +# Copyright Iris contributors +# +# This file is part of Iris and is released under the LGPL license. +# See COPYING and COPYING.LESSER in the root of the repository for full +# licensing details. +""" +Integration tests for the +:mod:`iris.fileformats.ugrid_cf_reader` package. + +""" diff --git a/lib/iris/tests/integration/ugrid_cf_reader/test_ugrid_load.py b/lib/iris/tests/integration/ugrid_cf_reader/test_ugrid_load.py new file mode 100644 index 0000000000..bc5d5c3b9b --- /dev/null +++ b/lib/iris/tests/integration/ugrid_cf_reader/test_ugrid_load.py @@ -0,0 +1,67 @@ +# Copyright Iris contributors +# +# This file is part of Iris and is released under the LGPL license. +# See COPYING and COPYING.LESSER in the root of the repository for full +# licensing details. +""" +Integration tests for the +:mod:`iris.fileformats.ugrid_cf_reader.UGridCFReader` class. + +""" + +# Import iris.tests first so that some things can be initialised before +# importing anything else. +import iris.tests as tests + +from gridded.pyugrid.ugrid import UGrid + +from iris.cube import CubeList +from iris.fileformats.ugrid_cf_reader import CubeUgrid +from iris import Constraint +from iris.fileformats.netcdf import load_cubes + + +@tests.skip_data +class TestUgrid(tests.IrisTest): + def test_basic_load(self): + file_path = tests.get_data_path( + ("NetCDF", "unstructured_grid", "theta_nodal_xios.nc") + ) + + # cube = iris.load_cube(file_path, "theta") + # Note: cannot use iris.load, as merge does not yet preserve + # the cube 'ugrid' properties. + + # Here's a thing that at least works. + loaded_cubes = CubeList(load_cubes(file_path, temp_xios_fix=True)) + + # Just check some expected details. + self.assertEqual(len(loaded_cubes), 2) + + (cube_0,) = loaded_cubes.extract(Constraint("theta")) + (cube_1,) = loaded_cubes.extract(Constraint("radius")) + + # Check the primary cube. + self.assertEqual(cube_0.var_name, "theta") + self.assertEqual(cube_0.long_name, "Potential Temperature") + self.assertEqual(cube_0.shape, (1, 6, 866)) + self.assertEqual( + cube_0.coord_dims(cube_0.coord("time", dim_coords=True)), (0,) + ) + self.assertEqual(cube_0.coord_dims("levels"), (1,)) + self.assertEqual(cube_0.coords(dimensions=2), []) + + # Check the cube.ugrid object. + cubegrid = cube_0.ugrid + self.assertIsInstance(cubegrid, CubeUgrid) + self.assertEqual(cubegrid.cube_dim, 2) + self.assertEqual(cubegrid.mesh_location, "node") + + # Check cube.ugrid.grid : a gridded Grid type. + ugrid = cubegrid.grid + self.assertIsInstance(ugrid, UGrid) + self.assertEqual(ugrid.mesh_name, "Mesh0") + + +if __name__ == "__main__": + tests.main() diff --git a/lib/iris/tests/test_cf.py b/lib/iris/tests/test_cf.py index 1bfd84ce7d..51ac0c6dae 100644 --- a/lib/iris/tests/test_cf.py +++ b/lib/iris/tests/test_cf.py @@ -13,6 +13,8 @@ from unittest import mock +import netCDF4 as nc + import iris import iris.fileformats.cf as cf @@ -54,7 +56,8 @@ def setUp(self): filename = tests.get_data_path( ("NetCDF", "rotated", "xyt", "small_rotPole_precipitation.nc") ) - self.cfr = cf.CFReader(filename) + ds = nc.Dataset(filename) + self.cfr = cf.CFReader(ds) def test_ancillary_variables_pass_0(self): self.assertEqual(self.cfr.cf_group.ancillary_variables, {}) @@ -328,7 +331,8 @@ def setUp(self): "A1B-99999a-river-sep-2070-2099.nc", ) ) - self.cfr = cf.CFReader(filename) + ds = nc.Dataset(filename) + self.cfr = cf.CFReader(ds) def test_bounds(self): time = self.cfr.cf_group["temp_dmax_tmean_abs"].cf_group.coordinates[ @@ -353,12 +357,14 @@ def setUp(self): "A1B-99999a-river-sep-2070-2099.nc", ) ) - self.cfr_start = cf.CFReader(filename) + ds = nc.Dataset(filename) + self.cfr_start = cf.CFReader(ds) filename = tests.get_data_path( ("NetCDF", "label_and_climate", "small_FC_167_mon_19601101.nc") ) - self.cfr_end = cf.CFReader(filename) + ds = nc.Dataset(filename) + self.cfr_end = cf.CFReader(ds) def test_label_dim_start(self): cf_data_var = self.cfr_start.cf_group["temp_dmax_tmean_abs"] diff --git a/lib/iris/tests/unit/fileformats/cf/test_CFReader.py b/lib/iris/tests/unit/fileformats/cf/test_CFReader.py index 7d4bf232fc..3d64136f23 100644 --- a/lib/iris/tests/unit/fileformats/cf/test_CFReader.py +++ b/lib/iris/tests/unit/fileformats/cf/test_CFReader.py @@ -70,11 +70,10 @@ def setUp(self): ) def test_create_global_attributes(self): - with mock.patch("netCDF4.Dataset", return_value=self.dataset): - global_attrs = CFReader("dummy").cf_group.global_attributes - self.assertEqual( - global_attrs["dimensions"], "something something_else" - ) + global_attrs = CFReader(self.dataset).cf_group.global_attributes + self.assertEqual( + global_attrs["dimensions"], "something something_else" + ) class Test_translate__formula_terms(tests.IrisTest): @@ -142,38 +141,37 @@ def setUp(self): self.addCleanup(reset_patch.stop) def test_create_formula_terms(self): - with mock.patch("netCDF4.Dataset", return_value=self.dataset): - cf_group = CFReader("dummy").cf_group - self.assertEqual(len(cf_group), len(self.variables)) - # Check there is a singular data variable. - group = cf_group.data_variables - self.assertEqual(len(group), 1) - self.assertEqual(list(group.keys()), ["temp"]) - self.assertIs(group["temp"].cf_data, self.temp) - # Check there are three coordinates. - group = cf_group.coordinates - self.assertEqual(len(group), 3) - coordinates = ["height", "lat", "lon"] - self.assertEqual(set(group.keys()), set(coordinates)) - for name in coordinates: - self.assertIs(group[name].cf_data, getattr(self, name)) - # Check there are three auxiliary coordinates. - group = cf_group.auxiliary_coordinates - self.assertEqual(len(group), 3) - aux_coordinates = ["delta", "sigma", "orography"] - self.assertEqual(set(group.keys()), set(aux_coordinates)) - for name in aux_coordinates: - self.assertIs(group[name].cf_data, getattr(self, name)) - # Check all the auxiliary coordinates are formula terms. - formula_terms = cf_group.formula_terms - self.assertEqual(set(group.items()), set(formula_terms.items())) - # Check there are three bounds. - group = cf_group.bounds - self.assertEqual(len(group), 3) - bounds = ["height_bnds", "delta_bnds", "sigma_bnds"] - self.assertEqual(set(group.keys()), set(bounds)) - for name in bounds: - self.assertEqual(group[name].cf_data, getattr(self, name)) + cf_group = CFReader(self.dataset).cf_group + self.assertEqual(len(cf_group), len(self.variables)) + # Check there is a singular data variable. + group = cf_group.data_variables + self.assertEqual(len(group), 1) + self.assertEqual(list(group.keys()), ["temp"]) + self.assertIs(group["temp"].cf_data, self.temp) + # Check there are three coordinates. + group = cf_group.coordinates + self.assertEqual(len(group), 3) + coordinates = ["height", "lat", "lon"] + self.assertEqual(set(group.keys()), set(coordinates)) + for name in coordinates: + self.assertIs(group[name].cf_data, getattr(self, name)) + # Check there are three auxiliary coordinates. + group = cf_group.auxiliary_coordinates + self.assertEqual(len(group), 3) + aux_coordinates = ["delta", "sigma", "orography"] + self.assertEqual(set(group.keys()), set(aux_coordinates)) + for name in aux_coordinates: + self.assertIs(group[name].cf_data, getattr(self, name)) + # Check all the auxiliary coordinates are formula terms. + formula_terms = cf_group.formula_terms + self.assertEqual(set(group.items()), set(formula_terms.items())) + # Check there are three bounds. + group = cf_group.bounds + self.assertEqual(len(group), 3) + bounds = ["height_bnds", "delta_bnds", "sigma_bnds"] + self.assertEqual(set(group.keys()), set(bounds)) + for name in bounds: + self.assertEqual(group[name].cf_data, getattr(self, name)) class Test_build_cf_groups__formula_terms(tests.IrisTest): @@ -241,78 +239,73 @@ def setUp(self): self.addCleanup(patcher.stop) def test_associate_formula_terms_with_data_variable(self): - with mock.patch("netCDF4.Dataset", return_value=self.dataset): - cf_group = CFReader("dummy").cf_group - self.assertEqual(len(cf_group), len(self.variables)) - # Check the cf-group associated with the data variable. - temp_cf_group = cf_group["temp"].cf_group - # Check the data variable is associated with eight variables. - self.assertEqual(len(temp_cf_group), 8) - # Check there are three coordinates. - group = temp_cf_group.coordinates - self.assertEqual(len(group), 3) - coordinates = ["height", "lat", "lon"] - self.assertEqual(set(group.keys()), set(coordinates)) - for name in coordinates: - self.assertIs(group[name].cf_data, getattr(self, name)) - # Check the height coordinate is bounded. - group = group["height"].cf_group - self.assertEqual(len(group.bounds), 1) - self.assertIn("height_bnds", group.bounds) - self.assertIs(group["height_bnds"].cf_data, self.height_bnds) - # Check there are five auxiliary coordinates. - group = temp_cf_group.auxiliary_coordinates - self.assertEqual(len(group), 5) - aux_coordinates = ["delta", "sigma", "orography", "x", "y"] - self.assertEqual(set(group.keys()), set(aux_coordinates)) - for name in aux_coordinates: - self.assertIs(group[name].cf_data, getattr(self, name)) - # Check all the auxiliary coordinates are formula terms. - formula_terms = cf_group.formula_terms - self.assertTrue( - set(formula_terms.items()).issubset(list(group.items())) + cf_group = CFReader(self.dataset).cf_group + self.assertEqual(len(cf_group), len(self.variables)) + # Check the cf-group associated with the data variable. + temp_cf_group = cf_group["temp"].cf_group + # Check the data variable is associated with eight variables. + self.assertEqual(len(temp_cf_group), 8) + # Check there are three coordinates. + group = temp_cf_group.coordinates + self.assertEqual(len(group), 3) + coordinates = ["height", "lat", "lon"] + self.assertEqual(set(group.keys()), set(coordinates)) + for name in coordinates: + self.assertIs(group[name].cf_data, getattr(self, name)) + # Check the height coordinate is bounded. + group = group["height"].cf_group + self.assertEqual(len(group.bounds), 1) + self.assertIn("height_bnds", group.bounds) + self.assertIs(group["height_bnds"].cf_data, self.height_bnds) + # Check there are five auxiliary coordinates. + group = temp_cf_group.auxiliary_coordinates + self.assertEqual(len(group), 5) + aux_coordinates = ["delta", "sigma", "orography", "x", "y"] + self.assertEqual(set(group.keys()), set(aux_coordinates)) + for name in aux_coordinates: + self.assertIs(group[name].cf_data, getattr(self, name)) + # Check all the auxiliary coordinates are formula terms. + formula_terms = cf_group.formula_terms + self.assertTrue( + set(formula_terms.items()).issubset(list(group.items())) + ) + # Check the terms by root. + for name, term in zip(aux_coordinates, ["a", "b", "orog"]): + self.assertEqual( + formula_terms[name].cf_terms_by_root, dict(height=term) + ) + # Check the bounded auxiliary coordinates. + for name, name_bnds in zip( + ["delta", "sigma"], ["delta_bnds", "sigma_bnds"] + ): + aux_coord_group = group[name].cf_group + self.assertEqual(len(aux_coord_group.bounds), 1) + self.assertIn(name_bnds, aux_coord_group.bounds) + self.assertIs( + aux_coord_group[name_bnds].cf_data, getattr(self, name_bnds), ) - # Check the terms by root. - for name, term in zip(aux_coordinates, ["a", "b", "orog"]): - self.assertEqual( - formula_terms[name].cf_terms_by_root, dict(height=term) - ) - # Check the bounded auxiliary coordinates. - for name, name_bnds in zip( - ["delta", "sigma"], ["delta_bnds", "sigma_bnds"] - ): - aux_coord_group = group[name].cf_group - self.assertEqual(len(aux_coord_group.bounds), 1) - self.assertIn(name_bnds, aux_coord_group.bounds) - self.assertIs( - aux_coord_group[name_bnds].cf_data, - getattr(self, name_bnds), - ) def test_promote_reference(self): - with mock.patch("netCDF4.Dataset", return_value=self.dataset): - cf_group = CFReader("dummy").cf_group - self.assertEqual(len(cf_group), len(self.variables)) - # Check the number of data variables. - self.assertEqual(len(cf_group.data_variables), 1) - self.assertEqual(list(cf_group.data_variables.keys()), ["temp"]) - # Check the number of promoted variables. - self.assertEqual(len(cf_group.promoted), 1) - self.assertEqual(list(cf_group.promoted.keys()), ["orography"]) - # Check the promoted variable dependencies. - group = cf_group.promoted["orography"].cf_group.coordinates - self.assertEqual(len(group), 2) - coordinates = ("lat", "lon") - self.assertEqual(set(group.keys()), set(coordinates)) - for name in coordinates: - self.assertIs(group[name].cf_data, getattr(self, name)) + cf_group = CFReader(self.dataset).cf_group + self.assertEqual(len(cf_group), len(self.variables)) + # Check the number of data variables. + self.assertEqual(len(cf_group.data_variables), 1) + self.assertEqual(list(cf_group.data_variables.keys()), ["temp"]) + # Check the number of promoted variables. + self.assertEqual(len(cf_group.promoted), 1) + self.assertEqual(list(cf_group.promoted.keys()), ["orography"]) + # Check the promoted variable dependencies. + group = cf_group.promoted["orography"].cf_group.coordinates + self.assertEqual(len(group), 2) + coordinates = ("lat", "lon") + self.assertEqual(set(group.keys()), set(coordinates)) + for name in coordinates: + self.assertIs(group[name].cf_data, getattr(self, name)) def test_formula_terms_ignore(self): self.orography.dimensions = ["lat", "wibble"] - with mock.patch( - "netCDF4.Dataset", return_value=self.dataset - ), mock.patch("warnings.warn") as warn: - cf_group = CFReader("dummy").cf_group + with mock.patch("warnings.warn") as warn: + cf_group = CFReader(self.dataset).cf_group group = cf_group.promoted self.assertEqual(list(group.keys()), ["orography"]) self.assertIs(group["orography"].cf_data, self.orography) @@ -320,10 +313,8 @@ def test_formula_terms_ignore(self): def test_auxiliary_ignore(self): self.x.dimensions = ["lat", "wibble"] - with mock.patch( - "netCDF4.Dataset", return_value=self.dataset - ), mock.patch("warnings.warn") as warn: - cf_group = CFReader("dummy").cf_group + with mock.patch("warnings.warn") as warn: + cf_group = CFReader(self.dataset).cf_group promoted = ["x", "orography"] group = cf_group.promoted self.assertEqual(set(group.keys()), set(promoted)) @@ -335,10 +326,8 @@ def test_promoted_auxiliary_ignore(self): self.wibble = netcdf_variable("wibble", "lat wibble", np.float) self.variables["wibble"] = self.wibble self.orography.coordinates = "wibble" - with mock.patch( - "netCDF4.Dataset", return_value=self.dataset - ), mock.patch("warnings.warn") as warn: - cf_group = CFReader("dummy").cf_group.promoted + with mock.patch("warnings.warn") as warn: + cf_group = CFReader(self.dataset).cf_group.promoted promoted = ["wibble", "orography"] self.assertEqual(set(cf_group.keys()), set(promoted)) for name in promoted: diff --git a/requirements/core.txt b/requirements/core.txt index c3f5775d7e..ddf4f4702b 100644 --- a/requirements/core.txt +++ b/requirements/core.txt @@ -10,5 +10,6 @@ cftime dask[array]>=2 #conda: dask>=2 matplotlib netcdf4 +gridded numpy>=1.14 scipy