Add fieldlist support (#15)

* Add fieldlist support

docs/conf.py

@@ -33,6 +33,7 @@
 extensions = [
     "sphinx_rtd_theme",
     "nbsphinx",
+    "earthkit.regrid.sphinxext.xref",
     "earthkit.regrid.sphinxext.module_output",
 ]
 
@@ -61,3 +62,16 @@
 html_css_files = ["style.css"]
 
 # html_logo = "_static/earthkit-regrid.png"
+
+
+xref_links = {
+    "earthkit": ("earthkit", "https://earthkit.readthedocs.io/en/latest/"),
+    "earthkit-data": (
+        "earthkit-data",
+        "https://earthkit-data.readthedocs.io/en/latest/",
+    ),
+    "fieldlist": (
+        "fieldlist",
+        "https://earthkit-data.readthedocs.io/en/latest/guide/data_format/grib.html",
+    ),
+}

docs/gridspec.rst

@@ -45,9 +45,9 @@ global regular latitude-longitude grid
 
 The ``grid`` format is::
 
-    [DX, DY]
+    [DLON, DLAT]
 
-where *DX* and *DY* are the increments in degrees in longitudes and latitudes, respectively. This grid definition uses the default scanning mode used at ECMWF: values start at North-West and follow in consecutive rows from West to East
+where *DLON* and *DLAT* are the increments in degrees in longitudes and latitudes, respectively. This grid definition uses the default scanning mode used at ECMWF: values start at North-West and follow in consecutive rows from West to East, where West is always the 0° meridian.
 
 Example:
 

docs/index.rst

@@ -6,7 +6,9 @@ Welcome to earthkit-regrids's documentation
     This project is **BETA** and will be **Experimental** for the foreseeable future. Interfaces and functionality are likely to change, and the project itself may be scrapped. **DO NOT** use this software in any project/software that is operational.
 
 
-**earthkit-regrid** is a Python package for regridding. It features the :func:`interpolate` function to regrid values stored in an ndarray. At the moment, regridding is only available for a **pre-defined** set of source and target global grid combinations.
+**earthkit-regrid** is a Python package for regridding. It is one of the components of :xref:`earthkit`.
+
+The **earthkit-regrid** API features the :func:`interpolate` function to regrid values stored in an ndarray or earthkit-data GRIB :xref:`fieldlist`. At the moment, regridding is only available for a **pre-defined** set of source and target global grid combinations.
 
 
 .. toctree::
@@ -29,10 +31,18 @@ Welcome to earthkit-regrids's documentation
    :caption: Installation
 
    install
+   release_notes/index
    development
    licence
 
 
+.. toctree::
+   :maxdepth: 1
+   :caption: Projects
+
+   earthkit <https://earthkit.readthedocs.io/en/latest>
+
+
 Indices and tables
 ==================
 

docs/interpolate.rst

@@ -1,25 +1,26 @@
 interpolate
 ==============
 
-.. py:function:: interpolate(values, source_gridspec, target_gridspec, matrix_source=None, method='linear')
+.. py:function:: interpolate(values, in_grid=None, out_grid=None, matrix_source=None, method='linear')
 
-    Interpolate the ``values`` from the ``source_gridspec`` onto the ``target_gridspec``.
+    Interpolate the ``values`` from the ``in_grid`` onto the ``out_grid``.
 
-    :param values: the input values
-    :type values: ndarray
-    :param source_gridspec: the :ref:`gridspec <gridspec>` describing the grid that ``values`` are defined on
-    :type source_gridspec: dict
-    :param target_gridspec: the :ref:`gridspec <gridspec>` describing the target grid that ``values`` will be interpolated onto
-    :type target_gridspec: dict
+    :param values: the input data. It can be an ndarray or an earthkit-data GRIB :xref:`fieldlist`. ndarrays are assumed to be defined on the ``in_grid``.
+    :type values: ndarray, :xref:`fieldlist`
+    :param in_grid: the :ref:`gridspec <gridspec>` describing the grid that ``values`` are defined on. When ``values`` is a :xref:`fieldlist` the input grid is automatically detected and ``in_grid`` cannot be specified.
+    :type in_grid: dict
+    :param out_grid: the :ref:`gridspec <gridspec>` describing the target grid that ``values`` will be interpolated onto
+    :type out_grid: dict
     :param method: the interpolation method. Possible values are ``linear`` and ``nearest-neighbour``. For ``nearest-neighbour`` the following aliases are also supported: ``nn``, ``nearest-neighbor``.
     :type method: str
     :param matrix_source: (experimental) the location of a user specified pre-generated matrix inventory. When it is None the default matrix inventory hosted on an ECMWF download server is used.
     :type matrix_source: str, None
-    :return: the interpolated values
-    :rtype: ndarray
+    :return: The same type of data as ``values`` containing the interpolated values.
+    :rtype: ndarray, :xref:`fieldlist`
     :raises ValueError: if a pre-generated interpolation matrix is not available
+    :raises ValueError: if ``in_grid`` is specified for a :xref:`fieldlist` input
 
-    The interpolation only works when a pre-generated interpolation matrix is available for the given ``source_gridspec``, ``target_gridspec`` and ``method`` combination.
+    The interpolation only works when a pre-generated interpolation matrix is available for the given ``in_grid``, ``out_grid`` and ``method`` combination.
 
     When ``matrix_source`` is None (default) the interpolation matrix is automatically downloaded and stored in a local cache and when it is needed again the cached version is used.
 

docs/release_notes/index.rst

@@ -0,0 +1,7 @@
+Release notes
+==============
+
+.. toctree::
+    :maxdepth: 1
+
+    version_0.3_updates

docs/release_notes/version_0.3_updates.rst

@@ -0,0 +1,17 @@
+Version 0.3 Updates
+/////////////////////////
+
+
+Version 0.3.0
+===============
+
+
+New features
+++++++++++++++++
+- restructured and regenerated matrix inventory
+- allow using the ``method`` keyword in :func:`interpolate` to specify the interpolation method
+- allow using earthkit-data GRIB :xref:`fieldlist` in :func:`interpolate` as input data. This only works when  the output grid is regular latitude-longitude grid.
+- added notebook examples:
+
+   - :ref:`/examples/healpix_fieldlist.ipynb`
+   - :ref:`/examples/octahedral_fieldlist.ipynb`

earthkit/regrid/interpolate.py

@@ -10,17 +10,85 @@
 from earthkit.regrid.db import find
 
 
-def interpolate(values, source_gridspec, target_gridspec, method="linear", **kwargs):
-    z, shape = find(source_gridspec, target_gridspec, method, **kwargs)
+def interpolate(values, in_grid=None, out_grid=None, method="linear", **kwargs):
+    interpolator = _find_interpolator(values)
+    if interpolator is None:
+        raise ValueError(f"Cannot interpolate data with type={type(values)}")
+
+    return interpolator(
+        values, in_grid=in_grid, out_grid=out_grid, method=method, **kwargs
+    )
+
+
+def _find_interpolator(values):
+    for interpolator in INTERPOLATORS:
+        if interpolator.match(values):
+            return interpolator
+    return None
+
+
+def _interpolate(values, in_grid, out_grid, method, **kwargs):
+    z, shape = find(in_grid, out_grid, method, **kwargs)
 
     if z is None:
-        raise ValueError(
-            f"No matrix found! {source_gridspec=} {target_gridspec=} {method=}"
-        )
+        raise ValueError(f"No matrix found! {in_grid=} {out_grid=} {method=}")
 
     # This should check for 1D (GG) and 2D (LL) matrices
     values = values.reshape(-1, 1)
 
     values = z @ values
 
     return values.reshape(shape)
+
+
+class NumpyInterpolator:
+    @staticmethod
+    def match(values):
+        import numpy as np
+
+        return isinstance(values, np.ndarray)
+
+    def __call__(self, values, **kwargs):
+        in_grid = kwargs.pop("in_grid")
+        out_grid = kwargs.pop("out_grid")
+        method = kwargs.pop("method")
+        return _interpolate(values, in_grid, out_grid, method, **kwargs)
+
+
+class FieldListInterpolator:
+    @staticmethod
+    def match(values):
+        try:
+            import earthkit.data
+
+            return isinstance(values, earthkit.data.FieldList)
+        except ImportError:
+            return False
+
+    def __call__(self, values, **kwargs):
+        import earthkit.data
+
+        ds = values
+        in_grid = kwargs.pop("in_grid")
+        if in_grid is not None:
+            raise ValueError("in_grid cannot be used for FieldList interpolation")
+        out_grid = kwargs.pop("out_grid")
+        method = kwargs.pop("method")
+
+        r = earthkit.data.FieldList()
+        for f in ds:
+            vv = f.to_numpy(flatten=True)
+            v_res = _interpolate(
+                vv,
+                f.metadata().gridspec,
+                out_grid,
+                method,
+                **kwargs,
+            )
+            md_res = f.metadata().override(gridspec=out_grid)
+            r += ds.from_numpy(v_res, md_res)
+
+        return r
+
+
+INTERPOLATORS = [NumpyInterpolator(), FieldListInterpolator()]

earthkit/regrid/sphinxext/xref.py

@@ -0,0 +1,48 @@
+# (C) Copyright 2022 ECMWF.
+#
+# This software is licensed under the terms of the Apache Licence Version 2.0
+# which can be obtained at http://www.apache.org/licenses/LICENSE-2.0.
+# In applying this licence, ECMWF does not waive the privileges and immunities
+# granted to it by virtue of its status as an intergovernmental organisation
+# nor does it submit to any jurisdiction.
+#
+
+# Taken from: https://github.com/michaeljones/sphinx-xref
+
+from docutils import nodes
+from sphinx.util import caption_ref_re
+
+
+def xref(typ, rawtext, text, lineno, inliner, options={}, content=[]):
+    title = target = text
+
+    # look if explicit title and target are given with `foo <bar>` syntax
+    brace = text.find("<")
+    if brace != -1:
+        m = caption_ref_re.match(text)
+        if m:
+            target = m.group(2)
+            title = m.group(1)
+        else:
+            # fallback: everything after '<' is the target
+            target = text[brace + 1 :]
+            title = text[:brace]
+
+    link = xref.links[target]
+
+    if brace != -1:
+        pnode = nodes.reference(target, title, refuri=link[1])
+    else:
+        pnode = nodes.reference(target, link[0], refuri=link[1])
+
+    return [pnode], []
+
+
+def get_refs(app):
+    xref.links = app.config.xref_links
+
+
+def setup(app):
+    app.add_config_value("xref_links", {}, True)
+    app.add_role("xref", xref)
+    app.connect("builder-inited", get_refs)