Feature/irregular grid mapper

.gitignore

@@ -13,7 +13,8 @@ site
 *.grib
 *.gif
 *.html
+*.req
+.mypy_cache
 example_eo
 example_mri
-.mypy_cache
-*.req
+polytope_python.egg-info

polytope/datacube/backends/datacube.py

@@ -127,11 +127,17 @@ def remap_path(self, path: DatacubePath):
         return path
 
     @staticmethod
-    def create(datacube, axis_options: dict, datacube_options={}):
+    def create(datacube, axis_options: dict, datacube_options={}, point_cloud_options=None):
         if isinstance(datacube, (xr.core.dataarray.DataArray, xr.core.dataset.Dataset)):
             from .xarray import XArrayDatacube
 
-            xadatacube = XArrayDatacube(datacube, axis_options, datacube_options)
+            xadatacube = XArrayDatacube(
+                datacube, axis_options, datacube_options, point_cloud_options=point_cloud_options
+            )
             return xadatacube
         else:
             return datacube
+
+    @abstractmethod
+    def find_point_cloud(self):
+        pass

polytope/datacube/backends/fdb.py

@@ -8,7 +8,7 @@
 
 
 class FDBDatacube(Datacube):
-    def __init__(self, config=None, axis_options=None, datacube_options=None):
+    def __init__(self, config=None, axis_options=None, datacube_options=None, point_cloud_options=None):
         if config is None:
             config = {}
 
@@ -17,6 +17,7 @@ def __init__(self, config=None, axis_options=None, datacube_options=None):
         logging.info("Created an FDB datacube with options: " + str(axis_options))
 
         self.unwanted_path = {}
+        self.has_point_cloud = point_cloud_options  # NOTE: here, will be True/False
 
         partial_request = config
         # Find values in the level 3 FDB datacube
@@ -43,6 +44,11 @@ def __init__(self, config=None, axis_options=None, datacube_options=None):
 
         logging.info("Polytope created axes for: " + str(self._axes.keys()))
 
+    def find_point_cloud(self):
+        # TODO: somehow, find the point cloud of irregular grid if it exists
+        if self.has_point_cloud:
+            return self.has_point_cloud
+
     def get(self, requests: IndexTree):
         fdb_requests = []
         fdb_requests_decoding_info = []
@@ -161,13 +167,15 @@ def get_2nd_last_values(self, requests, leaf_path=None):
 
         leaf_path_copy = deepcopy(leaf_path)
         leaf_path_copy.pop("values")
+        leaf_path_copy.pop("result")
         return (leaf_path_copy, range_lengths, current_start_idxs, fdb_node_ranges, lat_length)
 
     def get_last_layer_before_leaf(self, requests, leaf_path, range_l, current_idx, fdb_range_n):
         i = 0
         for c in requests.children:
             # now c are the leaves of the initial tree
             key_value_path = {c.axis.name: c.value}
+            leaf_path["result"] = c.result
             ax = c.axis
             (key_value_path, leaf_path, self.unwanted_path) = ax.unmap_path_key(
                 key_value_path, leaf_path, self.unwanted_path

polytope/datacube/backends/mock.py

@@ -67,3 +67,6 @@ def ax_vals(self, name):
 
     def _find_indexes_between(self, axis, indexes, low, up):
         pass
+
+    def find_point_cloud(self):
+        pass

polytope/datacube/backends/xarray.py

@@ -9,9 +9,10 @@
 class XArrayDatacube(Datacube):
     """Xarray arrays are labelled, axes can be defined as strings or integers (e.g. "time" or 0)."""
 
-    def __init__(self, dataarray: xr.DataArray, axis_options=None, datacube_options=None):
+    def __init__(self, dataarray: xr.DataArray, axis_options=None, datacube_options=None, point_cloud_options=None):
         super().__init__(axis_options, datacube_options)
         self.dataarray = dataarray
+        self.has_point_cloud = point_cloud_options
 
         for name, values in dataarray.coords.variables.items():
             if name in dataarray.dims:
@@ -37,24 +38,28 @@ def __init__(self, dataarray: xr.DataArray, axis_options=None, datacube_options=
                 val = self._axes[name].type
                 self._check_and_add_axes(options, name, val)
 
+    def find_point_cloud(self):
+        # TODO: somehow, find the point cloud of irregular grid if it exists
+        if self.has_point_cloud:
+            return self.has_point_cloud
+
     def get(self, requests: IndexTree):
         for r in requests.leaves:
-            path = r.flatten()
-            if len(path.items()) == self.axis_counter:
+            path = r.flatten_with_result()
+            if len(path.items()) == self.axis_counter + 1:
                 # first, find the grid mapper transform
                 unmapped_path = {}
                 path_copy = deepcopy(path)
                 for key in path_copy:
-                    axis = self._axes[key]
-                    key_value_path = {key: path_copy[key]}
-                    # (path, unmapped_path) = axis.unmap_to_datacube(path, unmapped_path)
-                    (key_value_path, path, unmapped_path) = axis.unmap_path_key(key_value_path, path, unmapped_path)
+                    if key != "result":
+                        axis = self._axes[key]
+                        key_value_path = {key: path_copy[key]}
+                        (key_value_path, path, unmapped_path) = axis.unmap_path_key(key_value_path, path, unmapped_path)
                 path.update(key_value_path)
                 path.update(unmapped_path)
 
                 unmapped_path = {}
                 self.refit_path(path, unmapped_path, path)
-
                 subxarray = self.dataarray.sel(path, method="nearest")
                 subxarray = subxarray.sel(unmapped_path)
                 value = subxarray.item()
@@ -83,9 +88,9 @@ def refit_path(self, path_copy, unmapped_path, path):
         for key in path.keys():
             if key not in self.dataarray.dims:
                 path_copy.pop(key)
-            if key not in self.dataarray.coords.dtypes:
+            elif key not in self.dataarray.coords.dtypes:
                 unmapped_path.update({key: path[key]})
-                path_copy.pop(key)
+                path_copy.pop(key, None)
             for key in self.dataarray.coords.dtypes:
                 key_dtype = self.dataarray.coords.dtypes[key]
                 if key_dtype.type is np.str_ and key in path.keys():

polytope/datacube/index_tree.py

@@ -54,13 +54,6 @@ def copy_children_from_other(self, other):
             c.copy_children_from_other(o)
         return
 
-    def pprint_2(self, level=0):
-        if self.axis.name == "root":
-            print("\n")
-        print("\t" * level + "\u21b3" + str(self))
-        for child in self.children:
-            child.pprint_2(level + 1)
-
     def _collect_leaf_nodes_old(self, leaves):
         if len(self.children) == 0:
             leaves.append(self)
@@ -191,6 +184,13 @@ def pprint(self, level=0):
         if len(self.children) == 0:
             logging.debug("\t" * (level + 1) + "\u21b3" + str(self.result))
 
+    def pprint_2(self, level=0):
+        if self.axis.name == "root":
+            print("\n")
+        print("\t" * level + "\u21b3" + str(self))
+        for child in self.children:
+            child.pprint_2(level + 1)
+
     def remove_branch(self):
         if not self.is_root():
             old_parent = self._parent
@@ -204,6 +204,20 @@ def flatten(self):
         ancestors = self.get_ancestors()
         for ancestor in ancestors:
             path[ancestor.axis.name] = ancestor.value
+        # add the result to the path
+        # NOTE: this is useful for the quadtree engine when we stash the point idx in the result
+        # if len(ancestors) != 0:
+        #     if ancestors[len(ancestors)-1].result is not None:
+        #         path["result"] = ancestors[len(ancestors)-1].result
+        return path
+
+    def flatten_with_result(self):
+        path = DatacubePath()
+        ancestors = self.get_ancestors()
+        for ancestor in ancestors:
+            path[ancestor.axis.name] = ancestor.value
+        if len(ancestors) != 0:
+            path["result"] = ancestors[-1].result
         return path
 
     def flatten_with_ancestors(self):