Add support for get_NOBSperc flag for calculating observation percentage and refactor job options and preprocessing functions for enhanced flexibility.

Author: mbuchhorn

README.md

@@ -81,12 +81,20 @@ SLC_masking_algo: str = 'mask_scl_dilation'
         Masking method for Sentinel-2 optical data ('satio', 'mask_scl_dilation', None)
         Note: if set to None, no masking is applied and the S2 L2A data is used as is.
 
-apply_cloud_mask : bool = True
+apply_cloud_mask: bool = True
         if True, the Sentinel-2 or PlanetScope data is masked for clouds (based on Sentinel-2 QA band). 
-        If False, no masking is applied but the mask band is still created and added to the cube.
-        Note: no effect when 'mask_scl_dilation' parameter is set to None.     
-
-append : bool = True
+        If False, no masking is applied but the mask band is still created and added to the cube. (Note: the cloud 
+        mask band is removed from the feature cube as soon the time domain is aggregated.)
+        Note: no effect when 'mask_scl_dilation' parameter is set to None.
+
+get_NOBSperc: bool = False 
+        if True, the number of observations per pixel is calculated and added to the feature cube.
+        Note: this flag is only relevant for Sentinel-2 data as well as the 'band' NOBSperc is only added 
+        to cubes with time domain aggregation. 
+        Note: 'apply_cloud_mask' has to set to 'True' AND 'SLC_masking_algo' has to be 'mask_scl_dilation' 
+              for this flag to work.
+
+append: bool = True
         if the VI's are appended to the reflectance/radar time series cube OR replace them
 
 S2_scaling: list = [0, 10000, 0, 1.0]

src/eo_processing/config/settings.py

@@ -20,6 +20,7 @@
 S2_MAX_CLOUD_COVER: int = 95
 MASKING_ALGO: str = 'mask_scl_dilation'
 APPLY_CLOUD_MASK: bool = True
+GET_NOBS: bool = False
 S2_TILEID_LIST: Optional[List[str]] = None
 SKIP_CHECK_S1: bool = False
 SKIP_CHECK_S2: bool = False
@@ -68,7 +69,7 @@
 
 # ---------------------------------------------------
 # Job options for OpenEO
-OPENEO_EXTRACT_JOB_OPTIONS: Dict[str, str] = {
+OPENEO_EXTRACT_JOB_OPTIONS: Dict[str, Union[str, int, float, List]] = {
     "driver-memory": "4G",
     "driver-memoryOverhead": "8G",
     "driver-cores": 2,
@@ -80,7 +81,7 @@
     "stac-version":"1.1"
 }
 
-OPENEO_EXTRACT_CREO_JOB_OPTIONS: Dict[str, str] = {
+OPENEO_EXTRACT_CREO_JOB_OPTIONS: Dict[str, Union[str, int, float, List]] = {
     "driver-memory": "4G",
     "driver-memoryOverhead": "2G",
     "driver-cores": 1,
@@ -91,7 +92,7 @@
     "max-executors": 200
 }
 
-OPENEO_EXTRACT_CDSE_JOB_OPTIONS: Dict[str, str] = {
+OPENEO_EXTRACT_CDSE_JOB_OPTIONS: Dict[str, Union[str, int, float, List]] = {
     "driver-memory": "8G",
     "driver-memoryOverhead": "5G",
     "driver-cores": 1,
@@ -103,7 +104,7 @@
     "logging-threshold": "info"
 }
 
-OPENEO_INFERENCE_CDSE_JOB_OPTIONS: Dict[str, str] = {
+OPENEO_INFERENCE_CDSE_JOB_OPTIONS: Dict[str, Union[str, int, float, List]] = {
     "driver-memory": "1000m",
     "driver-memoryOverhead": "1000m",
     "driver-cores": 1,
@@ -118,7 +119,7 @@
     ]
 }
 
-OPENEO_POINTEXTRACTION_CDSE_JOB_OPTIONS: Dict[str, str] = {
+OPENEO_POINTEXTRACTION_CDSE_JOB_OPTIONS: Dict[str, Union[str, int, float, List]] = {
     "driver-memory": "4G",
     "driver-memoryOverhead": "2G",
     "driver-cores": 1,
@@ -130,7 +131,7 @@
     "logging-threshold": "info"
 }
 
-OPENEO_CUBEEXTRACTION_CDSE_JOB_OPTIONS: Dict[str, str] = {
+OPENEO_CUBEEXTRACTION_CDSE_JOB_OPTIONS: Dict[str, Union[str, int, float, List]] = {
     "driver-memory": "4G",
     "driver-memoryOverhead": "4G",
     "driver-cores": 1,
@@ -168,7 +169,7 @@
     'S1_collection': "SENTINEL1_GRD"
 }
 
-def _get_default_job_options() -> Dict[str, str]:
+def _get_default_job_options() -> Dict[str, Union[str, int, float, List]]:
     """
     Retrieves the default job options for OpenEO extract operations.
 
@@ -181,7 +182,7 @@ def _get_default_job_options() -> Dict[str, str]:
     """
     return OPENEO_EXTRACT_JOB_OPTIONS.copy()
 
-def get_job_options(provider: str = None, task: str = 'raw_extraction') -> Dict[str, str]:
+def get_job_options(provider: str = None, task: str = 'raw_extraction') -> Dict[str, Union[str, int, float, List]]:
     """
     Retrieve job options based on the specified provider and task.
 
@@ -288,6 +289,7 @@ def get_standard_processing_options(provider: str, task: str = 'raw_extraction')
             "skip_check_S1": SKIP_CHECK_S1,
             "skip_check_S2": SKIP_CHECK_S2,
             "apply_cloud_mask": APPLY_CLOUD_MASK,
+            "get_NOBSperc": GET_NOBS,
             "openeo_chunk_size": CHUNK_SIZE,
         }
     elif (task == 'feature_generation') or (task == 'vi_generation'):
@@ -307,6 +309,7 @@ def get_standard_processing_options(provider: str, task: str = 'raw_extraction')
             "skip_check_S1": SKIP_CHECK_S1,
             "skip_check_S2": SKIP_CHECK_S2,
             "apply_cloud_mask": APPLY_CLOUD_MASK,
+            "get_NOBSperc": GET_NOBS,
             "optical_vi_list" : VI_LIST,
             "radar_vi_list" : RADAR_LIST,
             "S2_scaling" : S2_SCALING,
@@ -325,62 +328,57 @@ def get_advanced_options(provider: str, s1_orbitdirection: Optional[str] = S1_OR
                          S1_temporal_reducer: str = S1_TEMPORAL_REDUCER, slc_masking: Optional[str] = MASKING_ALGO,
                          S2_bands: List[str] = S2_BANDS, s2_tileid_list: Optional[List[str]] = S2_TILEID_LIST,
                          skip_check_S1: bool = SKIP_CHECK_S1, skip_check_S2: bool = SKIP_CHECK_S2,
-                         apply_cloud_mask: bool = APPLY_CLOUD_MASK, S2_max_cloud_cover: int = S2_MAX_CLOUD_COVER,
+                         apply_cloud_mask: bool = APPLY_CLOUD_MASK, get_NOBSperc: bool = GET_NOBS,
+                         S2_max_cloud_cover: int = S2_MAX_CLOUD_COVER,
                          optical_vi_list: List[str] = VI_LIST, radar_vi_list: List[str] = RADAR_LIST,
                          S2_scaling: List[int | float] = S2_SCALING, S1_db_rescale: bool = True,
                          append: bool = True) -> Dict[str, Union[str, bool, int | float, List[str], List[int | float]]]:
     """
-    Generate a dictionary of advanced options for processing satellite imagery.
-
-    This function validates input parameters and prepares a dictionary of options
-    required for processing satellite imagery data. Validation is applied for
-    various parameters to ensure that the provided values conform to the expected
-    formats or constraints. The function checks for valid values for key options
-    like satellite orbit direction, target coordinate reference system (CRS),
-    resolution, and others before constructing the dictionary.
-
-    :param provider: The provider of satellite data.
-    :param s1_orbitdirection: Direction of Sentinel-1 pass, can be 'ASCENDING',
-                              'DESCENDING' or None.
-    :param target_crs: An integer representing the target coordinate reference
-                       system (CRS).
-    :param resolution: The spatial resolution of the output in integer or float
-                       format.
-    :param ts_interpolation: Boolean indicating whether to apply linear time-series
-                             interpolation.
-    :param ts_interval: Interval for time series aggregation (temporal binning). Accepted values
-                        are 'day', 'week', 'dekad', 'month', 'season', 'year', or
-                        None.
-    :param S2_temporal_reducer: Method for reducing temporal data for Sentinel-2.
-                                Valid options include 'median', 'mean', 'max',
-                                'min', 'first', 'last', 'product', 'sd', 'sum',
-                                or 'variance'.
-    :param S1_temporal_reducer: Method for reducing temporal data for Sentinel-1.
-                                Valid options include 'median', 'mean', 'max',
-                                'min', 'first', 'last', 'product', 'sd', 'sum',
-                                or 'variance'.
-    :param slc_masking: Masking approach to be applied. Valid options are
-                        'mask_scl_dilation', 'satio', or None.
-    :param S2_bands: List of selected Sentinel-2 reflectance bands.
-    :param s2_tileid_list: Optional list of Sentinel-2 tiles for processing. Can
-                           be None if tiles are not specified.
-    :param skip_check_S1: Boolean indicating whether to skip validation checks
-                          for Sentinel-1 data.
-    :param skip_check_S2: Boolean indicating whether to skip validation checks
-                          for Sentinel-2 data.
-    :param apply_cloud_mask: Boolean indicating whether to apply cloud masking or to add it as own band.
-    :param S2_max_cloud_cover: Maximum allowable cloud cover percentage for
-                               Sentinel-2 data. Acceptable values are integers
-                               between 0 and 100.
-    :param optical_vi_list: List of selected vegetation indices for optical data.
-    :param radar_vi_list: List of selected vegetation indices for radar data.
-    :param S2_scaling: List of scaling factors to be applied to Sentinel-2 data.
-    :param S1_db_rescale: Boolean controlling whether Sentinel-1 data must be
-                          rescaled in decibels.
-    :param append: Boolean indicating whether to append the VI's to the reflectance cube or to replace them.
-    :param openeo_chunk_size: Chunk size for OpenEO processing.
-
-    :return: A dictionary containing all validated options as key-value pairs.
+    Constructs a dictionary of advanced processing options for remote sensing data.
+
+    This function validates input parameters and generates an options dictionary
+    tailored for processing remote sensing datasets, including Sentinel-1 and
+    Sentinel-2 data. It handles settings like temporal reducers, orbit directions,
+    spatial resolutions, and others. Input validation ensures proper parameter
+    types and value ranges.
+
+    :param provider: (str) The data provider identifier.
+    :param s1_orbitdirection: (Optional[str]) The orbit direction for Sentinel-1
+        data, can be 'ASCENDING', 'DESCENDING', or None.
+    :param target_crs: (Optional[int]) Target coordinate reference system code.
+    :param resolution: (int or float) Spatial resolution for the output data.
+    :param ts_interpolation: (bool) Whether to enable time-series interpolation.
+    :param ts_interval: (Optional[str]) Time-series interval, allowed values:
+        'day', 'week', 'dekad', 'month', 'season', 'year', or None.
+    :param S2_temporal_reducer: (str) Temporal reducer applied to Sentinel-2 data,
+        permitted values: 'median', 'mean', 'max', 'min', 'first', 'last',
+        'product', 'sd', 'sum', 'variance'.
+    :param openeo_chunk_size: (int) Chunk size for OpenEO processing.
+    :param S1_temporal_reducer: (str) Temporal reducer applied to Sentinel-1 data,
+        permitted values: 'median', 'mean', 'max', 'min', 'first', 'last',
+        'product', 'sd', 'sum', 'variance'.
+    :param slc_masking: (Optional[str]) Masking algorithm for SLC, either
+        'mask_scl_dilation', 'satio', or None.
+    :param S2_bands: (List[str]) List of Sentinel-2 reflectance bands to include.
+    :param s2_tileid_list: (Optional[List[str]]) List of Sentinel-2 tile IDs to
+        process, or None.
+    :param skip_check_S1: (bool) Whether to skip checking for Sentinel-1 data.
+    :param skip_check_S2: (bool) Whether to skip checking for Sentinel-2 data.
+    :param apply_cloud_mask: (bool) Whether to apply a cloud mask to Sentinel-2
+        data.
+    :param get_NOBSperc: (bool) Whether to calculate the percentage of observations.
+    :param S2_max_cloud_cover: (int) Maximum cloud cover percentage allowed
+        for Sentinel-2 scenes (0-100 inclusive).
+    :param optical_vi_list: (List[str]) List of vegetation indices to calculate
+        for optical data.
+    :param radar_vi_list: (List[str]) List of indices to calculate for radar data.
+    :param S2_scaling: (List[int or float]) List of scaling factors for Sentinel-2 data.
+    :param S1_db_rescale: (bool) Whether to rescale Sentinel-1 data to decibels.
+    :param append: (bool) Whether to append new processing options to an existing
+        configuration dictionary.
+
+    :return: (Dict[str, Union[str, bool, int or float, List[str], List[int or float]]])
+        A dictionary containing all validated and formatted processing options.
     """
 
     if s1_orbitdirection not in ['ASCENDING', 'DESCENDING', None]:
@@ -448,6 +446,13 @@ def get_advanced_options(provider: str, s1_orbitdirection: Optional[str] = S1_OR
     if type(apply_cloud_mask) != bool:
         raise ValueError(f'parameter for apply_cloud_mask must be an boolean.')
 
+    if type(get_NOBSperc) != bool:
+        raise ValueError(f'parameter for get_NOBSperc must be an boolean.')
+
+    if (get_NOBSperc == True) and ((apply_cloud_mask == False) or (slc_masking != 'mask_scl_dilation')):
+        raise ValueError(f'parameter for get_NOBSperc can be only set to True if apply_cloud_mask is set to True '
+                         f'AND slc_masking is set to mask_scl_dilation.')
+
     proc_opt = {
         "provider": provider,
         "s1_orbitdirection": s1_orbitdirection,
@@ -464,6 +469,7 @@ def get_advanced_options(provider: str, s1_orbitdirection: Optional[str] = S1_OR
         "skip_check_S1": skip_check_S1,
         "skip_check_S2": skip_check_S2,
         "apply_cloud_mask": apply_cloud_mask,
+        "get_NOBSperc": get_NOBSperc,
         "optical_vi_list": optical_vi_list,
         "radar_vi_list": radar_vi_list,
         "S2_scaling": S2_scaling,

src/eo_processing/openeo/preprocessing.py

@@ -1,5 +1,6 @@
 from __future__ import annotations
 import openeo
+from openeo import DataCube
 from openeo.processes import array_create, if_, is_nodata, power, array_contains
 from openeo.rest.datacube import DataCube
 
@@ -8,7 +9,7 @@
                                                  catalogue_check_CDSE_S1, catalogue_check_CDSE_S2)
 from eo_processing.config.settings import S2_BANDS, PLANET_BANDS, CHUNK_SIZE
 
-from typing import Optional, Dict, Union, List, TYPE_CHECKING
+from typing import Optional, Dict, Union, List, TYPE_CHECKING, Tuple
 if TYPE_CHECKING:
     from eo_processing.config.data_formats import openEO_bbox_format
 
@@ -29,11 +30,15 @@ def ts_datacube_extraction(
             resolution, ts_interval, time_interpolation, SLC_masking_algo, s1_orbitdirection, S2_bands)
     :return: DataCube
     """
-    # get the Sentinel-2 datacube as starting point
-    bands = extract_S2_datacube(connection, bbox, start, end,
-                                S2_collection=S2_collection,
-                                **processing_options)
+    # Note: for cubes with a time dimension, no nobs_perc band can be attached
+    get_NOBSperc: bool = processing_options.get("get_NOBSperc", False)
+    if get_NOBSperc:
+        processing_options["get_NOBSperc"] = False
 
+    # get the Sentinel-2 datacube as starting point
+    bands: DataCube = extract_S2_datacube(connection, bbox, start, end,
+                                          S2_collection=S2_collection,
+                                          **processing_options)
     # add the Sentinel-1 data
     if S1_collection is not None:
         bands = bands.merge_cubes(extract_S1_datacube(connection, bbox, start, end,
@@ -170,7 +175,8 @@ def extract_S1_datacube(
 def extract_S2_datacube(
         connection: openeo.Connection, bbox: Optional[openEO_bbox_format], start: str, end: str,
         S2_collection: str='SENTINEL2_L2A',
-        **processing_options: Dict[str, Union[str, bool, int | float, List[str], List[int | float]]]) -> DataCube:
+        **processing_options: Dict[str, Union[str, bool, int | float, List[str], List[int | float]]]) \
+        -> tuple[DataCube, DataCube | None] | DataCube:
     """ extract the Sentinel-2 data for requested time period and preprocess the data
 
     :param connection: active openEO connection object
@@ -179,7 +185,7 @@ def extract_S2_datacube(
     :param end: str, End date for requested input data (yyyy-mm-dd)
     :param S2_collection: (str, optional): Collection name for S2 data
     :param processing_options: (dict, optional), processing options for preprocessing routine (provider, target_crs,
-            resolution, ts_interval, time_interpolation, SLC_masking_algo, S2_bands)
+            resolution, ts_interval, time_interpolation, SLC_masking_algo, S2_bands, apply_cloud_mask, getNOBSperc)
     :return: DataCube
     """
     # evaluate additional processing_options
@@ -203,6 +209,7 @@ def extract_S2_datacube(
     max_cloud_max: int = processing_options.get("S2_max_cloud_cover", 95)
     chunk_size: int = processing_options.get("openeo_chunk_size", CHUNK_SIZE)
     s2_tileid_list: Optional[List[str]] = processing_options.get("s2_tileid_list", None)
+    get_NOBSperc: bool = processing_options.get("get_NOBSperc", False)
 
     # check if the masking parameter is valid
     if masking not in ['satio', 'mask_scl_dilation', None]:
@@ -240,6 +247,9 @@ def extract_S2_datacube(
     else:
         bands = bands.resample_spatial(resolution=target_res)
 
+    # for the NOBSperc to work we need to init the cube as None
+    nobs_perc_band = None
+
     # apply cloud masking
     if masking == 'mask_scl_dilation':
         # we have to load the SCL mask as an extra cube to get it correctly working
@@ -271,8 +281,22 @@ def extract_S2_datacube(
         ).rename_labels("bands", ["S2-CLOUD-MASK"])
 
         if apply_mask:
+            if get_NOBSperc:
+                # get temporal count of valid pixels before masking
+                tobs = sub_collection.count_time().rename_labels("bands", ["tobs"]).convert_data_type('float32')
+
+            # apply the cloud mask to the reference data
             bands = bands.mask(scl_dilated_mask) # here I do trust the automatic resampling of the mask
+
+            if get_NOBSperc:
+                # get temporal count after applying cloud mask = nobs
+                nobs = bands.count_time().rename_labels("bands", ["nobs"]).convert_data_type('float32')
+                # make sure TOBS is in right projection/resolution
+                tobs = tobs.resample_cube_spatial(target=nobs, method="near")
+                # calculate nobs_perc
+                nobs_perc_band = nobs.divide(tobs).multiply(100.).rename_labels("bands", ["nobs_perc"])
         else:
+            scl_dilated_mask = scl_dilated_mask.resample_cube_spatial(target=bands, method="near")
             bands = bands.merge_cubes(scl_dilated_mask)
     elif masking == 'satio':
         # Apply satio-based mask
@@ -287,6 +311,7 @@ def extract_S2_datacube(
         if apply_mask:
             bands = bands.mask(mask) # masks are automatically resampled/warped
         else:
+            mask = mask.resample_cube_spatial(target=bands, method="near")
             bands = bands.merge_cubes(mask)
 
     # time aggregation if wished
@@ -299,8 +324,12 @@ def extract_S2_datacube(
                                       process="array_interpolate_linear")
     # forcing 16bit --> UInt16
     bands = bands.linear_scale_range(0, 65534, 0, 65534)
-
-    return bands
+    
+    # the return is tricky since we have to be backwards compatible and still get nobs_perc working
+    if get_NOBSperc:
+        return bands, nobs_perc_band
+    else:
+        return bands
 
 def extract_planet_datacube( 
         connection: openeo.Connection, bbox: Optional[openEO_bbox_format], start: str, end: str,