name: rondonia-deforestation
version: "0.2"
category: change-detection
themes: [vegetation]
family: deforestation-forest-loss-detection
cites_papers: ["souza-2020"]   # reproduces MapBiomas Landsat LULC method (Souza et al. 2020)
description: >
  Two-epoch NDVI loss over the Rondônia fishbone clearing pattern in the
  Brazilian Amazon, along the BR-364 corridor. Landsat surface reflectance is
  composited over an early-2000s baseline epoch and a recent epoch; the drop in
  NDVI between them is the clearing proxy.

# WHAT THIS MEASURES, AND WHAT IT DOES NOT. READ BEFORE CITING A NUMBER.
#   dNDVI < -0.2 is a CLEARING PROXY, not a forest-loss measurement. It flags
#   pixels whose greenness fell between two composites. Selective logging,
#   drought, burn scars, pasture senescence, a crop cycle caught at a different
#   phenological phase and cloud-shadow residue all move NDVI the same way, and
#   the proxy cannot separate them from a clearing. Calling the output "forest
#   loss" would need a forest baseline mask, a per-year time axis and a
#   validated threshold — none of which this recipe carries.
#   MapBiomas (cited above) is the METHOD FAMILY this belongs to; this is not a
#   reproduction of MapBiomas, whose classifier uses a much larger feature set.
#
#   No per-year axis. Each epoch is ONE median composite over its window, so the
#   answer is "cleared between epoch A and epoch B", not an annual series. A
#   genuine annual series needs a per-year grouping reduction; no shipped op
#   provides one (tracker defect 63), and `raster_transition_matrix` — the
#   nearest categorical crosstab — consumes pre-binned H3 transition tables plus
#   a supply-shed assignment, not rasters (its `mode: raster` arm is declared
#   and not implemented).
#
#   UPSTREAM BUG (recipe-library tracker, defects 61/62). Some Landsat 7 (LE07)
#   assets on the Planetary Computer answer an HTML error page where a COG
#   should be, and the read dies as "Unknown byte order marker: [3C, 21]" (`<!`).
#   Confirmed 2026-08-17 on LE07_L2SP_036035_20221027_02_T1. There is no
#   platform filter to exclude L7 with: the source params carry none, and the
#   `query.platform` block the previous version of this recipe declared was
#   silently ignored. The only lever an author has is the calendar. Windows that
#   overlap 1999-2023 can pull L7 and can therefore hit this.
#
#   No month filter. The engine's source params have no month selector, so the
#   dry-season windows below are expressed as calendar spans and wet-season
#   scenes do enter the median. A median over many scenes is robust to that;
#   a mean would not be.

settings:
  dndvi_clearing_threshold: -0.2   # dNDVI below this → flagged as a clearing candidate
  cloud_cover_max_pct: 30          # eo:cloud_cover upper bound on the STAC search
  pixel_area_m2: 900               # 30 m Landsat pixel — multiply `count` by this for area

inputs:
  aoi:
    bbox: [-63.9, -11.8, -62.1, -10.2]
    crs_target: EPSG:32720
    notes: BR-364 highway corridor — fishbone roads cut E-W, farms clear N-S.

  baseline_epoch:
    start: "2000-07-01"
    end: "2002-11-01"
    notes: Three Rondônia dry seasons at the start of the Landsat 7 era.

  recent_epoch:
    start: "2023-07-01"
    end: "2025-11-01"
    notes: Three recent dry seasons, Landsat 8/9 OLI.

  collection:
    uri: "@usgs/landsat/c2-l2"

layers:

  # ── SOURCE ──────────────────────────────────────────────────────────────────
  #
  # Bands are CANONICAL concept codes from registry/datasets/usgs/landsat/c2-l2/
  # _dataset.yaml (N = nir, R = red, QA = qa-pixel). Each distribution maps them
  # onto its own producer asset keys, so one band list is correct across TM,
  # ETM+ and OLI — which the previous SR_B5/SR_B4 spelling was not.
  #
  # The FUSED per-scene form (`scenes: true` + `expression:` + `mask_band:`) is
  # what makes the index correct: the distribution's declared scale/offset is
  # applied per band before the expression, fill DNs become NaN, and QA_PIXEL is
  # the platform mask band, excluded from the b1..bN binding. `temporal_reduce`
  # then takes the median across the surviving scenes.
  #
  # Do NOT collapse this to a bare `bands: [N, R, QA]` source and difference the
  # result. A source read emits the DECLARED BANDS (measured 2026-08-17: three
  # bands out for a three-band request), so `b2 - b1` there is red-minus-nir on
  # raw reflectance, not dNDVI. It plans green and answers plausibly.

  source/ndvi-baseline-scenes:
    type: raster
    description: Per-scene NDVI, QA-masked, early-2000s baseline epoch.
    uri: "@usgs/landsat/c2-l2"
    params:
      scenes: true
      bands: [N, R, QA]
      mask_band: QA
      expression: "(b1 - b2) / (b1 + b2)"
      bbox: [-63.9, -11.8, -62.1, -10.2]
      temporal: ["2000-07-01", "2002-11-01"]
      cloud_lt: 30
      resolution: 30
      limit: 60

  source/ndvi-recent-scenes:
    type: raster
    description: Per-scene NDVI, QA-masked, recent epoch.
    uri: "@usgs/landsat/c2-l2"
    params:
      scenes: true
      bands: [N, R, QA]
      mask_band: QA
      expression: "(b1 - b2) / (b1 + b2)"
      bbox: [-63.9, -11.8, -62.1, -10.2]
      temporal: ["2023-07-01", "2025-11-01"]
      cloud_lt: 30
      resolution: 30
      limit: 60

  # ── COMPUTE ─────────────────────────────────────────────────────────────────

  compute/ndvi-baseline:
    type: raster
    description: Median NDVI over the baseline epoch's scenes.
    compute:
      op: temporal_reduce
      inputs:
        time_series: { layer: source/ndvi-baseline-scenes }
      params:
        method: median

  compute/ndvi-recent:
    type: raster
    description: Median NDVI over the recent epoch's scenes.
    compute:
      op: temporal_reduce
      inputs:
        time_series: { layer: source/ndvi-recent-scenes }
      params:
        method: median

  compute/ndvi-pair:
    type: raster
    description: The two NDVI composites stacked into one 2-band raster (b1 = baseline, b2 = recent).
    compute:
      op: stack_concat
      inputs:
        layers:
          - { layer: compute/ndvi-baseline }
          - { layer: compute/ndvi-recent }

  compute/ndvi-change:
    type: raster
    description: >
      dNDVI = recent - baseline. Negative = the surface lost greenness between
      epochs. This is the inspectable layer behind the answer.
    compute:
      op: raster_calc
      inputs:
        raster: { layer: compute/ndvi-pair }
      params:
        expression: "b2 - b1"
        output_band: dndvi

  compute/clearing-proxy:
    type: raster
    description: >
      Clearing-candidate footprint. Pixels below the dNDVI threshold keep the
      value 1; everything else becomes NaN (`b1 / 0.0`), which the stats kernel
      skips — so the answer layer's `count` is the flagged-pixel count.
    compute:
      op: raster_calc
      inputs:
        raster: { layer: compute/ndvi-change }
      params:
        expression: "where(b1 < ${settings.dndvi_clearing_threshold}, 1.0, b1 / 0.0)"
        output_band: clearing_proxy

  # ── ANSWER ──────────────────────────────────────────────────────────────────

  answer/cleared-extent:
    type: answer
    description: >
      Count of clearing-candidate pixels between the two epochs. Area =
      count x 900 m² (0.09 ha per pixel). A proxy extent, not a forest-loss area.
    compute:
      op: summary
      inputs:
        raster: { layer: compute/clearing-proxy }
      params:
        stats: [count]

outputs:
  cleared_extent:
    layer: answer/cleared-extent
    value_kind: answer
    description: Clearing-candidate pixel count between the baseline and recent epochs.
  ndvi_change:
    layer: compute/ndvi-change
    value_kind: raster
    description: dNDVI raster, recent minus baseline (show-your-work).
  clearing_proxy:
    layer: compute/clearing-proxy
    value_kind: raster
    description: The thresholded clearing-candidate mask (show-your-work).

license:
  recipe: Apache-2.0
  output_dataset: CC-BY-4.0
  upstream_data:
    - { source: "@usgs/landsat/c2-l2", license: usgs-public-domain }