diff --git a/contrib/backtest_forecast.py b/contrib/backtest_forecast.py new file mode 100644 index 0000000..f14b8d1 --- /dev/null +++ b/contrib/backtest_forecast.py @@ -0,0 +1,412 @@ +#!/usr/bin/env python3 +"""Score the thermal forecast against what the handle actually did, over +every recorded session. + +Three model changes in one evening were each checked by replaying one or +two hand-picked moments, and each found its counter-example within the +hour. This replays the whole history and scores every way the forecast can +be made against the plateau the handle actually reached, per scenario, so +a change to the model is judged on all of it at once. + +Two questions, scored separately: + +1. **In-run forecast.** While the current holds steady, how far off is the + projected plateau — the live forecast — as a function of how long it has + held? Trajectory basis (the run's own samples so far, the install's tau) + against model basis (an ambient plus the current law). + +2. **Cross-current prediction.** At every current change, and at every + session start, what would the plateau at the *new* current have been + predicted as from each ambient on offer — the LAN sensor, the ambient the + previous run's trajectory implies, the warmer of the two, the idle handle + before the session — under each current law: the I^2 prior, the + install's fitted exponent, and the fitted exponent with the ambient + term, each fitted with the scored session left out. This is the number + a restore, or the end of a calibration probe, is decided on. + +Ground truth is the plateau observed in a run that held its current for at +least --observe-tau time constants: an exponential fitted to the whole run +with tau free, the same fit the model's own parameters come from. Shorter +runs are still inputs (they imply an ambient) but never truth. + +Errors are predicted minus actual. Negative is optimistic — the handle ran +hotter than promised — which is the direction that trips the charger. + +Read-only. Point it at a copy of wallmonitor.db, or at the live file: SQLite +WAL readers do not block the writer. + +Example: + uv run python contrib/backtest_forecast.py --db ../wallmonitor.db + uv run python contrib/backtest_forecast.py --db ../wallmonitor.db --json out.json +""" + +from __future__ import annotations + +import argparse +import json +import sys +import time +from dataclasses import asdict, dataclass, field +from statistics import median + +from wallmonitor import thermal +from wallmonitor.db import Database + +MIN_RUN_S = 120.0 +TICK_S = 30.0 +FIRST_TICK_S = 120.0 +LAST_TICK_S = 1200.0 +OPTIMISTIC_C = 2.0 # a plateau under-read by this much is the dangerous kind of wrong +PROBE_CURRENT_FRAC = 0.75 +PROBE_MIN_S = 1800.0 +HOT_AMBIENT_C = 29.0 +WARM_START_C = 3.0 # handle this far above its idle level at a session's first run +BUCKETS_MIN = ((2, 5), (5, 10), (10, 20)) + + +@dataclass +class Run: + session_id: int + index: int + start_ts: float + end_ts: float + current_a: float + samples: list[tuple[float, float]] = field(repr=False) + handle_start_c: float = 0.0 + plateau_c: float | None = None # truth: see --truth + plateau_rmse_c: float | None = None + fit_c: float | None = None # exponential over the whole run, tau free + fit_tau_min: float | None = None + window_c: float | None = None # the same fit over the first PREFIX_SPAN_MIN_S — what the model's fitter sees + last_c: float | None = None # median handle over the run's last 3 min — model-free + max_c: float | None = None + sensor_ambient_c: float | None = None # LAN/car sensor at the run's start + idle_ambient_c: float | None = None # from the idle handle before the session (first run) + implied_ambient_c: float | None = None # this run's own trajectory, through the current law + kind: str = "" # cold_start | warm_start | step_down | step_up + probe: bool = False + hot: bool = False + + @property + def span_s(self) -> float: + return self.end_ts - self.start_ts + + +@dataclass +class TickScore: + session_id: int + kind: str + probe: bool + hot: bool + minutes: float + basis: str + error_c: float + + +@dataclass +class BoundaryScore: + session_id: int + kind: str + probe: bool + hot: bool + from_a: float | None + to_a: float + actual_c: float + predictions: dict[str, float] # "/" -> predicted plateau + + +# --------------------------------------------------------------------------- + + +def split_runs(rows: list[dict]) -> list[list[dict]]: + """Steady-current runs, the way the live forecast sees them: a run ends + when current leaves the band around its first sample, when the + contactor opens, or after a charging gap. Ramp-up samples fall into + runs too short to keep.""" + runs: list[list[dict]] = [] + current: list[dict] = [] + for row in rows: + amps = row.get("vehicle_current_a") or 0.0 + if not row.get("contactor_closed") or amps < 6.0 or row.get("handle_temp_c") is None: + if current: + runs.append(current) + current = [] + continue + if current: + ref = current[0]["vehicle_current_a"] + if abs(amps - ref) > max(2.0, 0.1 * ref) or row["ts"] - current[-1]["ts"] > thermal.SEGMENT_SPLIT_GAP_S: + runs.append(current) + current = [] + current.append(row) + if current: + runs.append(current) + return [run for run in runs if run[-1]["ts"] - run[0]["ts"] >= MIN_RUN_S] + + +def observe_plateau(run: Run, tau_min: float, observe_tau: float, truth: str) -> None: + """Fill in what the run actually reached. `truth` picks which reading + becomes plateau_c: "fit" — an exponential over the whole run with tau + free (needs >= observe_tau time constants; rejects a poor fit); "last" — + the handle's own median over the last three minutes, model-free (needs + >= observe_tau + 2 time constants, so the lag has worked itself out).""" + samples = run.samples + span = run.span_s + run.max_c = max(temp for _, temp in samples) + tail = [temp for ts, temp in samples if samples[-1][0] - ts <= 180.0] + run.last_c = median(tail) if tail else None + fit = thermal._fit_exponential(samples) + if fit is not None: + run.fit_tau_min, run.fit_c, run.plateau_rmse_c = fit[0] / 60.0, fit[1], fit[2] + head = [(ts, temp) for ts, temp in samples if ts - samples[0][0] <= thermal.PREFIX_SPAN_MIN_S] + if len(head) >= thermal.MIN_SEGMENT_SAMPLES: + window_fit = thermal._fit_exponential(head) + run.window_c = window_fit[1] if window_fit is not None else None + if truth == "fit": + if span >= observe_tau * tau_min * 60.0 and fit is not None and fit[2] <= thermal.MAX_FIT_RMSE_C: + run.plateau_c = fit[1] + elif span >= (observe_tau + 2.0) * tau_min * 60.0 and run.last_c is not None: + run.plateau_c = run.last_c + + +LAWS = { + # name -> (exponent pinned?, ambient_coef pinned?) — None means fitted + "I2": dict(exponent=thermal.DEFAULT_CURRENT_EXP, ambient_coef=thermal.DEFAULT_AMBIENT_COEF), + "n": dict(ambient_coef=thermal.DEFAULT_AMBIENT_COEF), + "n+k": dict(), +} + + +def params_without(fits: list[dict], sid: int, **law) -> thermal.ThermalParams: + """Model parameters from every fit but the scored session's own, under a + current law with the given terms pinned (see LAWS).""" + others = [dict(fit) for fit in fits if fit["session_id"] != sid] + return thermal.params_from_fits(others, **law) + + +def build_runs(db: Database, sess: dict, params: thermal.ThermalParams, observe_tau: float, + idle_model: thermal.IdleOffset, truth: str = "fit") -> list[Run]: + rows = db.vitals_range(sess["start_ts"] - 1, sess["end_ts"] + 1, 500_000) + runs: list[Run] = [] + for index, raw in enumerate(split_runs(rows)): + samples = [(row["ts"], row["handle_temp_c"]) for row in raw] + run = Run( + session_id=sess["id"], index=index, start_ts=samples[0][0], end_ts=samples[-1][0], + current_a=median(row["vehicle_current_a"] for row in raw), samples=samples, + handle_start_c=samples[0][1], + ) + observe_plateau(run, params.tau_min, observe_tau, truth) + measured = thermal._measured_ambient(db, run.start_ts - thermal.MEASURED_AMBIENT_WINDOW_S, run.start_ts + 60) + run.sensor_ambient_c = measured[0] if measured is not None else None + if index == 0: + run.idle_ambient_c = thermal._ambient_before(db, sess["start_ts"], idle_model) + if len(samples) >= thermal.TRAJECTORY_MIN_SAMPLES: + t_inf, _se = thermal._project_t_inf(samples, params.tau_min) + run.implied_ambient_c = params.ambient_from_plateau(t_inf, run.current_a) + run.probe = run.current_a <= PROBE_CURRENT_FRAC * thermal.REF_CURRENT_A and run.span_s >= PROBE_MIN_S + ambient_for_hot = run.sensor_ambient_c if run.sensor_ambient_c is not None else run.idle_ambient_c + run.hot = ambient_for_hot is not None and ambient_for_hot >= HOT_AMBIENT_C + if index == 0: + baseline = ambient_for_hot + idle_level = thermal.idle_handle_c(baseline, idle_model) if baseline is not None else None + run.kind = ( + "warm_start" if idle_level is not None and run.handle_start_c > idle_level + WARM_START_C + else "cold_start" + ) + else: + run.kind = "step_down" if run.current_a < runs[-1].current_a else "step_up" + runs.append(run) + return runs + + +def score_ticks(run: Run, params: thermal.ThermalParams) -> list[TickScore]: + """The live forecast, tick by tick, against the plateau this run reached.""" + if run.plateau_c is None: + return [] + scores: list[TickScore] = [] + t0 = run.start_ts + model_sensor = ( + params.plateau_at(run.current_a, run.sensor_ambient_c) if run.sensor_ambient_c is not None else None + ) + tick = FIRST_TICK_S + while tick <= min(LAST_TICK_S, run.span_s): + window = [(ts, temp) for ts, temp in run.samples if ts - t0 <= tick] + minutes = tick / 60.0 + if len(window) >= thermal.TRAJECTORY_MIN_SAMPLES: + t_inf, _se = thermal._project_t_inf(window, params.tau_min) + scores.append(TickScore(run.session_id, run.kind, run.probe, run.hot, minutes, "trajectory", t_inf - run.plateau_c)) + if model_sensor is not None: + scores.append(TickScore(run.session_id, run.kind, run.probe, run.hot, minutes, "model/sensor", model_sensor - run.plateau_c)) + tick += TICK_S + return scores + + +def score_boundary(run: Run, prev: Run | None, laws: dict[str, thermal.ThermalParams]) -> BoundaryScore | None: + """What each ambient, under each current law, would have predicted for + the plateau at this run's current, from what was known when it started.""" + if run.plateau_c is None: + return None + ambients: dict[str, float] = {} + if run.sensor_ambient_c is not None: + ambients["sensor"] = run.sensor_ambient_c + if prev is None: + if run.idle_ambient_c is not None: + ambients["idle"] = run.idle_ambient_c + else: + if prev.implied_ambient_c is not None: + ambients["implied"] = prev.implied_ambient_c + if run.sensor_ambient_c is not None: + ambients["warmer"] = max(run.sensor_ambient_c, prev.implied_ambient_c) + if not ambients: + return None + predictions = { + f"{name}/{law}": params.plateau_at(run.current_a, ambient) + for name, ambient in ambients.items() + for law, params in laws.items() + } + return BoundaryScore( + run.session_id, run.kind, run.probe, run.hot, + prev.current_a if prev is not None else None, run.current_a, run.plateau_c, predictions, + ) + + +# --------------------------------------------------------------------------- + + +def _pct(values: list[float], p: float) -> float: + ordered = sorted(values) + return ordered[min(len(ordered) - 1, int(round(p * (len(ordered) - 1))))] + + +def _stats(errors: list[float]) -> str: + if not errors: + return f"{'-':>5} {'-':>6} {'-':>6} {'-':>6} {'-':>5}" + absolute = [abs(e) for e in errors] + optimistic = sum(1 for e in errors if e < -OPTIMISTIC_C) / len(errors) + return f"{len(errors):>5} {median(errors):>+6.2f} {median(absolute):>6.2f} {_pct(absolute, 0.9):>6.2f} {optimistic:>5.0%}" + + +STATS_HEADER = f"{'n':>5} {'bias':>6} {'|med|':>6} {'|p90|':>6} {'opt':>5}" + + +def report(runs: list[Run], ticks: list[TickScore], boundaries: list[BoundaryScore], + params: thermal.ThermalParams, observe_tau: float) -> None: + observed = [run for run in runs if run.plateau_c is not None] + print( + f"model: tau {params.tau_min:.2f} min, rise {params.rise_ref_c:.1f} C at {thermal.REF_CURRENT_A:g} A " + f"and {thermal.AMBIENT_REF_C:g} C, n = {params.current_exp:.2f}, k = {params.ambient_coef:+.3f} C/C " + f"({params.current_exp_fits} fits)" + ) + print( + f"runs: {len(runs)} steady-current runs in {len({run.session_id for run in runs})} sessions; " + f"{len(observed)} held >= {observe_tau:g} tau and are scored as truth" + ) + kinds = sorted({run.kind for run in observed}) + print(" by kind:", ", ".join(f"{kind} {sum(1 for r in observed if r.kind == kind)}" for kind in kinds), + f"| probe {sum(1 for r in observed if r.probe)}, hot {sum(1 for r in observed if r.hot)}") + print() + print("== In-run forecast: projected plateau vs observed, by minutes at steady current") + print(" (bias = median signed error, predicted - actual; opt = share optimistic by > " + f"{OPTIMISTIC_C:g} C)") + bases = sorted({tick.basis for tick in ticks}) + print(f"{'basis':<14}" + "".join(f" {lo:>2}-{hi:<2} min {STATS_HEADER}" for lo, hi in BUCKETS_MIN)) + for basis in bases: + row = f"{basis:<14}" + for lo, hi in BUCKETS_MIN: + errors = [t.error_c for t in ticks if t.basis == basis and lo <= t.minutes < hi] + row += f" {'':>10}{_stats(errors)}" + print(row) + print() + print(" trajectory basis by scenario, 5-10 min:") + for kind in kinds: + errors = [t.error_c for t in ticks if t.basis == "trajectory" and t.kind == kind and 5 <= t.minutes < 10] + print(f" {kind:<12} {_stats(errors)}") + for flag in ("probe", "hot"): + errors = [t.error_c for t in ticks if t.basis == "trajectory" and getattr(t, flag) and 5 <= t.minutes < 10] + print(f" {flag:<12} {_stats(errors)}") + print() + print("== Cross-current prediction: plateau at the new current, from what was known at the change") + print(" ambient: sensor = LAN/car sensor at the change; implied = previous run's trajectory through the") + print(" current law; warmer = max(sensor, implied); idle = handle before the session (session start only)") + print(" law: I2 = the I^2 prior; n = fitted exponent; n+k = fitted exponent and ambient term;") + print(" each fitted with the scored session left out") + methods = sorted({m for b in boundaries for m in b.predictions}) + groups: list[tuple[str, list[BoundaryScore]]] = [(kind, [b for b in boundaries if b.kind == kind]) for kind in kinds] + groups += [("probe", [b for b in boundaries if b.probe]), ("hot", [b for b in boundaries if b.hot]), + ("all", boundaries)] + for label, group in groups: + if not group: + continue + print(f" {label} ({len(group)} changes) {STATS_HEADER}") + for method in methods: + errors = [b.predictions[method] - b.actual_c for b in group if method in b.predictions] + if errors: + print(f" {method:<14} {_stats(errors)}") + print() + + +def main(argv: list[str] | None = None) -> int: + parser = argparse.ArgumentParser(description=__doc__.splitlines()[0]) + parser.add_argument("--db", required=True, help="path to a wallmonitor.db (read-only)") + parser.add_argument("--lookback-days", type=float, default=180.0) + parser.add_argument("--observe-tau", type=float, default=3.0, + help="a run must hold this many time constants for its plateau to count as truth (default %(default)s)") + parser.add_argument("--truth", choices=("fit", "last"), default="fit", + help="what counts as the plateau a run reached: an exponential fitted over the whole run " + "(fit), or the handle's own last three minutes, needing two more time constants (last)") + parser.add_argument("--json", help="also write every run, tick score and boundary score here") + parser.add_argument("--verbose", action="store_true", help="list every run") + args = parser.parse_args(argv) + + db = Database(args.db) + now = time.time() + fits = thermal.fit_sessions(db, now, lookback_days=args.lookback_days) + params = thermal.fit_history(db, now, fits=fits) + idle_model = thermal.load_idle_offset(db) + sessions = [ + sess for sess in db.sessions_range(now - args.lookback_days * 86400, now) + if sess.get("end_ts") and (sess.get("charging_s") or 0) >= thermal.MIN_SEGMENT_S + ] + sessions.sort(key=lambda sess: sess["start_ts"]) + + runs: list[Run] = [] + ticks: list[TickScore] = [] + boundaries: list[BoundaryScore] = [] + for sess in sessions: + laws = {name: params_without(fits, sess["id"], **pins) for name, pins in LAWS.items()} + session_runs = build_runs(db, sess, laws["n+k"], args.observe_tau, idle_model, args.truth) + for index, run in enumerate(session_runs): + ticks.extend(score_ticks(run, laws["n+k"])) + boundary = score_boundary(run, session_runs[index - 1] if index else None, laws) + if boundary is not None: + boundaries.append(boundary) + if args.verbose: + fmt = lambda value: f"{value:5.1f}" if value is not None else " -" # noqa: E731 + print( + f"s{run.session_id:<4} run{run.index} {run.kind:<10} {run.current_a:5.1f}A " + f"{run.span_s / 60:5.0f}min handle {run.handle_start_c:5.1f} -> last {fmt(run.last_c)} " + f"max {fmt(run.max_c)} | fit {fmt(run.fit_c)} (tau {fmt(run.fit_tau_min)}) " + f"30min-fit {fmt(run.window_c)} | truth {fmt(run.plateau_c)} | " + f"sensor {fmt(run.sensor_ambient_c)} implied {fmt(run.implied_ambient_c)}" + f"{' probe' if run.probe else ''}{' hot' if run.hot else ''}" + ) + runs.extend(session_runs) + + report(runs, ticks, boundaries, params, args.observe_tau) + if args.json: + with open(args.json, "w") as fh: + json.dump( + { + "model": params.as_dict(), + "runs": [{k: v for k, v in asdict(run).items() if k != "samples"} for run in runs], + "ticks": [asdict(tick) for tick in ticks], + "boundaries": [asdict(b) for b in boundaries], + }, + fh, + indent=1, + ) + print(f"wrote {args.json}", file=sys.stderr) + return 0 + + +if __name__ == "__main__": + sys.exit(main()) diff --git a/docs/thermal-model.md b/docs/thermal-model.md index f9d02d5..7c229aa 100644 --- a/docs/thermal-model.md +++ b/docs/thermal-model.md @@ -31,9 +31,22 @@ the trip happens. forecast at an off-reference current inherited the error — including the current the amp controller was told to restore to. Once the free-running fits span ≥ 6 A of current, the exponent *n* in - rise = rise₄₈ · (I/48)ⁿ is fitted by log-log regression (`current_exp` - in `/api/thermal`, with its standard error), every fit's `rise_ref_c` is - re-normalized with it, and the model note says so. + rise = rise₄₈ · (I/48)ⁿ is fitted (`current_exp` in `/api/thermal`, with + its standard error), every fit's `rise_ref_c` is re-normalized with it, + and the model note says so. +- **And so is how it scales with ambient.** The sensor reads the garage + air; the handle's environment runs hotter than the air by an amount that + grows with the heat — sun on the wall, a heat-soaked structure and cable. + [Backtested](#measuring-the-forecast) over 74 sessions, a model without + that term was optimistic by ~1 °C on mild days and 3–5 °C above 30 °C. + So the rise carries an ambient term, k · (ambient − 25 °C), fitted + jointly with *n* once the free-running fits span ≥ 4 °C of ambient + (`ambient_coef`, with its standard error; `rise_ref_c` is the rise at + 48 A *and* 25 °C). On the install above k = 0.25 °C/°C. The + [degradation watch](#the-confounder-the-fits-cant-remove) had been + reporting the same slope as an unexplained "ambient coefficient" for + weeks; with the term in the model, that coefficient reads ~0 and the + watch is checking the model rather than doing its job for it. - **The charger is its own thermometer.** Idle, the handle sits ~1–2 °C above ambient (an ambient-dependent offset), so ambient can be read without any extra sensor. The offset model ships as a seed from one install and is @@ -160,6 +173,54 @@ alert-40 raise to within seconds. `/api/thermal` returns the fitted model, the live forecast, every per-segment fit, and the drift verdict. +## Measuring the forecast + +Every change to the model is judged against the whole recorded history, not +against the incident that prompted it: + +```bash +uv run python contrib/backtest_forecast.py --db /path/to/wallmonitor.db +``` + +The tool splits every session into steady-current runs, takes the plateau +each long-enough run actually reached as truth (an exponential fitted over +the whole run, or with `--truth last` the handle's own final minutes — the +two bracket the answer), and scores two things against it: the live +forecast tick by tick as the current holds, and — at every current change +and session start — the plateau that would have been predicted at the new +current from each ambient on offer under each current law, with the scored +session left out of the fit. Errors are predicted minus actual; negative is +optimistic, the direction that trips the charger. + +What it showed on its first run, over 74 sessions: the model-basis +forecast was optimistic by **3–4 °C** on median and by more than 2 °C four +times in five; every cross-current method by 2.5–5 °C at a step-down. The +error grew with ambient — about 1 °C on mild days, 3–5 °C above 30 °C — +which was the same ambient coefficient the degradation watch's regression +had been reporting (0.35 °C/°C) and that its confidence interval alone had +not made convincing. And the whole-run τ ran 1–2 min longer than the +30-minute fit windows' 11.25, which is where the trajectory's residual +optimism came from. + +Two model changes followed, each scored by the tool before it shipped: the +ambient term above, and fit windows of 4 τ instead of 30 min. Before → +after, as (whole-run-fit truth / model-free truth): + +| | bias, °C | optimistic by > 2 °C | +|---|---|---| +| model basis, in-run | −3.3 / −3.7 → **−2.1 / −1.2** | 79 % / 80 % → **52 % / 20 %** | +| trajectory, 10–20 min in | −1.0 / −0.6 → −0.6 / +0.4 | 31 % / 31 % → 21 % / 18 % | +| step-down, from the sensor | −4.6 / −4.0 → **−2.2 / −1.3** | 92 % / 100 % → **64 % / 20 %** | +| probe end, from the sensor | −3.8 / −3.9 → **−0.6 / −0.8** | 67 % / 100 % → **0 % / 0 %** | + +The remaining optimism at a step-down (1–2 °C) is a history effect no +static ambient carries: the cable is still warm from the higher current +that preceded it. The `SUGGEST_MARGIN_C` of 2 °C covers the median of it. +The same run settled a rule the anecdotes could not: at a step-down the +sensor's ambient beats the one the previous run's trajectory implies, and +"the warmer of the two" ties the sensor within 0.1 °C — so the sensor is +used, and nothing cleverer. + ## Degradation watch The same per-segment fits feed a trend. Rising heat at unchanged current diff --git a/tests/test_backtest_forecast.py b/tests/test_backtest_forecast.py new file mode 100644 index 0000000..9d45cce --- /dev/null +++ b/tests/test_backtest_forecast.py @@ -0,0 +1,116 @@ +"""The forecast backtest on a synthetic history whose truth is known: a +cold-start run at full rate followed by a step-down, both long enough to +plateau. The tool must find both runs, read both plateaus, score the live +forecast against them, and score the cross-current prediction at the step +from what was known at the change.""" + +import importlib.util +import math +import pathlib +import sys +import time + +import pytest + +from wallmonitor import thermal +from wallmonitor.db import Database + +spec = importlib.util.spec_from_file_location( + "backtest_forecast", pathlib.Path(__file__).parent.parent / "contrib" / "backtest_forecast.py" +) +bt = importlib.util.module_from_spec(spec) +sys.modules[spec.name] = bt +spec.loader.exec_module(bt) + +TAU_S = 720.0 +RISE_REF = 36.0 + + +@pytest.fixture +def db(tmp_path): + database = Database(str(tmp_path / "test.db")) + yield database + database.close() + + +def _seed(db, start_ts, ambient_c, steps, dt=10.0): + """Idle lead-in, then charging segments [(amps, seconds), ...] that + follow the first-order model with the I^2 law, integrated so a + step-down cools toward its lower plateau. Sensor samples every minute.""" + ts = start_ts - 1800 + while ts < start_ts: + db.insert_vitals(ts, { + "vehicle_connected": 0, "contactor_closed": 0, "vehicle_current_a": 0.0, + "handle_temp_c": round(thermal.idle_handle_c(ambient_c), 2), "pcba_temp_c": 38.0, "mcu_temp_c": 46.0, + }, None, 0.0) + ts += dt + for minute in range(-30, int(sum(s for _, s in steps) / 60) + 1): + db.insert_ambient(start_ts + 60 * minute, ambient_c) + sid = db.start_session(start_ts) + temp = thermal.idle_handle_c(ambient_c) + ts = start_ts + for amps, seconds in steps: + end = ts + seconds + while ts < end: + db.insert_vitals(ts, { + "vehicle_connected": 1, "contactor_closed": 1, "vehicle_current_a": amps, + "handle_temp_c": round(temp, 3), "pcba_temp_c": 55.0, "mcu_temp_c": 50.0, + }, sid, amps * 233.0) + t_inf = ambient_c + RISE_REF * (amps / thermal.REF_CURRENT_A) ** 2 + temp += dt * (t_inf - temp) / TAU_S + ts += dt + db.close_session(sid, ts, "vehicle_disconnected") + return sid + + +def test_backtest_scores_a_cold_start_and_a_step_down(db, tmp_path, capsys): + now = time.time() + # History for the model to fit from, then the session under test. + for i in range(4): + _seed(db, now - (8 - i) * 7200, ambient_c=25.0, steps=[(48.6, 2400)]) + sid = _seed(db, now - 7200, ambient_c=25.0, steps=[(48.6, 2700), (40.0, 2700)]) + + out = tmp_path / "bt.json" + assert bt.main(["--db", str(tmp_path / "test.db"), "--json", str(out)]) == 0 + text = capsys.readouterr().out + assert "cold_start" in text and "step_down" in text + + import json + data = json.loads(out.read_text()) + runs = [r for r in data["runs"] if r["session_id"] == sid] + assert [r["kind"] for r in runs] == ["cold_start", "step_down"] + assert [round(r["current_a"], 1) for r in runs] == [48.6, 40.0] + # Both runs held >= 3 tau, so both plateaus are observed, and they are + # the seeded ones. + expected = [25.0 + RISE_REF * (a / 48.0) ** 2 for a in (48.6, 40.0)] + for run, plateau in zip(runs, expected): + assert run["plateau_c"] is not None and abs(run["plateau_c"] - plateau) < 1.0 + # The live forecast converges onto the plateau as the run matures. + late = [t["error_c"] for t in data["ticks"] + if t["session_id"] == sid and t["basis"] == "trajectory" and 10 <= t["minutes"] < 20] + assert late and max(abs(e) for e in late) < 1.5 + # At the step-down, every ambient on offer predicts the lower plateau — + # the history is I^2 and isothermal, so sensor, implied and warmer agree. + step = [b for b in data["boundaries"] if b["session_id"] == sid and b["kind"] == "step_down"] + assert len(step) == 1 and step[0]["from_a"] == pytest.approx(48.6, abs=0.1) + for method, predicted in step[0]["predictions"].items(): + assert abs(predicted - step[0]["actual_c"]) < 1.5, (method, predicted, step[0]["actual_c"]) + assert {m.split("/")[0] for m in step[0]["predictions"]} == {"sensor", "implied", "warmer"} + assert {m.split("/")[1] for m in step[0]["predictions"]} == {"I2", "n", "n+k"} + # ...and the session start is scored from the sensor and the idle handle. + start = [b for b in data["boundaries"] if b["session_id"] == sid and b["kind"] == "cold_start"] + assert len(start) == 1 and start[0]["from_a"] is None + assert {m.split("/")[0] for m in start[0]["predictions"]} == {"sensor", "idle"} + + +def test_split_runs_drops_ramp_samples_and_splits_on_a_current_change(): + def row(i, amps, closed=1): + return {"ts": 1000.0 + 2.0 * i, "contactor_closed": closed, "vehicle_current_a": amps, "handle_temp_c": 30.0} + + rows = [row(i, a) for i, a in enumerate([6.1, 14.6, 23.9, 31.1, 32.6, 32.7, 32.9, 37.6, 43.3] + [32.5] * 200 + [44.0] * 100)] + runs = bt.split_runs(rows) + assert [len(r) for r in runs] == [200, 100] + assert runs[0][0]["vehicle_current_a"] == 32.5 and runs[1][0]["vehicle_current_a"] == 44.0 + # A contactor drop ends a run; a short remainder is dropped. + rows = [row(i, 48.0) for i in range(100)] + [row(100, 0.0, closed=0)] + [row(101 + i, 48.0) for i in range(10)] + assert [len(r) for r in bt.split_runs(rows)] == [100] diff --git a/tests/test_wallmonitor.py b/tests/test_wallmonitor.py index 761c5b7..2bc249e 100644 --- a/tests/test_wallmonitor.py +++ b/tests/test_wallmonitor.py @@ -336,7 +336,7 @@ def _seed_idle(db, t_from, t_to, ambient_c, dt=10.0): def _seed_thermal_session(db, start_ts, ambient_c, tau_s=720.0, rise_ref_c=36.0, amps=48.6, charge_s=1500.0, dt=10.0, ambient_end_c=None, cooldown_s=0.0, sag_to_a=None, - current_exp=2.0): + current_exp=2.0, ambient_coef=0.0): """Idle lead-in plus a charging ramp that follows the first-order model. With ambient_end_c set, ambient drifts linearly across the charge (the @@ -354,7 +354,7 @@ def _seed_thermal_session(db, start_ts, ambient_c, tau_s=720.0, rise_ref_c=36.0, _seed_idle(db, start_ts - 1800, start_ts, ambient_c, dt) sid = db.start_session(start_ts) t0_temp = thermal.idle_handle_c(ambient_c) - rise_at = rise_ref_c * (amps / thermal.REF_CURRENT_A) ** current_exp + rise_at = rise_ref_c * (amps / thermal.REF_CURRENT_A) ** current_exp + ambient_coef * (ambient_c - 25.0) integrate = ambient_end_c is not None or sag_to_a is not None temp = t0_temp ts = start_ts @@ -371,7 +371,7 @@ def _seed_thermal_session(db, start_ts, ambient_c, tau_s=720.0, rise_ref_c=36.0, if integrate: ambient_now = ambient_c if ambient_end_c is None else \ ambient_c + (ambient_end_c - ambient_c) * elapsed - rise_now = rise_ref_c * (amps_now / thermal.REF_CURRENT_A) ** current_exp + rise_now = rise_ref_c * (amps_now / thermal.REF_CURRENT_A) ** current_exp + ambient_coef * (ambient_now - 25.0) temp += dt * ((ambient_now + rise_now - temp) / tau_s) ts += dt db.close_session(sid, start_ts + charge_s, "vehicle_disconnected") @@ -410,7 +410,8 @@ async def test_thermal_fits_the_current_exponent_from_a_current_spread(db): params = thermal.fit_history(db, now, fits=fits) assert params.current_exp_fits == 6 assert abs(params.current_exp - 1.5) < 0.15 - assert params.current_exp_se is not None and params.current_exp_se < 0.15 + assert params.current_exp_se is not None and params.current_exp_se < 0.2 + assert params.ambient_coef == 0.0 and params.ambient_coef_se is None # one ambient: not identifiable assert abs(params.rise_ref_c - 36.0) < 2.0 for fit in fits: assert abs(fit["rise_ref_c"] - 36.0) < 2.5, fit @@ -419,6 +420,38 @@ async def test_thermal_fits_the_current_exponent_from_a_current_spread(db): assert all(fit["rise_c"] * (48.0 / fit["current_a"]) ** 2 > 40.0 for fit in low) +async def test_thermal_fits_the_ambient_term_and_predicts_hot_days(db): + # Backtested over 74 real sessions, the forecast was ~1 C optimistic on + # mild days and 3-5 C above 30 C: the handle's environment runs hotter + # than the garage air by an amount that grows with the heat. Seeded here + # with k = 0.3 C/C across a spread of ambients and two currents; the + # fitter must recover k and n together, normalize every fit to 48 A and + # 25 C, and the forecast must then land on a hot day it has never seen. + now = time.time() + plan = [(22.0, 48.6), (25.0, 48.6), (30.0, 48.6), (34.0, 48.6), (27.0, 40.0), (32.0, 40.0), (24.0, 40.0)] + for i, (ambient, amps) in enumerate(plan): + _seed_thermal_session(db, now - (len(plan) + 1 - i) * 7200, ambient_c=ambient, amps=amps, + current_exp=1.5, ambient_coef=0.3) + fits = thermal.fit_sessions(db, now) + params = thermal.fit_history(db, now, fits=fits) + assert params.current_exp_fits == len(plan) + assert abs(params.ambient_coef - 0.3) < 0.08, params.ambient_coef + assert params.ambient_coef_se is not None and params.ambient_coef_se < 0.08 + assert abs(params.current_exp - 1.5) < 0.2 + assert abs(params.rise_ref_c - 36.0) < 2.0 + for fit in fits: + assert abs(fit["rise_ref_c"] - 36.0) < 2.5, fit + # A 36 C day at full rate: the seeded truth is 36 + 36.4 + 3.3 = 75.7. + truth = 36.0 + 36.0 * (48.6 / 48.0) ** 1.5 + 0.3 * 11.0 + assert abs(params.plateau_at(48.6, 36.0) - truth) < 1.5 + assert abs(params.ambient_from_plateau(truth, 48.6) - 36.0) < 1.0 + # ...and the sustainable current shrinks accordingly: without k, 36 C + # leaves 27 C of headroom; with k = 0.3 it leaves 23.7. + without = thermal.ThermalParams(rise_ref_c=params.rise_ref_c, current_exp=params.current_exp) + assert thermal.sustainable_max_current(36.0, params) < thermal.sustainable_max_current(36.0, without) + assert params.safe_ambient_max_c() < without.safe_ambient_max_c() + + def test_steady_prefix_restarts_after_a_ramp_up_overshoot(): # 2026-09-14, first live calibration probe: the car overshot toward 48 A # for two samples while the 32 A cap was taking effect, then held 32.5 A diff --git a/wallmonitor/static/app.js b/wallmonitor/static/app.js index a5e4bf6..1dd0618 100644 --- a/wallmonitor/static/app.js +++ b/wallmonitor/static/app.js @@ -1145,8 +1145,10 @@ async function viewLive(root) { // The current law is a prior (I²) until this install's own fits span // enough current to measure it; say which one every forecast rests on. const expNote = model.current_exp_fits > 0 - ? ` Heat rise scales as I^${fmtNum(model.current_exp, 2)} here (fitted across ${model.current_exp_fits} free-running ` + - "sessions; the I² prior over-reads how much lower currents cool the handle)." + ? ` Heat rise scales as I^${fmtNum(model.current_exp, 2)} here` + + (model.ambient_coef ? ` and grows ${fmtNum(model.ambient_coef, 2)} °C per °C of ambient above ${fmtNum(model.ambient_ref_c, 0)} °C` : "") + + ` (fitted across ${model.current_exp_fits} free-running sessions; the I² prior over-reads how much lower currents ` + + "cool the handle, and the garage air under-reads how hot the handle's surroundings run on a hot day)." : ""; const modelNote = `Model: τ ≈ ${fmtNum(model.tau_min, 1)} min, +${fmtNum(model.rise_ref_c, 0)} °C at ${fmtNum(model.ref_current_a, 0)} A — ` + (model.fitted ? `fitted from ${model.tau_fits} recorded session ramp${model.tau_fits === 1 ? "" : "s"}.` + priorNote + expNote diff --git a/wallmonitor/thermal.py b/wallmonitor/thermal.py index 09e7699..55bc9a5 100644 --- a/wallmonitor/thermal.py +++ b/wallmonitor/thermal.py @@ -199,6 +199,21 @@ def ambient_from_idle_handle(handle_c: float, model: IdleOffset = BUILTIN_IDLE_O CURRENT_EXP_MIN_FITS = 4 CURRENT_EXP_MIN_SPAN_A = 6.0 CURRENT_EXP_RANGE = (1.0, 2.5) +CURRENT_EXP_GRID_STEP = 0.02 + +# The rise also grows with ambient: rise(I, a) = rise_ref * (I/48)^n + +# k * (a - 25). The sensor reads the garage air, and the handle's +# environment runs hotter than the air by an amount that scales with the +# heat — sun on the wall, a heat-soaked structure and cable. Backtested over +# 74 sessions with k = 0, the model-basis forecast was optimistic by ~1 C on +# mild days and 3-5 C above 30 C, in 79% of ticks by more than 2 C, and the +# degradation watch's regression had been reporting the same slope +# (0.35 C/C) as an "ambient coefficient" for weeks. So it is a model term, +# fitted from the same free-running fits once they span enough ambient, +# and rise_ref_c means the rise at 48 A and 25 C. +AMBIENT_REF_C = 25.0 +DEFAULT_AMBIENT_COEF = 0.0 +AMBIENT_COEF_MIN_SPREAD_C = 4.0 # Fit acceptance gates: a segment must actually contain a thermal ramp and # the exponential must describe it well, or it teaches the model nothing. @@ -220,7 +235,12 @@ def ambient_from_idle_handle(handle_c: float, model: IdleOffset = BUILTIN_IDLE_O # scales with the same tau estimate (PREFIX_SPAN_TAU) so a slow-tau install # is not starved of fits by a fixed cap it can never clear. MIN_SPAN_TAU = 1.8 -PREFIX_SPAN_TAU = 2.5 +# The fit window runs to 4 tau (was 2.5, which the 30 min floor made moot at +# a typical tau). Whole-run fits over the recorded history read tau 1-2 min +# longer than the 30 min windows did, and that gap is where the trajectory +# forecast's residual optimism came from: a window that ends at 2.7 tau has +# seen 93% of the rise and still trades a little tau for a little rise. +PREFIX_SPAN_TAU = 4.0 PREFIX_SPAN_MIN_S = 1800.0 # A steady run that ends within its first minute did not end; the ramp-up # was still wobbling. Vehicles overshoot on the way to a cap: one session @@ -294,19 +314,36 @@ class ThermalParams: current_exp: float = DEFAULT_CURRENT_EXP current_exp_fits: int = 0 current_exp_se: float | None = None + ambient_coef: float = DEFAULT_AMBIENT_COEF + ambient_coef_se: float | None = None @property def fitted(self) -> bool: return self.tau_fits > 0 and self.rise_fits > 0 - def rise_at(self, current_a: float) -> float: - """Steady-state rise above ambient at a charge current.""" - return self.rise_ref_c * (current_a / REF_CURRENT_A) ** self.current_exp + def rise_at(self, current_a: float, ambient_c: float = AMBIENT_REF_C) -> float: + """Steady-state rise above ambient at a charge current and ambient.""" + return ( + self.rise_ref_c * (current_a / REF_CURRENT_A) ** self.current_exp + + self.ambient_coef * (ambient_c - AMBIENT_REF_C) + ) + + def plateau_at(self, current_a: float, ambient_c: float) -> float: + return ambient_c + self.rise_at(current_a, ambient_c) + + def ambient_from_plateau(self, plateau_c: float, current_a: float) -> float: + """plateau_at inverted: the ambient at which this current settles here.""" + base = self.rise_ref_c * (current_a / REF_CURRENT_A) ** self.current_exp + return (plateau_c - base + self.ambient_coef * AMBIENT_REF_C) / (1.0 + self.ambient_coef) def current_for_rise(self, rise_c: float) -> float: - """The charge current whose steady-state rise is rise_c — rise_at inverted.""" + """The charge current whose current-dependent rise is rise_c.""" return REF_CURRENT_A * (rise_c / self.rise_ref_c) ** (1.0 / self.current_exp) + def safe_ambient_max_c(self) -> float: + """The ambient above which a full-rate charge settles at the trip point.""" + return (TRIP_HANDLE_C - self.rise_ref_c + self.ambient_coef * AMBIENT_REF_C) / (1.0 + self.ambient_coef) + # How far a fitted value may sit from the default before the dashboard # says the priors were a poor fit for this install. A heuristic, not a # statistic: 30% is roughly where the default-driven forecast's plateau @@ -352,6 +389,9 @@ def as_dict(self) -> dict: "current_exp": round(self.current_exp, 2), "current_exp_fits": self.current_exp_fits, "current_exp_se": round(self.current_exp_se, 2) if self.current_exp_se is not None else None, + "ambient_coef": round(self.ambient_coef, 3), + "ambient_coef_se": round(self.ambient_coef_se, 3) if self.ambient_coef_se is not None else None, + "ambient_ref_c": AMBIENT_REF_C, } @@ -810,51 +850,105 @@ def fit_sessions(db: Database, now: float, lookback_days: float = 120.0) -> list return fits -def _fit_current_exponent(fits: list[dict]) -> tuple[float, int, float | None]: - """The exponent n in rise = rise_ref * (I/48)^n, from free-running fits: - (n, fits used, standard error). The prior with no fits used when the - history cannot identify it — too few free-running fits, or all at one - current, where any n explains the data equally well. - - Log-log least squares: ln(rise) = ln(rise_ref) + n * ln(I/48). Only - free-running windows count — a regulated window's plateau is lower for a - reason that has nothing to do with current, and those windows cluster at - high current, which would bend n downward for the wrong reason.""" +@dataclass(frozen=True) +class CurrentLaw: + """How rise depends on current and ambient: rise = rise_ref * (I/48)^n + + k * (ambient - 25), with how many fits identified it and the standard + error of each fitted term (None for a term held at its prior).""" + + exponent: float = DEFAULT_CURRENT_EXP + ambient_coef: float = DEFAULT_AMBIENT_COEF + fits: int = 0 + exponent_se: float | None = None + ambient_coef_se: float | None = None + + +def _fit_current_law(fits: list[dict], exponent: float | None = None, + ambient_coef: float | None = None) -> CurrentLaw: + """Fit n and k from free-running fits; either can be pinned instead. + + n is profiled over a grid: for each n the rise is linear in rise_ref + and k, so ordinary least squares gives both and the n with the smallest + residual wins. Its standard error comes from the curvature of the + residual sum of squares at the optimum. A term is only fitted when the + history moved enough to identify it — >= CURRENT_EXP_MIN_SPAN_A of + current for n, >= AMBIENT_COEF_MIN_SPREAD_C of ambient for k — and + stays at its prior otherwise. Only free-running windows count: a + regulated window's plateau is lower for a reason that has nothing to do + with current or ambient, and those windows cluster at high current on + hot days, which would bend both terms for the wrong reason.""" points = [ - (math.log(fit["current_a"] / REF_CURRENT_A), math.log(fit["rise_c"])) + (fit["current_a"], fit["ambient_c"], fit["rise_c"]) for fit in fits - if fit.get("rise_c") is not None and fit["rise_c"] > 0 - and fit.get("current_a") and fit.get("free_plateau", True) + if fit.get("rise_c") is not None and fit["rise_c"] > 0 and fit.get("current_a") + and fit.get("ambient_c") is not None and fit.get("free_plateau", True) ] + pinned = CurrentLaw( + exponent=DEFAULT_CURRENT_EXP if exponent is None else exponent, + ambient_coef=DEFAULT_AMBIENT_COEF if ambient_coef is None else ambient_coef, + ) if len(points) < CURRENT_EXP_MIN_FITS: - return DEFAULT_CURRENT_EXP, 0, None - currents = [REF_CURRENT_A * math.exp(x) for x, _ in points] - if max(currents) - min(currents) < CURRENT_EXP_MIN_SPAN_A: - return DEFAULT_CURRENT_EXP, 0, None - result = _ols([y for _, y in points], [[1.0, x] for x, _ in points]) - if result is None: - return DEFAULT_CURRENT_EXP, 0, None - beta, se, _resid, _dof = result - exponent = min(CURRENT_EXP_RANGE[1], max(CURRENT_EXP_RANGE[0], beta[1])) - return exponent, len(points), se[1] - - -def fit_history(db: Database, now: float, lookback_days: float = 120.0, - fits: list[dict] | None = None) -> ThermalParams: - """Aggregate per-session fits into model parameters; defaults where thin. - - Fits the install's current exponent first and re-normalizes every fit's - rise_ref_c with it, in place — so the fits handed on to the API and the - degradation watch, and the rise_ref_c median taken here, all share one - current law. With the I^2 prior, off-reference fits carry a bias that the - watch's current term then has to absorb; with a fitted exponent that - term has nothing left to explain.""" - if fits is None: - fits = fit_sessions(db, now, lookback_days) - exponent, exp_fits, exp_se = _fit_current_exponent(fits) + return pinned + currents = [current for current, _, _ in points] + ambients = [ambient for _, ambient, _ in points] + rises = [rise for _, _, rise in points] + fit_n = exponent is None and max(currents) - min(currents) >= CURRENT_EXP_MIN_SPAN_A + fit_k = ambient_coef is None and max(ambients) - min(ambients) >= AMBIENT_COEF_MIN_SPREAD_C + if not fit_n and not fit_k: + return pinned + k_pinned = pinned.ambient_coef + + def solve(n: float) -> tuple[list[float], list[float], float, int] | None: + # With k pinned its contribution is moved to the left-hand side. + y = rises if fit_k else [rise - k_pinned * (ambient - AMBIENT_REF_C) for ambient, rise in zip(ambients, rises)] + design = [ + [(current / REF_CURRENT_A) ** n] + ([ambient - AMBIENT_REF_C] if fit_k else []) + for current, ambient in zip(currents, ambients) + ] + return _ols(y, design) + + if not fit_n: + result = solve(pinned.exponent) + if result is None: + return pinned + beta, se, _resid, _dof = result + return CurrentLaw(pinned.exponent, beta[1], len(points), None, se[1]) + + steps = int(round((CURRENT_EXP_RANGE[1] - CURRENT_EXP_RANGE[0]) / CURRENT_EXP_GRID_STEP)) + grid = [CURRENT_EXP_RANGE[0] + i * CURRENT_EXP_GRID_STEP for i in range(steps + 1)] + solved = [(n, solve(n)) for n in grid] + candidates = [(n, result) for n, result in solved if result is not None] + if not candidates: + return pinned + best_index = min(range(len(candidates)), key=lambda i: candidates[i][1][2]) + n, (beta, se, resid_sd, dof) = candidates[best_index] + n_se = None + if 0 < best_index < len(candidates) - 1: + # SSE(n) is locally quadratic; Var(n) = 2 sigma^2 / SSE''(n). + sse = [candidates[i][1][2] ** 2 * dof for i in (best_index - 1, best_index, best_index + 1)] + curvature = (sse[0] - 2 * sse[1] + sse[2]) / CURRENT_EXP_GRID_STEP**2 + if curvature > 0: + n_se = math.sqrt(2 * resid_sd**2 / curvature) + return CurrentLaw(n, beta[1] if fit_k else k_pinned, len(points), n_se, se[1] if fit_k else None) + + +def params_from_fits(fits: list[dict], exponent: float | None = None, + ambient_coef: float | None = None) -> ThermalParams: + """Model parameters from per-segment fits; defaults where thin. + + Fits the current law first (either term can be pinned, for comparing + laws) and re-normalizes every fit's rise_ref_c under it, in place — so + the fits handed on to the API and the degradation watch, and the + rise_ref_c median taken here, all share one law. Under the I^2, k = 0 + prior, off-reference and hot-day fits carry a bias that the watch's + current and ambient terms then have to absorb; under the fitted law + those terms have nothing left to explain.""" + law = _fit_current_law(fits, exponent, ambient_coef) for fit in fits: if fit.get("rise_c") is not None and fit.get("current_a"): - fit["rise_ref_c"] = round(fit["rise_c"] * (REF_CURRENT_A / fit["current_a"]) ** exponent, 2) + ambient = fit.get("ambient_c") + adjusted = fit["rise_c"] - (law.ambient_coef * (ambient - AMBIENT_REF_C) if ambient is not None else 0.0) + fit["rise_ref_c"] = round(adjusted * (REF_CURRENT_A / fit["current_a"]) ** law.exponent, 2) taus = [fit["tau_min"] for fit in fits] rises = [fit["rise_ref_c"] for fit in fits if fit["rise_ref_c"] is not None] rmses = [fit["rmse_c"] for fit in fits] @@ -864,12 +958,22 @@ def fit_history(db: Database, now: float, lookback_days: float = 120.0, tau_fits=len(taus), rise_fits=len(rises), fit_rmse_c=median(rmses) if rmses else None, - current_exp=exponent, - current_exp_fits=exp_fits, - current_exp_se=exp_se, + current_exp=law.exponent, + current_exp_fits=law.fits, + current_exp_se=law.exponent_se, + ambient_coef=law.ambient_coef, + ambient_coef_se=law.ambient_coef_se, ) +def fit_history(db: Database, now: float, lookback_days: float = 120.0, + fits: list[dict] | None = None) -> ThermalParams: + """Aggregate per-session fits into model parameters; defaults where thin.""" + if fits is None: + fits = fit_sessions(db, now, lookback_days) + return params_from_fits(fits) + + # ---------------- degradation watch ---------------- # A loose lug or degrading contact shows up as extra resistance: more heat @@ -1255,7 +1359,7 @@ def sustainable_max_current(ambient_c: float, params: ThermalParams) -> float | trim the last amp or two. Worked from a measured ambient when a sensor is reporting — see predict() for why the trajectory-implied one is only the fallback.""" - headroom = TRIP_HANDLE_C - SUGGEST_MARGIN_C - ambient_c + headroom = TRIP_HANDLE_C - SUGGEST_MARGIN_C - ambient_c - params.ambient_coef * (ambient_c - AMBIENT_REF_C) if headroom <= 0: return None amps = math.floor(params.current_for_rise(headroom)) @@ -1346,7 +1450,7 @@ def _recent_steady_ambient(recent: list[dict], params: ThermalParams) -> float | window = [(sample["ts"], sample["handle_temp_c"]) for sample in run] t_inf, _se = _project_t_inf(window, params.tau_min) run_current = median(sample["vehicle_current_a"] for sample in run) - ambient = t_inf - params.rise_at(run_current) + ambient = params.ambient_from_plateau(t_inf, run_current) if -30.0 <= ambient <= TRIP_HANDLE_C: return ambient return None @@ -1427,7 +1531,7 @@ def predict(db: Database, now: float, params: ThermalParams) -> dict: gap_reason = "warming_up" out["forecast"] = {"basis": "insufficient", "will_trip": None, "reason": gap_reason} return out - t_inf = ambient + params.rise_at(current) + t_inf = params.plateau_at(current, ambient) forecast["basis"] = "model" forecast["ambient_source"] = source # A model-basis plateau is only as good as its ambient: a sensor @@ -1447,7 +1551,7 @@ def predict(db: Database, now: float, params: ThermalParams) -> dict: # reporting — interpolation inside the fitted data rather than a # round trip to 32 A and back (on one install, after a 32 A probe, # the implied route named 47 A; the measured one 44 A, which held). - implied_ambient = t_inf - params.rise_at(current) + implied_ambient = params.ambient_from_plateau(t_inf, current) measured = _latest_measured_ambient(db, now) sustain_ambient, sustain_source = measured if measured is not None else (implied_ambient, "implied") forecast.update( @@ -1485,7 +1589,7 @@ def predict(db: Database, now: float, params: ThermalParams) -> dict: out["ambient_se_c"] = 0.3 if measured is not None else idle_model.ambient_se_c out["ambient_stable"] = stable # Hypothetical: a full-rate session started right now. - t_inf = ambient + params.rise_ref_c + t_inf = params.plateau_at(REF_CURRENT_A, ambient) minutes = _minutes_to_trip(last["handle_temp_c"], t_inf, tau_min) out["forecast"] = { "basis": "hypothetical", @@ -1493,7 +1597,7 @@ def predict(db: Database, now: float, params: ThermalParams) -> dict: "will_trip": minutes is not None, "minutes_to_trip": round(minutes, 1) if minutes is not None else None, "trip_ts": None, - "safe_ambient_max_c": round(TRIP_HANDLE_C - params.rise_ref_c, 1), + "safe_ambient_max_c": round(params.safe_ambient_max_c(), 1), "suggested_max_a": suggest_max_current(ambient, params) if minutes is not None else None, "sustainable_max_a": sustainable_max_current(ambient, params), }