name: crop-phenophase version: "0.1" status: draft # compute graph VALIDATED on synthetic data via # examples/reference_impl.py; engine-wiring # pending. See PROVENANCE. category: time-series-trends themes: [agriculture, imagery] family: crop-phenology-monitoring cites_papers: ["jonsson-2004", "zhang-2003"] # jonsson-2004: TIMESAT — the threshold-crossing # method (greenup/senescence as % of seasonal # amplitude) this recipe reproduces in Python. # zhang-2003: piecewise logistic for SOS/EOS, # the complementary parametric approach that # validates the amplitude-threshold dates. description: > Per-pixel and per-field phenophase transition dates from a Sentinel-2 NDVI season: green-up (spring upward threshold crossing at 50 % of seasonal amplitude), peak greenness (day-of-year of maximum NDVI from the harmonic fit), senescence (autumn downward crossing at 50 % of amplitude), and season length (senescence DOY minus green-up DOY). Crop-calendar detail that raw harmonic phase alone does not give: the phase tells you when the season peaks, but phenophase dates tell you when it starts and ends — the information growers, insurers, and agronomists need for planting-date compliance, growing-degree-day audits, and yield-gap attribution. Builds directly on crop-monitoring (same S2 source + harmonic_fit + zonal_stats graph) by adding the amplitude-threshold extraction layer between the harmonic fit and the field aggregation. Swap the threshold fraction (default 0.50) or the index (NDVI -> EVI/SAVI) without touching the rest. # ───────────────────────────────────────────────────────────────────────────── # PROVENANCE (why this recipe exists) # Authored 2026-07-16 as the crop-calendar detail layer that crop-monitoring # lacks. crop-monitoring (family: crop-phenology-monitoring) yields amplitude # + phase + trend; this recipe extends the same compute chain with # TIMESAT-style threshold crossings (Jonsson & Eklundh 2004) to produce the # four agronomic dates growers and insurers actually request: green-up, # peak, senescence, season-length. The threshold extraction is a pure-Python # op over the harmonic_fit output, so no new engine primitive is required. # # Papers: jonsson-2004 (TIMESAT, the canonical threshold-crossing algorithm); # zhang-2003 (piecewise logistic SOS/EOS, the other mainstream phenology # method for validation cross-checks). # # STATUS: draft. The op composition follows each op's documented signature. # Validate with `folia run` on the example AOI, reconcile any param/shape # mismatches, add a validations/ block, then drop status:draft. # ───────────────────────────────────────────────────────────────────────────── # ───────────────────────────────────────────────────────────────────────────── # PARAMETERIZATION # settings.threshold_fraction: the TIMESAT amplitude fraction that defines # green-up and senescence. 0.50 = 50 % of seasonal amplitude, the # Jonsson & Eklundh default. 0.20 gives an earlier/later outer envelope. # settings.n_harmonics / fit_trend / period: the harmonic model (same as # crop-monitoring; 2 harmonics captures double-cropping). # inputs.fields: the field-boundary vector — zonal aggregation zones. # inputs.aoi / date_range / collection: the bound instance. # ───────────────────────────────────────────────────────────────────────────── settings: index: slug: ndvi expression: "((B08 - reflectance_offset) - (B04 - reflectance_offset)) / ((B08 - reflectance_offset) + (B04 - reflectance_offset))" bands: [B04, B08] # 10 m native, no grid alignment needed range: [-1, 1] colormap: greens threshold_fraction: 0.50 # TIMESAT default: 50 % of seasonal amplitude # for green-up (spring) and senescence (fall) n_harmonics: 2 # annual + semi-annual (handles double-cropping) fit_trend: true # multi-year greening/browning period: 365.25 # DOY time axis # doy_grid: dense DOY vector (1..365) used to evaluate the fitted curve for # threshold crossing. 1-day resolution matches TIMESAT. doy_resolution: 1 # days between evaluated curve points field_stats: [mean, std, min, max] # per-field aggregation of each phenophase band reflectance_offset: 1000 # >=04.00 BOA offset; scale cancels in ND max_bad_fraction: 0.30 # AOI-level cloud gate bad_scl_classes: [0, 1, 3, 8, 9, 10] inputs: aoi: name: San Joaquin Valley cropland block — Fresno County bbox: [-120.05, 36.55, -119.85, 36.72] crs_target: EPSG:32611 mgrs_tile: "11SKB" fields: # The zones phenophase dates are aggregated into. Bound to Fields of the # World (@ftw/global-field-boundaries) — global, CC-BY, fiboa schema. # Any polygon set works (parcels, CAP/FSA, Overture land_use). name: Fields of the World — field boundaries (AOI subset) uri: "@ftw/global-field-boundaries" confidence_attr: confidence # per-field model certainty (0..1) geometry: registry/recipes/crop-monitoring/examples/fresno-fields.geojson date_range: start: "2023-01-01" # >=2 full growing seasons for stable harmonic fit end: "2026-06-12" collection: uri: "@esa/sentinel-2/l2a" assets: [B04, B08, SCL] # NDVI bands + SCL cloud gate layers: # ── SOURCE ────────────────────────────────────────────────────────────────── source/s2: type: raster_timeseries uri: "@esa/sentinel-2/l2a" description: > Every clear S2 L2A pass over the AOI (SCL-gated). One item per acquisition date; shared with crop-monitoring. params: aoi: { bbox: "${inputs.aoi.bbox}" } datetime: "${inputs.date_range.start}/${inputs.date_range.end}" assets: "${inputs.collection.assets}" dedup_by: acquisition_date aoi_cloud_gate: scl_asset: SCL bad_classes: "${settings.bad_scl_classes}" max_bad_fraction: "${settings.max_bad_fraction}" # ── COMPUTE ────────────────────────────────────────────────────────────────── compute/ndvi: type: raster_timeseries description: > NDVI per scene (raster_calc over B04/B08). The (t,y,x) stack feeds the harmonic fit. Identical to crop-monitoring; swap expression for EVI/SAVI. compute: op: raster_calc inputs: data: { layer: source/s2 } params: expression: "${settings.index.expression}" reflectance_offset: "${settings.reflectance_offset}" output_band: ndvi compute/harmonic: type: raster description: > Per-pixel seasonal decomposition (harmonic_fit, n_harmonics=2). Returns amplitude (seasonal NDVI swing), phase (DOY of peak greenness), trend (multi-year), fitted (t,y,x reconstructed curve), and coefficients. The fitted surface is what the threshold-crossing layer evaluates. compute: op: harmonic_fit inputs: collection: { layer: compute/ndvi } params: n_harmonics: "${settings.n_harmonics}" fit_trend: "${settings.fit_trend}" period: "${settings.period}" compute/phenophase_dates: type: raster description: > Four-band raster of phenophase transition dates (all as DOY): band 1 — greenup_doy: first spring crossing upward through (baseline + threshold_fraction * amplitude) band 2 — peak_doy: DOY of maximum fitted NDVI (= phase from harmonic_fit, re-evaluated on the dense grid) band 3 — senescence_doy: last autumn crossing downward through the same threshold band 4 — season_length: senescence_doy - greenup_doy (days) Method follows Jonsson & Eklundh 2004 (TIMESAT): evaluate the fitted harmonic curve on a dense 1-day DOY grid, find the seasonal range (max - min of the fitted curve), set the threshold at baseline + fraction * range, and detect the first upward and last downward crossings. Pixels where amplitude < 0.05 (bare soil / water) are set to NaN — no season detected. compute: op: raster_calc # implemented as a custom expression in # the reference_impl; engine maps this to # a multi-band output op once wired. inputs: harmonic: { layer: compute/harmonic } params: threshold_fraction: "${settings.threshold_fraction}" doy_resolution: "${settings.doy_resolution}" min_amplitude: 0.05 # below this = no detectable season compute/field_phenophase: type: table description: > The crop-phenophase product: one row per field with mean/std/min/max of each phenophase band (greenup_doy, peak_doy, senescence_doy, season_length) within the field polygon. zonal_stats over the four-band phenophase_dates raster keyed by inputs.fields. Sort by greenup_doy to rank early vs. late planters; threshold season_length to flag short-season stress; compare peak_doy to the regional crop calendar to detect off-schedule management. compute: op: zonal_stats inputs: raster: { layer: compute/phenophase_dates } zones: { feature: "${inputs.fields.geometry}" } params: stats: "${settings.field_stats}" outputs: primary: - compute/field_phenophase # geojson — per-field phenophase dates (THE product) show_your_work: - compute/ndvi # zarr (t,y,x) — the NDVI stack behind the fit - compute/harmonic # cog (amplitude, phase, trend, fitted) - compute/phenophase_dates # cog (greenup, peak, senescence, season_length) provenance: - source_scene_ids # STAC item id per timestamp # ── DASHBOARD ─────────────────────────────────────────────────────────────── dashboard: patterns: - slug: time-scrub-availability # scrub NDVI over the season; clear/clouded strip slot: temporal layer: compute/ndvi - slug: pixel-trajectory-climatology # click a field -> NDVI season vs harmonic fit slot: temporal layer: compute/ndvi - slug: choropleth-ranked # fields shaded by greenup_doy or season_length slot: recipe-manipulation layer: compute/field_phenophase # Recommended default: greenup_doy, colormap RdYlGn reversed # (early green-up = green; late = red). Lets grower rank fields by # earliest spring activation. Switch to season_length for stress triage. default_band: greenup_doy colormap: RdYlGn_r - slug: legend-as-histogram # distribution of greenup DOY across the AOI slot: recipe-manipulation layer: compute/field_phenophase license: recipe: Apache-2.0 output_dataset: CC-BY-4.0 upstream_data: - { source: sentinel-2-l2a, license: usgs-public-domain } # Related: crop-monitoring (parent recipe — same S2 source + harmonic_fit + # zonal_stats graph; this recipe adds the threshold-crossing extraction). # Papers: jonsson-2004 (TIMESAT), zhang-2003 (piecewise logistic SOS/EOS).