Fix multidimensional coordinates (#1768)

* Add test case for 1763.

* Add whitespace and wrap long lines to pass flake8 tests.

* Fix test case - time should not be a coordinate!

* When encoding coordinates, change the criteria so that the dimensions of the coordinate must be a subset (or equal to) the dimensions of the variable.

* Add note about bug fix (#1763) to the docs.

Author: neishm

doc/whats-new.rst

@@ -56,6 +56,9 @@ Enhancements
 
 Bug fixes
 ~~~~~~~~~
+- Fixed encoding of multi-dimensional coordinates in
+  :py:meth:`~Dataset.to_netcdf` (:issue:`1763`).
+  By `Mike Neish <https://github.com/neishm>`_.
 
 - Bug fix in open_dataset(engine='pydap') (:issue:`1775`)
   By `Keisuke Fujii <https://github.com/fujiisoup>`_.

xarray/conventions.py

@@ -663,7 +663,7 @@ def _encode_coordinates(variables, attributes, non_dim_coord_names):
         target_dims = variables[coord_name].dims
         for k, v in variables.items():
             if (k not in non_dim_coord_names and k not in v.dims and
-                    any(d in target_dims for d in v.dims)):
+                    set(target_dims) <= set(v.dims)):
                 variable_coordinates[k].add(coord_name)
                 global_coordinates.discard(coord_name)
 

xarray/tests/test_conventions.py

@@ -200,6 +200,41 @@ def test_missing_fillvalue(self):
         with self.assertWarns('floating point data as an integer'):
             conventions.encode_cf_variable(v)
 
+    def test_multidimensional_coordinates(self):
+        # regression test for GH1763
+        # Set up test case with coordinates that have overlapping (but not
+        # identical) dimensions.
+        zeros1 = np.zeros((1, 5, 3))
+        zeros2 = np.zeros((1, 6, 3))
+        zeros3 = np.zeros((1, 5, 4))
+        orig = Dataset({
+            'lon1': (['x1', 'y1'], zeros1.squeeze(0), {}),
+            'lon2': (['x2', 'y1'], zeros2.squeeze(0), {}),
+            'lon3': (['x1', 'y2'], zeros3.squeeze(0), {}),
+            'lat1': (['x1', 'y1'], zeros1.squeeze(0), {}),
+            'lat2': (['x2', 'y1'], zeros2.squeeze(0), {}),
+            'lat3': (['x1', 'y2'], zeros3.squeeze(0), {}),
+            'foo1': (['time', 'x1', 'y1'], zeros1,
+                     {'coordinates': 'lon1 lat1'}),
+            'foo2': (['time', 'x2', 'y1'], zeros2,
+                     {'coordinates': 'lon2 lat2'}),
+            'foo3': (['time', 'x1', 'y2'], zeros3,
+                     {'coordinates': 'lon3 lat3'}),
+            'time': ('time', [0.], {'units': 'hours since 2017-01-01'}),
+        })
+        orig = conventions.decode_cf(orig)
+        # Encode the coordinates, as they would be in a netCDF output file.
+        enc, attrs = conventions.encode_dataset_coordinates(orig)
+        # Make sure we have the right coordinates for each variable.
+        foo1_coords = enc['foo1'].attrs.get('coordinates', '')
+        foo2_coords = enc['foo2'].attrs.get('coordinates', '')
+        foo3_coords = enc['foo3'].attrs.get('coordinates', '')
+        assert set(foo1_coords.split()) == set(['lat1', 'lon1'])
+        assert set(foo2_coords.split()) == set(['lat2', 'lon2'])
+        assert set(foo3_coords.split()) == set(['lat3', 'lon3'])
+        # Should not have any global coordinates.
+        assert 'coordinates' not in attrs
+
 
 @requires_netCDF4
 class TestDecodeCF(TestCase):