Release 1.0.3

.github/workflows/ci.yaml

@@ -1,25 +1,19 @@
 name: ci
-
 on:
   # Trigger the workflow on push to master or develop, except tag creation
   push:
     branches:
       - 'main'
       - 'develop'
-
   # Trigger the workflow on pull request
   pull_request: ~
-
   # Trigger the workflow manually
   workflow_dispatch: ~
-
   # Trigger after public PR approved for CI
   pull_request_target:
     types: [labeled]
-
   release:
     types: [created]
-
 jobs:
   qa:
     name: qa
@@ -48,7 +42,6 @@ jobs:
 
       - name: Check flake8
         run: flake8 .
-
   setup:
     name: setup
     runs-on: ubuntu-20.04
@@ -95,7 +88,6 @@ jobs:
         run: |
           echo inputs=$(echo "${{ inputs.build_package_inputs || '{}' }}" | yq eval '.' --output-format json --indent 0 -) >> $GITHUB_OUTPUT
           echo inputs-for-ubuntu=$(echo "${{ inputs.build_package_inputs || '{}' }}" | yq eval '. * {"os":"ubuntu-20.04","compiler":"gnu-10","compiler_cc":"gcc-10","compiler_cxx":"g++-10","compiler_fc":"gfortran-10"}' --output-format json --indent 0 -) >> $GITHUB_OUTPUT
-
   test:
     name: test
     needs:
@@ -106,20 +98,16 @@ jobs:
       matrix: ${{ fromJson(needs.setup.outputs.matrix) }}
     runs-on: ${{ matrix.os }}
     steps:
-      - uses: actions/checkout@v3
+      - uses: actions/checkout@v4
+      - name: Install eccodes and Dependencies
+        id: install-dependencies
+        uses: ecmwf-actions/build-package@v2
         with:
-          fetch-depth: 3 # enable full clone to download Git LFS files
-      - name: Change write permission
-        run: chmod -R 777 ./tests/data/
-      - name: Install Git LFS
-        run: |
-          sudo apt-get install git-lfs
-          git lfs install
-      - name: Download Git LFS file
-        run: |
-          git lfs pull
-      - name: Install eccodes
-        run: sudo apt-get install -y libeccodes-dev
+          self_build: false
+          dependencies: |
+            ecmwf/ecbuild@develop
+            MathisRosenhauer/libaec@master
+            ecmwf/eccodes@develop
 
       - name: Setup Python
         uses: actions/setup-python@v4
@@ -143,23 +131,18 @@ jobs:
           LD_LIBRARY_PATH: ${{ steps.install-dependencies.outputs.lib_path }}
         shell: bash -eux {0}
         run: |
-          DYLD_LIBRARY_PATH=${{ env.LD_LIBRARY_PATH }} python -m pytest tests --cov=./ --cov-report=xml
+          DYLD_LIBRARY_PATH=${{ env.LD_LIBRARY_PATH }} python -m pytest -m "not fdb" tests --cov=./ --cov-report=xml
           python -m coverage report
 
       - name: Upload coverage to Codecov
         uses: codecov/codecov-action@v3
         with:
           files: coverage.xml
-
   deploy:
     needs: test
-
     if: ${{ github.event_name == 'release' }}
-
     name: Upload to Pypi
-
     runs-on: ubuntu-latest
-
     steps:
       - uses: actions/checkout@v3
       - name: Set up Python
@@ -176,4 +159,4 @@ jobs:
           TWINE_PASSWORD: ${{ secrets.PYPI_PASSWORD }}
         run: |
           python setup.py sdist
-          twine upload dist/*
+          twine upload dist/*

.gitignore

@@ -5,4 +5,15 @@
 polytope.egg-info
 .pytest_cache
 *.prof
-*.idx
+*.idx
+*.grib
+*.xml
+site
+.coverage
+*.grib
+*.gif
+*.html
+example_eo
+example_mri
+.mypy_cache
+*.req

ACKNOWLEDGEMENTS.rst

@@ -0,0 +1,8 @@
+Acknowledgements
+================
+
+Destination Earth
+    This software is developed with co-funding by the European Union under the Destination Earth initiative.
+
+LEXIS 
+    This software is part of a project that has received funding from the European Union’s Horizon 2020 research and innovation programme under grant agreement No 825532.

CONTRIBUTING.rst

@@ -0,0 +1,20 @@
+Contributing to *polytope*
+========================================
+
+Please report bug_ reports or pull-requests_ on GitHub_
+
+.. _bug: https://github.com/ecmwf/polytope/issues
+
+.. _pull-requests: https://github.com/ecmwf/polytope/pulls
+
+.. _GitHub: https://github.com/ecmwf/polytope
+
+We want your feedback, please e-mail: user-services@ecmwf.int
+
+The package is installed from PyPI with::
+
+    $ pip install polytope-python
+
+How you could use polytope:
+ * plot outputs, examples are provided in the examples folder
+ * go further than our examples and let us know how it goes

MANIFEST.in

@@ -0,0 +1 @@
+include requirements.txt

examples/3D_shipping_route.py

@@ -8,14 +8,14 @@
 
 from polytope.engine.hullslicer import HullSlicer
 from polytope.polytope import Polytope, Request
-from polytope.shapes import Ellipsoid, Path
+from polytope.shapes import Ellipsoid, Path, Select
 
 
 class Test:
     def setup_method(self):
         ds = data.from_source("file", "./examples/data/winds.grib")
         array = ds.to_xarray()
-        array = array.isel(time=0).isel(surface=0).isel(number=0)
+        array = array.isel(time=0).isel(surface=0).isel(number=0).u10
         self.array = array
         self.slicer = HullSlicer()
         self.API = Polytope(datacube=array, engine=self.slicer)
@@ -58,29 +58,25 @@ def test_slice_shipping_route(self):
         # Then somehow make this list of points into just a sequence of points
 
         ship_route_polytope = Path(["latitude", "longitude", "step"], initial_shape, *new_points)
-        request = Request(ship_route_polytope)
+        request = Request(
+            ship_route_polytope, Select("number", [0]), Select("surface", [0]), Select("time", ["2022-09-30T12:00:00"])
+        )
         result = self.API.retrieve(request)
 
         # Associate the results to the lat/long points in an array
         lats = []
         longs = []
         parameter_values = []
-        winds_u = []
-        winds_v = []
         for i in range(len(result.leaves)):
             cubepath = result.leaves[i].flatten()
             lat = cubepath["latitude"]
             long = cubepath["longitude"]
             lats.append(lat)
             longs.append(long)
 
-            u10_idx = result.leaves[i].result["u10"]
+            u10_idx = result.leaves[i].result[1]
             wind_u = u10_idx
-            v10_idx = result.leaves[i].result["v10"]
-            wind_v = v10_idx
-            winds_u.append(wind_u)
-            winds_v.append(wind_v)
-            parameter_values.append(math.sqrt(wind_u**2 + wind_v**2))
+            parameter_values.append(wind_u)
 
         parameter_values = np.array(parameter_values)
         # Plot this last array according to different colors for the result on a world map

examples/3D_shipping_route_wave_model.py

@@ -0,0 +1,131 @@
+import geopandas as gpd
+import matplotlib.pyplot as plt
+import pandas as pd
+import pytest
+from eccodes import codes_grib_find_nearest, codes_grib_new_from_file
+
+from polytope.datacube.backends.fdb import FDBDatacube
+from polytope.engine.hullslicer import HullSlicer
+from polytope.polytope import Polytope, Request
+from polytope.shapes import Ellipsoid, Path, Select
+from tests.helper_functions import download_test_data
+
+
+class TestReducedLatLonGrid:
+    def setup_method(self, method):
+        nexus_url = "https://get.ecmwf.int/test-data/polytope/test-data/wave.grib"
+        download_test_data(nexus_url, "wave.grib")
+        self.options = {
+            "values": {
+                "transformation": {
+                    "mapper": {"type": "reduced_ll", "resolution": 1441, "axes": ["latitude", "longitude"]}
+                }
+            },
+            "date": {"transformation": {"merge": {"with": "time", "linkers": ["T", "00"]}}},
+            "step": {"transformation": {"type_change": "int"}},
+            "number": {"transformation": {"type_change": "int"}},
+            "longitude": {"transformation": {"cyclic": [0, 360]}},
+        }
+        self.config = {"class": "od", "stream": "wave"}
+        self.fdbdatacube = FDBDatacube(self.config, axis_options=self.options)
+        self.slicer = HullSlicer()
+        self.API = Polytope(datacube=self.fdbdatacube, engine=self.slicer, axis_options=self.options)
+
+    def find_nearest_latlon(self, grib_file, target_lat, target_lon):
+        messages = grib_file
+
+        # Find the nearest grid points
+        nearest_points = []
+        for message in [messages[0]]:
+            nearest_index = codes_grib_find_nearest(message, target_lat, target_lon)
+            nearest_points.append(nearest_index)
+
+        return nearest_points
+
+    @pytest.mark.internet
+    @pytest.mark.skip(reason="can't install fdb branch on CI")
+    def test_reduced_ll_grid(self):
+        shapefile = gpd.read_file("./examples/data/Shipping-Lanes-v1.shp")
+        geometry_multiline = shapefile.iloc[2]
+        geometry_object = geometry_multiline["geometry"]
+
+        lines = []
+        i = 0
+
+        for line in geometry_object[:7]:
+            for point in line.coords:
+                point_list = list(point)
+                if list(point)[0] < 0:
+                    point_list[0] = list(point)[0] + 360
+                lines.append(point_list)
+
+        # Append for each point a corresponding step
+
+        new_points = []
+        for point in lines[:7]:
+            new_points.append([point[1], point[0], 1])
+
+        # Pad the shipping route with an initial shape
+
+        padded_point_upper = [0.24, 0.24, 1]
+        padded_point_lower = [-0.24, -0.24, 1]
+        initial_shape = Ellipsoid(["latitude", "longitude", "step"], padded_point_lower, padded_point_upper)
+
+        # Then somehow make this list of points into just a sequence of points
+
+        ship_route_polytope = Path(["latitude", "longitude", "step"], initial_shape, *new_points)
+
+        request = Request(
+            ship_route_polytope,
+            Select("date", [pd.Timestamp("20231129T000000")]),
+            Select("domain", ["g"]),
+            Select("expver", ["0001"]),
+            Select("param", ["140251"]),
+            Select("direction", ["1"]),
+            Select("frequency", ["1"]),
+            Select("class", ["od"]),
+            Select("stream", ["wave"]),
+            Select("levtype", ["sfc"]),
+            Select("type", ["fc"]),
+        )
+        result = self.API.retrieve(request)
+        result.pprint()
+
+        lats = []
+        lons = []
+        eccodes_lats = []
+        eccodes_lons = []
+        tol = 1e-8
+        f = open("./tests/data/wave.grib", "rb")
+        messages = []
+        message = codes_grib_new_from_file(f)
+        messages.append(message)
+
+        leaves = result.leaves
+        for i in range(len(leaves)):
+            cubepath = leaves[i].flatten()
+            lat = cubepath["latitude"]
+            lon = cubepath["longitude"]
+            del cubepath
+            lats.append(lat)
+            lons.append(lon)
+            nearest_points = codes_grib_find_nearest(message, lat, lon)[0]
+            eccodes_lat = nearest_points.lat
+            eccodes_lon = nearest_points.lon
+            eccodes_lats.append(eccodes_lat)
+            eccodes_lons.append(eccodes_lon)
+            assert eccodes_lat - tol <= lat
+            assert lat <= eccodes_lat + tol
+            assert eccodes_lon - tol <= lon
+            assert lon <= eccodes_lon + tol
+            print(i)
+        f.close()
+
+        worldmap = gpd.read_file(gpd.datasets.get_path("naturalearth_lowres"))
+        fig, ax = plt.subplots(figsize=(12, 6))
+        worldmap.plot(color="darkgrey", ax=ax)
+
+        plt.scatter(lons, lats, s=18, c="red", cmap="YlOrRd")
+        plt.scatter(eccodes_lons, eccodes_lats, s=6, c="green")
+        plt.colorbar(label="Temperature")
+        plt.show()

examples/4D_flight_path.py

@@ -4,24 +4,24 @@
 from earthkit import data
 from PIL import Image
 
-from polytope.datacube.xarray import XArrayDatacube
+from polytope.datacube.backends.xarray import XArrayDatacube
 from polytope.engine.hullslicer import HullSlicer
 from polytope.polytope import Polytope, Request
-from polytope.shapes import Box, Path
+from polytope.shapes import Box, Path, Select
 
 
 class Test:
     def setup_method(self):
         ds = data.from_source("file", "./examples/data/temp_model_levels.grib")
         array = ds.to_xarray()
-        array = array.isel(time=0)
-        options = {"longitude": {"Cyclic": [0, 360.0]}}
+        array = array.isel(time=0).t
+        axis_options = {"longitude": {"cyclic": [0, 360.0]}}
         self.xarraydatacube = XArrayDatacube(array)
         for dim in array.dims:
             array = array.sortby(dim)
         self.array = array
         self.slicer = HullSlicer()
-        self.API = Polytope(datacube=array, engine=self.slicer, options=options)
+        self.API = Polytope(datacube=array, engine=self.slicer, axis_options=axis_options)
 
     def test_slice_shipping_route(self):
         colorscale = [
@@ -84,7 +84,7 @@ def sphere(size, texture):
 
         flight_route_polytope = Path(["latitude", "longitude", "step", "hybrid"], initial_shape, *route_point_CDG_LHR)
 
-        request = Request(flight_route_polytope)
+        request = Request(flight_route_polytope, Select("time", ["2022-12-02T12:00:00"]))
         result = self.API.retrieve(request)
 
         lats = []
@@ -99,7 +99,7 @@ def sphere(size, texture):
             lats.append(lat)
             longs.append(long)
             levels.append(level)
-            t_idx = result.leaves[i].result["t"]
+            t_idx = result.leaves[i].result[1]
             parameter_values.append(t_idx)
         parameter_values = np.array(parameter_values)