issue #884 Update for mask_scl_dilation function

Open-EO · Dec 11, 2024 · 1a8aa38 · 1a8aa38
1 parent f7ff0e3
commit 1a8aa38
Show file tree

Hide file tree

Showing 12 changed files with 162 additions and 48 deletions.
diff --git a/notebooks/AI4FOOD_Whittaker_OpenEO_NDVI_S2.ipynb b/notebooks/AI4FOOD_Whittaker_OpenEO_NDVI_S2.ipynb
@@ -332,4 +332,4 @@
  },
  "nbformat": 4,
  "nbformat_minor": 5
-}
+}
diff --git a/notebooks/Local/FuseTS - Peak Valley Detection.ipynb b/notebooks/Local/FuseTS - Peak Valley Detection.ipynb
@@ -57,7 +57,8 @@
     "More information on the usage of openEO's Python client can be found on [GitHub](https://github.com/Open-EO/openeo-python-client).\n",
     "\n",
     "The first step is to connect to an openEO compatible backend.\n"
-   ]
+   ],
+   "id": "ccc826d7e3e794ad"
   },
   {
    "cell_type": "code",
@@ -74,14 +75,16 @@
    ],
    "source": [
     "connection = openeo.connect(\"openeo.vito.be\").authenticate_oidc()"
-   ]
+   ],
+   "id": "e6ad5560fc9b30dc"
   },
   {
    "cell_type": "markdown",
    "metadata": {},
    "source": [
     "Next we define the area of interest, in this case an extent, for which we would like to fetch time series data.\n"
-   ]
+   ],
+   "id": "efbacd1f48fedc4a"
   },
   {
    "cell_type": "code",
@@ -98,14 +101,16 @@
     ")\n",
     "spat_ext = dict(west=minx, east=maxx, north=maxy, south=miny, crs=4326)\n",
     "temp_ext = [\"2021-01-01\", \"2021-12-31\"]"
-   ]
+   ],
+   "id": "10a0472a737b4773"
   },
   {
    "cell_type": "markdown",
    "metadata": {},
    "source": [
     "We will create an openEO process to calculate the NDVI time series for our area of interest. We'll begin by using the SENTINEL2_L2A_SENTINELHUB collection, and apply a cloud masking algorithm to remove any interfering clouds before calculating the NDVI values.\n"
-   ]
+   ],
+   "id": "df3ee7d1a70f77ee"
   },
   {
    "cell_type": "code",
@@ -122,14 +127,16 @@
     ")\n",
     "s2 = s2.process(\"mask_scl_dilation\", data=s2, scl_band_name=\"SCL\")\n",
     "ndvi_cube = s2.ndvi(red=\"B04\", nir=\"B08\", target_band=\"NDVI\")"
-   ]
+   ],
+   "id": "5d862b567ae6f299"
   },
   {
    "cell_type": "markdown",
    "metadata": {},
    "source": [
     "Now that we have calculated the NDVI time series for our area of interest, we can request openEO to download the result to our local storage. This will allow us to access the file and use it for further analysis in this notebook. However, if we have already downloaded the file, we can use the existing time series to continue our analysis without the need for a new download.\n"
-   ]
+   ],
+   "id": "d7c22a63e09838e8"
   },
   {
    "cell_type": "code",
@@ -174,7 +181,8 @@
    "metadata": {},
    "source": [
     "### Plot Image Chips Across Time\n"
-   ]
+   ],
+   "id": "63390c062b62fadc"
   },
   {
    "cell_type": "code",
@@ -209,14 +217,16 @@
     "    ax.set_title(img.t.dt.date.to_numpy())\n",
     "\n",
     "plt.tight_layout()"
-   ]
+   ],
+   "id": "efc69fe5bfb141e5"
   },
   {
    "cell_type": "markdown",
    "metadata": {},
    "source": [
     "### Plot Time Series for Specific Polygon Area\n"
-   ]
+   ],
+   "id": "a7ed4de53954b38f"
   },
   {
    "cell_type": "code",
@@ -258,7 +268,8 @@
     "ax.plot(*np.array(geometry[\"coordinates\"][0]).T, \"C3\", lw=2, ls=\"dashed\")\n",
     "ax.set_title(ndvi_filtered.isel(t=7).t.dt.date.to_numpy())\n",
     "ax.set_aspect(\"auto\")"
-   ]
+   ],
+   "id": "b09851b79009b1f3"
   },
   {
    "cell_type": "code",
@@ -299,7 +310,8 @@
     "ax.imshow(ndvi_masked.isel(t=7), extent=[minx, maxx, miny, maxy])\n",
     "ax.set_title(ndvi_masked.isel(t=7).t.dt.date.to_numpy())\n",
     "ax.set_aspect(\"auto\")"
-   ]
+   ],
+   "id": "26f11754efa62e6b"
   },
   {
    "cell_type": "code",
@@ -356,7 +368,8 @@
     "ax.plot([], [], \"C0\", label=\"NDVI series of random pixels\")\n",
     "ax.axvspan(None, None, color=\"C3\", alpha=0.5, label=\"Detected peak-valley events\")\n",
     "ax.legend()"
-   ]
+   ],
+   "id": "1546c3e816acc487"
   }
  ],
  "metadata": {

diff --git a/notebooks/OpenEO/FuseTS - MOGPR Multi Source Fusion.ipynb b/notebooks/OpenEO/FuseTS - MOGPR Multi Source Fusion.ipynb
@@ -67,8 +67,8 @@
       "Requirement already satisfied: packaging>=21.3 in /Users/bramjanssen/projects/vito/FuseTS/venv_clean_v2/lib/python3.8/site-packages (from xarray>=0.12.3->openeo) (23.1)\n",
       "Requirement already satisfied: six>=1.5 in /Users/bramjanssen/projects/vito/FuseTS/venv_clean_v2/lib/python3.8/site-packages (from python-dateutil>=2.8.2->pandas>0.20.0->openeo) (1.16.0)\n",
       "\n",
-      "\u001b[1m[\u001b[0m\u001b[34;49mnotice\u001b[0m\u001b[1;39;49m]\u001b[0m\u001b[39;49m A new release of pip available: \u001b[0m\u001b[31;49m22.3.1\u001b[0m\u001b[39;49m -> \u001b[0m\u001b[32;49m24.0\u001b[0m\n",
-      "\u001b[1m[\u001b[0m\u001b[34;49mnotice\u001b[0m\u001b[1;39;49m]\u001b[0m\u001b[39;49m To update, run: \u001b[0m\u001b[32;49mpip install --upgrade pip\u001b[0m\n"
+      "\u001B[1m[\u001B[0m\u001B[34;49mnotice\u001B[0m\u001B[1;39;49m]\u001B[0m\u001B[39;49m A new release of pip available: \u001B[0m\u001B[31;49m22.3.1\u001B[0m\u001B[39;49m -> \u001B[0m\u001B[32;49m24.0\u001B[0m\n",
+      "\u001B[1m[\u001B[0m\u001B[34;49mnotice\u001B[0m\u001B[1;39;49m]\u001B[0m\u001B[39;49m To update, run: \u001B[0m\u001B[32;49mpip install --upgrade pip\u001B[0m\n"
      ]
     }
    ],

diff --git a/src/fusets/openeo/__init__.py b/src/fusets/openeo/__init__.py
@@ -72,10 +72,20 @@ def load_xarray(collection_id, spatial_extent, temporal_extent, properties=None,
     if openeo_connection == None:
         openeo_connection = openeo.connect("openeo.cloud").authenticate_oidc()
     print(spatial_extent)
+    scl = openeo_connection.load_collection(
+        collection_id, temporal_extent=temporal_extent, bands=["SCL"]
+    ).filter_bbox(spatial_extent)
+    cloud_mask = scl.process(
+        "to_scl_dilation_mask",
+        data=scl,
+        kernel1_size=17, kernel2_size=77,
+        mask1_values=[2, 4, 5, 6, 7],
+        mask2_values=[3, 8, 9, 10, 11],
+        erosion_kernel_size=3)
     data = openeo_connection.load_collection(
-        collection_id, temporal_extent=temporal_extent, bands=["B02", "B03", "B04", "SCL"]
+        collection_id, temporal_extent=temporal_extent, bands=["B02", "B03", "B04"]
     ).filter_bbox(spatial_extent)
-    data = data.process("mask_scl_dilation", data=data, scl_band_name="SCL")
+    data = data.mask(cloud_mask)
     job = data.execute_batch(out_format="netCDF")
     results = job.get_results()
     base_path = Path(job.job_id)

diff --git a/src/fusets/openeo/services/descriptions/phenology.md b/src/fusets/openeo/services/descriptions/phenology.md
@@ -60,7 +60,7 @@ s2 = connection.load_collection('SENTINEL2_L2A_SENTINELHUB',
                                 spatial_extent=spat_ext,
                                 temporal_extent=temp_ext,
                                 bands=["B04", "B08", "SCL"])
-s2 = s2.process("mask_scl_dilation", data=s2, scl_band_name="SCL")
+s2 = s2.process("mask_scl_dilation", data=s2, scl_band_name="SCL") # do md file also need updating?
 s2 = s2.mask_polygon(spat_ext)
 base_ndvi = s2.ndvi(red="B04", nir="B08")
 

diff --git a/src/fusets/openeo/services/publish_mogpr.py b/src/fusets/openeo/services/publish_mogpr.py
@@ -29,10 +29,20 @@ def execute_udf():
     temp_ext = ["2023-01-01", "2023-03-31"]
 
     # Setup NDVI cube
+    scl = connection.load_collection(
+        "SENTINEL2_L2A", spatial_extent=spat_ext, temporal_extent=temp_ext, bands=["SCL"]
+    )
+    cloud_mask = scl.process(
+        "to_scl_dilation_mask",
+        data=scl,
+        kernel1_size=17, kernel2_size=77,
+        mask1_values=[2, 4, 5, 6, 7],
+        mask2_values=[3, 8, 9, 10, 11],
+        erosion_kernel_size=3)
     base_s2 = connection.load_collection(
-        "SENTINEL2_L2A", spatial_extent=spat_ext, temporal_extent=temp_ext, bands=["B04", "B08", "SCL"]
+        "SENTINEL2_L2A", spatial_extent=spat_ext, temporal_extent=temp_ext, bands=["B04", "B08"]
     )
-    base_s2 = base_s2.process("mask_scl_dilation", data=base_s2, scl_band_name="SCL")
+    base_s2 = base_s2.mask(cloud_mask)
     base_s2 = base_s2.ndvi(red="B04", nir="B08", target_band="NDVI")
     base_s2 = base_s2.filter_bands(bands=["NDVI"])
     base_s2 = base_s2.mask_polygon(spat_ext)

diff --git a/src/fusets/openeo/services/publish_mogpr_s1_s2.py b/src/fusets/openeo/services/publish_mogpr_s1_s2.py
@@ -151,10 +151,20 @@ def _load_ndvi(connection, polygon, date):
     :param date:
     :return:
     """
+    scl = connection.load_collection(
+        "SENTINEL2_L2A", spatial_extent=polygon, temporal_extent=date, bands=["SCL"]
+    )
+    cloud_mask = scl.process(
+        "to_scl_dilation_mask",
+        data=scl,
+        kernel1_size=17, kernel2_size=77,
+        mask1_values=[2, 4, 5, 6, 7],
+        mask2_values=[3, 8, 9, 10, 11],
+        erosion_kernel_size=3)
     base_s2 = connection.load_collection(
-        "SENTINEL2_L2A", spatial_extent=polygon, temporal_extent=date, bands=["B04", "B08", "SCL"]
+        "SENTINEL2_L2A", spatial_extent=polygon, temporal_extent=date, bands=["B04", "B08"]
     )
-    base_s2 = base_s2.process("mask_scl_dilation", data=base_s2, scl_band_name="SCL")
+    base_s2 = base_s2.mask(cloud_mask)
     ndvi = base_s2.ndvi(red="B04", nir="B08", target_band="NDVI")
     ndvi = ndvi.filter_bands(bands=["NDVI"])
     return ndvi.mask_polygon(polygon)
@@ -182,13 +192,26 @@ def _load_evi(connection, polygon, date):
     :param date: Time of interest
     :return:
     """
+    scl = connection.load_collection(
+        collection_id="SENTINEL2_L2A",
+        spatial_extent=polygon,
+        temporal_extent=date,
+        bands=["SCL"],
+    )
+    cloud_mask = scl.process(
+        "to_scl_dilation_mask",
+        data=scl,
+        kernel1_size=17, kernel2_size=77,
+        mask1_values=[2, 4, 5, 6, 7],
+        mask2_values=[3, 8, 9, 10, 11],
+        erosion_kernel_size=3)
     base_s2 = connection.load_collection(
         collection_id="SENTINEL2_L2A",
         spatial_extent=polygon,
         temporal_extent=date,
-        bands=["B02", "B04", "B08", "SCL"],
+        bands=["B02", "B04", "B08"],
     )
-    base_s2 = base_s2.process("mask_scl_dilation", data=base_s2, scl_band_name="SCL")
+    base_s2 = base_s2.mask(cloud_mask)
 
     B02 = base_s2.band("B04")
     B04 = base_s2.band("B04")

diff --git a/src/fusets/openeo/services/publish_phenology.py b/src/fusets/openeo/services/publish_phenology.py
@@ -50,10 +50,20 @@ def test_udf():
     }
     temp_ext = ["2022-01-01", "2022-12-31"]
     smoothing_lambda = 10000
+    scl = connection.load_collection(
+        "SENTINEL2_L2A_SENTINELHUB", spatial_extent=spat_ext, temporal_extent=temp_ext, bands=["SCL"]
+    )
+    cloud_mask = scl.process(
+        "to_scl_dilation_mask",
+        data=scl,
+        kernel1_size=17, kernel2_size=77,
+        mask1_values=[2, 4, 5, 6, 7],
+        mask2_values=[3, 8, 9, 10, 11],
+        erosion_kernel_size=3)
     base = connection.load_collection(
-        "SENTINEL2_L2A_SENTINELHUB", spatial_extent=spat_ext, temporal_extent=temp_ext, bands=["B04", "B08", "SCL"]
+        "SENTINEL2_L2A_SENTINELHUB", spatial_extent=spat_ext, temporal_extent=temp_ext, bands=["B04", "B08"]
     )
-    base_cloudmasked = base.process("mask_scl_dilation", data=base, scl_band_name="SCL")
+    base_cloudmasked = base.mask(cloud_mask)
     base_ndvi = base_cloudmasked.ndvi(red="B04", nir="B08")
     phenology = base_ndvi.apply_dimension(
         process=lambda x: run_udf(x, udf=load_phenology_udf(), runtime="Python"),

diff --git a/tests/conftest.py b/tests/conftest.py
@@ -112,11 +112,20 @@ def wetland_sentinel2_ndvi(areas):
     import openeo
 
     openeo_connection = openeo.connect("openeo-dev.vito.be").authenticate_oidc()
-
+    scl = openeo_connection.load_collection(
+        "SENTINEL2_L2A", temporal_extent=("2020-01-01", "2021-01-01"), bands=["SCL"]
+    ).filter_bbox(areas["wetland"])
+    cloud_mask = scl.process(
+        "to_scl_dilation_mask",
+        data=scl,
+        kernel1_size=17, kernel2_size=77,
+        mask1_values=[2, 4, 5, 6, 7],
+        mask2_values=[3, 8, 9, 10, 11],
+        erosion_kernel_size=3)
     data = openeo_connection.load_collection(
-        "SENTINEL2_L2A", temporal_extent=("2020-01-01", "2021-01-01"), bands=["B08", "B04", "SCL"]
+        "SENTINEL2_L2A", temporal_extent=("2020-01-01", "2021-01-01"), bands=["B08", "B04"]
     ).filter_bbox(areas["wetland"])
-    data = data.process("mask_scl_dilation", data=data, scl_band_name="SCL")
+    data = data.mask(cloud_mask)
 
     return data.ndvi(nir="B08", red="B04")
 

diff --git a/tests/fusets_openeo_tests/test_performance.py b/tests/fusets_openeo_tests/test_performance.py
@@ -41,13 +41,26 @@ def start_job(data, context: dict, **kwargs) -> openeo.BatchJob:
             **context["jobinfo"],
         )
     else:
+        scl = connection.load_collection(
+            "SENTINEL2_L2A",
+            spatial_extent=aoi,
+            temporal_extent=context["params"]["temp-ext"],
+            bands=["SCL"],
+        )
+        cloud_mask = scl.process(
+            "to_scl_dilation_mask",
+            data=scl,
+            kernel1_size=17, kernel2_size=77,
+            mask1_values=[2, 4, 5, 6, 7],
+            mask2_values=[3, 8, 9, 10, 11],
+            erosion_kernel_size=3)
         base = connection.load_collection(
             "SENTINEL2_L2A",
             spatial_extent=aoi,
             temporal_extent=context["params"]["temp-ext"],
-            bands=["B04", "B08", "SCL"],
+            bands=["B04", "B08"],
         )
-        base_cloudmasked = base.process("mask_scl_dilation", data=base, scl_band_name="SCL")
+        base_cloudmasked = base.mask(cloud_mask)
         base_ndvi = base_cloudmasked.ndvi(red="B04", nir="B08")
         service_dc = connection.datacube_from_process(data=base_ndvi, **context["jobinfo"])
-Original file line number
+Diff line change
@@ Expand Up / @@ -332,4 +332,4 @@ @@
      },
      "nbformat": 4,
      "nbformat_minor": 5
-    }
+    }