Fix path feature

ecmwf · Aug 8, 2024 · 86ba62e · 86ba62e
1 parent 2d17f14
commit 86ba62e
Show file tree

Hide file tree

Showing 2 changed files with 67 additions and 11 deletions.
diff --git a/covjsonkit/decoder/Path.py b/covjsonkit/decoder/Path.py
@@ -36,21 +36,64 @@ def to_geopandas(self):
         pass
 
     def to_xarray(self):
-        dims = ["points"]
+        dims = ["datetimes", "number", "steps", "points"]
         dataarraydict = {}
 
         # Get coordinates
         x = []
         y = []
         t = []
         for coord in self.get_coordinates()["composite"]["values"]:
-            x.append(float(coord[0]))
-            y.append(float(coord[1]))
+            x.append(float(coord[1]))
+            y.append(float(coord[2]))
             t.append(coord[0])
 
-        # Get values
+        values = {}
         for parameter in self.parameters:
-            dataarray = xr.DataArray(self.get_values()[parameter][0], dims=dims)
+            values[parameter] = {}
+
+        datetimes = []
+        numbers = []
+        steps = []
+        for coverage in self.coverages:
+            if "number" not in coverage["mars:metadata"]:
+                coverage["mars:metadata"]["number"] = 0
+            numbers.append(coverage["mars:metadata"]["number"])
+            if "step" not in coverage["mars:metadata"]:
+                coverage["mars:metadata"]["step"] = 0
+            steps.append(coverage["mars:metadata"]["step"])
+            datetimes.append(coverage["mars:metadata"]["Forecast date"])
+            for parameter in self.parameters:
+                # values[parameter].append(coverage["ranges"][parameter]["values"])
+                if coverage["mars:metadata"]["Forecast date"] not in values[parameter]:
+                    values[parameter][coverage["mars:metadata"]["Forecast date"]] = {}
+                if (
+                    coverage["mars:metadata"]["number"]
+                    not in values[parameter][coverage["mars:metadata"]["Forecast date"]]
+                ):
+                    values[parameter][coverage["mars:metadata"]["Forecast date"]][
+                        coverage["mars:metadata"]["number"]
+                    ] = {}
+                values[parameter][coverage["mars:metadata"]["Forecast date"]][coverage["mars:metadata"]["number"]][
+                    coverage["mars:metadata"]["step"]
+                ] = coverage["ranges"][parameter]["values"]
+
+        datetimes = list(set(datetimes))
+        numbers = list(set(numbers))
+        steps = list(set(steps))
+
+        new_values = {}
+        for parameter in values.keys():
+            new_values[parameter] = []
+            for i, datetime in enumerate(datetimes):
+                new_values[parameter].append([])
+                for j, number in enumerate(numbers):
+                    new_values[parameter][i].append([])
+                    for k, step in enumerate(steps):
+                        new_values[parameter][i][j].append(values[parameter][datetime][number][step])
+
+        for parameter in self.parameters:
+            dataarray = xr.DataArray(new_values[parameter], dims=dims)
             dataarray.attrs["type"] = self.get_parameter_metadata(parameter)["type"]
             dataarray.attrs["units"] = self.get_parameter_metadata(parameter)["unit"]["symbol"]
             dataarray.attrs["long_name"] = self.get_parameter_metadata(parameter)["observedProperty"]["id"]
@@ -59,10 +102,19 @@ def to_xarray(self):
         ds = xr.Dataset(
             dataarraydict,
             coords=dict(
-                points=(["points"], list(range(0, len(x)))), x=(["points"], x), y=(["points"], y), t=(["points"], t)
+                datetimes=(["datetimes"], datetimes),
+                number=(["number"], numbers),
+                steps=(["steps"], steps),
+                points=(["points"], list(range(0, len(x)))),
+                x=(["points"], x),
+                y=(["points"], y),
+                t=(["points"], t),
             ),
         )
         for mars_metadata in self.mars_metadata[0]:
             ds.attrs[mars_metadata] = self.mars_metadata[0][mars_metadata]
 
+        # Add date attribute
+        # ds.attrs["date"] = self.get_coordinates()["t"]["values"][0]
+
         return ds
diff --git a/covjsonkit/encoder/Path.py b/covjsonkit/encoder/Path.py
@@ -25,7 +25,7 @@ def add_domain(self, coverage, coords):
         coverage["domain"]["axes"] = {}
         coverage["domain"]["axes"]["composite"] = {}
         coverage["domain"]["axes"]["composite"]["dataType"] = "tuple"
-        coverage["domain"]["axes"]["composite"]["coordinates"] = self.covjson["referencing"]
+        coverage["domain"]["axes"]["composite"]["coordinates"] = self.covjson["referencing"][0]["coordinates"]
         coverage["domain"]["axes"]["composite"]["values"] = coords["composite"]
 
     def add_range(self, coverage, values):
@@ -50,7 +50,7 @@ def from_xarray(self, dataset):
 
         self.add_reference(
             {
-                "coordinates": ["x", "y", "z"],
+                "coordinates": ["t", "x", "y", "z"],
                 "system": {
                     "type": "GeographicCRS",
                     "id": "http://www.opengis.net/def/crs/OGC/1.3/CRS84",
@@ -89,17 +89,21 @@ def from_polytope(self, result):
 
         self.add_reference(
             {
-                "coordinates": ["x", "y", "z"],
+                "coordinates": ["t", "x", "y", "z"],
                 "system": {
                     "type": "GeographicCRS",
                     "id": "http://www.opengis.net/def/crs/OGC/1.3/CRS84",
                 },
             }
         )
-        print(coords)
-        print(range_dict)
 
         for date in range_dict.keys():
+            new_coords = []
+            for val in coords[date]["composite"]:
+                val.insert(0, date)
+                new_coords.append(val)
+            coords[date]["composite"] = new_coords
+
             for num in range_dict[date].keys():
                 val_dict = {}
                 for step in range_dict[date][num][self.parameters[0]].keys():