From 165805dfa5e06436ebdf1e418d0b25b20a2c37e2 Mon Sep 17 00:00:00 2001 From: Manuel Schlund Date: Fri, 5 Mar 2021 14:40:44 +0100 Subject: [PATCH 1/5] Added AtmosphereSigmaFactory --- lib/iris/aux_factory.py | 179 +++++++++++ lib/iris/fileformats/netcdf.py | 17 +- lib/iris/tests/integration/test_netcdf.py | 33 ++ .../netcdf/TestAtmosphereSigma/save.cdl | 62 ++++ .../test_AtmosphereSigmaFactory.py | 285 ++++++++++++++++++ 5 files changed, 575 insertions(+), 1 deletion(-) create mode 100644 lib/iris/tests/results/integration/netcdf/TestAtmosphereSigma/save.cdl create mode 100644 lib/iris/tests/unit/aux_factory/test_AtmosphereSigmaFactory.py diff --git a/lib/iris/aux_factory.py b/lib/iris/aux_factory.py index 962b46e9e2..e0445f36e7 100644 --- a/lib/iris/aux_factory.py +++ b/lib/iris/aux_factory.py @@ -365,6 +365,185 @@ def _remap_with_bounds(self, dependency_dims, derived_dims): return nd_values_by_key +class AtmosphereSigmaFactory(AuxCoordFactory): + """Defines an atmosphere sigma coordinate factory.""" + + def __init__( + self, pressure_at_top=None, sigma=None, surface_air_pressure=None + ): + """Create class instance. + + Creates an atmosphere sigma coordinate factory with the formula: + + p(n, k, j, i) = pressure_at_top + sigma(k) * + (surface_air_pressure(n, j, i) - pressure_at_top) + + """ + # Configure the metadata manager. + self._metadata_manager = metadata_manager_factory(CoordMetadata) + super().__init__() + + # Check that provided coordinates meet necessary conditions. + self._check_dependencies(pressure_at_top, sigma, surface_air_pressure) + + # Initialize instance attributes + self.units = pressure_at_top.units + self.pressure_at_top = pressure_at_top + self.sigma = sigma + self.surface_air_pressure = surface_air_pressure + self.standard_name = "air_pressure" + self.attributes = {"positive": "down"} + + @staticmethod + def _check_dependencies(pressure_at_top, sigma, surface_air_pressure): + """Check for sufficient coordinates.""" + if any( + [ + pressure_at_top is None, + sigma is None, + surface_air_pressure is None, + ] + ): + raise ValueError( + "Unable to construct atmosphere sigma coordinate factory due " + "to insufficient source coordinates" + ) + + # Check dimensions + if pressure_at_top.shape not in ((), (1,)): + raise ValueError( + f"Expected scalar 'pressure_at_top' coordinate, got shape " + f"{pressure_at_top.shape}" + ) + + # Check bounds + if sigma.nbounds not in (0, 2): + raise ValueError( + f"Invalid 'sigma' coordinate: must have either 0 or 2 bounds, " + f"got {sigma.nbounds:d}" + ) + for coord in (pressure_at_top, surface_air_pressure): + if coord.nbounds: + msg = ( + f"Coordinate '{coord.name()}' has bounds. These will " + "be disregarded" + ) + warnings.warn(msg, UserWarning, stacklevel=2) + + # Check units + if sigma.units.is_unknown(): + # Be graceful, and promote unknown to dimensionless units. + sigma.units = cf_units.Unit("1") + if not sigma.units.is_dimensionless(): + raise ValueError( + f"Invalid units: 'sigma' must be dimensionless, got " + f"'{sigma.units}'" + ) + if pressure_at_top.units != surface_air_pressure.units: + raise ValueError( + f"Incompatible units: 'pressure_at_top' and " + f"'surface_air_pressure' must have the same units, got " + f"'{pressure_at_top.units}' and " + f"'{surface_air_pressure.units}'" + ) + if not pressure_at_top.units.is_convertible("Pa"): + raise ValueError( + "Invalid units: 'pressure_at_top' and 'surface_air_pressure' " + "must have units of pressure" + ) + + @property + def dependencies(self): + """Return dependencies.""" + dependencies = { + "pressure_at_top": self.pressure_at_top, + "sigma": self.sigma, + "surface_air_pressure": self.surface_air_pressure, + } + return dependencies + + @staticmethod + def _derive(pressure_at_top, sigma, surface_air_pressure): + """Derive coordinate.""" + return pressure_at_top + sigma * ( + surface_air_pressure - pressure_at_top + ) + + def make_coord(self, coord_dims_func): + """Make new :class:`iris.coords.AuxCoord`.""" + # Which dimensions are relevant? + derived_dims = self.derived_dims(coord_dims_func) + dependency_dims = self._dependency_dims(coord_dims_func) + + # Build the points array + nd_points_by_key = self._remap(dependency_dims, derived_dims) + points = self._derive( + nd_points_by_key["pressure_at_top"], + nd_points_by_key["sigma"], + nd_points_by_key["surface_air_pressure"], + ) + + # Bounds + bounds = None + if self.sigma.nbounds: + nd_values_by_key = self._remap_with_bounds( + dependency_dims, derived_dims + ) + pressure_at_top = nd_values_by_key["pressure_at_top"] + sigma = nd_values_by_key["sigma"] + surface_air_pressure = nd_values_by_key["surface_air_pressure"] + ok_bound_shapes = [(), (1,), (2,)] + if sigma.shape[-1:] not in ok_bound_shapes: + raise ValueError("Invalid sigma coordinate bounds") + if pressure_at_top.shape[-1:] not in [(), (1,)]: + warnings.warn( + "Pressure at top coordinate has bounds. These are being " + "disregarded" + ) + pressure_at_top_pts = nd_points_by_key["pressure_at_top"] + bds_shape = list(pressure_at_top_pts.shape) + [1] + pressure_at_top = pressure_at_top_pts.reshape(bds_shape) + if surface_air_pressure.shape[-1:] not in [(), (1,)]: + warnings.warn( + "Surface pressure coordinate has bounds. These are being " + "disregarded" + ) + surface_air_pressure_pts = nd_points_by_key[ + "surface_air_pressure" + ] + bds_shape = list(surface_air_pressure_pts.shape) + [1] + surface_air_pressure = surface_air_pressure_pts.reshape( + bds_shape + ) + bounds = self._derive(pressure_at_top, sigma, surface_air_pressure) + + # Create coordinate + return iris.coords.AuxCoord( + points, + standard_name=self.standard_name, + long_name=self.long_name, + var_name=self.var_name, + units=self.units, + bounds=bounds, + attributes=self.attributes, + coord_system=self.coord_system, + ) + + def update(self, old_coord, new_coord=None): + """Notify the factory of the removal/replacement of a coordinate.""" + new_dependencies = self.dependencies + for (name, coord) in self.dependencies.items(): + if old_coord is coord: + new_dependencies[name] = new_coord + try: + self._check_dependencies(**new_dependencies) + except ValueError as exc: + raise ValueError(f"Failed to update dependencies: {exc}") + else: + setattr(self, name, new_coord) + break + + class HybridHeightFactory(AuxCoordFactory): """ Defines a hybrid-height coordinate factory with the formula: diff --git a/lib/iris/fileformats/netcdf.py b/lib/iris/fileformats/netcdf.py index 14dbab8054..73945f715c 100644 --- a/lib/iris/fileformats/netcdf.py +++ b/lib/iris/fileformats/netcdf.py @@ -28,6 +28,7 @@ from iris._lazy_data import as_lazy_data from iris.aux_factory import ( + AtmosphereSigmaFactory, HybridHeightFactory, HybridPressureFactory, OceanSFactory, @@ -111,6 +112,12 @@ "_FactoryDefn", ("primary", "std_name", "formula_terms_format") ) _FACTORY_DEFNS = { + AtmosphereSigmaFactory: _FactoryDefn( + primary="sigma", + std_name="atmosphere_sigma_coordinate", + formula_terms_format="ptop: {pressure_at_top} sigma: {sigma} " + "ps: {surface_air_pressure}", + ), HybridHeightFactory: _FactoryDefn( primary="delta", std_name="atmosphere_hybrid_height_coordinate", @@ -652,6 +659,7 @@ def _load_aux_factory(engine, cube): """ formula_type = engine.requires.get("formula_type") if formula_type in [ + "atmosphere_sigma_coordinate", "atmosphere_hybrid_height_coordinate", "atmosphere_hybrid_sigma_pressure_coordinate", "ocean_sigma_z_coordinate", @@ -673,7 +681,14 @@ def coord_from_term(term): "{!r}".format(name) ) - if formula_type == "atmosphere_hybrid_height_coordinate": + if formula_type == "atmosphere_sigma_coordinate": + pressure_at_top = coord_from_term("ptop") + sigma = coord_from_term("sigma") + surface_air_pressure = coord_from_term("ps") + factory = AtmosphereSigmaFactory( + pressure_at_top, sigma, surface_air_pressure + ) + elif formula_type == "atmosphere_hybrid_height_coordinate": delta = coord_from_term("a") sigma = coord_from_term("b") orography = coord_from_term("orog") diff --git a/lib/iris/tests/integration/test_netcdf.py b/lib/iris/tests/integration/test_netcdf.py index 3ff5bbb19d..9bbac189fa 100644 --- a/lib/iris/tests/integration/test_netcdf.py +++ b/lib/iris/tests/integration/test_netcdf.py @@ -34,6 +34,39 @@ import iris.tests.stock as stock +@tests.skip_data +class TestAtmosphereSigma(tests.IrisTest): + def setUp(self): + # Modify stock cube so it is suitable to have a atmosphere sigma + # factory added to it. + cube = stock.realistic_4d_no_derived() + cube.coord("surface_altitude").rename("surface_air_pressure") + cube.coord("surface_air_pressure").units = "Pa" + cube.coord("sigma").units = "1" + ptop_coord = iris.coords.AuxCoord(1000.0, var_name="ptop", units="Pa") + cube.add_aux_coord(ptop_coord, ()) + cube.remove_coord("level_height") + # Construct and add atmosphere sigma factory. + factory = iris.aux_factory.AtmosphereSigmaFactory( + cube.coord("ptop"), + cube.coord("sigma"), + cube.coord("surface_air_pressure"), + ) + cube.add_aux_factory(factory) + self.cube = cube + + def test_save(self): + with self.temp_filename(suffix=".nc") as filename: + iris.save(self.cube, filename) + self.assertCDL(filename) + + def test_save_load_loop(self): + with self.temp_filename(suffix=".nc") as filename: + iris.save(self.cube, filename) + cube = iris.load_cube(filename, "air_potential_temperature") + assert cube.coords("air_pressure") + + @tests.skip_data class TestHybridPressure(tests.IrisTest): def setUp(self): diff --git a/lib/iris/tests/results/integration/netcdf/TestAtmosphereSigma/save.cdl b/lib/iris/tests/results/integration/netcdf/TestAtmosphereSigma/save.cdl new file mode 100644 index 0000000000..cfb3143050 --- /dev/null +++ b/lib/iris/tests/results/integration/netcdf/TestAtmosphereSigma/save.cdl @@ -0,0 +1,62 @@ +dimensions: + bnds = 2 ; + grid_latitude = 100 ; + grid_longitude = 100 ; + model_level_number = 70 ; + time = 6 ; +variables: + float air_potential_temperature(time, model_level_number, grid_latitude, grid_longitude) ; + air_potential_temperature:standard_name = "air_potential_temperature" ; + air_potential_temperature:units = "K" ; + air_potential_temperature:grid_mapping = "rotated_latitude_longitude" ; + air_potential_temperature:coordinates = "forecast_period ptop sigma surface_air_pressure" ; + int rotated_latitude_longitude ; + rotated_latitude_longitude:grid_mapping_name = "rotated_latitude_longitude" ; + rotated_latitude_longitude:longitude_of_prime_meridian = 0. ; + rotated_latitude_longitude:earth_radius = 6371229. ; + rotated_latitude_longitude:grid_north_pole_latitude = 37.5 ; + rotated_latitude_longitude:grid_north_pole_longitude = 177.5 ; + rotated_latitude_longitude:north_pole_grid_longitude = 0. ; + double time(time) ; + time:axis = "T" ; + time:units = "hours since 1970-01-01 00:00:00" ; + time:standard_name = "time" ; + time:calendar = "gregorian" ; + int model_level_number(model_level_number) ; + model_level_number:axis = "Z" ; + model_level_number:units = "1" ; + model_level_number:standard_name = "model_level_number" ; + model_level_number:positive = "up" ; + float grid_latitude(grid_latitude) ; + grid_latitude:axis = "Y" ; + grid_latitude:bounds = "grid_latitude_bnds" ; + grid_latitude:units = "degrees" ; + grid_latitude:standard_name = "grid_latitude" ; + float grid_latitude_bnds(grid_latitude, bnds) ; + float grid_longitude(grid_longitude) ; + grid_longitude:axis = "X" ; + grid_longitude:bounds = "grid_longitude_bnds" ; + grid_longitude:units = "degrees" ; + grid_longitude:standard_name = "grid_longitude" ; + float grid_longitude_bnds(grid_longitude, bnds) ; + double forecast_period ; + forecast_period:units = "hours" ; + forecast_period:standard_name = "forecast_period" ; + double ptop ; + ptop:units = "Pa" ; + float sigma(model_level_number) ; + sigma:bounds = "sigma_bnds" ; + sigma:units = "1" ; + sigma:long_name = "sigma" ; + sigma:standard_name = "atmosphere_sigma_coordinate" ; + sigma:axis = "Z" ; + sigma:formula_terms = "ptop: ptop sigma: sigma ps: surface_air_pressure" ; + float sigma_bnds(model_level_number, bnds) ; + float surface_air_pressure(grid_latitude, grid_longitude) ; + surface_air_pressure:units = "Pa" ; + surface_air_pressure:standard_name = "surface_air_pressure" ; + +// global attributes: + :source = "Iris test case" ; + :Conventions = "CF-1.7" ; +} diff --git a/lib/iris/tests/unit/aux_factory/test_AtmosphereSigmaFactory.py b/lib/iris/tests/unit/aux_factory/test_AtmosphereSigmaFactory.py new file mode 100644 index 0000000000..c0bb733278 --- /dev/null +++ b/lib/iris/tests/unit/aux_factory/test_AtmosphereSigmaFactory.py @@ -0,0 +1,285 @@ +# 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. +""" +Unit tests for the +`iris.aux_factory.AtmosphereSigmaFactory` class. + +""" + +# Import iris.tests first so that some things can be initialised before +# importing anything else. +import iris.tests as tests + +from unittest import mock + +from cf_units import Unit +import numpy as np + +from iris.aux_factory import AtmosphereSigmaFactory +from iris.coords import AuxCoord, DimCoord + + +class Test___init__(tests.IrisTest): + def setUp(self): + self.pressure_at_top = mock.Mock(units=Unit("Pa"), nbounds=0, shape=()) + self.sigma = mock.Mock(units=Unit("1"), nbounds=0) + self.surface_air_pressure = mock.Mock(units=Unit("Pa"), nbounds=0) + self.kwargs = dict( + pressure_at_top=self.pressure_at_top, + sigma=self.sigma, + surface_air_pressure=self.surface_air_pressure, + ) + + def test_insufficient_coordinates(self): + with self.assertRaises(ValueError): + AtmosphereSigmaFactory() + with self.assertRaises(ValueError): + AtmosphereSigmaFactory( + pressure_at_top=None, + sigma=self.sigma, + surface_air_pressure=self.surface_air_pressure, + ) + with self.assertRaises(ValueError): + AtmosphereSigmaFactory( + pressure_at_top=self.pressure_at_top, + sigma=None, + surface_air_pressure=self.surface_air_pressure, + ) + with self.assertRaises(ValueError): + AtmosphereSigmaFactory( + pressure_at_top=self.pressure_at_top, + sigma=self.sigma, + surface_air_pressure=None, + ) + + def test_ptop_shapes(self): + for shape in [(), (1,)]: + self.pressure_at_top.shape = shape + AtmosphereSigmaFactory(**self.kwargs) + + for shape in [(2,), (1, 1)]: + self.pressure_at_top.shape = shape + with self.assertRaises(ValueError): + AtmosphereSigmaFactory(**self.kwargs) + + def test_sigma_bounds(self): + for n_bounds in [0, 2]: + self.sigma.nbounds = n_bounds + AtmosphereSigmaFactory(**self.kwargs) + + for n_bounds in [-1, 1, 3]: + self.sigma.nbounds = n_bounds + with self.assertRaises(ValueError): + AtmosphereSigmaFactory(**self.kwargs) + + def test_sigma_units(self): + for units in ["1", "unknown", None]: + self.sigma.units = Unit(units) + AtmosphereSigmaFactory(**self.kwargs) + + for units in ["Pa", "m"]: + self.sigma.units = Unit(units) + with self.assertRaises(ValueError): + AtmosphereSigmaFactory(**self.kwargs) + + def test_ptop_ps_units(self): + for units in [("Pa", "Pa")]: + self.pressure_at_top.units = Unit(units[0]) + self.surface_air_pressure.units = Unit(units[1]) + AtmosphereSigmaFactory(**self.kwargs) + + for units in [("Pa", "1"), ("1", "Pa"), ("bar", "Pa"), ("Pa", "hPa")]: + self.pressure_at_top.units = Unit(units[0]) + self.surface_air_pressure.units = Unit(units[1]) + with self.assertRaises(ValueError): + AtmosphereSigmaFactory(**self.kwargs) + + def test_ptop_units(self): + for units in ["Pa", "bar", "mbar", "hPa"]: + self.pressure_at_top.units = Unit(units) + self.surface_air_pressure.units = Unit(units) + AtmosphereSigmaFactory(**self.kwargs) + + for units in ["1", "m", "kg", None]: + self.pressure_at_top.units = Unit(units) + self.surface_air_pressure.units = Unit(units) + with self.assertRaises(ValueError): + AtmosphereSigmaFactory(**self.kwargs) + + +class Test_dependencies(tests.IrisTest): + def setUp(self): + self.pressure_at_top = mock.Mock(units=Unit("Pa"), nbounds=0, shape=()) + self.sigma = mock.Mock(units=Unit("1"), nbounds=0) + self.surface_air_pressure = mock.Mock(units=Unit("Pa"), nbounds=0) + self.kwargs = dict( + pressure_at_top=self.pressure_at_top, + sigma=self.sigma, + surface_air_pressure=self.surface_air_pressure, + ) + + def test_values(self): + factory = AtmosphereSigmaFactory(**self.kwargs) + self.assertEqual(factory.dependencies, self.kwargs) + + +class Test__derive(tests.IrisTest): + def test_function(self): + assert AtmosphereSigmaFactory._derive(0, 0, 0) == 0 + assert AtmosphereSigmaFactory._derive(3, 0, 0) == 3 + assert AtmosphereSigmaFactory._derive(0, 5, 0) == 0 + assert AtmosphereSigmaFactory._derive(0, 0, 7) == 0 + assert AtmosphereSigmaFactory._derive(3, 5, 0) == -12 + assert AtmosphereSigmaFactory._derive(3, 0, 7) == 3 + assert AtmosphereSigmaFactory._derive(0, 5, 7) == 35 + assert AtmosphereSigmaFactory._derive(3, 5, 7) == 23 + + ptop = 3 + sigma = np.array([2, 4]) + ps = np.arange(4).reshape(2, 2) + np.testing.assert_equal( + AtmosphereSigmaFactory._derive(ptop, sigma, ps), + [[-3, -5], [1, 3]], + ) + + +class Test_make_coord(tests.IrisTest): + @staticmethod + def coord_dims(coord): + mapping = dict( + pressure_at_top=(), + sigma=(0,), + surface_air_pressure=(1, 2), + ) + return mapping[coord.name()] + + @staticmethod + def derive(pressure_at_top, sigma, surface_air_pressure, coord=True): + result = pressure_at_top + sigma * ( + surface_air_pressure - pressure_at_top + ) + if coord: + name = "air_pressure" + result = AuxCoord( + result, + standard_name=name, + units="Pa", + attributes=dict(positive="down"), + ) + return result + + def setUp(self): + self.pressure_at_top = AuxCoord( + [3.0], + long_name="pressure_at_top", + units="Pa", + ) + self.sigma = DimCoord( + [1.0, 0.4, 0.1], + bounds=[[1.0, 0.6], [0.6, 0.2], [0.2, 0.0]], + long_name="sigma", + units="1", + ) + self.surface_air_pressure = AuxCoord( + [[-1.0, 2.0], [1.0, 3.0]], + long_name="surface_air_pressure", + units="Pa", + ) + self.kwargs = dict( + pressure_at_top=self.pressure_at_top, + sigma=self.sigma, + surface_air_pressure=self.surface_air_pressure, + ) + + def test_derived_coord(self): + # Broadcast expected points given the known dimensional mapping + pressure_at_top = self.pressure_at_top.points[0] + sigma = self.sigma.points[..., np.newaxis, np.newaxis] + surface_air_pressure = self.surface_air_pressure.points[ + np.newaxis, ... + ] + + # Calculate the expected result + + expected_coord = self.derive( + pressure_at_top, sigma, surface_air_pressure + ) + + # Calculate the actual result + factory = AtmosphereSigmaFactory(**self.kwargs) + coord = factory.make_coord(self.coord_dims) + + # Check bounds + expected_bounds = [ + [[[-1.0, 0.6], [2.0, 2.4]], [[1.0, 1.8], [3.0, 3.0]]], + [[[0.6, 2.2], [2.4, 2.8]], [[1.8, 2.6], [3.0, 3.0]]], + [[[2.2, 3.0], [2.8, 3.0]], [[2.6, 3.0], [3.0, 3.0]]], + ] + np.testing.assert_allclose(coord.bounds, expected_bounds) + coord.bounds = None + + # Check points and metadata + self.assertEqual(expected_coord, coord) + + +class Test_update(tests.IrisTest): + def setUp(self): + self.pressure_at_top = mock.Mock(units=Unit("Pa"), nbounds=0, shape=()) + self.sigma = mock.Mock(units=Unit("1"), nbounds=0) + self.surface_air_pressure = mock.Mock(units=Unit("Pa"), nbounds=0) + self.kwargs = dict( + pressure_at_top=self.pressure_at_top, + sigma=self.sigma, + surface_air_pressure=self.surface_air_pressure, + ) + self.factory = AtmosphereSigmaFactory(**self.kwargs) + + def test_pressure_at_top(self): + new_pressure_at_top = mock.Mock(units=Unit("Pa"), nbounds=0, shape=()) + self.factory.update(self.pressure_at_top, new_pressure_at_top) + self.assertIs(self.factory.pressure_at_top, new_pressure_at_top) + + def test_pressure_at_top_wrong_shape(self): + new_pressure_at_top = mock.Mock( + units=Unit("Pa"), nbounds=0, shape=(2,) + ) + with self.assertRaises(ValueError): + self.factory.update(self.pressure_at_top, new_pressure_at_top) + + def test_sigma(self): + new_sigma = mock.Mock(units=Unit("1"), nbounds=0) + self.factory.update(self.sigma, new_sigma) + self.assertIs(self.factory.sigma, new_sigma) + + def test_sigma_too_many_bounds(self): + new_sigma = mock.Mock(units=Unit("1"), nbounds=4) + with self.assertRaises(ValueError): + self.factory.update(self.sigma, new_sigma) + + def test_sigma_incompatible_units(self): + new_sigma = mock.Mock(units=Unit("Pa"), nbounds=0) + with self.assertRaises(ValueError): + self.factory.update(self.sigma, new_sigma) + + def test_surface_air_pressure(self): + new_surface_air_pressure = mock.Mock(units=Unit("Pa"), nbounds=0) + self.factory.update( + self.surface_air_pressure, new_surface_air_pressure + ) + self.assertIs( + self.factory.surface_air_pressure, new_surface_air_pressure + ) + + def test_surface_air_pressure_incompatible_units(self): + new_surface_air_pressure = mock.Mock(units=Unit("mbar"), nbounds=0) + with self.assertRaises(ValueError): + self.factory.update( + self.surface_air_pressure, new_surface_air_pressure + ) + + +if __name__ == "__main__": + tests.main() From 0eaefdf69749f805768ca5747b93197805b2089e Mon Sep 17 00:00:00 2001 From: Manuel Schlund Date: Thu, 5 Aug 2021 09:50:49 +0200 Subject: [PATCH 2/5] Expanded tests for AtmosphereSigmaFactory --- lib/iris/tests/integration/test_netcdf.py | 2 ++ .../fileformats/nc_load_rules/actions/test__hybrid_formulae.py | 1 + 2 files changed, 3 insertions(+) diff --git a/lib/iris/tests/integration/test_netcdf.py b/lib/iris/tests/integration/test_netcdf.py index 9bbac189fa..01033eb96a 100644 --- a/lib/iris/tests/integration/test_netcdf.py +++ b/lib/iris/tests/integration/test_netcdf.py @@ -61,6 +61,8 @@ def test_save(self): self.assertCDL(filename) def test_save_load_loop(self): + # Ensure that the AtmosphereSigmaFactory is automatically loaded + # when loading the file. with self.temp_filename(suffix=".nc") as filename: iris.save(self.cube, filename) cube = iris.load_cube(filename, "air_potential_temperature") diff --git a/lib/iris/tests/unit/fileformats/nc_load_rules/actions/test__hybrid_formulae.py b/lib/iris/tests/unit/fileformats/nc_load_rules/actions/test__hybrid_formulae.py index f9a11ba403..ca92926542 100644 --- a/lib/iris/tests/unit/fileformats/nc_load_rules/actions/test__hybrid_formulae.py +++ b/lib/iris/tests/unit/fileformats/nc_load_rules/actions/test__hybrid_formulae.py @@ -245,6 +245,7 @@ def test_two_formulae(self): # support all of them. _SUPPORTED_FORMULA_TYPES = ( # NOTE: omit "atmosphere_hybrid_height_coordinate" : our basic testcase + "atmosphere_sigma_coordinate", "atmosphere_hybrid_sigma_pressure_coordinate", "ocean_sigma_z_coordinate", "ocean_sigma_coordinate", From b475984954a94b5d92730f2dbc2090d3101b2eb1 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Thu, 5 Aug 2021 07:52:50 +0000 Subject: [PATCH 3/5] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- .../tests/unit/aux_factory/test_AtmosphereSigmaFactory.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/lib/iris/tests/unit/aux_factory/test_AtmosphereSigmaFactory.py b/lib/iris/tests/unit/aux_factory/test_AtmosphereSigmaFactory.py index c0bb733278..6a5c9b2025 100644 --- a/lib/iris/tests/unit/aux_factory/test_AtmosphereSigmaFactory.py +++ b/lib/iris/tests/unit/aux_factory/test_AtmosphereSigmaFactory.py @@ -9,10 +9,6 @@ """ -# Import iris.tests first so that some things can be initialised before -# importing anything else. -import iris.tests as tests - from unittest import mock from cf_units import Unit @@ -21,6 +17,10 @@ from iris.aux_factory import AtmosphereSigmaFactory from iris.coords import AuxCoord, DimCoord +# Import iris.tests first so that some things can be initialised before +# importing anything else. +import iris.tests as tests + class Test___init__(tests.IrisTest): def setUp(self): From dec0b5c4ab7441f735a1b9a8accdef22577aa370 Mon Sep 17 00:00:00 2001 From: Manuel Schlund Date: Fri, 13 Aug 2021 09:48:03 +0200 Subject: [PATCH 4/5] Addressed review comments --- lib/iris/aux_factory.py | 222 +++--------------- .../test_AtmosphereSigmaFactory.py | 25 +- 2 files changed, 53 insertions(+), 194 deletions(-) diff --git a/lib/iris/aux_factory.py b/lib/iris/aux_factory.py index e0445f36e7..08317ffdf3 100644 --- a/lib/iris/aux_factory.py +++ b/lib/iris/aux_factory.py @@ -97,20 +97,32 @@ def make_coord(self, coord_dims_func): """ - @abstractmethod def update(self, old_coord, new_coord=None): """ - Notifies the factory of a removal/replacement of a dependency. + Notifies the factory of the removal/replacement of a coordinate + which might be a dependency. Args: * old_coord: - The dependency coordinate to be removed/replaced. + The coordinate to be removed/replaced. * new_coord: - If None, the dependency using old_coord is removed, otherwise - the dependency is updated to use new_coord. + If None, any dependency using old_coord is removed, otherwise + any dependency using old_coord is updated to use new_coord. """ + new_dependencies = self.dependencies + for name, coord in self.dependencies.items(): + if old_coord is coord: + new_dependencies[name] = new_coord + try: + self._check_dependencies(**new_dependencies) + except ValueError as e: + msg = "Failed to update dependencies. " + str(e) + raise ValueError(msg) + else: + setattr(self, name, new_coord) + break def __repr__(self): def arg_text(item): @@ -366,14 +378,15 @@ def _remap_with_bounds(self, dependency_dims, derived_dims): class AtmosphereSigmaFactory(AuxCoordFactory): - """Defines an atmosphere sigma coordinate factory.""" + """Defines an atmosphere sigma coordinate factory with the formula: + p = ptop + sigma * (ps - ptop) + + """ def __init__( self, pressure_at_top=None, sigma=None, surface_air_pressure=None ): - """Create class instance. - - Creates an atmosphere sigma coordinate factory with the formula: + """Creates an atmosphere sigma coordinate factory with the formula: p(n, k, j, i) = pressure_at_top + sigma(k) * (surface_air_pressure(n, j, i) - pressure_at_top) @@ -392,7 +405,7 @@ def __init__( self.sigma = sigma self.surface_air_pressure = surface_air_pressure self.standard_name = "air_pressure" - self.attributes = {"positive": "down"} + self.attributes = {} @staticmethod def _check_dependencies(pressure_at_top, sigma, surface_air_pressure): @@ -470,7 +483,18 @@ def _derive(pressure_at_top, sigma, surface_air_pressure): ) def make_coord(self, coord_dims_func): - """Make new :class:`iris.coords.AuxCoord`.""" + """ + Returns a new :class:`iris.coords.AuxCoord` as defined by this + factory. + + Args: + + * coord_dims_func: + A callable which can return the list of dimensions relevant + to a given coordinate. + See :meth:`iris.cube.Cube.coord_dims()`. + + """ # Which dimensions are relevant? derived_dims = self.derived_dims(coord_dims_func) dependency_dims = self._dependency_dims(coord_dims_func) @@ -529,20 +553,6 @@ def make_coord(self, coord_dims_func): coord_system=self.coord_system, ) - def update(self, old_coord, new_coord=None): - """Notify the factory of the removal/replacement of a coordinate.""" - new_dependencies = self.dependencies - for (name, coord) in self.dependencies.items(): - if old_coord is coord: - new_dependencies[name] = new_coord - try: - self._check_dependencies(**new_dependencies) - except ValueError as exc: - raise ValueError(f"Failed to update dependencies: {exc}") - else: - setattr(self, name, new_coord) - break - class HybridHeightFactory(AuxCoordFactory): """ @@ -912,33 +922,6 @@ def make_coord(self, coord_dims_func): ) return hybrid_pressure - def update(self, old_coord, new_coord=None): - """ - Notifies the factory of the removal/replacement of a coordinate - which might be a dependency. - - Args: - - * old_coord: - The coordinate to be removed/replaced. - * new_coord: - If None, any dependency using old_coord is removed, otherwise - any dependency using old_coord is updated to use new_coord. - - """ - new_dependencies = self.dependencies - for name, coord in self.dependencies.items(): - if old_coord is coord: - new_dependencies[name] = new_coord - try: - self._check_dependencies(**new_dependencies) - except ValueError as e: - msg = "Failed to update dependencies. " + str(e) - raise ValueError(msg) - else: - setattr(self, name, new_coord) - break - class OceanSigmaZFactory(AuxCoordFactory): """Defines an ocean sigma over z coordinate factory.""" @@ -1232,33 +1215,6 @@ def make_coord(self, coord_dims_func): ) return coord - def update(self, old_coord, new_coord=None): - """ - Notifies the factory of the removal/replacement of a coordinate - which might be a dependency. - - Args: - - * old_coord: - The coordinate to be removed/replaced. - * new_coord: - If None, any dependency using old_coord is removed, otherwise - any dependency using old_coord is updated to use new_coord. - - """ - new_dependencies = self.dependencies - for name, coord in self.dependencies.items(): - if old_coord is coord: - new_dependencies[name] = new_coord - try: - self._check_dependencies(**new_dependencies) - except ValueError as e: - msg = "Failed to update dependencies. " + str(e) - raise ValueError(msg) - else: - setattr(self, name, new_coord) - break - class OceanSigmaFactory(AuxCoordFactory): """Defines an ocean sigma coordinate factory.""" @@ -1417,33 +1373,6 @@ def make_coord(self, coord_dims_func): ) return coord - def update(self, old_coord, new_coord=None): - """ - Notifies the factory of the removal/replacement of a coordinate - which might be a dependency. - - Args: - - * old_coord: - The coordinate to be removed/replaced. - * new_coord: - If None, any dependency using old_coord is removed, otherwise - any dependency using old_coord is updated to use new_coord. - - """ - new_dependencies = self.dependencies - for name, coord in self.dependencies.items(): - if old_coord is coord: - new_dependencies[name] = new_coord - try: - self._check_dependencies(**new_dependencies) - except ValueError as e: - msg = "Failed to update dependencies. " + str(e) - raise ValueError(msg) - else: - setattr(self, name, new_coord) - break - class OceanSg1Factory(AuxCoordFactory): """Defines an Ocean s-coordinate, generic form 1 factory.""" @@ -1640,33 +1569,6 @@ def make_coord(self, coord_dims_func): ) return coord - def update(self, old_coord, new_coord=None): - """ - Notifies the factory of the removal/replacement of a coordinate - which might be a dependency. - - Args: - - * old_coord: - The coordinate to be removed/replaced. - * new_coord: - If None, any dependency using old_coord is removed, otherwise - any dependency using old_coord is updated to use new_coord. - - """ - new_dependencies = self.dependencies - for name, coord in self.dependencies.items(): - if old_coord is coord: - new_dependencies[name] = new_coord - try: - self._check_dependencies(**new_dependencies) - except ValueError as e: - msg = "Failed to update dependencies. " + str(e) - raise ValueError(msg) - else: - setattr(self, name, new_coord) - break - class OceanSFactory(AuxCoordFactory): """Defines an Ocean s-coordinate factory.""" @@ -1865,33 +1767,6 @@ def make_coord(self, coord_dims_func): ) return coord - def update(self, old_coord, new_coord=None): - """ - Notifies the factory of the removal/replacement of a coordinate - which might be a dependency. - - Args: - - * old_coord: - The coordinate to be removed/replaced. - * new_coord: - If None, any dependency using old_coord is removed, otherwise - any dependency using old_coord is updated to use new_coord. - - """ - new_dependencies = self.dependencies - for name, coord in self.dependencies.items(): - if old_coord is coord: - new_dependencies[name] = new_coord - try: - self._check_dependencies(**new_dependencies) - except ValueError as e: - msg = "Failed to update dependencies. " + str(e) - raise ValueError(msg) - else: - setattr(self, name, new_coord) - break - class OceanSg2Factory(AuxCoordFactory): """Defines an Ocean s-coordinate, generic form 2 factory.""" @@ -2088,30 +1963,3 @@ def make_coord(self, coord_dims_func): coord_system=self.coord_system, ) return coord - - def update(self, old_coord, new_coord=None): - """ - Notifies the factory of the removal/replacement of a coordinate - which might be a dependency. - - Args: - - * old_coord: - The coordinate to be removed/replaced. - * new_coord: - If None, any dependency using old_coord is removed, otherwise - any dependency using old_coord is updated to use new_coord. - - """ - new_dependencies = self.dependencies - for name, coord in self.dependencies.items(): - if old_coord is coord: - new_dependencies[name] = new_coord - try: - self._check_dependencies(**new_dependencies) - except ValueError as e: - msg = "Failed to update dependencies. " + str(e) - raise ValueError(msg) - else: - setattr(self, name, new_coord) - break diff --git a/lib/iris/tests/unit/aux_factory/test_AtmosphereSigmaFactory.py b/lib/iris/tests/unit/aux_factory/test_AtmosphereSigmaFactory.py index 6a5c9b2025..6e417a3b38 100644 --- a/lib/iris/tests/unit/aux_factory/test_AtmosphereSigmaFactory.py +++ b/lib/iris/tests/unit/aux_factory/test_AtmosphereSigmaFactory.py @@ -9,6 +9,10 @@ """ +# Import iris.tests first so that some things can be initialised before +# importing anything else. +import iris.tests as tests # isort:skip + from unittest import mock from cf_units import Unit @@ -17,10 +21,6 @@ from iris.aux_factory import AtmosphereSigmaFactory from iris.coords import AuxCoord, DimCoord -# Import iris.tests first so that some things can be initialised before -# importing anything else. -import iris.tests as tests - class Test___init__(tests.IrisTest): def setUp(self): @@ -33,21 +33,27 @@ def setUp(self): surface_air_pressure=self.surface_air_pressure, ) - def test_insufficient_coordinates(self): + def test_insufficient_coordinates_no_args(self): with self.assertRaises(ValueError): AtmosphereSigmaFactory() + + def test_insufficient_coordinates_no_ptop(self): with self.assertRaises(ValueError): AtmosphereSigmaFactory( pressure_at_top=None, sigma=self.sigma, surface_air_pressure=self.surface_air_pressure, ) + + def test_insufficient_coordinates_no_sigma(self): with self.assertRaises(ValueError): AtmosphereSigmaFactory( pressure_at_top=self.pressure_at_top, sigma=None, surface_air_pressure=self.surface_air_pressure, ) + + def test_insufficient_coordinates_no_ps(self): with self.assertRaises(ValueError): AtmosphereSigmaFactory( pressure_at_top=self.pressure_at_top, @@ -60,6 +66,7 @@ def test_ptop_shapes(self): self.pressure_at_top.shape = shape AtmosphereSigmaFactory(**self.kwargs) + def test_ptop_invalid_shapes(self): for shape in [(2,), (1, 1)]: self.pressure_at_top.shape = shape with self.assertRaises(ValueError): @@ -70,6 +77,7 @@ def test_sigma_bounds(self): self.sigma.nbounds = n_bounds AtmosphereSigmaFactory(**self.kwargs) + def test_sigma_invalid_bounds(self): for n_bounds in [-1, 1, 3]: self.sigma.nbounds = n_bounds with self.assertRaises(ValueError): @@ -80,6 +88,7 @@ def test_sigma_units(self): self.sigma.units = Unit(units) AtmosphereSigmaFactory(**self.kwargs) + def test_sigma_invalid_units(self): for units in ["Pa", "m"]: self.sigma.units = Unit(units) with self.assertRaises(ValueError): @@ -91,6 +100,7 @@ def test_ptop_ps_units(self): self.surface_air_pressure.units = Unit(units[1]) AtmosphereSigmaFactory(**self.kwargs) + def test_ptop_ps_invalid_units(self): for units in [("Pa", "1"), ("1", "Pa"), ("bar", "Pa"), ("Pa", "hPa")]: self.pressure_at_top.units = Unit(units[0]) self.surface_air_pressure.units = Unit(units[1]) @@ -103,6 +113,7 @@ def test_ptop_units(self): self.surface_air_pressure.units = Unit(units) AtmosphereSigmaFactory(**self.kwargs) + def test_ptop_invalid_units(self): for units in ["1", "m", "kg", None]: self.pressure_at_top.units = Unit(units) self.surface_air_pressure.units = Unit(units) @@ -127,7 +138,7 @@ def test_values(self): class Test__derive(tests.IrisTest): - def test_function(self): + def test_function_scalar(self): assert AtmosphereSigmaFactory._derive(0, 0, 0) == 0 assert AtmosphereSigmaFactory._derive(3, 0, 0) == 3 assert AtmosphereSigmaFactory._derive(0, 5, 0) == 0 @@ -137,6 +148,7 @@ def test_function(self): assert AtmosphereSigmaFactory._derive(0, 5, 7) == 35 assert AtmosphereSigmaFactory._derive(3, 5, 7) == 23 + def test_function_array(self): ptop = 3 sigma = np.array([2, 4]) ps = np.arange(4).reshape(2, 2) @@ -167,7 +179,6 @@ def derive(pressure_at_top, sigma, surface_air_pressure, coord=True): result, standard_name=name, units="Pa", - attributes=dict(positive="down"), ) return result From a81b1cd67b43f49a7cee2d02ad88c65313097587 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Fri, 13 Aug 2021 07:49:51 +0000 Subject: [PATCH 5/5] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- lib/iris/aux_factory.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/lib/iris/aux_factory.py b/lib/iris/aux_factory.py index 08317ffdf3..f49de62b3f 100644 --- a/lib/iris/aux_factory.py +++ b/lib/iris/aux_factory.py @@ -379,7 +379,7 @@ def _remap_with_bounds(self, dependency_dims, derived_dims): class AtmosphereSigmaFactory(AuxCoordFactory): """Defines an atmosphere sigma coordinate factory with the formula: - p = ptop + sigma * (ps - ptop) + p = ptop + sigma * (ps - ptop) """