Skip to content

Commit

Permalink
ENH: Working zonal mean linear regridding for circular sources or wit…
Browse files Browse the repository at this point in the history
…h use of extrapolation (#3085)

* Working zonal mean from circular/extrapolated source to target

* MAINT: Refactor of zonal mean testing

* MAINT: Documentation changes from review
  • Loading branch information
Carwyn Pelley authored and DPeterK committed Jul 13, 2018
1 parent 22dd241 commit fcfea63
Show file tree
Hide file tree
Showing 2 changed files with 114 additions and 2 deletions.
5 changes: 3 additions & 2 deletions lib/iris/analysis/_regrid.py
Original file line number Diff line number Diff line change
Expand Up @@ -231,8 +231,9 @@ def _regrid(src_data, x_dim, y_dim,
data = np.empty(shape, dtype=dtype)

# The interpolation class requires monotonically increasing
# coordinates, so flip the coordinate(s) and data if the aren't.
reverse_x = src_x_coord.points[0] > src_x_coord.points[1]
# coordinates, so flip the coordinate(s) and data if they aren't.
reverse_x = (src_x_coord.points[0] > src_x_coord.points[1] if
src_x_coord.points.size > 1 else False)
reverse_y = src_y_coord.points[0] > src_y_coord.points[1]
flip_index = [slice(None)] * src_data.ndim
if reverse_x:
Expand Down
111 changes: 111 additions & 0 deletions lib/iris/tests/integration/test_regridding.py
Original file line number Diff line number Diff line change
Expand Up @@ -109,5 +109,116 @@ def test_nearest(self):
self.assertArrayShapeStats(res, (1, 6, 3, 4), 315.890808, 11.000724)


class TestZonalMean_global(tests.IrisTest):
def setUp(self):
np.random.seed(0)
self.src = iris.cube.Cube(np.random.random_integers(0, 10, (140, 1)))
s_crs = iris.coord_systems.GeogCS(6371229.0)
sy_coord = iris.coords.DimCoord(
np.linspace(-90, 90, 140), standard_name='latitude',
units='degrees', coord_system=s_crs)
sx_coord = iris.coords.DimCoord(
-180, bounds=[-180, 180], standard_name='longitude',
units='degrees', circular=True, coord_system=s_crs)
self.src.add_dim_coord(sy_coord, 0)
self.src.add_dim_coord(sx_coord, 1)

def test_linear_same_crs_global(self):
# Regrid the zonal mean onto an identical coordinate system target, but
# on a different set of longitudes - which should result in no change.
points = [-150, -90, -30, 30, 90, 150]
bounds = [[-180, -120], [-120, -60], [-60, 0], [0, 60], [60, 120],
[120, 180]]
sx_coord = self.src.coord(axis='x')
sy_coord = self.src.coord(axis='y')
x_coord = sx_coord.copy(points, bounds=bounds)
grid = iris.cube.Cube(
np.zeros([sy_coord.points.size, x_coord.points.size]))
grid.add_dim_coord(sy_coord, 0)
grid.add_dim_coord(x_coord, 1)

res = self.src.regrid(grid, iris.analysis.Linear())

# Ensure data remains unchanged.
# (the same along each column)
self.assertTrue(
np.array([(res.data[:, 0]-res.data[:, i]).max() for i in
range(1, res.shape[1])]).max() < 1e-10)
self.assertArrayAlmostEqual(res.data[:, 0], self.src.data.reshape(-1))


class TestZonalMean_regional(TestZonalMean_global, tests.IrisTest):
def setUp(self):
super(TestZonalMean_regional, self).setUp()

# Define a target grid and a target result (what we expect the
# regridder to return).
sx_coord = self.src.coord(axis='x')
sy_coord = self.src.coord(axis='y')
grid_crs = iris.coord_systems.RotatedGeogCS(
37.5, 177.5, ellipsoid=iris.coord_systems.GeogCS(6371229.0))
grid_x = sx_coord.copy(np.linspace(350, 370, 100))
grid_x.circular = False
grid_x.coord_system = grid_crs
grid_y = sy_coord.copy(np.linspace(-10, 10, 100))
grid_y.coord_system = grid_crs
grid = iris.cube.Cube(
np.zeros([grid_y.points.size, grid_x.points.size]))
grid.add_dim_coord(grid_y, 0)
grid.add_dim_coord(grid_x, 1)

# The target result is derived by regridding a multi-column version of
# the source to the target (i.e. turning a zonal mean regrid into a
# conventional regrid).
self.tar = self.zonal_mean_as_multi_column(self.src).regrid(
grid, iris.analysis.Linear())
self.grid = grid

def zonal_mean_as_multi_column(self, src_cube):
# Munge the source (duplicate source latitudes) so that we can
# utilise linear regridding as a conventional problem (that is, to
# duplicate columns so that it is no longer a zonal mean problem).
src_cube2 = src_cube.copy()
src_cube2.coord(axis='x').points = -90
src_cube2.coord(axis='x').bounds = [-180, 0]
src_cube.coord(axis='x').points = 90
src_cube.coord(axis='x').bounds = [0, 180]
src_cubes = iris.cube.CubeList([src_cube, src_cube2])
return src_cubes.concatenate_cube()

def test_linear_rotated_regional(self):
# Ensure that zonal mean source data is linearly interpolated onto a
# high resolution target.
regridder = iris.analysis.Linear()
res = self.src.regrid(self.grid, regridder)
self.assertArrayAlmostEqual(res.data, self.tar.data)

def test_linear_rotated_regional_no_extrapolation(self):
# Capture the case where our source remains circular but we don't use
# extrapolation.
regridder = iris.analysis.Linear(extrapolation_mode='nan')
res = self.src.regrid(self.grid, regridder)
self.assertArrayAlmostEqual(res.data, self.tar.data)

def test_linear_rotated_regional_not_circular(self):
# Capture the case where our source is not circular but we utilise
# extrapolation.
regridder = iris.analysis.Linear()
self.src.coord(axis='x').circular = False
res = self.src.regrid(self.grid, regridder)
self.assertArrayAlmostEqual(res.data, self.tar.data)

def test_linear_rotated_regional_no_extrapolation_not_circular(self):
# Confirm how zonal mean actually works in so far as, that
# extrapolation and circular source handling is the means by which
# these usecases are supported.
# In the case where the source is neither using extrapolation and is
# not circular, then 'nan' values will result (as we would expect).
regridder = iris.analysis.Linear(extrapolation_mode='nan')
self.src.coord(axis='x').circular = False
res = self.src.regrid(self.grid, regridder)
self.assertTrue(np.isnan(res.data).all())


if __name__ == "__main__":
tests.main()

0 comments on commit fcfea63

Please sign in to comment.