#!/usr/bin/env catnip
# Calibration d'une sonde thermique par ajustement non linéaire.
#
# La sonde plongée dans un bain répond avec une montée exponentielle :
# T(t) = T_bain + (T0 - T_bain) * exp(-t/tau). Trois paramètres sont ajustés
# (T_bain, amplitude, constante de temps tau) sur des mesures bruitées.
# Catnip porte le protocole (valeurs vraies, tolérances de fiche technique,
# verdict) ; SciPy ajuste et fournit la covariance des paramètres.
#
# Contrôles : chaque paramètre ajusté doit tomber à moins de 3 sigma de sa
# valeur vraie, et les résidus rester compatibles avec le bruit injecté.
#
# DEPS: numpy scipy matplotlib
# OFFLINE: données synthétiques déterministes

numpy = import('numpy')
import('scipy.optimize', 'curve_fit')
mpl = import('matplotlib')
mpl.use('Agg')
plt = import('matplotlib.pyplot')
import('pathlib', 'Path')

script_dir = Path(META.file).parent
output_dir = script_dir / 'output'
output_dir.mkdir(exist_ok=True)

# Vérité terrain : bain à 80 °C, sonde initialement à 20 °C, tau = 45 s.

T_BATH = 80.0
T0 = 20.0
TAU = 45.0
NOISE = 0.3

rng = numpy.random.default_rng(42)
t = numpy.linspace(0.0, 300.0, 60)
measured = T_BATH + (T0 - T_BATH) * numpy.exp(-t / TAU) + rng.normal(0.0, NOISE, len(t))

# Modèle transmis au solveur : lambda Catnip vectorisée via numpy.

model = (t, t_bath, delta, tau) => { t_bath + delta * numpy.exp(-t / tau) }

# Point de départ volontairement décalé : 70 °C, -50 °C, 30 s.

popt, pcov = curve_fit(model, t, measured, list(70.0, -50.0, 30.0))
perr = numpy.sqrt(numpy.diag(pcov))

union Verdict {
    nominal; suspect

    label(self): str => {
        match self {
            Verdict.nominal => { "nominal" }
            Verdict.suspect => { "suspect" }
        }
    }
}

struct Param {
    name: str; truth: float; fitted: float; sigma: float

    bias_sigmas(self): float => { abs(self.fitted - self.truth) / self.sigma }

    verdict(self): Verdict => {
        if self.bias_sigmas() < 3.0 { Verdict.nominal } else { Verdict.suspect }
    }

    display(self): str => {
        f"{self.name:<7} vrai {self.truth:>7.2f}  ajusté {self.fitted:>8.3f} ± {self.sigma:<6.3f} biais {self.bias_sigmas():>5.2f}σ  {self.verdict().label()}"
    }
}

params = list(
    Param('T_bain', T_BATH, popt[0], perr[0]),
    Param('delta', T0 - T_BATH, popt[1], perr[1]),
    Param('tau', TAU, popt[2], perr[2]),
)

print("⇒ Paramètres ajustés (point de départ : 70, -50, 30)")
for p in params {
    print(f"  {p.display()}")
}

# Fiche technique : constante de temps annoncée 45 s ± 5 s.

tau_ok = abs(popt[2] - TAU) <= 5.0
print()
print(f"⇒ Fiche technique (tau = 45 s ± 5 s) : {tau_ok}")

# Résidus : doivent rester du même ordre que le bruit injecté (0.3 °C).

fitted_curve = model(t, popt[0], popt[1], popt[2])
residuals = measured - fitted_curve
rms = float(numpy.sqrt(numpy.mean(residuals * residuals)))
print(f"⇒ Résidus RMS : {round(rms, 3)} °C (bruit injecté : {NOISE} °C)")

oracle_ok = tau_ok and rms < 0.5
for p in params {
    if p.verdict() != Verdict.nominal { oracle_ok = False }
}

print()
print(f"⇒ Oracle (3σ sur chaque paramètre + résidus) : {oracle_ok}")

# Planche : mesures et courbe ajustée en haut, résidus en bas.

fig, axes = plt.subplots(2, 1, figsize=tuple(10, 7), sharex=True)

axes[0].scatter(t, measured, s=12, color='#4c78a8', label='mesures')
axes[0].plot(t, fitted_curve, color='#f58518', label=f"fit (tau = {round(popt[2], 1)} s)")
axes[0].set_ylabel('température (°C)')
axes[0].legend(loc='lower right')
axes[0].set_title("Calibration d'une sonde thermique (échelon 20 °C → 80 °C)")

axes[1].axhline(0.0, color='#888888', linewidth=0.8)
axes[1].scatter(t, residuals, s=12, color='#54a24b')
axes[1].set_xlabel('temps (s)')
axes[1].set_ylabel('résidu (°C)')

out_path = output_dir / 'curve_fit_calibration.png'
fig.savefig(str(out_path), dpi=130, bbox_inches='tight')
print(f"⇒ Planche écrite : {out_path}")

# --no-browser garde un chemin headless (le PNG reste écrit ci-dessus).
if '--no-browser' not in import('sys').argv {
    http = import('http')
    b64 = import('base64').b64encode(out_path.read_bytes()).decode('ascii')
    http.serve(f'<!doctype html><meta charset="utf-8"><title>Calibration sonde (curve_fit)</title><body style="margin:0;background:#0d1117"><img style="max-width:100%;display:block;margin:0 auto" src="data:image/png;base64,{b64}">',
        0, 'text/html; charset=utf-8', True)
}