diff --git a/lib/iris/aux_factory.py b/lib/iris/aux_factory.py
index efd6607fc6..6d5427ba65 100644
--- a/lib/iris/aux_factory.py
+++ b/lib/iris/aux_factory.py
@@ -26,6 +26,7 @@
from abc import ABCMeta, abstractmethod, abstractproperty
import warnings
+import dask.array as da
import numpy as np
from iris._cube_coord_common import CFVariableMixin
@@ -110,7 +111,7 @@ def arg_text(item):
def derived_dims(self, coord_dims_func):
"""
- Returns the virtual dim-mapping for the derived coordinate.
+ Returns the cube dimensions for the derived coordinate.
Args:
@@ -119,6 +120,10 @@ def derived_dims(self, coord_dims_func):
to a given coordinate.
See :meth:`iris.cube.Cube.coord_dims()`.
+ Returns:
+
+ A sorted list of cube dimension numbers.
+
"""
# Which dimensions are relevant?
# e.g. If sigma -> [1] and orog -> [2, 3] then result = [1, 2, 3]
@@ -173,8 +178,15 @@ def _dependency_dims(self, coord_dims_func):
@staticmethod
def _nd_bounds(coord, dims, ndim):
"""
- Returns the coord's bounds in Cube-orientation and
- broadcastable to N dimensions.
+ Return a lazy bounds array for a dependency coordinate, 'coord'.
+
+ The result is aligned to the first 'ndim' cube dimensions, and
+ expanded to the full ('ndim'+1)-dimensional shape.
+
+ The value of 'ndim' must be >= the highest cube dimension of the
+ dependency coordinate.
+
+ The extra final result dimension ('ndim'-th) is the bounds dimension.
Example:
coord.shape == (70,)
@@ -188,8 +200,8 @@ def _nd_bounds(coord, dims, ndim):
# Transpose to be consistent with the Cube.
sorted_pairs = sorted(enumerate(dims), key=lambda pair: pair[1])
transpose_order = [pair[0] for pair in sorted_pairs] + [len(dims)]
- bounds = coord.core_bounds()
- if dims:
+ bounds = coord.lazy_bounds()
+ if dims and transpose_order != list(range(len(dims))):
bounds = bounds.transpose(transpose_order)
# Figure out the n-dimensional shape.
@@ -202,8 +214,13 @@ def _nd_bounds(coord, dims, ndim):
@staticmethod
def _nd_points(coord, dims, ndim):
"""
- Returns the coord's points in Cube-orientation and
- broadcastable to N dimensions.
+ Return a lazy points array for a dependency coordinate, 'coord'.
+
+ The result is aligned to the first 'ndim' cube dimensions, and
+ expanded to the full 'ndim'-dimensional shape.
+
+ The value of 'ndim' must be >= the highest cube dimension of the
+ dependency coordinate.
Example:
coord.shape == (4, 3)
@@ -216,7 +233,7 @@ def _nd_points(coord, dims, ndim):
# Transpose to be consistent with the Cube.
sorted_pairs = sorted(enumerate(dims), key=lambda pair: pair[1])
transpose_order = [pair[0] for pair in sorted_pairs]
- points = coord.core_points()
+ points = coord.lazy_points()
if dims and transpose_order != list(range(len(dims))):
points = points.transpose(transpose_order)
@@ -235,6 +252,16 @@ def _nd_points(coord, dims, ndim):
return points
def _remap(self, dependency_dims, derived_dims):
+ """
+ Return a mapping from dependency names to coordinate points arrays.
+
+ For dependencies that are present, the values are all expanded and
+ aligned to the same dimensions, which is the full set of all the
+ dependency dimensions.
+ These non-missing values are all lazy arrays.
+ Missing dependencies, however, are assigned a scalar value of 0.0.
+
+ """
if derived_dims:
ndim = max(derived_dims) + 1
else:
@@ -263,6 +290,20 @@ def _remap(self, dependency_dims, derived_dims):
return nd_points_by_key
def _remap_with_bounds(self, dependency_dims, derived_dims):
+ """
+ Return a mapping from dependency names to coordinate bounds arrays.
+
+ For dependencies that are present, the values are all expanded and
+ aligned to the same dimensions, which is the full set of all the
+ dependency dimensions, plus an extra bounds dimension.
+ These non-missing values are all lazy arrays.
+ Missing dependencies, however, are assigned a scalar value of 0.0.
+
+ Where a dependency coordinate has no bounds, then the associated value
+ is taken from its points array, but reshaped to have an extra bounds
+ dimension of length 1.
+
+ """
if derived_dims:
ndim = max(derived_dims) + 1
else:
@@ -307,31 +348,6 @@ def _remap_with_bounds(self, dependency_dims, derived_dims):
nd_values_by_key[key] = nd_values
return nd_values_by_key
- def _shape(self, nd_values_by_key):
- nd_values = sorted(nd_values_by_key.values(),
- key=lambda value: value.ndim)
- shape = list(nd_values.pop().shape)
- for array in nd_values:
- for i, size in enumerate(array.shape):
- if size > 1:
- # NB. If there's an inconsistency it can only come
- # from a mismatch in the number of bounds (the Cube
- # ensures the other dimensions must match).
- # But we can't afford to raise an error now - it'd
- # break Cube.derived_coords. Instead, we let the
- # error happen when the derived coordinate's bounds
- # are accessed.
- shape[i] = size
- return shape
-
- def _dtype(self, arrays_by_key, **other_args):
- dummy_args = {}
- for key, array in six.iteritems(arrays_by_key):
- dummy_args[key] = np.zeros(1, dtype=array.dtype)
- dummy_args.update(other_args)
- dummy_data = self._derive(**dummy_args)
- return dummy_data.dtype
-
class HybridHeightFactory(AuxCoordFactory):
"""
@@ -428,27 +444,23 @@ def make_coord(self, coord_dims_func):
# Build the bounds array.
nd_values_by_key = self._remap_with_bounds(dependency_dims,
derived_dims)
-
- # Define the function here to obtain a closure.
- def calc_bounds():
- delta = nd_values_by_key['delta']
- sigma = nd_values_by_key['sigma']
- orography = nd_values_by_key['orography']
- ok_bound_shapes = [(), (1,), (2,)]
- if delta.shape[-1:] not in ok_bound_shapes:
- raise ValueError('Invalid delta coordinate bounds.')
- if sigma.shape[-1:] not in ok_bound_shapes:
- raise ValueError('Invalid sigma coordinate bounds.')
- if orography.shape[-1:] not in [(), (1,)]:
- warnings.warn('Orography coordinate has bounds. '
- 'These are being disregarded.',
- UserWarning, stacklevel=2)
- orography_pts = nd_points_by_key['orography']
- orography_pts_shape = list(orography_pts.shape)
- orography = orography_pts.reshape(
- orography_pts_shape.append(1))
- return self._derive(delta, sigma, orography)
- bounds = calc_bounds()
+ delta = nd_values_by_key['delta']
+ sigma = nd_values_by_key['sigma']
+ orography = nd_values_by_key['orography']
+ ok_bound_shapes = [(), (1,), (2,)]
+ if delta.shape[-1:] not in ok_bound_shapes:
+ raise ValueError('Invalid delta coordinate bounds.')
+ if sigma.shape[-1:] not in ok_bound_shapes:
+ raise ValueError('Invalid sigma coordinate bounds.')
+ if orography.shape[-1:] not in [(), (1,)]:
+ warnings.warn('Orography coordinate has bounds. '
+ 'These are being disregarded.',
+ UserWarning, stacklevel=2)
+ orography_pts = nd_points_by_key['orography']
+ bds_shape = list(orography_pts.shape) + [1]
+ orography = orography_pts.reshape(bds_shape)
+
+ bounds = self._derive(delta, sigma, orography)
hybrid_height = iris.coords.AuxCoord(points,
standard_name=self.standard_name,
@@ -618,28 +630,24 @@ def make_coord(self, coord_dims_func):
# Build the bounds array.
nd_values_by_key = self._remap_with_bounds(dependency_dims,
derived_dims)
-
- # Define the function here to obtain a closure.
- def calc_bounds():
- delta = nd_values_by_key['delta']
- sigma = nd_values_by_key['sigma']
- surface_air_pressure = nd_values_by_key['surface_air_pressure']
- ok_bound_shapes = [(), (1,), (2,)]
- if delta.shape[-1:] not in ok_bound_shapes:
- raise ValueError('Invalid delta coordinate bounds.')
- if sigma.shape[-1:] not in ok_bound_shapes:
- raise ValueError('Invalid sigma coordinate bounds.')
- 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']
- surface_air_pressure_pts_shape = list(
- surface_air_pressure_pts.shape)
- surface_air_pressure = surface_air_pressure_pts.reshape(
- surface_air_pressure_pts_shape.append(1))
- return self._derive(delta, sigma, surface_air_pressure)
- bounds = calc_bounds()
+ delta = nd_values_by_key['delta']
+ sigma = nd_values_by_key['sigma']
+ surface_air_pressure = nd_values_by_key['surface_air_pressure']
+ ok_bound_shapes = [(), (1,), (2,)]
+ if delta.shape[-1:] not in ok_bound_shapes:
+ raise ValueError('Invalid delta coordinate bounds.')
+ if sigma.shape[-1:] not in ok_bound_shapes:
+ raise ValueError('Invalid sigma coordinate bounds.')
+ 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(delta, sigma, surface_air_pressure)
hybrid_pressure = iris.coords.AuxCoord(
points, standard_name=self.standard_name, long_name=self.long_name,
@@ -784,20 +792,65 @@ def dependencies(self):
depth_c=self.depth_c, nsigma=self.nsigma, zlev=self.zlev)
def _derive(self, sigma, eta, depth, depth_c,
- zlev, shape, nsigma_slice):
+ zlev, nsigma, coord_dims_func):
+ # Calculate the index of the 'z' dimension in the input arrays.
+ # First find the cube 'z' dimension ...
+ [cube_z_dim] = coord_dims_func(self.dependencies['zlev'])
+ # ... then calculate the corresponding dependency dimension.
+ derived_cubedims = self.derived_dims(coord_dims_func)
+ z_dim = derived_cubedims.index(cube_z_dim)
+
+ # Calculate the result shape as a combination of all the inputs.
+ # Note: all the inputs have the same number of dimensions >= 1, except
+ # for any missing dependencies, which have scalar values.
+ allshapes = np.array(
+ [el.shape
+ for el in (sigma, eta, depth, depth_c, zlev)
+ if el.ndim > 0])
+ result_shape = list(np.max(allshapes, axis=0))
+ ndims = len(result_shape)
+
+ # Make a slice tuple to index the first nsigma z-levels.
+ z_slices_nsigma = [slice(None)] * ndims
+ z_slices_nsigma[z_dim] = slice(0, int(nsigma))
+ z_slices_nsigma = tuple(z_slices_nsigma)
+ # Make a slice tuple to index the remaining z-levels.
+ z_slices_rest = [slice(None)] * ndims
+ z_slices_rest[z_dim] = slice(int(nsigma), None)
+ z_slices_rest = tuple(z_slices_rest)
+
# Perform the ocean sigma over z coordinate nsigma slice.
if eta.ndim:
- eta = eta[nsigma_slice]
+ eta = eta[z_slices_nsigma]
if sigma.ndim:
- sigma = sigma[nsigma_slice]
+ sigma = sigma[z_slices_nsigma]
if depth.ndim:
- depth = depth[nsigma_slice]
+ depth = depth[z_slices_nsigma]
# Note that, this performs a point-wise minimum.
- temp = eta + sigma * (np.minimum(depth_c, depth) + eta)
- # Calculate the final derived result.
- result = np.ones(shape, dtype=temp.dtype) * zlev
- result[nsigma_slice] = temp
-
+ nsigma_levs = eta + sigma * (da.minimum(depth_c, depth) + eta)
+
+ # Make a result-shaped lazy "ones" array for expanding partial results.
+ # Note: for the 'chunks' arg, we try to use [1, 1, ... ny, nx].
+ # This calculation could be assuming too much in some cases, as we
+ # don't actually check the dimensions of our dependencies anywhere.
+ result_chunks = result_shape
+ if len(result_shape) > 1:
+ result_chunks = [1] * len(result_shape)
+ result_chunks[-2:] = result_shape[-2:]
+ ones_full_result = da.ones(result_shape, chunks=result_chunks,
+ dtype=zlev.dtype)
+
+ # Expand nsigma_levs to its full required shape : needed as the
+ # calculated result may have a fixed size of 1 in some dimensions.
+ result_nsigma_levs = nsigma_levs * ones_full_result[z_slices_nsigma]
+
+ # Likewise, expand zlev to its full required shape.
+ result_rest_levs = (zlev[z_slices_rest] *
+ ones_full_result[z_slices_rest])
+
+ # Combine nsigma and 'rest' levels for the final result.
+ result = da.concatenate([result_nsigma_levs, result_rest_levs],
+ axis=z_dim)
return result
def make_coord(self, coord_dims_func):
@@ -817,59 +870,47 @@ def make_coord(self, coord_dims_func):
# Build the points array.
nd_points_by_key = self._remap(dependency_dims, derived_dims)
- points_shape = self._shape(nd_points_by_key)
-
- # Calculate the nsigma slice.
- nsigma_slice = [slice(None)] * len(derived_dims)
- dim, = dependency_dims['zlev']
- index = derived_dims.index(dim)
- nsigma_slice[index] = slice(0, int(nd_points_by_key['nsigma']))
- nsigma_slice = tuple(nsigma_slice)
+ [nsigma] = nd_points_by_key['nsigma']
points = self._derive(nd_points_by_key['sigma'],
nd_points_by_key['eta'],
nd_points_by_key['depth'],
nd_points_by_key['depth_c'],
nd_points_by_key['zlev'],
- points_shape,
- nsigma_slice)
+ nsigma,
+ coord_dims_func)
bounds = None
if self.zlev.nbounds or (self.sigma and self.sigma.nbounds):
# Build the bounds array.
nd_values_by_key = self._remap_with_bounds(dependency_dims,
derived_dims)
- bounds_shape = self._shape(nd_values_by_key)
- nsigma_slice_bounds = nsigma_slice + (slice(None),)
-
- # Define the function here to obtain a closure.
- def calc_bounds():
- valid_shapes = [(), (1,), (2,)]
- for key in ('sigma', 'zlev'):
- if nd_values_by_key[key].shape[-1:] not in valid_shapes:
- name = self.dependencies[key].name()
- msg = 'Invalid bounds for {} ' \
- 'coordinate {!r}.'.format(key, name)
- raise ValueError(msg)
- valid_shapes.pop()
- for key in ('eta', 'depth', 'depth_c', 'nsigma'):
- if nd_values_by_key[key].shape[-1:] not in valid_shapes:
- name = self.dependencies[key].name()
- msg = 'The {} coordinate {!r} has bounds. ' \
- 'These are being disregarded.'.format(key, name)
- warnings.warn(msg, UserWarning, stacklevel=2)
- # Swap bounds with points.
- shape = list(nd_points_by_key[key].shape)
- bounds = nd_points_by_key[key].reshape(shape.append(1))
- nd_values_by_key[key] = bounds
- return self._derive(nd_values_by_key['sigma'],
- nd_values_by_key['eta'],
- nd_values_by_key['depth'],
- nd_values_by_key['depth_c'],
- nd_values_by_key['zlev'],
- bounds_shape,
- nsigma_slice_bounds)
- bounds = calc_bounds()
+ valid_shapes = [(), (1,), (2,)]
+ for key in ('sigma', 'zlev'):
+ if nd_values_by_key[key].shape[-1:] not in valid_shapes:
+ name = self.dependencies[key].name()
+ msg = 'Invalid bounds for {} ' \
+ 'coordinate {!r}.'.format(key, name)
+ raise ValueError(msg)
+ valid_shapes.pop()
+ for key in ('eta', 'depth', 'depth_c', 'nsigma'):
+ if nd_values_by_key[key].shape[-1:] not in valid_shapes:
+ name = self.dependencies[key].name()
+ msg = 'The {} coordinate {!r} has bounds. ' \
+ 'These are being disregarded.'.format(key, name)
+ warnings.warn(msg, UserWarning, stacklevel=2)
+ # Swap bounds with points.
+ bds_shape = list(nd_points_by_key[key].shape) + [1]
+ bounds = nd_points_by_key[key].reshape(bds_shape)
+ nd_values_by_key[key] = bounds
+
+ bounds = self._derive(nd_values_by_key['sigma'],
+ nd_values_by_key['eta'],
+ nd_values_by_key['depth'],
+ nd_values_by_key['depth_c'],
+ nd_values_by_key['zlev'],
+ nsigma,
+ coord_dims_func)
coord = iris.coords.AuxCoord(points,
standard_name=self.standard_name,
@@ -1010,32 +1051,28 @@ def make_coord(self, coord_dims_func):
# Build the bounds array.
nd_values_by_key = self._remap_with_bounds(dependency_dims,
derived_dims)
-
- # Define the function here to obtain a closure.
- def calc_bounds():
- valid_shapes = [(), (1,), (2,)]
- key = 'sigma'
+ valid_shapes = [(), (1,), (2,)]
+ key = 'sigma'
+ if nd_values_by_key[key].shape[-1:] not in valid_shapes:
+ name = self.dependencies[key].name()
+ msg = 'Invalid bounds for {} ' \
+ 'coordinate {!r}.'.format(key, name)
+ raise ValueError(msg)
+ valid_shapes.pop()
+ for key in ('eta', 'depth'):
if nd_values_by_key[key].shape[-1:] not in valid_shapes:
name = self.dependencies[key].name()
- msg = 'Invalid bounds for {} ' \
- 'coordinate {!r}.'.format(key, name)
- raise ValueError(msg)
- valid_shapes.pop()
- for key in ('eta', 'depth'):
- if nd_values_by_key[key].shape[-1:] not in valid_shapes:
- name = self.dependencies[key].name()
- msg = 'The {} coordinate {!r} has bounds. ' \
- 'These are being disregarded.'.format(key, name)
- warnings.warn(msg, UserWarning, stacklevel=2)
- # Swap bounds with points.
- shape = list(nd_points_by_key[key].shape)
- bounds = nd_points_by_key[key].reshape(shape.append(1))
- nd_values_by_key[key] = bounds
- return self._derive(nd_values_by_key['sigma'],
- nd_values_by_key['eta'],
- nd_values_by_key['depth'])
-
- bounds = calc_bounds()
+ msg = 'The {} coordinate {!r} has bounds. ' \
+ 'These are being disregarded.'.format(key, name)
+ warnings.warn(msg, UserWarning, stacklevel=2)
+ # Swap bounds with points.
+ bds_shape = list(nd_points_by_key[key].shape) + [1]
+ bounds = nd_points_by_key[key].reshape(bds_shape)
+ nd_values_by_key[key] = bounds
+
+ bounds = self._derive(nd_values_by_key['sigma'],
+ nd_values_by_key['eta'],
+ nd_values_by_key['depth'])
coord = iris.coords.AuxCoord(points,
standard_name=self.standard_name,
@@ -1198,33 +1235,30 @@ def make_coord(self, coord_dims_func):
# Build the bounds array.
nd_values_by_key = self._remap_with_bounds(dependency_dims,
derived_dims)
-
- # Define the function here to obtain a closure.
- def calc_bounds():
- valid_shapes = [(), (1,), (2,)]
- key = 's'
+ valid_shapes = [(), (1,), (2,)]
+ key = 's'
+ if nd_values_by_key[key].shape[-1:] not in valid_shapes:
+ name = self.dependencies[key].name()
+ msg = 'Invalid bounds for {} ' \
+ 'coordinate {!r}.'.format(key, name)
+ raise ValueError(msg)
+ valid_shapes.pop()
+ for key in ('eta', 'depth', 'depth_c'):
if nd_values_by_key[key].shape[-1:] not in valid_shapes:
name = self.dependencies[key].name()
- msg = 'Invalid bounds for {} ' \
- 'coordinate {!r}.'.format(key, name)
- raise ValueError(msg)
- valid_shapes.pop()
- for key in ('eta', 'depth', 'depth_c'):
- if nd_values_by_key[key].shape[-1:] not in valid_shapes:
- name = self.dependencies[key].name()
- msg = 'The {} coordinate {!r} has bounds. ' \
- 'These are being disregarded.'.format(key, name)
- warnings.warn(msg, UserWarning, stacklevel=2)
- # Swap bounds with points.
- shape = list(nd_points_by_key[key].shape)
- bounds = nd_points_by_key[key].reshape(shape.append(1))
- nd_values_by_key[key] = bounds
- return self._derive(nd_values_by_key['s'],
- nd_values_by_key['c'],
- nd_values_by_key['eta'],
- nd_values_by_key['depth'],
- nd_values_by_key['depth_c'])
- bounds = calc_bounds()
+ msg = 'The {} coordinate {!r} has bounds. ' \
+ 'These are being disregarded.'.format(key, name)
+ warnings.warn(msg, UserWarning, stacklevel=2)
+ # Swap bounds with points.
+ bds_shape = list(nd_points_by_key[key].shape) + [1]
+ bounds = nd_points_by_key[key].reshape(bds_shape)
+ nd_values_by_key[key] = bounds
+
+ bounds = self._derive(nd_values_by_key['s'],
+ nd_values_by_key['c'],
+ nd_values_by_key['eta'],
+ nd_values_by_key['depth'],
+ nd_values_by_key['depth_c'])
coord = iris.coords.AuxCoord(points,
standard_name=self.standard_name,
@@ -1353,8 +1387,8 @@ def dependencies(self):
b=self.b, depth_c=self.depth_c)
def _derive(self, s, eta, depth, a, b, depth_c):
- c = ((1 - b) * np.sinh(a * s) / np.sinh(a) + b *
- (np.tanh(a * (s + 0.5)) / (2 * np.tanh(0.5 * a)) - 0.5))
+ c = ((1 - b) * da.sinh(a * s) / da.sinh(a) + b *
+ (da.tanh(a * (s + 0.5)) / (2 * da.tanh(0.5 * a)) - 0.5))
return eta * (1 + s) + depth_c * s + (depth - depth_c) * c
def make_coord(self, coord_dims_func):
@@ -1386,34 +1420,31 @@ def make_coord(self, coord_dims_func):
# Build the bounds array.
nd_values_by_key = self._remap_with_bounds(dependency_dims,
derived_dims)
-
- # Define the function here to obtain a closure.
- def calc_bounds():
- valid_shapes = [(), (1,), (2,)]
- key = 's'
+ valid_shapes = [(), (1,), (2,)]
+ key = 's'
+ if nd_values_by_key[key].shape[-1:] not in valid_shapes:
+ name = self.dependencies[key].name()
+ msg = 'Invalid bounds for {} ' \
+ 'coordinate {!r}.'.format(key, name)
+ raise ValueError(msg)
+ valid_shapes.pop()
+ for key in ('eta', 'depth', 'a', 'b', 'depth_c'):
if nd_values_by_key[key].shape[-1:] not in valid_shapes:
name = self.dependencies[key].name()
- msg = 'Invalid bounds for {} ' \
- 'coordinate {!r}.'.format(key, name)
- raise ValueError(msg)
- valid_shapes.pop()
- for key in ('eta', 'depth', 'a', 'b', 'depth_c'):
- if nd_values_by_key[key].shape[-1:] not in valid_shapes:
- name = self.dependencies[key].name()
- msg = 'The {} coordinate {!r} has bounds. ' \
- 'These are being disregarded.'.format(key, name)
- warnings.warn(msg, UserWarning, stacklevel=2)
- # Swap bounds with points.
- shape = list(nd_points_by_key[key].shape)
- bounds = nd_points_by_key[key].reshape(shape.append(1))
- nd_values_by_key[key] = bounds
- return self._derive(nd_values_by_key['s'],
- nd_values_by_key['eta'],
- nd_values_by_key['depth'],
- nd_values_by_key['a'],
- nd_values_by_key['b'],
- nd_values_by_key['depth_c'])
- bounds = calc_bounds()
+ msg = 'The {} coordinate {!r} has bounds. ' \
+ 'These are being disregarded.'.format(key, name)
+ warnings.warn(msg, UserWarning, stacklevel=2)
+ # Swap bounds with points.
+ bds_shape = list(nd_points_by_key[key].shape) + [1]
+ bounds = nd_points_by_key[key].reshape(bds_shape)
+ nd_values_by_key[key] = bounds
+
+ bounds = self._derive(nd_values_by_key['s'],
+ nd_values_by_key['eta'],
+ nd_values_by_key['depth'],
+ nd_values_by_key['a'],
+ nd_values_by_key['b'],
+ nd_values_by_key['depth_c'])
coord = iris.coords.AuxCoord(points,
standard_name=self.standard_name,
@@ -1577,33 +1608,30 @@ def make_coord(self, coord_dims_func):
# Build the bounds array.
nd_values_by_key = self._remap_with_bounds(dependency_dims,
derived_dims)
-
- # Define the function here to obtain a closure.
- def calc_bounds():
- valid_shapes = [(), (1,), (2,)]
- key = 's'
+ valid_shapes = [(), (1,), (2,)]
+ key = 's'
+ if nd_values_by_key[key].shape[-1:] not in valid_shapes:
+ name = self.dependencies[key].name()
+ msg = 'Invalid bounds for {} ' \
+ 'coordinate {!r}.'.format(key, name)
+ raise ValueError(msg)
+ valid_shapes.pop()
+ for key in ('eta', 'depth', 'depth_c'):
if nd_values_by_key[key].shape[-1:] not in valid_shapes:
name = self.dependencies[key].name()
- msg = 'Invalid bounds for {} ' \
- 'coordinate {!r}.'.format(key, name)
- raise ValueError(msg)
- valid_shapes.pop()
- for key in ('eta', 'depth', 'depth_c'):
- if nd_values_by_key[key].shape[-1:] not in valid_shapes:
- name = self.dependencies[key].name()
- msg = 'The {} coordinate {!r} has bounds. ' \
- 'These are being disregarded.'.format(key, name)
- warnings.warn(msg, UserWarning, stacklevel=2)
- # Swap bounds with points.
- shape = list(nd_points_by_key[key].shape)
- bounds = nd_points_by_key[key].reshape(shape.append(1))
- nd_values_by_key[key] = bounds
- return self._derive(nd_values_by_key['s'],
- nd_values_by_key['c'],
- nd_values_by_key['eta'],
- nd_values_by_key['depth'],
- nd_values_by_key['depth_c'])
- bounds = calc_bounds()
+ msg = 'The {} coordinate {!r} has bounds. ' \
+ 'These are being disregarded.'.format(key, name)
+ warnings.warn(msg, UserWarning, stacklevel=2)
+ # Swap bounds with points.
+ bds_shape = list(nd_points_by_key[key].shape) + [1]
+ bounds = nd_points_by_key[key].reshape(bds_shape)
+ nd_values_by_key[key] = bounds
+
+ bounds = self._derive(nd_values_by_key['s'],
+ nd_values_by_key['c'],
+ nd_values_by_key['eta'],
+ nd_values_by_key['depth'],
+ nd_values_by_key['depth_c'])
coord = iris.coords.AuxCoord(points,
standard_name=self.standard_name,
diff --git a/lib/iris/tests/__init__.py b/lib/iris/tests/__init__.py
index d59b2a1678..b420bba961 100644
--- a/lib/iris/tests/__init__.py
+++ b/lib/iris/tests/__init__.py
@@ -512,6 +512,20 @@ def assertXMLElement(self, obj, reference_filename):
def assertArrayEqual(self, a, b, err_msg=''):
np.testing.assert_array_equal(a, b, err_msg=err_msg)
+ def assertRaisesRegexp(self, *args, **kwargs):
+ """
+ Emulate the old :meth:`unittest.TestCase.assertRaisesRegexp`.
+
+ Because the original function is now deprecated in Python 3.
+ Now calls :meth:`six.assertRaisesRegex()` (no final "p") instead.
+ It is the same, except for providing an additional 'msg' argument.
+
+ """
+ # Note: invoke via parent class to avoid recursion as, in Python 2,
+ # "six.assertRaisesRegex" calls getattr(self, 'assertRaisesRegexp').
+ return six.assertRaisesRegex(super(IrisTest_nometa, self),
+ *args, **kwargs)
+
def _assertMaskedArray(self, assertion, a, b, strict, **kwargs):
# Define helper function to extract unmasked values as a 1d
# array.
diff --git a/lib/iris/tests/integration/aux_factory/__init__.py b/lib/iris/tests/integration/aux_factory/__init__.py
new file mode 100644
index 0000000000..f49ea330c1
--- /dev/null
+++ b/lib/iris/tests/integration/aux_factory/__init__.py
@@ -0,0 +1,20 @@
+# (C) British Crown Copyright 2017, Met Office
+#
+# This file is part of Iris.
+#
+# Iris is free software: you can redistribute it and/or modify it under
+# the terms of the GNU Lesser General Public License as published by the
+# Free Software Foundation, either version 3 of the License, or
+# (at your option) any later version.
+#
+# Iris is distributed in the hope that it will be useful,
+# but WITHOUT ANY WARRANTY; without even the implied warranty of
+# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+# GNU Lesser General Public License for more details.
+#
+# You should have received a copy of the GNU Lesser General Public License
+# along with Iris. If not, see .
+"""Integration tests for the :mod:`iris.aux_factory` package."""
+
+from __future__ import (absolute_import, division, print_function)
+from six.moves import (filter, input, map, range, zip) # noqa
diff --git a/lib/iris/tests/integration/aux_factory/test_OceanSigmaZFactory.py b/lib/iris/tests/integration/aux_factory/test_OceanSigmaZFactory.py
new file mode 100644
index 0000000000..a9245e0333
--- /dev/null
+++ b/lib/iris/tests/integration/aux_factory/test_OceanSigmaZFactory.py
@@ -0,0 +1,189 @@
+# (C) British Crown Copyright 2017, Met Office
+#
+# This file is part of Iris.
+#
+# Iris is free software: you can redistribute it and/or modify it under
+# the terms of the GNU Lesser General Public License as published by the
+# Free Software Foundation, either version 3 of the License, or
+# (at your option) any later version.
+#
+# Iris is distributed in the hope that it will be useful,
+# but WITHOUT ANY WARRANTY; without even the implied warranty of
+# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+# GNU Lesser General Public License for more details.
+#
+# You should have received a copy of the GNU Lesser General Public License
+# along with Iris. If not, see .
+"""
+Integratation tests for the
+`iris.aux_factory.OceanSigmaZFactory` class.
+
+"""
+
+from __future__ import (absolute_import, division, print_function)
+from six.moves import (filter, input, map, range, zip) # noqa
+
+# Import iris.tests first so that some things can be initialised before
+# importing anything else.
+import iris.tests as tests
+
+import itertools
+
+import numpy as np
+
+from iris._lazy_data import as_lazy_data
+from iris.tests.stock import ocean_sigma_z as stock_sample_osz
+import iris.util
+
+
+class Test_sample(tests.IrisTest):
+ def setUp(self):
+ self.cube = stock_sample_osz()
+ # Snapshot result, printed with ...
+ # >>> np.set_printoptions(linewidth=180,
+ # formatter={'float':lambda x:'{:-09.3f}'.format(x)})
+ # >>> print repr(coord.points)
+ self.basic_derived_result = np.array(
+ [[[[-0000.632, -0000.526, -0000.421, -0000.316],
+ [-0000.789, -0000.684, -0000.579, -0000.474],
+ [-0000.947, -0000.842, -0000.737, -0000.632]],
+ [[-0014.358, -0014.264, -0014.169, -0014.074],
+ [-0014.501, -0014.406, -0014.311, -0014.216],
+ [-0014.643, -0014.548, -0014.453, -0014.358]],
+ [[-0082.993, -0082.951, -0082.908, -0082.866],
+ [-0083.056, -0083.014, -0082.972, -0082.929],
+ [-0083.119, -0083.077, -0083.035, -0082.993]],
+ [[-0368.400, -0368.400, -0368.400, -0368.400],
+ [-0368.400, -0368.400, -0368.400, -0368.400],
+ [-0368.400, -0368.400, -0368.400, -0368.400]],
+ [[-1495.600, -1495.600, -1495.600, -1495.600],
+ [-1495.600, -1495.600, -1495.600, -1495.600],
+ [-1495.600, -1495.600, -1495.600, -1495.600]]],
+
+ [[[-0000.842, -0000.737, -0000.632, -0000.526],
+ [-0001.000, -0000.895, -0000.789, -0000.684],
+ [-0001.158, -0001.053, -0000.947, -0000.842]],
+ [[-0014.548, -0014.453, -0014.358, -0014.264],
+ [-0014.690, -0014.595, -0014.501, -0014.406],
+ [-0014.832, -0014.737, -0014.643, -0014.548]],
+ [[-0083.077, -0083.035, -0082.993, -0082.951],
+ [-0083.140, -0083.098, -0083.056, -0083.014],
+ [-0083.203, -0083.161, -0083.119, -0083.077]],
+ [[-0368.400, -0368.400, -0368.400, -0368.400],
+ [-0368.400, -0368.400, -0368.400, -0368.400],
+ [-0368.400, -0368.400, -0368.400, -0368.400]],
+ [[-1495.600, -1495.600, -1495.600, -1495.600],
+ [-1495.600, -1495.600, -1495.600, -1495.600],
+ [-1495.600, -1495.600, -1495.600, -1495.600]]]])
+
+ self.derived_coord_name = \
+ 'sea_surface_height_above_reference_ellipsoid'
+
+ def _check_result(self, cube, expected_result=None, **kwargs):
+ if expected_result is None:
+ expected_result = self.basic_derived_result
+ coord = cube.coord(self.derived_coord_name)
+ result = coord.points
+ self.assertArrayAllClose(result, expected_result, atol=0.005, **kwargs)
+
+ def test_basic(self):
+ self._check_result(self.cube)
+
+ def _lazy_testcube(self):
+ cube = self.cube
+ for dep_name in ('depth', 'layer_depth', 'ocean_sigma_z_coordinate'):
+ coord = cube.coord(dep_name)
+ coord.points = as_lazy_data(coord.points, coord.shape)
+ return cube
+
+ def test_nonlazy_cube_has_lazy_derived(self):
+ # Check same results when key coords are made lazy.
+ cube = self.cube
+ self.assertEqual(
+ cube.coord('depth').has_lazy_points(),
+ False)
+ self.assertEqual(
+ cube.coord(self.derived_coord_name).has_lazy_points(),
+ True)
+
+ def test_lazy_cube_same_result(self):
+ cube = self._lazy_testcube()
+ self.assertEqual(
+ cube.coord('depth').has_lazy_points(),
+ True)
+ self.assertEqual(
+ cube.coord(self.derived_coord_name).has_lazy_points(),
+ True)
+ self._check_result(cube)
+
+ def test_transpose(self):
+ # Check it works with all possible dimension orders.
+ for dims_list in itertools.permutations(range(self.cube.ndim)):
+ cube = self.cube.copy()
+ cube.transpose(dims_list)
+ expected = self.basic_derived_result.transpose(dims_list)
+ msg = 'Unexpected result when cube transposed by {}'
+ msg = msg.format(dims_list)
+ self._check_result(cube, expected,
+ err_msg=msg)
+
+ def test_lazy_transpose(self):
+ # Check lazy calc works with all possible dimension orders.
+ for dims_list in itertools.permutations(range(self.cube.ndim)):
+ cube = self._lazy_testcube().copy()
+ cube.transpose(dims_list)
+ expected = self.basic_derived_result.transpose(dims_list)
+ msg = 'Unexpected result when cube transposed by {}'
+ msg = msg.format(dims_list)
+ self._check_result(cube, expected,
+ err_msg=msg)
+
+ def test_extra_dims(self):
+ # Insert some extra cube dimensions + check it still works.
+ cube = self.cube
+ cube = iris.util.new_axis(cube)
+ cube = iris.util.new_axis(cube)
+ cube = iris.util.new_axis(cube)
+ # N.B. shape is now (1, 1, 1, t, z, y, x)
+ cube.transpose((0, 3, 1, 4, 5, 2, 6))
+ # N.B. shape is now (1, t, 1, z, y, 1, x)
+ # Should get same original result, as derived dims are the same.
+ self._check_result(cube)
+
+ def test_no_sigma(self):
+ # Check it still works when 'sigma' is removed.
+ # NOTE: the unit test for this does not cover all cases because it
+ # doesn't provide a time dimension.
+
+ # Set all sigma points to zero + snapshot the resulting derived points.
+ trial_cube = self.cube.copy()
+ trial_cube.coord('ocean_sigma_z_coordinate').points[:] = 0.0
+ expected = trial_cube.coord(self.derived_coord_name).points
+
+ # Remove sigma altogether + check the result is the same.
+ cube = self.cube
+ cube.remove_coord('ocean_sigma_z_coordinate')
+ self._check_result(cube, expected)
+
+ def test_no_eta(self):
+ # Check it still works when 'eta' is removed.
+ # NOTE: the unit test for this does not cover all cases because it
+ # doesn't provide a time dimension.
+
+ # Set all sigma points to zero + snapshot the resulting derived points.
+ trial_cube = self.cube.copy()
+ trial_cube.coord('sea_surface_height').points[:] = 0.0
+ expected = trial_cube.coord(self.derived_coord_name).points
+ # Check this has no variation between the two timepoints.
+ self.assertArrayAllClose(expected[0], expected[1])
+ # Take first time, as no sigma --> result *has* no time dimension.
+ expected = expected[0]
+
+ # Remove eta altogether + check the result is the same.
+ cube = self.cube
+ cube.remove_coord('sea_surface_height')
+ self._check_result(cube, expected)
+
+
+if __name__ == "__main__":
+ tests.main()
diff --git a/lib/iris/tests/results/analysis/regrid/linear_circular_grid.cml b/lib/iris/tests/results/analysis/regrid/linear_circular_grid.cml
index 09e68c9f4c..08c5a81e1f 100644
--- a/lib/iris/tests/results/analysis/regrid/linear_circular_grid.cml
+++ b/lib/iris/tests/results/analysis/regrid/linear_circular_grid.cml
@@ -7,10 +7,10 @@
-
+
diff --git a/lib/iris/tests/results/analysis/regrid/linear_masked_altitude.cml b/lib/iris/tests/results/analysis/regrid/linear_masked_altitude.cml
index f4e9914aa6..e9b0aab378 100644
--- a/lib/iris/tests/results/analysis/regrid/linear_masked_altitude.cml
+++ b/lib/iris/tests/results/analysis/regrid/linear_masked_altitude.cml
@@ -8,60 +8,60 @@
+ [[424.423, 398.043, --, --, --],
+ [368.688, 343.878, --, --, --],
+ [375.091, 347.861, --, --, --],
+ [446.161, 414.22, --, --, --]]]" shape="(2, 4, 5)" standard_name="altitude" units="Unit('m')" value_type="float32">
diff --git a/lib/iris/tests/results/analysis/regrid/linear_non_circular.cml b/lib/iris/tests/results/analysis/regrid/linear_non_circular.cml
index 32a80bcce5..ded7ffd106 100644
--- a/lib/iris/tests/results/analysis/regrid/linear_non_circular.cml
+++ b/lib/iris/tests/results/analysis/regrid/linear_non_circular.cml
@@ -7,10 +7,10 @@
-
+
diff --git a/lib/iris/tests/results/analysis/regrid/linear_partial_overlap.cml b/lib/iris/tests/results/analysis/regrid/linear_partial_overlap.cml
index e55e259bc8..147f5c69a1 100644
--- a/lib/iris/tests/results/analysis/regrid/linear_partial_overlap.cml
+++ b/lib/iris/tests/results/analysis/regrid/linear_partial_overlap.cml
@@ -6,54 +6,54 @@
-
+ [[--, --, 367.726, 355.63],
+ [--, --, 340.443, 407.574],
+ [--, --, 336.602, 419.093],
+ [--, --, 376.39, 341.021]]]" shape="(2, 4, 4)" standard_name="altitude" units="Unit('m')" value_type="float32">
diff --git a/lib/iris/tests/results/analysis/regrid/nearest_circular_grid.cml b/lib/iris/tests/results/analysis/regrid/nearest_circular_grid.cml
index 7f040a6ff7..d73639f66d 100644
--- a/lib/iris/tests/results/analysis/regrid/nearest_circular_grid.cml
+++ b/lib/iris/tests/results/analysis/regrid/nearest_circular_grid.cml
@@ -7,10 +7,10 @@
-
+
diff --git a/lib/iris/tests/results/analysis/regrid/nearest_masked_altitude.cml b/lib/iris/tests/results/analysis/regrid/nearest_masked_altitude.cml
index 83ceaea0a9..cfe47b3f8b 100644
--- a/lib/iris/tests/results/analysis/regrid/nearest_masked_altitude.cml
+++ b/lib/iris/tests/results/analysis/regrid/nearest_masked_altitude.cml
@@ -14,14 +14,14 @@
[[325.116, 337.95],
[289.584, 302.472],
- [nan, nan],
- [nan, nan],
+ [--, --],
+ [--, --],
[424.168, 436.85]],
[[325.116, 337.95],
[289.584, 302.472],
- [nan, nan],
- [nan, nan],
+ [--, --],
+ [--, --],
[424.168, 436.85]],
[[440.66, 453.316],
@@ -39,14 +39,14 @@
[[337.95, 357.202],
[302.472, 321.806],
- [nan, nan],
- [nan, nan],
+ [--, --],
+ [--, --],
[436.85, 455.873]],
[[337.95, 357.202],
[302.472, 321.806],
- [nan, nan],
- [nan, nan],
+ [--, --],
+ [--, --],
[436.85, 455.873]],
[[453.316, 472.302],
@@ -54,13 +54,13 @@
[315.637, 334.94],
[315.637, 334.94],
[272.726, 292.128]]]]" id="9041e969" points="[[[418.698, 379.592, 202.986, 202.986, 333.61],
- [329.929, 294.417, nan, nan, 428.923],
- [329.929, 294.417, nan, nan, 428.923],
+ [329.929, 294.417, --, --, 428.923],
+ [329.929, 294.417, --, --, 428.923],
[445.406, 398.973, 307.595, 307.595, 264.642]],
[[434.57, 395.539, 219.272, 219.272, 349.646],
- [345.971, 310.528, nan, nan, 444.776],
- [345.971, 310.528, nan, nan, 444.776],
+ [345.971, 310.528, --, --, 444.776],
+ [345.971, 310.528, --, --, 444.776],
[461.227, 414.883, 323.68, 323.68, 280.81]]]" shape="(2, 4, 5)" standard_name="altitude" units="Unit('m')" value_type="float32">
diff --git a/lib/iris/tests/results/analysis/regrid/nearest_non_circular.cml b/lib/iris/tests/results/analysis/regrid/nearest_non_circular.cml
index 1d8e87c7e4..5b842c8811 100644
--- a/lib/iris/tests/results/analysis/regrid/nearest_non_circular.cml
+++ b/lib/iris/tests/results/analysis/regrid/nearest_non_circular.cml
@@ -7,10 +7,10 @@
-
+
diff --git a/lib/iris/tests/results/analysis/regrid/nearest_partial_overlap.cml b/lib/iris/tests/results/analysis/regrid/nearest_partial_overlap.cml
index 351e86a2ac..3763379f50 100644
--- a/lib/iris/tests/results/analysis/regrid/nearest_partial_overlap.cml
+++ b/lib/iris/tests/results/analysis/regrid/nearest_partial_overlap.cml
@@ -6,54 +6,54 @@
-
+ [[--, --, 395.539, 349.646],
+ [--, --, 310.528, 444.776],
+ [--, --, 310.528, 444.776],
+ [--, --, 414.883, 280.81]]]" shape="(2, 4, 4)" standard_name="altitude" units="Unit('m')" value_type="float32">
diff --git a/lib/iris/tests/results/analysis/regrid/no_overlap.cml b/lib/iris/tests/results/analysis/regrid/no_overlap.cml
index 5168c88ca5..a4503e2532 100644
--- a/lib/iris/tests/results/analysis/regrid/no_overlap.cml
+++ b/lib/iris/tests/results/analysis/regrid/no_overlap.cml
@@ -6,54 +6,54 @@
-
+ [[--, --, --, --],
+ [--, --, --, --],
+ [--, --, --, --],
+ [--, --, --, --]]]" shape="(2, 4, 4)" standard_name="altitude" units="Unit('m')" value_type="float32">
diff --git a/lib/iris/tests/stock.py b/lib/iris/tests/stock.py
index df63472690..25df8eda33 100644
--- a/lib/iris/tests/stock.py
+++ b/lib/iris/tests/stock.py
@@ -1,4 +1,4 @@
-# (C) British Crown Copyright 2010 - 2016, Met Office
+# (C) British Crown Copyright 2010 - 2017, Met Office
#
# This file is part of Iris.
#
@@ -28,10 +28,12 @@
import numpy as np
import numpy.ma as ma
+from cf_units import Unit
from iris.cube import Cube
import iris.aux_factory
import iris.coords
import iris.coords as icoords
+from iris.coords import DimCoord, AuxCoord
import iris.tests as tests
from iris.coord_systems import GeogCS, RotatedGeogCS
@@ -482,7 +484,6 @@ def realistic_4d():
level_height_bnds, model_level_pts, sigma_pts, sigma_bnds, time_pts, \
_source_pts, forecast_period_pts, data, orography = arrays
-
ll_cs = RotatedGeogCS(37.5, 177.5, ellipsoid=GeogCS(6371229.0))
lat = icoords.DimCoord(lat_pts, standard_name='grid_latitude', units='degrees',
@@ -558,3 +559,58 @@ def global_grib2():
path = tests.get_data_path(('GRIB', 'global_t', 'global.grib2'))
cube = iris.load_cube(path)
return cube
+
+
+def ocean_sigma_z():
+ """
+ Return a sample cube with an
+ :class:`iris.aux_factory.OceanSigmaZFactory` vertical coordinate.
+
+ This is a fairly small cube with real coordinate arrays. The coordinate
+ values are derived from the sample data linked at
+ https://github.com/SciTools/iris/pull/509#issuecomment-23565381.
+
+ """
+ co_time = DimCoord([0.0, 1.0], standard_name='time', units='')
+ co_lats = DimCoord([-58.1, -52.7, -46.9],
+ standard_name='latitude', units=Unit('degrees'))
+ co_lons = DimCoord([65.1, 72.9, 83.7, 96.5],
+ standard_name='longitude', units=Unit('degrees'))
+ co_ssh = AuxCoord([[[-0.63157895, -0.52631579, -0.42105263, -0.31578947],
+ [-0.78947368, -0.68421053, -0.57894737, -0.47368421],
+ [-0.94736842, -0.84210526, -0.73684211, -0.63157895]],
+ [[-0.84210526, -0.73684211, -0.63157895, -0.52631579],
+ [-1.00000000, -0.89473684, -0.78947368, -0.68421053],
+ [-1.15789474, -1.05263158, -0.94736842, -0.84210526]]],
+ standard_name=u'sea_surface_height', units=Unit('m'))
+
+ co_sigma = AuxCoord([0., -0.1, -0.6, -1., -1.],
+ standard_name=u'ocean_sigma_z_coordinate',
+ units=Unit('1'),
+ attributes={'positive': 'up'})
+
+ co_zlay = AuxCoord([-137.2, -137.3, -137.4, -368.4, -1495.6],
+ long_name='layer_depth', units=Unit('m'))
+ co_depth = AuxCoord([[1625.7, 3921.2, 4106.4, 5243.5],
+ [3615.4, 4942.6, 3883.6, 4823.1],
+ [3263.2, 2816.3, 2741.8, 3883.6]],
+ standard_name=u'depth', units=Unit('m'))
+ co_depthc = DimCoord(137.9, long_name='depth_c', units=Unit('m'))
+ co_nsigma = DimCoord(3, long_name='nsigma')
+
+ cube = Cube(np.zeros((2, 5, 3, 4)))
+ cube.add_dim_coord(co_time, 0)
+ cube.add_dim_coord(co_lats, 2)
+ cube.add_dim_coord(co_lons, 3)
+ cube.add_aux_coord(co_zlay, 1)
+ cube.add_aux_coord(co_sigma, 1)
+ cube.add_aux_coord(co_ssh, (0, 2, 3))
+ cube.add_aux_coord(co_depth, (2, 3))
+ cube.add_aux_coord(co_depthc)
+ cube.add_aux_coord(co_nsigma)
+
+ fact = iris.aux_factory.OceanSigmaZFactory(
+ depth=co_depth, eta=co_ssh, depth_c=co_depthc, zlev=co_zlay,
+ sigma=co_sigma, nsigma=co_nsigma)
+ cube.add_aux_factory(fact)
+ return cube
diff --git a/lib/iris/tests/test_hybrid.py b/lib/iris/tests/test_hybrid.py
index a3240d8578..6c83df9502 100644
--- a/lib/iris/tests/test_hybrid.py
+++ b/lib/iris/tests/test_hybrid.py
@@ -151,14 +151,15 @@ def test_bounded_orography(self):
# Make sure altitude still works OK if orography was messed
# with *after* altitude was created.
- altitude = self.cube.coord('altitude')
orog.bounds = np.zeros(orog.shape + (4,))
- self.assertIsInstance(altitude.bounds, np.ndarray)
- # Make sure altitude.bounds now raises an error.
- exp_emsg = 'operands could not be broadcast together'
- with self.assertRaisesRegexp(ValueError, exp_emsg):
- self.cube.coord('altitude')
+ # Check that altitude derivation now issues a warning.
+ msg = 'Orography.* bounds.* being disregarded'
+ with warnings.catch_warnings():
+ # Cause all warnings to raise Exceptions
+ warnings.simplefilter("error")
+ with self.assertRaisesRegexp(UserWarning, msg):
+ self.cube.coord('altitude')
@tests.skip_data
@@ -230,15 +231,15 @@ def test_bounded_surface_pressure(self):
# Make sure pressure still works OK if surface pressure was messed
# with *after* pressure was created.
- pressure = self.cube.coord('air_pressure')
surface_pressure.bounds = np.zeros(surface_pressure.shape + (4,))
- self.assertIsInstance(pressure.bounds, np.ndarray)
-
- # Make sure pressure.bounds now raises an error.
- exp_emsg = 'operands could not be broadcast together'
- with self.assertRaisesRegexp(ValueError, exp_emsg):
- self.cube.coord('air_pressure')
+ # Check that air_pressure derivation now issues a warning.
+ msg = 'Surface pressure.* bounds.* being disregarded'
+ with warnings.catch_warnings():
+ # Cause all warnings to raise Exceptions
+ warnings.simplefilter("error")
+ with self.assertRaisesRegexp(UserWarning, msg):
+ self.cube.coord('air_pressure')
if __name__ == "__main__":
tests.main()