From 6fc68fe1aa8bda617f20506ba9288c3035f56208 Mon Sep 17 00:00:00 2001
From: Lars Bogner <upyxo@student.kit.edu>
Date: Mon, 22 Jun 2026 14:41:45 +0200
Subject: [PATCH] Export plots for knowledge base

---
 analysis/export_photon_edep.py      | 40 +++++++++++++++
 analysis/export_validation_plots.py | 80 +++++++++++++++++++++++++++++
 2 files changed, 120 insertions(+)
 create mode 100644 analysis/export_photon_edep.py
 create mode 100644 analysis/export_validation_plots.py

diff --git a/analysis/export_photon_edep.py b/analysis/export_photon_edep.py
new file mode 100644
index 0000000..c940c85
--- /dev/null
+++ b/analysis/export_photon_edep.py
@@ -0,0 +1,40 @@
+"""Zoomed-in real-vs-generated edep histogram for photons only (pdg=22),
+to characterize the KL spike flagged by plot_kl_bars_pl(group_by='pdg')."""
+
+from pathlib import Path
+
+import matplotlib.pyplot as plt
+import numpy as np
+
+import giant.analysis as a
+
+FILE = "/home/lars/Programming/giant/pbwo4_10k_9_predicted_local.parquet"
+OUT = Path("/home/lars/knowledge-base/meta/attachments")
+
+samples = a.load_predicted_local(FILE, sample_frac=0.15)
+mask = samples.pdg == 22
+edep_idx = a.RAW_TARGET_NAMES.index("edep")
+real = samples.real_raw[mask, edep_idx]
+gen = samples.gen_raw[mask, edep_idx]
+print("n photon rows:", mask.sum())
+print("real: mean", real.mean(), "std", real.std(), "max", real.max(), "frac==0", (real == 0).mean())
+print("gen:  mean", gen.mean(), "std", gen.std(), "max", gen.max(), "frac==0", (gen == 0).mean())
+for q in [0.5, 0.9, 0.99, 0.999]:
+    print(f"q={q}: real={np.quantile(real, q):.4f} gen={np.quantile(gen, q):.4f}")
+
+fig, axes = plt.subplots(1, 2, figsize=(10, 4))
+bins = np.linspace(0, np.quantile(real, 0.999), 80)
+axes[0].hist(real, bins=bins, alpha=0.6, label="real", density=True)
+axes[0].hist(gen, bins=bins, alpha=0.6, label="gen", density=True)
+axes[0].set_yscale("log")
+axes[0].set_xlabel("edep (photons, pdg=22)")
+axes[0].legend()
+
+axes[1].hist(real, bins=bins, alpha=0.6, label="real", density=True, cumulative=True, histtype="step")
+axes[1].hist(gen, bins=bins, alpha=0.6, label="gen", density=True, cumulative=True, histtype="step")
+axes[1].set_xlabel("edep (photons, pdg=22) - CDF")
+axes[1].legend()
+
+fig.tight_layout()
+fig.savefig(OUT / "giant-h1024n8d0.1lr3e-4-photon-edep-zoom.png", dpi=150, bbox_inches="tight")
+print("saved photon-edep-zoom")
diff --git a/analysis/export_validation_plots.py b/analysis/export_validation_plots.py
new file mode 100644
index 0000000..1a08d15
--- /dev/null
+++ b/analysis/export_validation_plots.py
@@ -0,0 +1,80 @@
+"""One-off export of validation plots for checkpoints/scan/h1024_n8_d0.1_lr0.0003/best.pt
+into the knowledge-base attachments folder. Not part of the package; run manually."""
+
+from pathlib import Path
+
+import giant.analysis as a
+
+FILE = "/home/lars/Programming/giant/pbwo4_10k_9_predicted_local.parquet"
+OUT = Path("/home/lars/knowledge-base/meta/attachments")
+PREFIX = "giant-h1024n8d0.1lr3e-4"
+
+PDG_NAMES = {
+    11: "e-",
+    -11: "e+",
+    22: "gamma",
+    2112: "n",
+    2212: "p",
+}
+
+
+def pdg_label(code: int) -> str:
+    if code in PDG_NAMES:
+        return PDG_NAMES[code]
+    if code > 1000000000:
+        return f"ion{code}"
+    return str(code)
+
+
+print("=== KL bar plots (lazy, full dataset) ===")
+for grouping in [None, "energy", "pdg", "material"]:
+    fig = a.plot_kl_bars_pl(FILE, group_by=grouping)
+    name = f"{PREFIX}-kl-bars-{grouping or 'all'}.png"
+    fig.savefig(OUT / name, dpi=150, bbox_inches="tight")
+    print("saved", name)
+
+print("=== loading sampled SampleCollection ===")
+samples = a.load_predicted_local(FILE, sample_frac=0.15)
+print("n rows:", len(samples.gen_raw))
+
+print("=== marginals ===")
+fig = a.plot_marginals(samples)
+fig.savefig(OUT / f"{PREFIX}-marginals-all.png", dpi=150, bbox_inches="tight")
+
+fig = a.plot_marginals(samples, group_by="energy")
+fig.savefig(OUT / f"{PREFIX}-marginals-energy.png", dpi=150, bbox_inches="tight")
+
+fig = a.plot_marginals(samples, group_by="pdg")
+fig.savefig(OUT / f"{PREFIX}-marginals-pdg.png", dpi=150, bbox_inches="tight")
+print("saved marginals")
+
+print("=== correlation matrices ===")
+fig = a.plot_correlation_matrices(samples)
+fig.savefig(OUT / f"{PREFIX}-correlation.png", dpi=150, bbox_inches="tight")
+
+print("=== pairwise ===")
+fig = a.plot_pairwise(samples, n_sample=5000)
+fig.savefig(OUT / f"{PREFIX}-pairwise.png", dpi=150, bbox_inches="tight")
+
+print("=== direction alignment ===")
+fig = a.plot_direction_alignment(samples)
+fig.savefig(OUT / f"{PREFIX}-direction-alignment.png", dpi=150, bbox_inches="tight")
+
+print("=== constraint violations ===")
+fig = a.plot_constraint_violations(samples)
+fig.savefig(OUT / f"{PREFIX}-constraint-violations.png", dpi=150, bbox_inches="tight")
+
+print("=== marginal_table aggregate ===")
+agg = a.marginal_table(samples)
+print(agg.to_string())
+
+print("=== marginal_table by pdg (top rows incl. photon) ===")
+by_pdg = a.marginal_table(samples, group_by="pdg")
+by_pdg["particle"] = by_pdg["group"].astype(str)
+print(by_pdg.to_string())
+
+print("=== photon-only rows ===")
+photon_rows = by_pdg[by_pdg["group"].astype(str) == "pdg=22"]
+print(photon_rows.to_string())
+
+print("DONE")