Interpolator API

[VeckoTheGecko]

621a948 (#1946) introduced the interpolation template which allows users to pass in a function which performs the interpolation.

Conceptually, interpolation is the interface between the grid (the particles position in the grid), as well as the underlying dataset itself.

When the interpolator, the user needs to make assumptions about the underlying dataset (e.g., labels, dimension ordering etc), as well as the grid object itself. It would be good to have these assumptions being verified on Field initialisation so that things can fail early with an explainable error message, and so that people writing interpolators can document what is expected.

I suggest we implement the following API:

from typing import Any
from abc import ABC, abstractmethod


class Interpolator(ABC):
    """
    Base class for interpolators.
    """

    def __init__(self, *args, **kwargs):
        pass

    @abstractmethod
    def assert_compatible(data: Any, grid: Any) -> None: ...

    @abstractmethod
    def interpolate(self, data: Any, grid: Any, ti: int, ei: int, tau: float, t: float, z: float, y: float, x: float):
        # particle_position = grid.search(z, y, x, ei=ei)
        # xi, xsi = particle_position["X"]
        # yi, eta = particle_position["Y"]
        # zi, zeta = particle_position["Z"]
        # data[...]
        ...

Then, in field.py, roughly the following changes will take place:

diff --git a/parcels/field.py b/parcels/field.py
index 51ade1e2..cf5d78ed 100644
--- a/parcels/field.py
+++ b/parcels/field.py
@@ -101,52 +101,13 @@ class Field:
 
     """
 
-    @staticmethod
-    def _interp_template(
-        self,
-        ti: int,
-        ei: int,
-        bcoords: np.ndarray,
-        tau: np.float32 | np.float64,
-        t: np.float32 | np.float64,
-        z: np.float32 | np.float64,
-        y: np.float32 | np.float64,
-        x: np.float32 | np.float64,
-    ) -> np.float32 | np.float64:
-        """Template function used for the signature check of the lateral interpolation methods."""
-        return 0.0
-
-    def _validate_interp_function(self, func: Callable) -> bool:
-        """Ensures that the function has the correct signature."""
-        template_sig = inspect.signature(self._interp_template)
-        func_sig = inspect.signature(func)
-
-        if len(template_sig.parameters) != len(func_sig.parameters):
-            return False
-
-        for (_name1, param1), (_name2, param2) in zip(
-            template_sig.parameters.items(), func_sig.parameters.items(), strict=False
-        ):
-            if param1.kind != param2.kind:
-                return False
-            if param1.annotation != param2.annotation:
-                return False
-
-        return_annotation = func_sig.return_annotation
-        template_return = template_sig.return_annotation
-
-        if return_annotation != template_return:
-            return False
-
-        return True
-
     def __init__(
         self,
         name: str,
         data: xr.DataArray | ux.UxDataArray,
         grid: UxGrid | XGrid,
         mesh_type: Mesh = "flat",
-        interp_method: Callable | None = None,
+        interp_method: Interpolator | None = None,
     ):
         if not isinstance(data, (ux.UxDataArray, xr.DataArray)):
             raise ValueError(
@@ -189,7 +150,7 @@ class Field:
         if interp_method is None:
             self._interp_method = self._interp_template  # Default to method that returns 0 always
         else:
-            self._validate_interp_function(interp_method)
+            interp_method.assert_compatible(data, grid)
             self._interp_method = interp_method
 
         self.igrid = -1  # Default the grid index to -1

(_validate_interp_function wouldn't really be needed as this would be the purpose of the base class)

Thoughts @fluidnumerics-joe and @erikvansebille ? I think that this approach would be compatible with unstructured grids as well?

[fluidnumerics-joe]

I think this is a good way to clean things up. Though this may be implied already, want to make sure to point out the VectorField could also have an Interpolator.

Just making sure I understand that this would push the search method completely underneath the interpolate method ? Given our discussion in #2048, do we want to push time search under the grid and remove ti and tau ?

For users to provide their own interpolator, this would require that someone defines a type extension of the Interpolator class that ovverides the Interpolator.interpolate method. This is a little different than what we have currently, but I don't really see an issue.

[erikvansebille]

I agree that this seems a good way to move forward. So please go ahead and implement; perhaps considering moving ti and tau into the search, as @fluidnumerics-joe also suggests, to make the interpolate signature shorter/simpler

[VeckoTheGecko]

Just had a meeting with Erik, and we have decided to keep time out of the grid (i.e., also out of the search). If in the future we want to include it, we will scope out the implications that would have in a separate issue.