diff --git a/docs/iris/src/whatsnew/contributions_1.11/newfeature_2015-Aug-06_lambert_azimuthal_equal_area.txt b/docs/iris/src/whatsnew/contributions_1.11/newfeature_2015-Aug-06_lambert_azimuthal_equal_area.txt
new file mode 100644
index 0000000000..700b599b65
--- /dev/null
+++ b/docs/iris/src/whatsnew/contributions_1.11/newfeature_2015-Aug-06_lambert_azimuthal_equal_area.txt
@@ -0,0 +1 @@
+* The coordinate system :class:`iris.coord_systems.LambertAzimuthalEqualArea` has been added with NetCDF saving support.
diff --git a/lib/iris/coord_systems.py b/lib/iris/coord_systems.py
index 6808f717ce..66a93c5ad5 100644
--- a/lib/iris/coord_systems.py
+++ b/lib/iris/coord_systems.py
@@ -494,14 +494,14 @@ def __init__(self, latitude_of_projection_origin,
* longitude_of_projection_origin:
True longitude of planar origin in degrees.
+ Kwargs:
+
* false_easting
X offset from planar origin in metres. Defaults to 0.
* false_northing
Y offset from planar origin in metres. Defaults to 0.
- Kwargs:
-
* ellipsoid
:class:`GeogCS` defining the ellipsoid.
@@ -577,14 +577,14 @@ def __init__(self, latitude_of_projection_origin,
Altitude of satellite in metres above the surface of the
ellipsoid.
+ Kwargs:
+
* false_easting
X offset from planar origin in metres. Defaults to 0.
* false_northing
Y offset from planar origin in metres. Defaults to 0.
- Kwargs:
-
* ellipsoid
:class:`GeogCS` defining the ellipsoid.
@@ -667,14 +667,14 @@ def __init__(self, central_lat, central_lon,
* central_lon
The central longitude, which aligns with the y axis.
+ Kwargs:
+
* false_easting
X offset from planar origin in metres. Defaults to 0.
* false_northing
Y offset from planar origin in metres. Defaults to 0.
- Kwargs:
-
* true_scale_lat
Latitude of true scale.
@@ -739,7 +739,7 @@ def __init__(self, central_lat=39.0, central_lon=-96.0,
"""
Constructs a LambertConformal coord system.
- Args:
+ Kwargs:
* central_lat
The latitude of "unitary scale".
@@ -753,8 +753,6 @@ def __init__(self, central_lat=39.0, central_lon=-96.0,
* false_northing
Y offset from planar origin in metres.
- Kwargs:
-
* secant_latitudes
Latitudes of secant intersection.
@@ -837,13 +835,9 @@ def __init__(self, longitude_of_projection_origin=0, ellipsoid=None):
"""
Constructs a Mercator coord system.
- Args:
-
+ Kwargs:
* longitude_of_projection_origin
True longitude of planar origin in degrees.
-
- Kwargs:
-
* ellipsoid
:class:`GeogCS` defining the ellipsoid.
@@ -870,3 +864,73 @@ def as_cartopy_crs(self):
def as_cartopy_projection(self):
return self.as_cartopy_crs()
+
+
+class LambertAzimuthalEqualArea(CoordSystem):
+ """
+ A coordinate system in the Lambert Azimuthal Equal Area projection.
+
+ """
+
+ grid_mapping_name = "lambert_azimuthal_equal_area"
+
+ def __init__(self, latitude_of_projection_origin=0.0,
+ longitude_of_projection_origin=0.0,
+ false_easting=0.0, false_northing=0.0,
+ ellipsoid=None):
+ """
+ Constructs a Lambert Azimuthal Equal Area coord system.
+
+ Kwargs:
+
+ * latitude_of_projection_origin
+ True latitude of planar origin in degrees. Defaults to 0.
+
+ * longitude_of_projection_origin
+ True longitude of planar origin in degrees. Defaults to 0.
+
+ * false_easting
+ X offset from planar origin in metres. Defaults to 0.
+
+ * false_northing
+ Y offset from planar origin in metres. Defaults to 0.
+
+ * ellipsoid
+ :class:`GeogCS` defining the ellipsoid.
+
+ """
+ #: True latitude of planar origin in degrees.
+ self.latitude_of_projection_origin = latitude_of_projection_origin
+ #: True longitude of planar origin in degrees.
+ self.longitude_of_projection_origin = longitude_of_projection_origin
+ #: X offset from planar origin in metres.
+ self.false_easting = false_easting
+ #: Y offset from planar origin in metres.
+ self.false_northing = false_northing
+ #: Ellipsoid definition.
+ self.ellipsoid = ellipsoid
+
+ def __repr__(self):
+ return "LambertAzimuthalEqualArea(latitude_of_projection_origin={!r}, "\
+ "longitude_of_projection_origin={!r}, false_easting={!r}, "\
+ "false_northing={!r}, ellipsoid={!r})".format(
+ self.latitude_of_projection_origin,
+ self.longitude_of_projection_origin,
+ self.false_easting,
+ self.false_northing,
+ self.ellipsoid)
+
+ def as_cartopy_crs(self):
+ if self.ellipsoid is not None:
+ globe = self.ellipsoid.as_cartopy_globe()
+ else:
+ globe = ccrs.Globe()
+ return ccrs.LambertAzimuthalEqualArea(
+ central_longitude=self.longitude_of_projection_origin,
+ central_latitude=self.latitude_of_projection_origin,
+ false_easting=self.false_easting,
+ false_northing=self.false_northing,
+ globe=globe)
+
+ def as_cartopy_projection(self):
+ return self.as_cartopy_crs()
diff --git a/lib/iris/fileformats/_pyke_rules/fc_rules_cf.krb b/lib/iris/fileformats/_pyke_rules/fc_rules_cf.krb
index 28da8cd836..cef4d35de7 100644
--- a/lib/iris/fileformats/_pyke_rules/fc_rules_cf.krb
+++ b/lib/iris/fileformats/_pyke_rules/fc_rules_cf.krb
@@ -157,6 +157,26 @@ fc_provides_grid_mapping_lambert_conformal
facts_cf.provides(coordinate_system, lambert_conformal)
python engine.rule_triggered.add(rule.name)
+#
+# Context:
+# This rule will trigger iff a grid_mapping() case specific fact
+# has been asserted that refers to a lambert azimuthal equal area.
+#
+# Purpose:
+# Creates the lambert azimuthal equal area coordinate system.
+#
+fc_provides_grid_mapping_lambert_azimuthal_equal_area
+ foreach
+ facts_cf.grid_mapping($grid_mapping)
+ check is_grid_mapping(engine, $grid_mapping, CF_GRID_MAPPING_LAMBERT_AZIMUTHAL)
+ assert
+ python cf_grid_var = engine.cf_var.cf_group.grid_mappings[$grid_mapping]
+ python coordinate_system = build_lambert_azimuthal_equal_area_coordinate_system(engine, cf_grid_var)
+ python engine.provides['coordinate_system'] = coordinate_system
+ facts_cf.provides(coordinate_system, lambert_azimuthal_equal_area)
+ python engine.rule_triggered.add(rule.name)
+
+
#
# Context:
# This rule will trigger iff a coordinate() case specific fact
@@ -715,6 +735,45 @@ fc_build_coordinate_projection_y_stereographic
python engine.rule_triggered.add(rule.name)
+#
+# Context:
+# This rule will trigger iff a projection_x_coordinate coordinate exists and
+# a lambert azimuthal equal area coordinate system exists.
+#
+# Purpose:
+# Add the projection_x_coordinate coordinate into the cube.
+#
+fc_build_coordinate_projection_x_lambert_azimuthal_equal_area
+ foreach
+ facts_cf.provides(coordinate, projection_x_coordinate, $coordinate)
+ facts_cf.provides(coordinate_system, lambert_azimuthal_equal_area)
+ assert
+ python cf_coord_var = engine.cf_var.cf_group.coordinates[$coordinate]
+ python build_dimension_coordinate(engine, cf_coord_var,
+ coord_name=CF_VALUE_STD_NAME_PROJ_X,
+ coord_system=engine.provides['coordinate_system'])
+ python engine.rule_triggered.add(rule.name)
+
+
+#
+# Context:
+# This rule will trigger iff a projection_y_coordinate coordinate exists and
+# a lambert azimuthal equal area coordinate system exists.
+#
+# Purpose:
+# Add the projection_y_coordinate coordinate into the cube.
+#
+fc_build_coordinate_projection_y_lambert_azimuthal_equal_area
+ foreach
+ facts_cf.provides(coordinate, projection_y_coordinate, $coordinate)
+ facts_cf.provides(coordinate_system, lambert_azimuthal_equal_area)
+ assert
+ python cf_coord_var = engine.cf_var.cf_group.coordinates[$coordinate]
+ python build_dimension_coordinate(engine, cf_coord_var,
+ coord_name=CF_VALUE_STD_NAME_PROJ_Y,
+ coord_system=engine.provides['coordinate_system'])
+ python engine.rule_triggered.add(rule.name)
+
#
# Context:
# This rule will trigger iff a CF time coordinate exists.
@@ -1298,6 +1357,37 @@ fc_extras
return cs
+ ################################################################################
+ def build_lambert_azimuthal_equal_area_coordinate_system(engine, cf_grid_var):
+ """
+ Create a lambert azimuthal equal area coordinate system from the CF-netCDF
+ grid mapping variable.
+
+ """
+ major, minor, inverse_flattening = _get_ellipsoid(cf_grid_var)
+
+ latitude_of_projection_origin = getattr(
+ cf_grid_var, CF_ATTR_GRID_LAT_OF_PROJ_ORIGIN, None)
+ longitude_of_projection_origin = getattr(
+ cf_grid_var, CF_ATTR_GRID_LON_OF_PROJ_ORIGIN, None)
+ false_easting = getattr(
+ cf_grid_var, CF_ATTR_GRID_FALSE_EASTING, None)
+ false_northing = getattr(
+ cf_grid_var, CF_ATTR_GRID_FALSE_NORTHING, None)
+
+ ellipsoid = None
+ if major is not None or minor is not None or \
+ inverse_flattening is not None:
+ ellipsoid = iris.coord_systems.GeogCS(major, minor,
+ inverse_flattening)
+
+ cs = iris.coord_systems.LambertAzimuthalEqualArea(
+ latitude_of_projection_origin, longitude_of_projection_origin,
+ false_easting, false_northing, ellipsoid)
+
+ return cs
+
+
################################################################################
def get_attr_units(cf_var, attributes):
attr_units = getattr(cf_var, CF_ATTR_UNITS, cf_units._UNIT_DIMENSIONLESS)
diff --git a/lib/iris/fileformats/netcdf.py b/lib/iris/fileformats/netcdf.py
index 62279ba20d..65af4eea85 100644
--- a/lib/iris/fileformats/netcdf.py
+++ b/lib/iris/fileformats/netcdf.py
@@ -1300,7 +1300,8 @@ def _get_dim_names(self, cube):
dimension_names.append(dim_name)
return dimension_names
- def _cf_coord_identity(self, coord):
+ @staticmethod
+ def _cf_coord_identity(coord):
"""
Determine a suitable units from a given coordinate.
@@ -1327,7 +1328,7 @@ def _cf_coord_identity(self, coord):
elif coord.standard_name == "longitude":
units = 'degrees_east'
- return (coord.standard_name, coord.long_name, units)
+ return coord.standard_name, coord.long_name, units
def _ensure_valid_dtype(self, values, src_name, src_object):
# NetCDF3 does not support int64 or unsigned ints, so we check
@@ -1779,6 +1780,18 @@ def add_ellipsoid(ellipsoid):
elif isinstance(cs, iris.coord_systems.OSGB):
warnings.warn('OSGB coordinate system not yet handled')
+ # lambert azimuthal equal area
+ elif isinstance(cs,
+ iris.coord_systems.LambertAzimuthalEqualArea):
+ if cs.ellipsoid:
+ add_ellipsoid(cs.ellipsoid)
+ cf_var_grid.longitude_of_projection_origin = (
+ cs.longitude_of_projection_origin)
+ cf_var_grid.latitude_of_projection_origin = (
+ cs.latitude_of_projection_origin)
+ cf_var_grid.false_easting = cs.false_easting
+ cf_var_grid.false_northing = cs.false_northing
+
# other
else:
warnings.warn('Unable to represent the horizontal '
diff --git a/lib/iris/tests/integration/test_netcdf.py b/lib/iris/tests/integration/test_netcdf.py
index ce58043686..16b73745d1 100644
--- a/lib/iris/tests/integration/test_netcdf.py
+++ b/lib/iris/tests/integration/test_netcdf.py
@@ -301,5 +301,13 @@ def test_unknown_method(self):
shutil.rmtree(temp_dirpath)
+class TestCoordSystem(tests.IrisTest):
+ def test_load_laea_grid(self):
+ cube = iris.load_cube(
+ tests.get_data_path(('NetCDF', 'lambert_azimuthal_equal_area',
+ 'euro_air_temp.nc')))
+ self.assertCML(cube, ('netcdf', 'netcdf_laea.cml'))
+
+
if __name__ == "__main__":
tests.main()
diff --git a/lib/iris/tests/results/netcdf/netcdf_laea.cml b/lib/iris/tests/results/netcdf/netcdf_laea.cml
new file mode 100644
index 0000000000..3dc66dc0e3
--- /dev/null
+++ b/lib/iris/tests/results/netcdf/netcdf_laea.cml
@@ -0,0 +1,72 @@
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
diff --git a/lib/iris/tests/unit/coord_systems/test_LambertAzimuthalEqualArea.py b/lib/iris/tests/unit/coord_systems/test_LambertAzimuthalEqualArea.py
new file mode 100644
index 0000000000..fe16156607
--- /dev/null
+++ b/lib/iris/tests/unit/coord_systems/test_LambertAzimuthalEqualArea.py
@@ -0,0 +1,94 @@
+# (C) British Crown Copyright 2015 - 2016, 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 .
+"""
+Unit tests for the :class:`iris.coord_systems.LambertAzimuthalEqualArea` 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 cartopy.crs as ccrs
+from iris.coord_systems import GeogCS, LambertAzimuthalEqualArea
+
+
+class Test_as_cartopy_crs(tests.IrisTest):
+ def setUp(self):
+ self.latitude_of_projection_origin = 90.0
+ self.longitude_of_projection_origin = 0.0
+ self.semi_major_axis = 6377563.396
+ self.semi_minor_axis = 6356256.909
+ self.false_easting = 0.0
+ self.false_northing = 0.0
+ self.ellipsoid = GeogCS(self.semi_major_axis, self.semi_minor_axis)
+ self.laea_cs = LambertAzimuthalEqualArea(
+ self.latitude_of_projection_origin,
+ self.longitude_of_projection_origin,
+ self.false_easting,
+ self.false_northing,
+ ellipsoid=self.ellipsoid)
+
+ def test_crs_creation(self):
+ res = self.laea_cs.as_cartopy_crs()
+ globe = ccrs.Globe(semimajor_axis=self.semi_major_axis,
+ semiminor_axis=self.semi_minor_axis,
+ ellipse=None)
+ expected = ccrs.LambertAzimuthalEqualArea(
+ self.longitude_of_projection_origin,
+ self.latitude_of_projection_origin,
+ self.false_easting,
+ self.false_northing,
+ globe=globe)
+ self.assertEqual(res, expected)
+
+
+class Test_as_cartopy_projection(tests.IrisTest):
+ def setUp(self):
+ self.latitude_of_projection_origin = 0.0
+ self.longitude_of_projection_origin = 0.0
+ self.semi_major_axis = 6377563.396
+ self.semi_minor_axis = 6356256.909
+ self.false_easting = 0.0
+ self.false_northing = 0.0
+ self.ellipsoid = GeogCS(self.semi_major_axis, self.semi_minor_axis)
+ self.laea_cs = LambertAzimuthalEqualArea(
+ self.latitude_of_projection_origin,
+ self.longitude_of_projection_origin,
+ self.false_easting,
+ self.false_northing,
+ ellipsoid=self.ellipsoid)
+
+ def test_projection_creation(self):
+ res = self.laea_cs.as_cartopy_projection()
+ globe = ccrs.Globe(semimajor_axis=self.semi_major_axis,
+ semiminor_axis=self.semi_minor_axis,
+ ellipse=None)
+ expected = ccrs.LambertAzimuthalEqualArea(
+ self.latitude_of_projection_origin,
+ self.longitude_of_projection_origin,
+ self.false_easting,
+ self.false_northing,
+ globe=globe)
+ self.assertEqual(res, expected)
+
+
+if __name__ == '__main__':
+ tests.main()
diff --git a/lib/iris/tests/unit/fileformats/netcdf/test_Saver.py b/lib/iris/tests/unit/fileformats/netcdf/test_Saver.py
index c58871e785..a06ee870f0 100644
--- a/lib/iris/tests/unit/fileformats/netcdf/test_Saver.py
+++ b/lib/iris/tests/unit/fileformats/netcdf/test_Saver.py
@@ -29,8 +29,8 @@
import iris
from iris.coord_systems import (GeogCS, TransverseMercator, RotatedGeogCS,
- LambertConformal, Mercator, Stereographic)
-
+ LambertConformal, Mercator, Stereographic,
+ LambertAzimuthalEqualArea)
from iris.coords import DimCoord
from iris.cube import Cube
from iris.fileformats.netcdf import Saver
@@ -458,6 +458,38 @@ def test_valid_range_and_valid_min_valid_max_provided(self):
saver.check_attribute_compliance(self.container, self.data)
+class Test__cf_coord_identity(tests.IrisTest):
+ def check_call(self, coord_name, coord_system, units, expected_units):
+ coord = iris.coords.DimCoord([30, 45], coord_name, units=units,
+ coord_system=coord_system)
+ result = Saver._cf_coord_identity(coord)
+ self.assertEqual(result, (coord.standard_name, coord.long_name,
+ expected_units))
+
+ def test_geogcs_latitude(self):
+ crs = iris.coord_systems.GeogCS(60, 0)
+ self.check_call('latitude', coord_system=crs, units='degrees',
+ expected_units='degrees_north')
+
+ def test_geogcs_longitude(self):
+ crs = iris.coord_systems.GeogCS(60, 0)
+ self.check_call('longitude', coord_system=crs, units='degrees',
+ expected_units='degrees_east')
+
+ def test_no_coord_system_latitude(self):
+ self.check_call('latitude', coord_system=None, units='degrees',
+ expected_units='degrees_north')
+
+ def test_no_coord_system_longitude(self):
+ self.check_call('longitude', coord_system=None, units='degrees',
+ expected_units='degrees_east')
+
+ def test_passthrough_units(self):
+ crs = iris.coord_systems.LambertConformal(0, 20)
+ self.check_call('projection_x_coordinate', coord_system=crs,
+ units='km', expected_units='km')
+
+
class Test__create_cf_grid_mapping(tests.IrisTest):
def _cube_with_cs(self, coord_system):
"""Return a simple 2D cube that uses the given coordinate system."""
@@ -556,6 +588,24 @@ def test_lambert_conformal(self):
}
self._test(coord_system, expected)
+ def test_laea_cs(self):
+ coord_system = LambertAzimuthalEqualArea(
+ latitude_of_projection_origin=52,
+ longitude_of_projection_origin=10,
+ false_easting=100,
+ false_northing=200,
+ ellipsoid=GeogCS(6377563.396, 6356256.909))
+ expected = {'grid_mapping_name': 'lambert_azimuthal_equal_area',
+ 'latitude_of_projection_origin': 52,
+ 'longitude_of_projection_origin': 10,
+ 'false_easting': 100,
+ 'false_northing': 200,
+ 'semi_major_axis': 6377563.396,
+ 'semi_minor_axis': 6356256.909,
+ 'longitude_of_prime_meridian': 0,
+ }
+ self._test(coord_system, expected)
+
if __name__ == "__main__":
tests.main()