gen_data/genera.py

# /// script
# requires-python = ">=3.11"
# dependencies = ["numpy"]
# ///
"""Genera un CSV con T e ΔT = T − T_amb con rumore realistico (vedi docs/SPEC.md)."""

import argparse
import re
import secrets
import sys
import tomllib
from pathlib import Path

import numpy as np

CARTELLA_OUTPUT = Path("output")
GRIGLIA_MS = 5.0  # passo della griglia fine per il filtro del sensore

DESCRIZIONE = """\
Genera un CSV con la temperatura T di un sensore IR e la differenza dT = T - T_amb,
con rumore, quantizzazione e jitter dei tempi simili a un sensore reale.

T parte dalla temperatura ambiente, sale con uno o più gradini di riscaldamento e
ridiscende con raffreddamento esponenziale fino a T_amb. Ogni corsa (seed) è diversa:
i parametri principali sono intervalli [min, max] sorteggiati a ogni corsa."""

EPILOGO = """\
esempi:
  genera.py config.toml                    scrive output/gen_<seed>.csv, seed casuale
  genera.py config.toml --seed 42          riproducibile: output/gen_42.csv
  genera.py config.toml -o dati/prova.csv  percorso a scelta
  genera.py config.toml -o - > out.csv     su stdout

formato dell'output (separatore ',', fine riga '\\n'):
  t_ms,T_C,dT_C
  151,22.90,0.20
  t_ms è intero; T_C e dT_C hanno due decimali. dT è calcolato dopo rumore e
  quantizzazione di entrambi i canali (T e ambiente).

cosa fa (in ordine):
  1. tempi         t_i = t_(i-1) + dt_ms +/- jitter_ms, fino a durata_s
  2. T vera        ambiente vero + gradini di riscaldamento (numero, istanti, ampiezze
                   sorteggiati) con raffreddamento esponenziale + disturbi a salto
  3. sensore       filtro del 1° ordine (tau_sensore_s) su griglia fine da 5 ms
  4. rumore di T   bianco (più forte a caldo) + deriva lenta, poi quantizzazione
  5. ambiente, dT  ambiente con rumore e quantizzazione propri; dT = T - ambiente

gruppi di chiavi di config.toml (unica fonte dei valori, unità nel suffisso:
_s secondi, _ms millisecondi, _c gradi C, _hz 1/s; ogni chiave è obbligatoria).
Le chiavi marcate [min, max] sono intervalli: a ogni corsa ne viene sorteggiato un
valore uniforme (n_gradini: intero, estremi inclusi):
  tempi            durata_s, dt_ms, jitter_ms
  ambiente         amb_inf_c [min, max], amb_A_c [min, max], amb_tau_s [min, max],
                   amb_sigma_c, amb_passo_c
  rumore di T      sigma_c [min, max], sigma_caldo_c [min, max], x_rif_c,
                   wander_sigma_c [min, max], wander_tau_s [min, max], passo_T_c
  gradini          n_gradini [min, max], t_primo_gradino_s [min, max],
                   distanza_gradini_s [min, max] (tra gradini consecutivi),
                   dT_gradino_c [min, max] (ampiezza di ogni gradino),
                   k_raffreddamento_hz [min, max], tau_sensore_s [min, max]
  salti            salti_rate_hz [min, max], salti_media_c, salti_sigma_c [min, max],
                   salti_tau_s, salti_x_min_c, salti_max_c (raffreddamento);
                   salti_risc_rate_hz [min, max], salti_risc_media_c,
                   salti_risc_sigma_c [min, max], salti_risc_tau_s (riscaldamento)

codici di uscita:
  0  ok
  1  errore, con messaggio in italiano su stderr: file di configurazione illeggibile
     o TOML non valido, chiave mancante, intervallo [min, max] non valido,
     file di uscita non scrivibile
  2  argomenti della riga di comando non validi

altri strumenti: grafico.py (disegna un CSV), pulisci.sh (svuota output/).
documentazione: docs/SPEC.md (modello e controlli), docs/FUNZIONALITA.md (funzionalità)."""

# chiavi che sono intervalli [min, max], sorteggiati uniformi a ogni corsa
INTERVALLI = [
    "amb_inf_c", "amb_A_c", "amb_tau_s",
    "sigma_c", "sigma_caldo_c", "wander_sigma_c", "wander_tau_s",
    "k_raffreddamento_hz", "tau_sensore_s",
    "salti_rate_hz", "salti_sigma_c", "salti_risc_rate_hz", "salti_risc_sigma_c",
]
# intervalli che descrivono i gradini (istanti e ampiezze derivati da questi)
INTERVALLI_GRADINI = ["n_gradini", "t_primo_gradino_s", "distanza_gradini_s", "dT_gradino_c"]

CHIAVI = [
    "durata_s", "dt_ms", "jitter_ms",
    "amb_inf_c", "amb_A_c", "amb_tau_s", "amb_sigma_c", "amb_passo_c",
    "sigma_c", "sigma_caldo_c", "x_rif_c", "wander_sigma_c", "wander_tau_s", "passo_T_c",
    "n_gradini", "t_primo_gradino_s", "distanza_gradini_s", "dT_gradino_c",
    "k_raffreddamento_hz", "tau_sensore_s",
    "salti_rate_hz", "salti_media_c", "salti_sigma_c", "salti_tau_s", "salti_x_min_c",
    "salti_max_c", "salti_risc_rate_hz", "salti_risc_media_c", "salti_risc_sigma_c", "salti_risc_tau_s",
]


class ParserIt(argparse.ArgumentParser):
    """ArgumentParser con i messaggi d'errore più comuni in italiano (exit 2)."""

    TRADUZIONI = [
        (r"the following arguments are required: (.*)", r"manca l'argomento obbligatorio: \1"),
        (r"unrecognized arguments: (.*)", r"argomenti non riconosciuti: \1"),
        (r"argument (.*): invalid int value: (.*)", r"argomento \1: \2 non è un numero intero"),
        (r"argument (.*): expected one argument", r"argomento \1: manca il valore"),
    ]

    def error(self, message):
        for pattern, sostituzione in self.TRADUZIONI:
            message = re.sub(pattern, sostituzione, message)
        self.print_usage(sys.stderr)
        print(f"{self.prog}: errore: {message}", file=sys.stderr)
        sys.exit(2)


def errore(msg):
    print(f"Errore: {msg}", file=sys.stderr)
    sys.exit(1)


def leggi_config(percorso):
    try:
        with open(percorso, "rb") as f:
            cfg = tomllib.load(f)
    except OSError as e:
        errore(f"impossibile leggere '{percorso}': {e.strerror}")
    except tomllib.TOMLDecodeError as e:
        errore(f"TOML non valido in '{percorso}': {e}")
    mancanti = [k for k in CHIAVI if k not in cfg]
    if mancanti:
        errore(f"chiavi mancanti in '{percorso}': {', '.join(mancanti)}")
    for k in INTERVALLI + INTERVALLI_GRADINI:
        v = cfg[k]
        if (not isinstance(v, list) or len(v) != 2
                or not all(isinstance(x, (int, float)) and not isinstance(x, bool) for x in v)
                or v[0] > v[1]):
            errore(f"'{k}' deve essere un intervallo [min, max] con min <= max")
    n_min, n_max = cfg["n_gradini"]
    if not (isinstance(n_min, int) and isinstance(n_max, int) and n_min >= 1):
        errore("'n_gradini' deve essere un intervallo [min, max] di interi con min >= 1")
    return cfg


def eccesso(t_s, cfg):
    """Eccesso termico x(t): un gradino per impulso, poi raffreddamento esponenziale."""
    x = np.zeros_like(t_s)
    for tj, dj in zip(cfg["t_impulsi_s"], cfg["dT_impulsi_c"]):
        dopo = t_s >= tj
        x[dopo] += dj * np.exp(-cfg["k_raffreddamento_hz"] * (t_s[dopo] - tj))
    return x


def sorteggia(cfg, rng):
    """Copia di cfg con ogni intervallo sostituito da un valore sorteggiato e i gradini generati (ogni corsa è diversa)."""
    valori = {k: rng.uniform(*cfg[k]) for k in INTERVALLI}
    n = int(rng.integers(cfg["n_gradini"][0], cfg["n_gradini"][1] + 1))
    t_primo = rng.uniform(*cfg["t_primo_gradino_s"])
    t = t_primo + np.concatenate(([0.0], np.cumsum(rng.uniform(*cfg["distanza_gradini_s"], n - 1))))
    d = rng.uniform(*cfg["dT_gradino_c"], n)
    return {**cfg, **valori, "t_impulsi_s": t, "dT_impulsi_c": d}


def eventi(g_s, cfg, rng, t_ini, t_fin, rate, media, sigma, tau, x_min):
    """Scarti casuali (Poisson) in [t_ini, t_fin], con rientro esponenziale, se x ≥ x_min."""
    s = np.zeros_like(g_s)
    n = rng.poisson(rate * (t_fin - t_ini))
    ts = np.sort(rng.uniform(t_ini, t_fin, n))
    ampiezze = np.clip(rng.normal(media, sigma, n), -cfg["salti_max_c"], cfg["salti_max_c"])
    for tj, aj in zip(ts, ampiezze):
        if eccesso(np.array([tj]), cfg)[0] < x_min:
            continue
        dopo = g_s >= tj
        s[dopo] += aj * np.exp(-(g_s[dopo] - tj) / tau)
    return s


def salti(g_s, cfg, rng):
    """Disturbi a salto: durante il riscaldamento (tra primo e ultimo gradino) e il raffreddamento."""
    if len(cfg["t_impulsi_s"]) == 0:
        return np.zeros_like(g_s)
    t_primo, t_ultimo = cfg["t_impulsi_s"][0], cfg["t_impulsi_s"][-1]
    risc = eventi(g_s, cfg, rng, t_primo, t_ultimo, cfg["salti_risc_rate_hz"],
                  cfg["salti_risc_media_c"], cfg["salti_risc_sigma_c"], cfg["salti_risc_tau_s"], 0.0)
    raff = eventi(g_s, cfg, rng, t_ultimo, cfg["durata_s"], cfg["salti_rate_hz"],
                  cfg["salti_media_c"], cfg["salti_sigma_c"], cfg["salti_tau_s"], cfg["salti_x_min_c"])
    return risc + raff


def ambiente_vero(t_s, cfg):
    return cfg["amb_inf_c"] + cfg["amb_A_c"] * np.exp(-t_s / cfg["amb_tau_s"])


def quantizza(v, passo):
    return np.round(v / passo) * passo


def genera(cfg, rng):
    cfg = sorteggia(cfg, rng)
    # 1. tempi
    n_max = int(cfg["durata_s"] * 1000 / (cfg["dt_ms"] - cfg["jitter_ms"])) + 2
    dt = cfg["dt_ms"] + rng.uniform(-cfg["jitter_ms"], cfg["jitter_ms"], n_max)
    t_ms = np.cumsum(dt)
    t_ms = t_ms[t_ms <= cfg["durata_s"] * 1000]
    dt = np.diff(t_ms, prepend=0.0)
    t_s = t_ms / 1000

    # 2. T vera su griglia fine = ambiente vero + eccesso termico + salti
    g_s = np.arange(0.0, cfg["durata_s"] + GRIGLIA_MS / 1000, GRIGLIA_MS / 1000)
    t_vera = ambiente_vero(g_s, cfg) + eccesso(g_s, cfg) + salti(g_s, cfg, rng)

    # 3. risposta del sensore: filtro del 1° ordine, interpolato ai t_i
    a = np.exp(-GRIGLIA_MS / 1000 / cfg["tau_sensore_s"])
    y = np.empty_like(t_vera)
    y[0] = t_vera[0]
    for n in range(1, len(y)):
        y[n] = a * y[n - 1] + (1 - a) * t_vera[n]
    y_i = np.interp(t_s, g_s, y)

    # 4. rumore di T: bianco con σ(x) + deriva Ornstein-Uhlenbeck, poi quantizzazione
    x_i = eccesso(t_s, cfg)
    sigma = cfg["sigma_c"] + (cfg["sigma_caldo_c"] - cfg["sigma_c"]) * np.minimum(1.0, x_i / cfg["x_rif_c"])
    phi = np.exp(-(dt / 1000) / cfg["wander_tau_s"])
    inn = rng.standard_normal(len(t_s))
    w = np.empty_like(t_s)
    w[0] = cfg["wander_sigma_c"] * inn[0]
    for i in range(1, len(w)):
        w[i] = phi[i] * w[i - 1] + cfg["wander_sigma_c"] * np.sqrt(1 - phi[i] ** 2) * inn[i]
    T = quantizza(y_i + w + sigma * rng.standard_normal(len(t_s)), cfg["passo_T_c"])

    # 5. ambiente misurato e ΔT (dopo rumore e quantizzazione di entrambi i canali)
    amb = quantizza(
        ambiente_vero(t_s, cfg) + cfg["amb_sigma_c"] * rng.standard_normal(len(t_s)),
        cfg["amb_passo_c"],
    )
    dT = quantizza(T - amb, cfg["passo_T_c"])
    return np.round(t_ms).astype(int), T, dT


def scrivi(f, t_ms, T, dT):
    f.write("t_ms,T_C,dT_C\n")
    for t, a, b in zip(t_ms, T, dT):
        # "+ 0.0" evita "-0.00"
        f.write(f"{t},{round(a, 2) + 0.0:.2f},{round(b, 2) + 0.0:.2f}\n")


def main():
    ap = ParserIt(
        prog="genera.py",
        usage="genera.py [-h] [-o FILE] [--seed N] config",
        description=DESCRIZIONE,
        epilog=EPILOGO,
        formatter_class=argparse.RawDescriptionHelpFormatter,
        add_help=False,
    )
    argomenti = ap.add_argument_group("argomenti")
    opzioni = ap.add_argument_group("opzioni")
    opzioni.add_argument("-h", "--help", action="help", help="mostra questo aiuto ed esce")
    argomenti.add_argument("config", help="file di configurazione TOML (es. config.toml); tutte le chiavi sono obbligatorie")
    opzioni.add_argument(
        "-o", "--output", metavar="FILE",
        help="file CSV di uscita; le cartelle mancanti vengono create. '-' scrive su stdout. "
             "Default: output/gen_<seed>.csv",
    )
    opzioni.add_argument(
        "--seed", type=int, metavar="N",
        help="seme del generatore casuale: stesso seed e stessa configurazione danno lo stesso file. "
             "Default: casuale, stampato su stderr",
    )
    args = ap.parse_args()

    cfg = leggi_config(args.config)
    seed = args.seed if args.seed is not None else secrets.randbelow(2**32)
    if args.seed is None:
        print(f"seed: {seed}", file=sys.stderr)

    t_ms, T, dT = genera(cfg, np.random.default_rng(seed))

    if args.output == "-":
        scrivi(sys.stdout, t_ms, T, dT)
        return
    percorso = Path(args.output) if args.output else CARTELLA_OUTPUT / f"gen_{seed}.csv"
    try:
        percorso.parent.mkdir(parents=True, exist_ok=True)
        with open(percorso, "w", newline="\n") as f:
            scrivi(f, t_ms, T, dT)
    except OSError as e:
        errore(f"impossibile scrivere '{percorso}': {e.strerror}")
    print(f"scritto: {percorso}", file=sys.stderr)


if __name__ == "__main__":
    main()