Repository: ioos/compliance-checker
Pull request #697 (2019-11-14 20:54:27)

Verify grid mapping coordinate vars under CF 1.7
Verifies that grid mapping coordinate variables that exist as part of a grid_mapping attribute are well formed and exist as variables within the dataset.  Implements https://cf-trac.llnl.gov/trac/ticket/70.


Patch:
diff --git a/compliance_checker/cf/cf.py b/compliance_checker/cf/cf.py
index d8c5907..a518070 100644
--- a/compliance_checker/cf/cf.py
+++ b/compliance_checker/cf/cf.py
@@ -139,6 +139,142 @@ class CFBaseCheck(BaseCheck):
         self._find_cf_standard_name_table(ds)
         self._find_geophysical_vars(ds)
 
+    def check_grid_mapping(self, ds):
+        """
+        5.6 When the coordinate variables for a horizontal grid are not
+        longitude and latitude, it is required that the true latitude and
+        longitude coordinates be supplied via the coordinates attribute. If in
+        addition it is desired to describe the mapping between the given
+        coordinate variables and the true latitude and longitude coordinates,
+        the attribute grid_mapping may be used to supply this description.
+
+        This attribute is attached to data variables so that variables with
+        different mappings may be present in a single file. The attribute takes
+        a string value which is the name of another variable in the file that
+        provides the description of the mapping via a collection of attached
+        attributes. This variable is called a grid mapping variable and is of
+        arbitrary type since it contains no data. Its purpose is to act as a
+        container for the attributes that define the mapping.
+
+        The one attribute that all grid mapping variables must have is
+        grid_mapping_name which takes a string value that contains the mapping's
+        name. The other attributes that define a specific mapping depend on the
+        value of grid_mapping_name. The valid values of grid_mapping_name along
+        with the attributes that provide specific map parameter values are
+        described in Appendix F, Grid Mappings.
+
+        When the coordinate variables for a horizontal grid are longitude and
+        latitude, a grid mapping variable with grid_mapping_name of
+        latitude_longitude may be used to specify the ellipsoid and prime
+        meridian.
+
+
+        In order to make use of a grid mapping to directly calculate latitude
+        and longitude values it is necessary to associate the coordinate
+        variables with the independent variables of the mapping. This is done by
+        assigning a standard_name to the coordinate variable. The appropriate
+        values of the standard_name depend on the grid mapping and are given in
+        Appendix F, Grid Mappings.
+
+        :param netCDF4.Dataset ds: An open netCDF dataset
+        :rtype: list
+        :return: List of results
+        """
+
+        ret_val = OrderedDict()
+        grid_mapping_variables = cfutil.get_grid_mapping_variables(ds)
+
+        # Check the grid_mapping attribute to be a non-empty string and that its reference exists
+        for variable in ds.get_variables_by_attributes(grid_mapping=lambda x: x is not None):
+            grid_mapping = getattr(variable, 'grid_mapping', None)
+            defines_grid_mapping = self.get_test_ctx(BaseCheck.HIGH,
+                                                     self.section_titles["5.6"],
+                                                     variable.name)
+            defines_grid_mapping.assert_true((isinstance(grid_mapping, basestring) and grid_mapping),
+                                             "{}'s grid_mapping attribute must be a "+\
+                                             "space-separated non-empty string".format(variable.name))
+            if isinstance(grid_mapping, basestring):
+                # TODO (badams): refactor functionality to split functionality
+                #                into requisite classes
+                if ':' in grid_mapping and self._cc_spec_version >= '1.7':
+                    colon_count = grid_mapping.count(':')
+                    re_all = regex.findall(r"(\w+):\s*((?:\w+\s+)*(?:\w+)(?![\w:]))",
+                                        grid_mapping)
+                    if colon_count != len(re_all):
+                        defines_grid_mapping.out_of += 1
+                        defines_grid_mapping.messages.append(
+                         "Could not consume entire grid_mapping expression, please check for well-formedness")
+                    else:
+                        for grid_var_name, coord_var_str in re_all:
+                            defines_grid_mapping.assert_true(grid_var_name in ds.variables,
+                                                        "grid mapping variable {} must exist in this dataset".format(grid_var_name))
+                            for ref_var in coord_var_str.split():
+                                defines_grid_mapping.assert_true(ref_var in ds.variables,
+                                                            "Coordinate-related variable {} referenced by grid_mapping variable {} must exist in this dataset".format(ref_var, grid_var_name))
+
+                else:
+                    for grid_var_name in grid_mapping.split():
+                        defines_grid_mapping.assert_true(grid_var_name in ds.variables,
+                                                    "grid mapping variable {} must exist in this dataset".format(grid_var_name))
+            ret_val[variable.name] = defines_grid_mapping.to_result()
+
+        # Check the grid mapping variables themselves
+        for grid_var_name in grid_mapping_variables:
+            valid_grid_mapping = self.get_test_ctx(BaseCheck.HIGH,
+                                                   self.section_titles["5.6"],
+                                                   grid_var_name)
+            grid_var = ds.variables[grid_var_name]
+
+            grid_mapping_name = getattr(grid_var, 'grid_mapping_name', None)
+
+            # Grid mapping name must be in appendix F
+            valid_grid_mapping.assert_true(grid_mapping_name in self.grid_mapping_dict,
+                                           "{} is not a valid grid_mapping_name.".format(grid_mapping_name)+\
+                                           " See Appendix F for valid grid mappings")
+
+            # The self.grid_mapping_dict has a values of:
+            # - required attributes
+            # - optional attributes (can't check)
+            # - required standard_names defined
+            # - at least one of these attributes must be defined
+
+            # We can't do any of the other grid mapping checks if it's not a valid grid mapping name
+            if grid_mapping_name not in self.grid_mapping_dict:
+                ret_val[grid_mapping_name] = valid_grid_mapping.to_result()
+                continue
+
+            grid_mapping = self.grid_mapping_dict[grid_mapping_name]
+            required_attrs = grid_mapping[0]
+            # Make sure all the required attributes are defined
+            for req in required_attrs:
+                valid_grid_mapping.assert_true(hasattr(grid_var, req),
+                                               "{} is a required attribute for grid mapping {}".format(req, grid_mapping_name))
+
+            # Make sure that exactly one of the exclusive attributes exist
+            if len(grid_mapping) == 4:
+                at_least_attr = grid_mapping[3]
+                number_found = 0
+                for attr in at_least_attr:
+                    if hasattr(grid_var, attr):
+                        number_found += 1
+                valid_grid_mapping.assert_true(number_found == 1,
+                                               "grid mapping {}".format(grid_mapping_name) +\
+                                               "must define exactly one of these attributes: "+\
+                                               "{}".format(' or '.join(at_least_attr)))
+
+            # Make sure that exactly one variable is defined for each of the required standard_names
+            expected_std_names = grid_mapping[2]
+            for expected_std_name in expected_std_names:
+                found_vars = ds.get_variables_by_attributes(standard_name=expected_std_name)
+                valid_grid_mapping.assert_true(len(found_vars) == 1,
+                                               "grid mapping {} requires exactly".format(grid_mapping_name)+\
+                                               "one variable with standard_name "+\
+                                               "{} to be defined".format(expected_std_name))
+
+            ret_val[grid_var_name] = valid_grid_mapping.to_result()
+
+        return ret_val
+
     def check_conventions_version(self, ds):
         '''
         CF §2.6.1 the NUG defined global attribute Conventions to the string
@@ -2414,123 +2550,6 @@ class CF1_6Check(CFNCCheck):
 
     # grid mapping dictionary, appendix F
 
-    def check_grid_mapping(self, ds):
-        """
-        5.6 When the coordinate variables for a horizontal grid are not
-        longitude and latitude, it is required that the true latitude and
-        longitude coordinates be supplied via the coordinates attribute. If in
-        addition it is desired to describe the mapping between the given
-        coordinate variables and the true latitude and longitude coordinates,
-        the attribute grid_mapping may be used to supply this description.
-
-        This attribute is attached to data variables so that variables with
-        different mappings may be present in a single file. The attribute takes
-        a string value which is the name of another variable in the file that
-        provides the description of the mapping via a collection of attached
-        attributes. This variable is called a grid mapping variable and is of
-        arbitrary type since it contains no data. Its purpose is to act as a
-        container for the attributes that define the mapping.
-
-        The one attribute that all grid mapping variables must have is
-        grid_mapping_name which takes a string value that contains the mapping's
-        name. The other attributes that define a specific mapping depend on the
-        value of grid_mapping_name. The valid values of grid_mapping_name along
-        with the attributes that provide specific map parameter values are
-        described in Appendix F, Grid Mappings.
-
-        When the coordinate variables for a horizontal grid are longitude and
-        latitude, a grid mapping variable with grid_mapping_name of
-        latitude_longitude may be used to specify the ellipsoid and prime
-        meridian.
-
-
-        In order to make use of a grid mapping to directly calculate latitude
-        and longitude values it is necessary to associate the coordinate
-        variables with the independent variables of the mapping. This is done by
-        assigning a standard_name to the coordinate variable. The appropriate
-        values of the standard_name depend on the grid mapping and are given in
-        Appendix F, Grid Mappings.
-
-        :param netCDF4.Dataset ds: An open netCDF dataset
-        :rtype: list
-        :return: List of results
-        """
-
-        ret_val = OrderedDict()
-        grid_mapping_variables = cfutil.get_grid_mapping_variables(ds)
-
-        # Check the grid_mapping attribute to be a non-empty string and that its reference exists
-        for variable in ds.get_variables_by_attributes(grid_mapping=lambda x: x is not None):
-            grid_mapping = getattr(variable, 'grid_mapping', None)
-            defines_grid_mapping = self.get_test_ctx(BaseCheck.HIGH,
-                                                     self.section_titles["5.6"],
-                                                     variable.name)
-            defines_grid_mapping.assert_true((isinstance(grid_mapping, basestring) and grid_mapping),
-                                             "{}'s grid_mapping attribute must be a "+\
-                                             "space-separated non-empty string".format(variable.name))
-
-            if isinstance(grid_mapping, basestring):
-                for grid_var_name in grid_mapping.split():
-                    defines_grid_mapping.assert_true(grid_var_name in ds.variables,
-                                                  "grid mapping variable {} must exist in this dataset".format(variable.name))
-            ret_val[variable.name] = defines_grid_mapping.to_result()
-
-        # Check the grid mapping variables themselves
-        for grid_var_name in grid_mapping_variables:
-            valid_grid_mapping = self.get_test_ctx(BaseCheck.HIGH,
-                                                   self.section_titles["5.6"],
-                                                   grid_var_name)
-            grid_var = ds.variables[grid_var_name]
-
-            grid_mapping_name = getattr(grid_var, 'grid_mapping_name', None)
-
-            # Grid mapping name must be in appendix F
-            valid_grid_mapping.assert_true(grid_mapping_name in self.grid_mapping_dict,
-                                           "{} is not a valid grid_mapping_name.".format(grid_mapping_name)+\
-                                           " See Appendix F for valid grid mappings")
-
-            # The self.grid_mapping_dict has a values of:
-            # - required attributes
-            # - optional attributes (can't check)
-            # - required standard_names defined
-            # - at least one of these attributes must be defined
-
-            # We can't do any of the other grid mapping checks if it's not a valid grid mapping name
-            if grid_mapping_name not in self.grid_mapping_dict:
-                ret_val[grid_mapping_name] = valid_grid_mapping.to_result()
-                continue
-
-            grid_mapping = self.grid_mapping_dict[grid_mapping_name]
-            required_attrs = grid_mapping[0]
-            # Make sure all the required attributes are defined
-            for req in required_attrs:
-                valid_grid_mapping.assert_true(hasattr(grid_var, req),
-                                               "{} is a required attribute for grid mapping {}".format(req, grid_mapping_name))
-
-            # Make sure that exactly one of the exclusive attributes exist
-            if len(grid_mapping) == 4:
-                at_least_attr = grid_mapping[3]
-                number_found = 0
-                for attr in at_least_attr:
-                    if hasattr(grid_var, attr):
-                        number_found += 1
-                valid_grid_mapping.assert_true(number_found == 1,
-                                               "grid mapping {}".format(grid_mapping_name) +\
-                                               "must define exactly one of these attributes: "+\
-                                               "{}".format(' or '.join(at_least_attr)))
-
-            # Make sure that exactly one variable is defined for each of the required standard_names
-            expected_std_names = grid_mapping[2]
-            for expected_std_name in expected_std_names:
-                found_vars = ds.get_variables_by_attributes(standard_name=expected_std_name)
-                valid_grid_mapping.assert_true(len(found_vars) == 1,
-                                               "grid mapping {} requires exactly".format(grid_mapping_name)+\
-                                               "one variable with standard_name "+\
-                                               "{} to be defined".format(expected_std_name))
-
-            ret_val[grid_var_name] = valid_grid_mapping.to_result()
-
-        return ret_val
 
     ###############################################################################
     # Chapter 6: Labels and Alternative Coordinates
@@ -3846,17 +3865,17 @@ class CF1_7Check(CF1_6Check):
             deprecated_units,
             self._check_dimensionless_vertical_coordinate_1_6,
             dimless_vertical_coordinates_1_7)
-        ) 
+        )
 
         ret_val.extend(self._check_dimensionless_vertical_coordinates(
             ds,
             deprecated_units,
             self._check_dimensionless_vertical_coordinate_1_7,
             dimless_vertical_coordinates_1_7)
-        ) 
+        )
 
         return ret_val
-     
+
 
 class CFNCCheck(BaseNCCheck, CFBaseCheck):
 


Test patch:
diff --git a/compliance_checker/tests/data/grid_mapping_coordinates.cdl b/compliance_checker/tests/data/grid_mapping_coordinates.cdl
new file mode 100644
index 0000000..82bb211
--- /dev/null
+++ b/compliance_checker/tests/data/grid_mapping_coordinates.cdl
@@ -0,0 +1,51 @@
+netcdf example_grid_mapping {
+  dimensions:
+    z = 1 ;
+    y = 1 ;
+    x = 1 ;
+  variables:
+    double x(x) ;
+      x:standard_name = "projection_x_coordinate" ;
+      x:long_name = "Easting" ;
+      x:units = "m" ;
+    double y(y) ;
+      y:standard_name = "projection_y_coordinate" ;
+      y:long_name = "Northing" ;
+      y:units = "m" ;
+    double z(z) ;
+      z:standard_name = "height_above_reference_ellipsoid" ;
+      z:long_name = "height_above_osgb_newlyn_datum_masl" ;
+      z:units = "m" ;
+    double lat(y, x) ;
+      lat:standard_name = "latitude" ;
+      lat:units = "degrees_north" ;
+    double lon(y, x) ;
+      lon:standard_name = "longitude" ;
+      lon:units = "degrees_east" ;
+    float temp(z, y, x) ;
+      temp:standard_name = "air_temperature" ;
+      temp:units = "K" ;
+      temp:coordinates = "lat lon" ;
+      temp:grid_mapping = "crsOSGB: x y crsWGS84: lat lon" ;
+    float pres(z, y, x) ;
+      temp:standard_name = "air_pressure" ;
+      temp:units = "Pa" ;
+      temp:coordinates = "lat lon" ;
+      temp:grid_mapping = "crsOSGB: x y crsWGS84: lat lon" ;
+    int crsOSGB ;
+      crsOSGB:grid_mapping_name = "transverse_mercator";
+      crsOSGB:semi_major_axis = 6377563.396 ;
+      crsOSGB:inverse_flattening = 299.3249646 ;
+      crsOSGB:longitude_of_prime_meridian = 0.0 ;
+      crsOSGB:latitude_of_projection_origin = 49.0 ;
+      crsOSGB:longitude_of_central_meridian = -2.0 ;
+      crsOSGB:scale_factor_at_central_meridian = 0.9996012717 ;
+      crsOSGB:false_easting = 400000.0 ;
+      crsOSGB:false_northing = -100000.0 ;
+      crsOSGB:unit = "metre" ;
+    int crsWGS84 ;
+      crsWGS84:grid_mapping_name = "latitude_longitude";
+      crsWGS84:longitude_of_prime_meridian = 0.0 ;
+      crsWGS84:semi_major_axis = 6378137.0 ;
+      crsWGS84:inverse_flattening = 298.257223563 ;
+}
diff --git a/compliance_checker/tests/resources.py b/compliance_checker/tests/resources.py
index 71714a7..d125ff6 100644
--- a/compliance_checker/tests/resources.py
+++ b/compliance_checker/tests/resources.py
@@ -52,6 +52,7 @@ STATIC_FILES = {
     'ghrsst'                               : get_filename('tests/data/20160919092000-ABOM-L3S_GHRSST-SSTfnd-AVHRR_D-1d_dn_truncate.cdl'),
     'glcfs'                                : get_filename('tests/data/examples/glcfs.cdl'),
     'grid-boundaries'                      : get_filename('tests/data/grid-boundaries.cdl'),
+    'grid_mapping_coordinates'             : get_filename('tests/data/grid_mapping_coordinates.cdl'),
     'hycom_global'                         : get_filename('tests/data/examples/hycom_global.cdl'),
     'h_point'                              : get_filename('tests/data/appendix_h/point.cdl'),
     'h_timeseries-incomplete'              : get_filename('tests/data/appendix_h/timeseries-incomplete.cdl'),
diff --git a/compliance_checker/tests/test_cf.py b/compliance_checker/tests/test_cf.py
index 9a3143f..510c7c9 100644
--- a/compliance_checker/tests/test_cf.py
+++ b/compliance_checker/tests/test_cf.py
@@ -1395,6 +1395,31 @@ class TestCF1_7(BaseTestCase):
          'false_easting is a required attribute for grid mapping stereographic'])
         self.assertEqual(msg_len, 0)
 
+    def test_check_grid_mapping_coordinates(self):
+        """
+        Checks that coordinates variables referred to by a grid mapping
+        are well-formed and exist.
+        """
+        dataset = self.load_dataset(STATIC_FILES['grid_mapping_coordinates'])
+        valid_grid_mapping = copy.deepcopy(self.cf)
+        valid_grid_mapping_2 = copy.deepcopy(self.cf)
+        dataset.variables['temp'] = MockVariable(dataset.variables['temp'])
+        results = self.cf.check_grid_mapping(dataset)
+        self.assertEqual(results['temp'].value[0], results['temp'].value[1])
+        malformed_sep = "crsOSGB: x y : lat lon"
+        dataset.variables['temp'].grid_mapping = malformed_sep
+        results = valid_grid_mapping.check_grid_mapping(dataset)
+        self.assertIn("Could not consume entire grid_mapping expression, please check for well-formedness",
+                      results['temp'].msgs)
+        self.assertLess(*results['temp'].value)
+        malformed_var = "crsOSGB: x y_null z_null"
+        dataset.variables['temp'].grid_mapping = malformed_var
+        results = valid_grid_mapping_2.check_grid_mapping(dataset)
+        self.assertEqual(['Coordinate-related variable y_null referenced by grid_mapping variable crsOSGB must exist in this dataset',
+                          'Coordinate-related variable z_null referenced by grid_mapping variable crsOSGB must exist in this dataset'],
+                         results['temp'].msgs)
+        self.assertLess(*results['temp'].value)
+
     def test_check_grid_mapping_vert_datum_geoid_name(self):
         """Checks that geoid_name works proerly"""
         dataset = self.load_dataset(STATIC_FILES['mapping'])
@@ -1473,7 +1498,7 @@ class TestCF1_7(BaseTestCase):
             dataset.createVariable('lev', 'd') # dtype=double, dims=1
             dataset.variables['lev'].setncattr('standard_name', 'atmosphere_sigma_coordinate')
             dataset.variables['lev'].setncattr('formula_terms', 'sigma: lev ps: PS ptop: PTOP')
-            
+
             dataset.createVariable('PS', 'd', ('time',)) # dtype=double, dims=time
             dataset.createVariable('PTOP', 'd', ('time',)) # dtype=double, dims=time
 
@@ -1500,7 +1525,7 @@ class TestCF1_7(BaseTestCase):
             # computed_standard_name is assigned, should pass
             score, out_of, messages = get_results(ret_val)
             assert score == out_of
-        
+
     def test_dimensionless_vertical(self):
         '''
         Section 4.3.2 check, but for CF-1.7 implementation. With the refactor in