# RUNS (validated 2026-08-11 against live STAC; ADR-0086 closed the capability gap this
# recipe was blocked on). The source declares `bands_mode: preserve` — every declared band
# per scene, a [time, band] axis-product tensor — the mask+NDVI layer is ONE per-scene
# `raster_calc` binding b1..b3, and `temporal_reduce` reduces the surviving time axis.
# First live run: 4 scenes x 3 bands (12 planes), 17.7 MiB / 34 requests / 7.4 s; composite
# NDVI in [-0.08, 0.88] with ~30% of pixels clear across the 4 early-March scenes.
#
# The registry's worked example of PER-PIXEL SCL cloud masking.
#
# Two things share the name "cloud mask" and they are not the same operation:
#
#   WHOLE-SCENE GATE   `aoi_cloud_gate:` in a source layer's params — drop an entire
#                      acquisition when too much of the AOI is bad. One tiny SCL read
#                      decides the scene. See registry/recipes/index-timeseries.
#   PER-PIXEL MASK     this recipe — keep the scene, void the bad pixels inside it.
#
# They compose (gate cheaply, then mask exactly); they do not substitute for each other. A
# gate alone leaves cloud in every surviving scene. A mask alone pays to read scenes that are
# 90% cloud. This file deliberately uses only the mask, so the pattern is visible on its own.
#
# A dedicated op DOES exist: `cloud_mask_sentinel2`. It is one layer instead of two, and it is
# the better choice when you want the conventional mask and nothing more. It is not used here
# because its class set is compiled into the Rust kernel (it NaN-s SCL 3, 8, 9, 10 and keeps
# everything else) and `scl_classes:` is Python-backend only, so it cannot express this file's
# class set at all. See docs/RECIPE_AUTHORING.md "Two routes, and how to choose".
#
# The composition is what you fall back to when the class set is a decision of the analysis
# rather than a default. It puts both decisions on the page: WHICH SCL classes you call bad
# (a `settings:` value here) and whether you are gating or masking.

name: openeo-scl-cloud-mask
category: time-series-trends
themes: [imagery, vegetation]
version: "0.2.0"
description: >
  Per-pixel Sentinel-2 SCL cloud masking, then an NDVI temporal-mean composite. For every
  L2A pass over the AOI, evaluate the scene classification layer per pixel, void the NDVI
  pixels it rejects, and average what survives across time. The load -> mask half reproduces
  the openEO cookbook's SCL masking example ("Mask Out Specific Values", keep SCL {4, 5});
  the NDVI and the temporal mean are this recipe's own completion of it. Each step is its
  own inspectable layer.

settings:
  # SCL classes to KEEP. Two different class sets get called "the SCL cloud mask", and this
  # file implements the FIRST of them. Both are labelled here because picking the wrong one
  # is silent: either produces a plausible composite of the wrong pixels.
  #
  #   0 nodata  1 saturated  2 dark-area  3 cloud-shadow  4 vegetation  5 bare-soil
  #   6 water   7 unclassified  8 cloud-medium  9 cloud-high  10 cirrus  11 snow/ice
  #
  # (A) THE COOKBOOK'S EXAMPLE — what this recipe reproduces. openEO's cookbook, under
  #     "Mask Out Specific Values", writes exactly:
  #
  #         SCL = cube_s2.band("SCL")
  #         cube_s2_mask = ~ ((SCL == 4) | (SCL == 5))
  #         cube_s2_masked = cube_s2.mask(cube_s2_mask)
  #
  #     Keep {4, 5}. This is narrower than a cloud mask and worth knowing you are choosing:
  #     it also drops water (6), unclassified (7) and snow (11), so a lake or a snowfield
  #     reads as NaN, not as itself. Faithful to the tutorial; often not what you want.
  #
  # (B) THE CLOUD MASK — the other common reading, which drops only the cloud-related
  #     classes and keeps everything else. Swap the expression in compute/ndvi-masked for:
  #
  #         expression: "where((b3 == 3) | (b3 == 8) | (b3 == 9) | (b3 == 10), 0/0, (b2 - b1) / (b2 + b1))"
  #
  #     That drop set is exactly what the built-in `cloud_mask_sentinel2` op applies, so if
  #     you want (B) you do not need this composition at all — use the op (one layer).
  #
  # Named here rather than buried in the expression so a reader can see, and change, exactly
  # what "clear" means for this recipe.
  keep_scl_classes: [4, 5]

inputs:
  aoi:
    name: Pellworm, Schleswig-Holstein, DE
    bbox: [8.5464, 54.4473, 9.0724, 54.5685]
  date_range:
    start: "2021-03-01"
    end: "2021-03-10"

layers:

  # ── SOURCE ────────────────────────────────────────────────────────────────────────────
  # `scenes: true` selects the per-scene bridge; `bands_mode: preserve` (ADR-0086) keeps
  # EVERY declared band per scene instead of collapsing each scene to one index band — the
  # output is the [time, band] axis-product tensor the two consumers below are keyed to.
  # Bands bind POSITIONALLY inside the per-scene expression below, in the order declared
  # here: b1 = R, b2 = N, b3 = SCL. Preserve mode passes RAW band values through (no
  # reflectance scaling — scaling would destroy the SCL class codes); NDVI is a pure ratio,
  # so the missing scale cancels.
  source/s2-scenes:
    type: tensor
    uri: "@esa/sentinel-2/l2a"
    description: >
      Every Sentinel-2 L2A pass over the AOI in the window, red + NIR + the scene
      classification layer, one band-set per acquisition date.
    params:
      scenes: true
      bands_mode: preserve
      bands: [R, N, SCL]
      bbox: "${inputs.aoi.bbox}"
      temporal: ["${inputs.date_range.start}", "${inputs.date_range.end}"]

  # ── MASK + INDEX, PER SCENE ───────────────────────────────────────────────────────────
  # One per-scene evaluation (ADR-0086: raster_calc maps over the time axis, binding the
  # band axis coordinates to b1..b3 per scene): where this scene's OWN SCL says keep,
  # compute this scene's NDVI; everywhere else void the pixel (0/0 evaluates to NaN, the
  # engine's nodata). POLARITY, which is where openEO ports go wrong: openEO's `.mask()` is
  # TRUE-where-REMOVED, so the cookbook's ~((SCL == 4) | (SCL == 5)) becomes the UNNEGATED
  # keep-condition here. `where(cond, then, else)` is the dialect's conditional; there is no
  # ternary. (`folia/domains/geo/operations/_functions.yaml` is the whole function list.)
  # Voiding to NaN is what makes the reduce below correct — the reducer is NaN-aware and
  # SKIPS rejected pixels rather than averaging in a sentinel.
  compute/ndvi-masked:
    type: tensor
    description: >
      The per-scene NDVI stack with every SCL-rejected pixel voided — one plane per
      acquisition, each judged by its own SCL. This is the show-your-work layer: scrub it
      to see which pixels each pass actually contributed.
    compute:
      op: raster_calc
      inputs: { data: { layer: source/s2-scenes } }
      params:
        expression: "where((b3 == 4) | (b3 == 5), (b2 - b1) / (b2 + b1), 0/0)"
        output_band: ndvi_masked

  # ── COMPOSITE ─────────────────────────────────────────────────────────────────────────
  # This step is OURS, not the cookbook's. The cookbook's masking example stops at the mask:
  # it reduces over the `bands` dimension to build the mask and never reduces over time. The
  # temporal reduce is added here so the recipe answers something, and `mean` over `median`
  # is this file's choice, not a ported one. `method: median` is a one-word swap and is the
  # more robust pick when a few passes survive masking with residual haze.
  #
  # The input edge must be named `data` (or `stack`/`raster`) — the op reads those and
  # nothing else.
  compute/ndvi-composite:
    type: raster
    description: >
      Mean NDVI per pixel over the surviving observations. A pixel cloudy in every pass has
      no valid samples and answers NaN rather than a fabricated value.
    compute:
      op: temporal_reduce
      inputs: { data: { layer: compute/ndvi-masked } }
      params:
        method: mean

outputs:
  ndvi_composite:
    layer: compute/ndvi-composite
    value_kind: raster
    description: Cloud-masked mean NDVI composite.
  ndvi_stack:
    layer: compute/ndvi-masked
    value_kind: tensor
    description: Show-your-work — the masked per-scene NDVI stack behind the composite.

license:
  recipe: Apache-2.0
  output_dataset: CC-BY-4.0
  upstream_data:
    - { source: sentinel-2-l2a, license: cc-by-sa-3.0-igo }