Repository: ecmwf/earthkit-data
Pull request #595 (2025-01-22 15:50:02)

Feature/merge from main into develop


Patch:
diff --git a/docs/release_notes/version_0.12_updates.rst b/docs/release_notes/version_0.12_updates.rst
index 12640c3..0b5bdf4 100644
--- a/docs/release_notes/version_0.12_updates.rst
+++ b/docs/release_notes/version_0.12_updates.rst
@@ -1,6 +1,20 @@
 Version 0.12 Updates
 /////////////////////////
 
+Version 0.12.1
+===============
+
+Changes
+++++++++
+
+- Increased covjsonkit version to 0.1.1
+
+Fixes
+++++++++
+
+- Fixed issue when the earthkit Xarray engine prevented opening geotiff with :py:meth:`xarray.open_dataset()` (:pr:`591`)
+- Fixed issue when no "units" attribute was added to the Xarray dataset generated from a single GRIB variable (:pr:`592`)
+
 
 Version 0.12.0
 ===============
diff --git a/environment.yml b/environment.yml
index b60bb40..634169f 100644
--- a/environment.yml
+++ b/environment.yml
@@ -27,7 +27,7 @@ dependencies:
   - ecmwf-opendata>=0.1.2
   - polytope-client>=0.7.4
   - earthkit-meteo>=0.0.1
-  - covjsonkit>=0.0.28
+  - covjsonkit>=0.1.1
   - earthkit-geo>=0.2.0
 - tqdm>=4.63.0
 - lru-dict
diff --git a/pyproject.toml b/pyproject.toml
index af9ab44..7a6ba0f 100644
--- a/pyproject.toml
+++ b/pyproject.toml
@@ -49,7 +49,7 @@ dependencies = [
 optional-dependencies.all = [
   "cartopy",
   "cdsapi>=0.7.2",
-  "covjsonkit>=0.0.28",
+  "covjsonkit>=0.1.1",
   "earthkit-geo>=0.2",
   "ecmwf-api-client>=1.6.1",
   "ecmwf-opendata>=0.3.3",
@@ -61,11 +61,11 @@ optional-dependencies.all = [
 ]
 optional-dependencies.cds = [ "cdsapi>=0.7.2" ]
 optional-dependencies.ci = [ "numpy" ]
-optional-dependencies.covjsonkit = [ "covjsonkit>=0.0.28" ]
+optional-dependencies.covjsonkit = [ "covjsonkit>=0.1.1" ]
 optional-dependencies.dev = [
   "cartopy",
   "cdsapi>=0.7.2",
-  "covjsonkit>=0.0.28",
+  "covjsonkit>=0.1.1",
   "earthkit-data-demo-source",
   "earthkit-geo>=0.2",
   "ecmwf-api-client>=1.6.1",
diff --git a/src/earthkit/data/readers/grib/__init__.py b/src/earthkit/data/readers/grib/__init__.py
index 8b6162e..131b8c7 100644
--- a/src/earthkit/data/readers/grib/__init__.py
+++ b/src/earthkit/data/readers/grib/__init__.py
@@ -25,6 +25,19 @@ def _is_default(magic, content_type):
     return (magic is None or len(magic) == 0) and (content_type is None or len(content_type) == 0)
 
 
+def is_grib_file(path):
+    from earthkit.data.core.settings import SETTINGS
+
+    n_bytes = SETTINGS.get("reader-type-check-bytes")
+    magic = None
+    with open(path, "rb") as f:
+        magic = f.read(n_bytes)
+
+    if _match_magic(magic, False):
+        return True
+    return _match_magic(magic, True)
+
+
 def reader(source, path, *, magic=None, deeper_check=False, **kwargs):
     if _match_magic(magic, deeper_check):
         from .file import GRIBReader
diff --git a/src/earthkit/data/readers/grib/xarray.py b/src/earthkit/data/readers/grib/xarray.py
index 9d6e8a2..2d71698 100644
--- a/src/earthkit/data/readers/grib/xarray.py
+++ b/src/earthkit/data/readers/grib/xarray.py
@@ -231,10 +231,12 @@ class XarrayMixIn:
                 - "fixed": Use the attributes defined in ``variable_attrs`` as variables
                   attributes and ``global_attrs`` as global attributes.
                 - "unique": Use all the attributes defined in ``attrs``, ``variable_attrs``
-                  and ``global_attrs``. When an attribute has unique a value for a dataset
+                  and ``global_attrs``.  When an attribute has unique value for a dataset
                   it will be a global attribute, otherwise it will be a variable attribute.
-                  However keys in ``variable_attrs`` are always used as variable attributes,
-                  while keys in ``global_attrs`` are always used as global attributes.
+                  However, this logic is only applied if a unique variable attribute can be
+                  a global attribute according to the CF conventions Appendix A. (e.g. "units" cannot
+                  be a global attribute). Additionally, keys in ``variable_attrs`` are always used as
+                  variable attributes, while keys in ``global_attrs`` are always used as global attributes.
             * attrs: str, number, callable, dict or list of these, None
                 Attribute or list of attributes. Only used when ``attrs_mode`` is ``unique``.
                 Its default value (None) expands to [] unless the ``profile`` overwrites it.
diff --git a/src/earthkit/data/utils/xarray/attrs.py b/src/earthkit/data/utils/xarray/attrs.py
index 2515016..6762fd0 100644
--- a/src/earthkit/data/utils/xarray/attrs.py
+++ b/src/earthkit/data/utils/xarray/attrs.py
@@ -8,9 +8,11 @@
 #
 
 import logging
+import os
 from abc import ABCMeta
 from abc import abstractmethod
 from collections import defaultdict
+from functools import cached_property
 
 from earthkit.data.utils import ensure_dict
 from earthkit.data.utils import ensure_iterable
@@ -18,6 +20,36 @@ from earthkit.data.utils import ensure_iterable
 LOG = logging.getLogger(__name__)
 
 
+class CFAttrs:
+    def _load(self):
+        here = os.path.dirname(__file__)
+        path = os.path.join(here, "cf_attrs.yaml")
+        if os.path.exists(path):
+            import yaml
+
+            try:
+                with open(path, "r") as f:
+                    return yaml.safe_load(f)
+            except Exception as e:
+                LOG.exception(f"Failed read CF attributes file {path}. {e}")
+                raise
+        else:
+            raise ValueError(f"CF attributes file not found! path={path}")
+
+    @cached_property
+    def attrs(self):
+        return self._load()
+
+    def can_be_global(self, name):
+        item = self.attrs.get(name, None)
+        if item:
+            return "G" in item["use"]
+        return True
+
+
+CF_ATTRS = CFAttrs()
+
+
 class Attr:
     """Generic attribute class.
 
@@ -238,7 +270,7 @@ class AttrsBuilder(metaclass=ABCMeta):
                 global_attrs[item.name] = item.value()
 
         # TODO: make it optional
-        global_attrs.pop("units", None)
+        # global_attrs.pop("units", None)
 
         return global_attrs
 
@@ -257,13 +289,24 @@ class UniqueAttrBuilder(AttrsBuilder):
             if len(v) == 1 and k not in self.attrs.variable_attrs:
                 global_attrs[k] = list(v)[0]
 
-        for var_obj in t_vars.values():
-            var_obj.adjust_attrs(drop_keys=global_attrs.keys(), rename=rename)
+        # Some attrs cannot be global according to the CF convention.
+        # These are removed from global attrs and kept as variable attrs.
+        global_attrs_keys = list(global_attrs.keys())
+        global_attrs_renamed_keys = global_attrs_keys
+        if rename:
+            global_attrs_renamed_keys = list(rename(global_attrs).keys())
+
+        for k1, k2 in zip(global_attrs_keys, global_attrs_renamed_keys):
+            if not CF_ATTRS.can_be_global(k1) or not CF_ATTRS.can_be_global(k2):
+                global_attrs.pop(k1)
 
         for k in self.attrs.variable_attrs:
             if k in global_attrs:
                 global_attrs.pop(k)
 
+        for var_obj in t_vars.values():
+            var_obj.adjust_attrs(drop_keys=global_attrs.keys(), rename=rename)
+
         global_attrs = {k: v for k, v in global_attrs.items() if v is not None}
 
         return global_attrs
diff --git a/src/earthkit/data/utils/xarray/cf_attrs.yaml b/src/earthkit/data/utils/xarray/cf_attrs.yaml
new file mode 100644
index 0000000..03767eb
--- /dev/null
+++ b/src/earthkit/data/utils/xarray/cf_attrs.yaml
@@ -0,0 +1,195 @@
+# Based on CF Conventions Appendix A
+# All CF attributes are listed here except for those that are used to describe grid mappings. See Appendix F for the grid mapping attributes.
+# The 'Type' values are S for string, N for numeric, and D for the type of the data variable.
+# The 'Use' values are G for global, C for variables containing coordinate data, and D for variables containing non-coordinate data
+Conventions:
+  type: S
+  use: G
+_FillValue:
+  type: D
+  use:
+  - C
+  - D
+actual_range:
+  type: N
+  use:
+  - C
+  - D
+add_offset:
+  type: N
+  use:
+  - C
+  - D
+ancillary_variables:
+  type: S
+  use: D
+axis:
+  type: S
+  use: C
+bounds:
+  type: S
+  use: C
+calendar:
+  type: S
+  use: C
+cell_measures:
+  type: S
+  use: D
+cell_methods:
+  type: S
+  use: D
+cf_role:
+  type: S
+  use: C
+climatology:
+  type: S
+  use: C
+comment:
+  type: S
+  use:
+  - G
+  - C
+  - D
+compress:
+  type: S
+  use: C
+computed_standard_name:
+  type: S
+  use: C
+coordinates:
+  type: S
+  use:
+  - D
+  - M
+external_variables:
+  type: S
+  use: G
+featureType:
+  type: S
+  use: G
+flag_masks:
+  type: D
+  use: D
+flag_meanings:
+  type: S
+  use: D
+flag_values:
+  type: D
+  use: D
+formula_terms:
+  type: S
+  use: C
+geometry:
+  type: S
+  use:
+  - C
+  - D
+geometry_type:
+  type: S
+  use: M
+grid_mapping:
+  type: S
+  use:
+  - D
+  - M
+history:
+  type: S
+  use:
+  - G
+  - Gr
+instance_dimension:
+  type: S
+  use: "-"
+institution:
+  type: S
+  use:
+  - G
+  - D
+interior_ring:
+  type: S
+  use: M
+leap_month:
+  type: N
+  use: C
+leap_year:
+  type: N
+  use: C
+long_name:
+  type: S
+  use:
+  - C
+  - D
+missing_value:
+  type: D
+  use:
+  - C
+  - D
+month_lengths:
+  type: N
+  use: C
+node_coordinates:
+  type: S
+  use: M
+node_count:
+  type: S
+  use: M
+nodes:
+  type: S
+  use: C
+part_node_count:
+  type: S
+  use: M
+positive:
+  type: S
+  use: C
+references:
+  type: S
+  use:
+  - G
+  - D
+sample_dimension:
+  type: S
+  use: "-"
+scale_factor:
+  type: N
+  use:
+  - C
+  - D
+source:
+  type: S
+  use:
+  - G
+  - D
+standard_error_multiplier:
+  type: N
+  use: D
+standard_name:
+  type: S
+  use:
+  - C
+  - D
+title:
+  type: S
+  use:
+  - G
+  - Gr
+units:
+  type: S
+  use:
+  - C
+  - D
+valid_max:
+  type: N
+  use:
+  - C
+  - D
+valid_min:
+  type: N
+  use:
+  - C
+  - D
+valid_range:
+  type: N
+  use:
+  - C
+  - D
diff --git a/src/earthkit/data/utils/xarray/engine.py b/src/earthkit/data/utils/xarray/engine.py
index 86a9f8e..1dfaa20 100644
--- a/src/earthkit/data/utils/xarray/engine.py
+++ b/src/earthkit/data/utils/xarray/engine.py
@@ -188,10 +188,12 @@ class EarthkitBackendEntrypoint(BackendEntrypoint):
             - "fixed": Use the attributes defined in ``variable_attrs`` as variables
               attributes and ``global_attrs`` as global attributes.
             - "unique": Use all the attributes defined in ``attrs``, ``variable_attrs``
-              and ``global_attrs``. When an attribute has unique a value for a dataset
+              and ``global_attrs``. When an attribute has unique value for a dataset
               it will be a global attribute, otherwise it will be a variable attribute.
-              However keys in ``variable_attrs`` are always used as variable attributes,
-              while keys in ``global_attrs`` are always used as global attributes.
+              However, this logic is only applied if a unique variable attribute can be
+              a global attribute according to the CF conventions Appendix A. (e.g. "units" cannot
+              be a global attribute). Additionally, keys in ``variable_attrs`` are always used as
+              variable attributes, while keys in ``global_attrs`` are always used as global attributes.
         attrs: str, number, callable, dict or list of these, None
             Attribute or list of attributes. Only used when ``attrs_mode`` is ``unique``.
             Its default value (None) expands to [] unless the ``profile`` overwrites it.
@@ -293,7 +295,21 @@ class EarthkitBackendEntrypoint(BackendEntrypoint):
 
     @classmethod
     def guess_can_open(cls, filename_or_obj):
-        return True  # filename_or_obj.endswith(".grib")
+        try:
+            from earthkit.data.core import Base
+
+            if isinstance(filename_or_obj, Base):
+                return True
+            elif isinstance(filename_or_obj, str):
+                from earthkit.data.readers.grib import is_grib_file
+
+                return is_grib_file(filename_or_obj)
+        except Exception:
+            LOG.debug(
+                "Failed to guess if %s can be opened by the earthkit backend", filename_or_obj, exc_info=True
+            )
+
+        return False
 
     @staticmethod
     def _fieldlist(filename_or_obj, source_type):
@@ -301,15 +317,10 @@ class EarthkitBackendEntrypoint(BackendEntrypoint):
 
         if isinstance(filename_or_obj, Base):
             ds = filename_or_obj
-        # TODO: Add Path? or handle with try statement
         elif isinstance(filename_or_obj, str):
             from earthkit.data import from_source
 
             ds = from_source(source_type, filename_or_obj)
-        else:
-            from earthkit.data import from_object
-
-            ds = from_object(filename_or_obj)
         return ds
 
 


Test patch:
diff --git a/tests/environment-unit-tests.yml b/tests/environment-unit-tests.yml
index 55e841e..bece444 100644
--- a/tests/environment-unit-tests.yml
+++ b/tests/environment-unit-tests.yml
@@ -30,7 +30,7 @@ dependencies:
   - polytope-client>=0.7.4
   - earthkit-meteo>=0.0.1
   - git+https://github.com/ecmwf/earthkit-data-demo-source
-  - covjsonkit>=0.0.28
+  - covjsonkit>=0.1.1
   - earthkit-geo>=0.2.0
 - tqdm>=4.63.0
 - lru-dict
diff --git a/tests/readers/test_geotiff_reader.py b/tests/readers/test_geotiff_reader.py
index b5fcaad..7726a39 100644
--- a/tests/readers/test_geotiff_reader.py
+++ b/tests/readers/test_geotiff_reader.py
@@ -32,6 +32,16 @@ def test_geotiff_reader_with_multiband():
     assert s[0].shape == (294, 315)
 
 
+@pytest.mark.skipif(NO_RIOXARRAY, reason="rioxarray not available")
+@pytest.mark.with_proj
+def test_geotiff_bypassing_xr_engine():
+    import xarray as xr
+
+    # the earthkit xarray engine should not try to handle this file but leave it to the rioxarray backend
+    da = xr.open_dataarray(earthkit_test_data_file("dgm50hs_col_32_368_5616_nw.tif"))
+    assert da.name == "band_data"
+
+
 if __name__ == "__main__":
     from earthkit.data.testing import main
 
diff --git a/tests/xr_engine/test_xr_engine.py b/tests/xr_engine/test_xr_engine.py
index 4aa3794..59a9e45 100644
--- a/tests/xr_engine/test_xr_engine.py
+++ b/tests/xr_engine/test_xr_engine.py
@@ -469,10 +469,14 @@ def test_xr_engine_single_field():
     lats = np.linspace(90, -90, 19)
     lons = np.linspace(0, 350, 36)
 
-    attrs_ref = {
-        "param": "t",
+    var_attrs_ref = {
         "standard_name": "air_temperature",
         "long_name": "Temperature",
+        "units": "K",
+    }
+
+    global_attrs_ref = {
+        "param": "t",
         "paramId": 130,
         "class": "od",
         "stream": "oper",
@@ -488,7 +492,7 @@ def test_xr_engine_single_field():
         "institution": "ECMWF",
     }
 
-    assert ds.attrs == attrs_ref
+    assert ds.attrs == global_attrs_ref
 
     data_vars = ["t"]
 
@@ -510,6 +514,9 @@ def test_xr_engine_single_field():
 
     da = ds["t"]
 
+    for k, v in var_attrs_ref.items():
+        assert da.attrs[k] == v
+
     r = da[:, :]
     r.shape == (19, 36)
     assert np.allclose(r.values, vals_ref)