From 77e148883041eb5de8976011da2b0dd6aff9d2ed Mon Sep 17 00:00:00 2001 From: Felix Delattre Date: Mon, 4 May 2026 09:53:12 +0200 Subject: [PATCH] foo --- metrics_stats.py | 67 ++++++++++++++ run.py | 91 ++++++++++--------- webapp/index.html | 57 ++++++++++-- webapp/metrics.html | 215 +++++++++++++++++++++++++++++++++++--------- 4 files changed, 343 insertions(+), 87 deletions(-) diff --git a/metrics_stats.py b/metrics_stats.py index e22bd57..865a702 100644 --- a/metrics_stats.py +++ b/metrics_stats.py @@ -110,6 +110,25 @@ def nse(y_true, y_pred): return float(1 - (numerator / denominator)) +def residual_vs_phenocam(fusion_ts, phenocam_ts): + """Stats of (fused_GCC − PhenoCam_GCC) on matched dates; None if too few points. + + Mean: positive → fusion systematically above PhenoCam; negative → below; ~0 → unbiased mean. + Compare BtI vs ItB means at same strategy/σ (``derived.bti_vs_itb_mean_residual``): closer to 0 → less mean bias vs PhenoCam. + """ + yf, yp, _dates = match_dates(fusion_ts, phenocam_ts) + if len(yf) < 2: + return None + r = yf - yp + return { + "mean": float(np.mean(r)), + "std": float(np.std(r)), + "mae": float(np.mean(np.abs(r))), + "rmse": float(np.sqrt(np.mean(r**2))), + "n_samples": int(len(r)), + } + + def calculate_temporal_metrics(fusion_ts, phenocam_ts): """Temporal metrics vs PhenoCam (nse_pc; nse is the same value).""" fusion_vals, phenocam_vals, dates = match_dates(fusion_ts, phenocam_ts) @@ -129,9 +148,54 @@ def calculate_temporal_metrics(fusion_ts, phenocam_ts): "n_samples": len(fusion_vals), "date_range": {"start": dates[0], "end": dates[-1]} if dates else None, } + rv = residual_vs_phenocam(fusion_ts, phenocam_ts) + if rv: + metrics["residual_vs_phenocam"] = rv return metrics +def derived_tier1(temporal: dict) -> dict: + """ΔNSE_PC (σ20 − σ30) and paired BtI vs ItB mean residual; needs temporal fusion keys. + + ΔNSE_PC > 0 → NSE_PC higher at σ=20 than σ=30 (tighter EFAST temporal kernel wins). + ΔNSE_PC < 0 → σ=30 wins (broader smoothing matches PhenoCam better). + """ + d_nse = {"bti": {}, "itb": {}} + for strategy in ("aggressive", "nonaggressive"): + for mode, suf in (("bti", ""), ("itb", "_itb")): + k20 = f"{strategy}_sigma20{suf}" + k30 = f"{strategy}_sigma30{suf}" + n20 = (temporal.get(k20) or {}).get("nse_pc") + n30 = (temporal.get(k30) or {}).get("nse_pc") + if isinstance(n20, (int, float)) and isinstance(n30, (int, float)): + d_nse[mode][strategy] = float(n20 - n30) + else: + d_nse[mode][strategy] = None + + paired = [] + for strategy in ("aggressive", "nonaggressive"): + for sig in (20, 30): + kb, ki = f"{strategy}_sigma{sig}", f"{strategy}_sigma{sig}_itb" + mb = (temporal.get(kb) or {}).get("residual_vs_phenocam", {}).get("mean") + mi = (temporal.get(ki) or {}).get("residual_vs_phenocam", {}).get("mean") + paired.append( + { + "strategy": strategy, + "sigma": sig, + "mean_residual_bti": float(mb) + if isinstance(mb, (int, float)) + else None, + "mean_residual_itb": float(mi) + if isinstance(mi, (int, float)) + else None, + } + ) + return { + "delta_nse_pc_sigma20_minus_sigma30": d_nse, + "bti_vs_itb_mean_residual": paired, + } + + def calculate_phenocam_stats(phenocam_ts): """Calculate phenocam summary statistics.""" values = [v for v in phenocam_ts.values() if v is not None] @@ -323,6 +387,9 @@ def calculate_all_metrics(season, site_name, site_position): if temporal_metrics: results["temporal"][scenario_name] = temporal_metrics + if results["temporal"]: + results["derived"] = derived_tier1(results["temporal"]) + # Save results output_path = Path(f"data/{site_name}/{season}/metrics.json") output_path.parent.mkdir(parents=True, exist_ok=True) diff --git a/run.py b/run.py index 73643a0..4e55259 100644 --- a/run.py +++ b/run.py @@ -1,64 +1,73 @@ -from fusion import run_all_efast_scenarios, run_all_efast_itb_scenarios -from postprocessing import ( - post_process_all_scenarios, - post_process_all_itb_scenarios, - post_process_timeseries, -) +"""Pipeline entry point. +Active snippet below only **regenerates metrics.json** (temporal, baseline, +`derived`, `residual_vs_phenocam`). Requires existing post-processed GCC +timeseries under `data/{site}/{season}/processed_*`. + +Un-comment imports and steps below for acquisition → fusion → post-process. +""" + +# from fusion import run_all_efast_scenarios, run_all_efast_itb_scenarios +# from postprocessing import ( +# post_process_all_scenarios, +# post_process_all_itb_scenarios, +# post_process_timeseries, +# ) # from acquisition_s2 import download_s2 # from acquisition_s3 import download_s3 # from acquisition_phenocam import download_phenocam -from preselection import create_timeseries -from preparation import ( - prepare_s2, - prepare_s3, - prepare_s2_gcc_for_itb, - prepare_s3_gcc_for_itb, -) -from metrics_indices import create_prepared_fusion_timeseries +# from preselection import create_timeseries +# from preparation import ( +# prepare_s2, +# prepare_s3, +# prepare_s2_gcc_for_itb, +# prepare_s3_gcc_for_itb, +# ) +# from metrics_indices import create_prepared_fusion_timeseries from metrics_stats import calculate_all_metrics -from phenology_timesat import write_phenocam_phenology_for_site + +# from phenology_timesat import write_phenocam_phenology_for_site def run_pipeline(season, site_position, site_name): - """Run pipeline.""" + """Run pipeline (metrics-only by default; see module docstring).""" try: # print(f"Downloading S2, S3, and PhenoCam: {site_name}, {season}") # download_s2(season, site_position, site_name) # download_s3(season, site_position, site_name) # download_phenocam(season, site_position, site_name) - print(f"PhenoCam phenology (50 % amplitude): {site_name}, {season}") - write_phenocam_phenology_for_site(site_name, season) + # print(f"PhenoCam phenology (50 % amplitude): {site_name}, {season}") + # write_phenocam_phenology_for_site(site_name, season) - print(f"Creating preselection timeseries: {site_name}, {season}") - create_timeseries(season, site_position, site_name) + # print(f"Creating preselection timeseries: {site_name}, {season}") + # create_timeseries(season, site_position, site_name) - print(f"Preparing S2 and S3 for fusion: {site_name}, {season}") - for strategy in ["aggressive", "nonaggressive"]: - prepare_s2(season, site_position, site_name, cleaning_strategy=strategy) - prepare_s3(season, site_position, site_name, cleaning_strategy=strategy) + # print(f"Preparing S2 and S3 for fusion: {site_name}, {season}") + # for strategy in ["aggressive", "nonaggressive"]: + # prepare_s2(season, site_position, site_name, cleaning_strategy=strategy) + # prepare_s3(season, site_position, site_name, cleaning_strategy=strategy) - print(f"Running EFAST fusion for all scenarios: {site_name}, {season}") - run_all_efast_scenarios(season, site_position, site_name) + # print(f"Running EFAST fusion for all scenarios: {site_name}, {season}") + # run_all_efast_scenarios(season, site_position, site_name) - print(f"Index-then-Blend (ItB): {site_name}, {season}") - for strategy in ["aggressive", "nonaggressive"]: - prepare_s2_gcc_for_itb( - season, site_position, site_name, cleaning_strategy=strategy - ) - prepare_s3_gcc_for_itb( - season, site_position, site_name, cleaning_strategy=strategy - ) - run_all_efast_itb_scenarios(season, site_position, site_name) - post_process_all_itb_scenarios(season, site_position, site_name) + # print(f"Index-then-Blend (ItB): {site_name}, {season}") + # for strategy in ["aggressive", "nonaggressive"]: + # prepare_s2_gcc_for_itb( + # season, site_position, site_name, cleaning_strategy=strategy + # ) + # prepare_s3_gcc_for_itb( + # season, site_position, site_name, cleaning_strategy=strategy + # ) + # run_all_efast_itb_scenarios(season, site_position, site_name) + # post_process_all_itb_scenarios(season, site_position, site_name) - print(f"Creating prepared/fusion timeseries: {site_name}, {season}") - create_prepared_fusion_timeseries(season, site_position, site_name) + # print(f"Creating prepared/fusion timeseries: {site_name}, {season}") + # create_prepared_fusion_timeseries(season, site_position, site_name) - print(f"Post-processing (crop): {site_name}, {season}") - post_process_all_scenarios(season, site_position, site_name) - post_process_timeseries(season, site_position, site_name) + # print(f"Post-processing (crop): {site_name}, {season}") + # post_process_all_scenarios(season, site_position, site_name) + # post_process_timeseries(season, site_position, site_name) print(f"Calculating metrics: {site_name}, {season}") calculate_all_metrics(season, site_name, site_position) diff --git a/webapp/index.html b/webapp/index.html index 8fad0bd..784e7a0 100644 --- a/webapp/index.html +++ b/webapp/index.html @@ -136,7 +136,7 @@
Metrics vs PhenoCam (fusion scenarios)

- R² (variance explained), nRMSE (RMSE normalised by PhenoCam σ), NSE_PC (Nash–Sutcliffe vs PhenoCam). Other metrics remain in metrics.json. + R² (variance explained), nRMSE (RMSE normalised by PhenoCam σ), NSE_PC (Nash–Sutcliffe vs PhenoCam). ΔNSE_PC = NSE_PC(σ20)−NSE_PC(σ30): positive → σ20 better; negative → σ30 better. Mean residual (fused−PhenoCam): positive → fusion high vs PhenoCam on average; negative → low; compare BtI vs ItB in the same row (closer to 0 = less mean bias). Tables at the top when metrics.json has derived (regenerate with metrics_stats.py / run.py).

@@ -605,13 +605,58 @@ container.innerHTML = "

Metrics not available. Run the pipeline (run.py) or metrics_stats.py to generate.

"; return; } - + + const fmtD = (v) => { + const n = Number(v); + return Number.isFinite(n) ? n.toFixed(3) : "—"; + }; + const d = metricsData.derived; + let html = ""; + if (d && d.delta_nse_pc_sigma20_minus_sigma30) { + const dn = d.delta_nse_pc_sigma20_minus_sigma30; + let anyDelta = false; + for (const mode of ["bti", "itb"]) { + for (const strat of ["aggressive", "nonaggressive"]) { + if (Number.isFinite(Number(dn[mode]?.[strat]))) anyDelta = true; + } + } + html += + "

ΔNSE_PC (σ20 − σ30)

"; + html += + "

NSE_PC(σ20) − NSE_PC(σ30): + → σ20 better, → σ30 better.

"; + html += + ""; + for (const mode of ["bti", "itb"]) { + for (const strat of ["aggressive", "nonaggressive"]) { + const v = dn[mode]?.[strat]; + html += ``; + } + } + html += "
ModeStrategyΔNSE_PC
${mode.toUpperCase()}${strat}${fmtD(v)}
"; + if (!anyDelta) { + html += + "

ΔNSE_PC needs both σ20 and σ30 rows in temporal (BtI and ItB) with nse_pc.

"; + } + } + if (d && d.bti_vs_itb_mean_residual && d.bti_vs_itb_mean_residual.length) { + html += + "

Mean residual (fused − PhenoCam): BtI vs ItB

"; + html += + "

Each cell: mean(fused−PhenoCam) on matched dates. + overestimates, underestimates; ~0 little mean bias (see R²/MAE for overall fit). Same row: column closer to 0 → less systematic offset vs PhenoCam (RQ1.1).

"; + html += + ""; + for (const row of d.bti_vs_itb_mean_residual) { + html += ``; + } + html += "
StrategyσBtIItB
${row.strategy}${row.sigma}${fmtD(row.mean_residual_bti)}${fmtD(row.mean_residual_itb)}
"; + } + const scenarios = ["aggressive_sigma20", "aggressive_sigma30", "nonaggressive_sigma20", "nonaggressive_sigma30"]; const scenarioNames = ["Aggressive σ20", "Aggressive σ30", "Non-aggressive σ20", "Non-aggressive σ30"]; const metrics = ["r_squared", "nrmse", "nse_pc"]; const metricLabels = { r_squared: "R²", nrmse: "nRMSE", nse_pc: "NSE_PC" }; - - let html = ""; + + html += "
"; html += ""; html += ""; metrics.forEach(m => html += ``); @@ -655,13 +700,13 @@ }); html += "
Scenario${metricLabels[m]}
"; - + // Add phenocam stats info if available if (metricsData.phenocam_stats) { const stats = metricsData.phenocam_stats; html += `

PhenoCam GCC samples (ground truth): n = ${stats.n_samples}

`; } - + container.innerHTML = html; } diff --git a/webapp/metrics.html b/webapp/metrics.html index 1dab574..e710b65 100644 --- a/webapp/metrics.html +++ b/webapp/metrics.html @@ -14,18 +14,36 @@ .selectors select { padding: 5px 10px; font-size: 14px; margin-right: 15px; } h1 { font-size: 22px; } h2 { font-size: 16px; margin-top: 24px; color: #333; } + h2:first-of-type { margin-top: 8px; } + h3 { font-size: 14px; margin: 14px 0 6px 0; color: #444; font-weight: 600; } table { border-collapse: collapse; width: 100%; font-size: 13px; margin-bottom: 12px; } th, td { border: 1px solid #ccc; padding: 6px 8px; text-align: left; } th { background: #f5f5f5; } td.num { text-align: right; font-variant-numeric: tabular-nums; } + .fusion-block table { margin-bottom: 4px; } + .fusion-block table + table { margin-top: 12px; } .section-note { font-size: 12px; color: #555; margin: -6px 0 8px 0; max-width: 720px; line-height: 1.45; } .section-note code { background: #f1f1f1; padding: 1px 4px; border-radius: 3px; font-size: 11px; } .intro { font-size: 13px; color: #333; background: #fafafa; border: 1px solid #e5e5e5; padding: 10px 12px; border-radius: 4px; margin-bottom: 18px; line-height: 1.5; } - .intro ul { margin: 6px 0 0 18px; padding: 0; } - .intro li { margin-bottom: 2px; } + .intro-short { margin-bottom: 0; } + details.definitions { margin-top: 28px; font-size: 13px; border: 1px solid #e5e5e5; border-radius: 4px; padding: 8px 12px; background: #fafafa; } + details.definitions summary { cursor: pointer; font-weight: 600; color: #333; } + details.definitions ul { margin: 8px 0 0 18px; padding: 0; } + details.definitions li { margin-bottom: 4px; } + .scenario-key { font-size: 11px; color: #666; font-weight: normal; } .empty { color: #666; font-style: italic; } .err { color: #a00; } + details.how-read { + font-size: 12px; color: #333; line-height: 1.5; max-width: 820px; margin: 0 0 18px 0; + padding: 8px 12px 10px; border: 1px solid #ccd; border-radius: 4px; background: #f8fafc; + } + details.how-read summary { + cursor: pointer; font-weight: 600; font-size: 13px; color: #111; margin-bottom: 0; + } + details.how-read ol { margin: 10px 0 0; padding-left: 1.35rem; } + details.how-read li { margin-bottom: 7px; } + details.how-read li:last-child { margin-bottom: 0; } @@ -56,6 +74,14 @@ nrmse: "nRMSE", nse_pc: "NSE_PC", }; + + const FUSION_BTI_ROWS = [ + ["aggressive_sigma20", "Aggressive", 20], + ["aggressive_sigma30", "Aggressive", 30], + ["nonaggressive_sigma20", "Non-aggressive", 20], + ["nonaggressive_sigma30", "Non-aggressive", 30], + ]; + function mv(m, c) { return c === "nse_pc" ? (m.nse_pc ?? m.nse) : m[c]; } @@ -81,21 +107,64 @@ return Number.isInteger(v) ? String(v) : v.toFixed(4); } - function tableSection(title, obj) { - const heading = title ? `

${title}

` : ""; - if (!obj || typeof obj !== "object" || !Object.keys(obj).length) { - return `${heading}

No data

`; - } - const keys = Object.keys(obj).sort(); - let head = `Scenario${DISPLAY_METRIC_COLS.map((c) => `${DISPLAY_METRIC_LABELS[c]}`).join("")}`; - const rows = keys.map((k) => { - const m = obj[k] || {}; - return `${k}${DISPLAY_METRIC_COLS.map((c) => `${fmtMetric(c, mv(m, c))}`).join("")}`; - }).join(""); - return `${heading}${head}${rows}
`; + function fusionMeanResidual(m) { + const x = m?.residual_vs_phenocam?.mean; + const n = Number(x); + return Number.isFinite(n) ? n : null; } - function baselineSection(b) { + function fusionSubTableRows(temporal, keysWithLabels, includeMeanResid) { + const parts = []; + for (const [key, stratLabel, sig] of keysWithLabels) { + const m = temporal[key]; + if (!m) continue; + const mr = fusionMeanResidual(m); + const meanCell = includeMeanResid + ? `${mr !== null ? mr.toFixed(3) : "—"}` + : ""; + parts.push( + `${stratLabel}, σ=${sig} (${key})${DISPLAY_METRIC_COLS.map((c) => `${fmtMetric(c, mv(m, c))}`).join("")}${meanCell}` + ); + } + return parts; + } + + function fusionTables(temporal) { + if (!temporal || typeof temporal !== "object") { + return `

No fusion temporal data

`; + } + const itbRows = FUSION_BTI_ROWS.map(([k, s, sig]) => [`${k}_itb`, s, sig]); + const allKeys = [...FUSION_BTI_ROWS.map((r) => r[0]), ...itbRows.map((r) => r[0])]; + let showMean = false; + for (const k of allKeys) { + if (fusionMeanResidual(temporal[k]) !== null) { + showMean = true; + break; + } + } + const btiBody = fusionSubTableRows(temporal, FUSION_BTI_ROWS, showMean); + const itbBody = fusionSubTableRows(temporal, itbRows, showMean); + if (!btiBody.length && !itbBody.length) { + return `

No fusion scenarios in temporal

`; + } + const meanTh = showMean ? `Mean resid.` : ""; + const head = `Setting${DISPLAY_METRIC_COLS.map((c) => `${DISPLAY_METRIC_LABELS[c]}`).join("")}${meanTh}`; + + let h = `
`; + if (btiBody.length) { + h += `

Bands-then-Index (BtI)

`; + h += `${head}${btiBody.join("")}
`; + } + if (itbBody.length) { + h += `

Index-then-Bands (ItB)

`; + h += `${head}${itbBody.join("")}
`; + } + h += `
`; + return h; + } + + /** Returns only <table>…</table> or empty string (no heading). */ + function baselineTable(b) { if (!b || typeof b !== "object") return ""; const rows = []; const pushRow = (label, m) => { @@ -111,8 +180,74 @@ pushRow(`S2 Whittaker λ=400 (${strat})`, b.s2_whittaker_lambda400?.[strat]); } if (!rows.length) return ""; - const head = `Baseline${DISPLAY_METRIC_COLS.map((c) => `${DISPLAY_METRIC_LABELS[c]}`).join("")}`; - return `

Baselines (temporal vs PhenoCam)

${head}${rows.join("")}
`; + const head = `Baseline${DISPLAY_METRIC_COLS.map((c) => `${DISPLAY_METRIC_LABELS[c]}`).join("")}`; + return `${head}${rows.join("")}
`; + } + + function fmtFixed3(v) { + const n = Number(v); + return Number.isFinite(n) ? n.toFixed(3) : "—"; + } + + function derivedSection(d) { + if (!d) return ""; + const dn = d.delta_nse_pc_sigma20_minus_sigma30; + const paired = d.bti_vs_itb_mean_residual || []; + if (!dn && !paired.length) return ""; + + let h = `

Summaries

`; + h += `

Same numbers as Fusion, condensed. First table: which σ fits PhenoCam better (NSE_PC only). Second: mean bias BtI vs ItB.

`; + if (dn) { + h += `

ΔNSE_PC = NSE_PC(σ20) − NSE_PC(σ30). + → σ20 better. → σ30 better.

`; + h += ``; + let anyDelta = false; + for (const mode of ["bti", "itb"]) { + for (const strat of ["aggressive", "nonaggressive"]) { + const v = dn[mode]?.[strat]; + if (Number.isFinite(Number(v))) anyDelta = true; + h += ``; + } + } + h += `
ModeStrategyΔNSE_PC
${mode.toUpperCase()}${strat}${fmtFixed3(v)}
`; + if (!anyDelta) { + h += `

ΔNSE_PC needs both σ20 and σ30 fusion rows in temporal (BtI and ItB). Re-run metrics_stats.

`; + } + } + if (paired.length) { + h += `

Mean(fused − PhenoCam) per row. + / = average over / under PhenoCam. Closer to 0 in a column = less bias for that workflow.

`; + h += ``; + for (const row of paired) { + h += ``; + } + h += `
StrategyσMean residual BtIMean residual ItB
${row.strategy}${row.sigma}${fmtFixed3(row.mean_residual_bti)}${fmtFixed3(row.mean_residual_itb)}
`; + } + return h; + } + + function howToReadBlock() { + return `
+ How to read +
    +
  1. All scores are satellite or fusion GCC vs PhenoCam GCC at the site 3×3 window, same calendar days only. Extra stats: metrics.json.
  2. +
  3. , NSE_PC: higher = better. nRMSE: lower = better.
  4. +
  5. Fusion: same row number in BtI and in ItB = same screening + same σ — compare left/right. Down one block = change screening or σ.
  6. +
  7. Mean resid. (if present): mean(fused − PhenoCam). Sign = average bias; use R² / nRMSE / NSE_PC for overall fit.
  8. +
  9. Summaries: ΔNSE_PC = NSE at σ20 minus NSE at σ30 (+ means σ20 wins). Paired table: closer to 0 = less mean bias.
  10. +
+
`; + } + + function definitionsDetails() { + return `
+ Definitions + +
`; } function render(data) { @@ -122,35 +257,35 @@ return; } let html = ""; - html += ` -
- All metrics compare a greenness index (GCC) from satellite products against PhenoCam - ground-truth GCC at the site's 3×3 pixel window. - -
`; + html += `
+ GCC at the 3×3 site window vs PhenoCam. Sections: PhenoCam → baselines → fusion (BtI, then ItB) → summaries. + data/${siteName}/${season}/metrics.json +
`; + html += howToReadBlock(); + if (data.phenocam_stats) { - html += `

PhenoCam

`; - html += `

Summary statistics of the PhenoCam GCC timeseries used as ground truth for this site and season.

`; + html += `

PhenoCam (ground truth)

`; + html += `

Camera ROI GCC (not compared to itself). Dates / SOS–EOS: Phenology.

`; html += ``; const p = data.phenocam_stats; html += `
meanstdminmaxn
${fmt(p.mean)}${fmt(p.std)}${fmt(p.min)}${fmt(p.max)}${fmt(p.n_samples)}
`; } - const baselineHtml = baselineSection(data.baseline); - if (baselineHtml) { - html += `

Baselines (temporal vs PhenoCam)

`; - html += `

Reference GCC series before any fusion: raw S2 (all dates and cloud-screened per strategy), S3 composite per strategy, and a Whittaker-smoothed S2 series (λ=400). Useful to see what fusion has to beat.

`; - html += baselineHtml.replace(/^

[^<]*<\/h2>/, ""); + + const baselineTbl = baselineTable(data.baseline); + if (baselineTbl) { + html += `

Baselines (vs PhenoCam)

`; + html += `

Same columns as fusion (vs PhenoCam). Higher R² / NSE_PC, lower nRMSE = better. S3 = coarse-only; Whittaker = smoothed S2-only.

`; + html += baselineTbl; } - html += `

Temporal (vs PhenoCam)

`; - html += `

Per-scenario agreement between the fusion GCC timeseries at the site 3×3 window and the PhenoCam GCC timeseries, across all matched dates. Scenarios ending in _itb are Index-then-Bands; the others are Bands-then-Index.

`; - html += tableSection("", data.temporal); + + html += `

Fusion (vs PhenoCam)

`; + html += `

BtI block vs ItB block: same row = same screening + σ. Within a block: four EFAST combinations.

`; + html += fusionTables(data.temporal || {}); + + html += derivedSection(data.derived); + + html += definitionsDetails(); + el.innerHTML = html || `

Empty metrics file

`; }