HormoneTrack-web/js/pk/pk-engine.js

/**
 * ─────────────────────────────────────────────────────────────────────────────
 * PharmacokineticEngine (partie 1/4 : noyau) — portage web de
 * `pk/PharmacokineticEngine.kt`. Le reste est découpé en modules :
 *  - pk-calibration.js : auto/manuelle par ester/période/modèle (#60/#61) ;
 *  - pk-reminders.js   : prochain déclenchement (grille Posologie, #52) ;
 *  - pk-extrema.js     : détection pics/creux (triangles du graphique).
 *
 * Chaque injection contribue `dose_mg × profil(dt)` où `profil` est la
 * réponse normalisée (pg/mL par mg) ; les contributions se superposent.
 * Trois modèles superposables :
 *  - **Estrannaise (ESE)** = tables horaires de l'asset (ODS) ;
 *  - **Transfem Science (TFS)** = forme close V3C de la méta-analyse (7 esters) ;
 *  - **WHSAH (WHS)** = fit « license-free » de Mona (6 esters, PEP non couvert).
 * Gel/patch/oral → modèle de Bateman paramétrable.
 *
 * ⚠️ INVARIANTS à préserver (doc Android §14 — bugs à ne pas réintroduire) :
 *  - `computeKa` : bisection DANS LE BON SENS (#19) ;
 *  - tables ESE : modèle strict + extrapolation terminale ≥ 1 % du pic (#20/#21) ;
 *  - `convertTToNgMl` appliquée à la calibration ET au rendu (#23) ;
 *  - `generateForecastDoses` : créneaux passés SAUTÉS (#35).
 *
 * PUR : aucune dépendance au DOM ni au localStorage → testable en Node.
 * ─────────────────────────────────────────────────────────────────────────────
 */

import * as TransfemScienceModels from './transfem-science-models.js';
import * as WhsahModels from './whsah-models.js';
import * as PKProfileStore from './pk-profile-store.js';

/** Une heure en ms. */
export const HOUR_MS = 3600000;

/** Un point de courbe : E2 en pg/mL, T en ng/mL. */
export class LevelPoint {
  /** @param {number} timestamp Epoch ms · @param {number} e2 pg/mL · @param {number} t ng/mL */
  constructor(timestamp, e2, t) {
    this.timestamp = timestamp;
    this.e2 = e2;
    this.t = t;
  }
}

/** Réglages du modèle T empirique (mêmes défauts que TConfig.kt). */
export class TConfig {
  /** @param {number} [base=6.0] T de base (ng/mL) à E2 ≈ 0 · @param {number} [floor=0.2] plancher · @param {number} [k=0.19] constante d'inhibition */
  constructor(base = 6.0, floor = 0.2, k = 0.19) {
    this.base = base;
    this.floor = floor;
    this.k = k;
  }
}

/** Constantes d'ester (miroir de `Esters.kt` — mêmes clés que les backups JSON). */
export const Esters = {
  NONE: 'NONE', EV: 'EV', EU: 'EU', EEN: 'EEN',
  EB: 'EB', EC: 'EC', ECS: 'ECS', PEP: 'PEP',
};

/** Constantes de modèles PK (miroir de `PKModels.kt`). */
export const PKModels = { ESTRANNAISE: 'ESE', TRANSFEM_SCIENCE: 'TFS', WHSAH: 'WHS' };

/**
 * Propriété DÉRIVÉE (miroir du getter Kotlin `usesProfileModel`) : le
 * traitement a un PROFIL PK s'il est injectable (IM/SC) ET avec un ester.
 * En JS ce n'est pas une propriété stockée (objets JSON plats, mêmes champs
 * que Gson) — le moteur la recalcule à la demande.
 *
 * @param {object} treatment
 * @returns {boolean}
 */
export function usesProfileModel(treatment) {
  return (treatment.route === 'INJECTION_IM' || treatment.route === 'INJECTION_SUBCUT')
    && treatment.esterType !== Esters.NONE;
}

// ---------------------------------------------------------------------------
// Bateman (gel / patch / oral)
// ---------------------------------------------------------------------------

/**
 * Résout `ln(ka/ke) = (ka−ke)·Tmax` par bisection (50 itérations).
 *
 * ⚠️ BUG #19 (épinglé depuis la session 1) : eq(mid) DÉCROÎT en mid et
 * s'annule en ka > ke → `eq > 0` signifie que la racine est AU-DESSUS de
 * mid ⇒ `lo = mid` (l'ancienne bisection inversée convergeait vers un ka
 * énorme, pic à ~0 h au lieu de Tmax).
 *
 * @param {number} tHalfHours  Demi-vie d'élimination (h)
 * @param {number} tMaxHours   Temps au pic souhaité (h)
 * @returns {number} ka (h⁻¹)
 */
export function computeKa(tHalfHours, tMaxHours) {
  if (tMaxHours <= 0 || tHalfHours <= 0) return 1.0;
  const ke = Math.log(2.0) / tHalfHours;
  if (tMaxHours < 0.01) return ke * 100.0;
  let lo = ke * 1.001;
  let hi = ke * 1000.0;
  for (let i = 0; i < 50; i++) {
    const mid = (lo + hi) / 2.0;
    const eq = Math.log(mid / ke) - (mid - ke) * tMaxHours;
    if (eq > 0) lo = mid; else hi = mid; // direction corrigée (bug #19)
  }
  return (lo + hi) / 2.0;
}

/**
 * Paramètres (ke, ka) d'un traitement Bateman, bornes défensives incluses.
 * @param {object} treatment  Treatment (data/models.js)
 * @returns {{ke:number, ka:number}}
 */
export function batemanParams(treatment) {
  const ke = Math.log(2.0) / Math.max(treatment.eliminationHalfLifeHours, 0.01);
  const ka = computeKa(
    Math.max(treatment.eliminationHalfLifeHours, 0.01),
    Math.max(treatment.absorptionHours, 0.01),
  );
  return { ke, ka };
}

// ---------------------------------------------------------------------------
// Contribution d'une dose
// ---------------------------------------------------------------------------

/**
 * Ester EFFECTIF d'une dose : override par injection (`dose.esterType`) ou
 * ester par défaut du traitement — permet de switcher d'ester d'une
 * injection à l'autre (comme dans le tableur d'origine).
 *
 * @param {object} treatment
 * @param {object} dose  DoseLog
 * @returns {string} clé ester ("EV"…, "NONE" si aucun)
 */
export function doseEster(treatment, dose) {
  return dose.esterType ?? treatment.esterType;
}

/**
 * Contribution d'UNE dose au niveau d'E2 à queryTimeMs (pg/mL), tous modèles :
 *  1. traitement à PROFIL PK (injection + ester ≠ NONE) :
 *     - TFS + ester couvert par le V3C → TransfemScienceModels ;
 *     - WHSAH + ester couvert par le fit Mona → WhsahModels ;
 *     - sinon → tables ODS via PKProfileStore (ester sans modèle → 0) ;
 *  2. sinon → Bateman (gel / patch / oral / custom).
 *
 * `modelOverride` force ESE/TFS/WHS pour CE calcul (le graphique dessine les
 * trois modèles côte à côte depuis le même traitement) ; Bateman n'est pas
 * concerné (les séries y sont identiques).
 *
 * @param {object} treatment
 * @param {object} dose
 * @param {number} queryTimeMs
 * @param {{ke:number,ka:number}|null} [bateman=null]  Params précalculés (cache)
 * @param {string|null} [modelOverride=null]  "ESE" | "TFS" | "WHS"
 * @returns {number} pg/mL
 */
export function concentrationOfDose(treatment, dose, queryTimeMs, bateman = null, modelOverride = null) {
  const dtH = (queryTimeMs - dose.timestamp) / 3600000.0;
  if (dtH <= 0.0) return 0.0;

  const mg = dose.doseAmount;
  if (mg <= 0.0) return 0.0;

  if (usesProfileModel(treatment)) {
    const ester = doseEster(treatment, dose);
    if (ester !== Esters.NONE) {
      const mod = modelOverride || treatment.pkModel;
      // Dispatch parallèle des 3 modèles (aucun ne remplace un autre) :
      // TFS (V3C méta-analyse) puis WHSAH (fit Mona), fallback tables ODS.
      if (mod === PKModels.TRANSFEM_SCIENCE && TransfemScienceModels.hasModel(ester)) {
        return TransfemScienceModels.sample(ester, dtH) * mg;
      }
      if (mod === PKModels.WHSAH && WhsahModels.hasModel(ester)) {
        return WhsahModels.sample(ester, dtH) * mg;
      }
      return PKProfileStore.sample(ester, mod, dtH) * mg;
    }
  }

  // Fallback Bateman (gel / patch / oral / custom)
  const p = bateman || batemanParams(treatment);
  const diff = p.ka - p.ke;
  const a = mg * treatment.bioavailabilityFraction;
  // Cas dégénéré ka≈ke : limite analytique (sinon division par ~0)
  const c = Math.abs(diff) < 1e-3
    ? a * p.ke * dtH * Math.exp(-p.ke * dtH)
    : (a * p.ka / diff) * (Math.exp(-p.ke * dtH) - Math.exp(-p.ka * dtH));
  return Math.max(0.0, c);
}

/**
 * Fenêtre de contribution d'un traitement en heures (au-delà : la dose est
 * ignorée par les boucles d'agrégation — performance).
 *  - TFS/WHS : 10 demi-vies TERMINALES du fit (ex. EV 30 j, PEP 284 j) ;
 *  - tables ODS : longueur de table (8001 h) ;
 *  - Bateman : 30 × t½. Toujours ≥ 24 h.
 *
 * @param {object} treatment
 * @returns {number} heures
 */
export function cutoffHours(treatment) {
  let profileH = 0.0;
  if (usesProfileModel(treatment)) {
    const mod = treatment.pkModel;
    if (mod === PKModels.WHSAH && WhsahModels.hasModel(treatment.esterType)) {
      profileH = WhsahModels.model(treatment.esterType).terminalHalfLifeDays * 24.0 * 10.0;
    } else if (mod === PKModels.TRANSFEM_SCIENCE && TransfemScienceModels.hasModel(treatment.esterType)) {
      profileH = TransfemScienceModels.model(treatment.esterType).terminalHalfLifeDays * 24.0 * 10.0;
    } else {
      profileH = PKProfileStore.profileLength(treatment.esterType, mod);
    }
  }
  const batemanH = usesProfileModel(treatment) ? 0.0 : 30.0 * treatment.eliminationHalfLifeHours;
  const cut = Math.max(profileH, batemanH);
  return cut >= 24.0 ? cut : 24.0;
}

// ---------------------------------------------------------------------------
// Niveaux agrégés
// ---------------------------------------------------------------------------

/**
 * Niveau E2 total à tMs : somme des contributions de toutes les doses E2
 * (traitements ESTRADIOL uniquement), chacune multipliée par SON facteur
 * d'échelle — `scalePerEster[ester]` si présent (calibration par période
 * d'ester, v1.2.1), sinon `scaleFactor` stocké du traitement.
 *
 * ⚠️ TOUS les traitements sont passés ici, actifs comme INACTIFS :
 * « isActive » est un drapeau administratif (plus de nouvelles doses),
 * JAMAIS un filtre de données (bug v1.2.4 — l'historique EV d'un traitement
 * inactivé doit rester simulé, cf régression #3).
 *
 * @param {object[]} treatments
 * @param {object[]} doseLogs
 * @param {number} tMs  Instant d'évaluation (epoch ms)
 * @param {string|null} [modelOverride=null]
 * @param {Object<string,number>|null} [scalePerEster=null]
 * @returns {number} pg/mL
 */
export function e2At(treatments, doseLogs, tMs, modelOverride = null, scalePerEster = null) {
  const batemanCache = new Map();
  let total = 0.0;
  for (const treatment of treatments) {
    if (treatment.type !== 'ESTRADIOL') continue; // AA/progestatifs → 0 en E2
    if (!batemanCache.has(treatment.id)) {
      batemanCache.set(treatment.id, batemanParams(treatment));
    }
    for (const dose of doseLogs) {
      if (dose.treatmentId !== treatment.id || dose.timestamp > tMs) continue;
      const dtH = (tMs - dose.timestamp) / 3600000.0;
      if (dtH > cutoffHours(treatment)) continue;
      const c = concentrationOfDose(treatment, dose, tMs, batemanCache.get(treatment.id), modelOverride);
      if (c > 0.0) {
        const scale = (scalePerEster && scalePerEster[doseEster(treatment, dose)]) ?? treatment.scaleFactor;
        total += c * scale;
      }
    }
  }
  return total;
}

/**
 * Modèle T EMPIRIQUE (non publié — étiqueté « estimation » partout) :
 * `T = floor + (base − floor) / (1 + k·E2)` en ng/mL.
 *
 * @param {number} e2Level  E2 calibrée (pg/mL)
 * @param {TConfig} config
 * @returns {number} ng/mL
 */
export function testosteroneAt(e2Level, config) {
  if (e2Level <= 0.0) return config.base;
  return config.floor + (config.base - config.floor) / (1.0 + config.k * e2Level);
}

/**
 * Normalisation des labs T vers ng/mL — appliquée à la calibration ET au
 * rendu du graphique (sinon l'axe T est faux d'un facteur 100, bug #23 :
 * labs 33/44 ng/dL). Branche défensive pg/mL ÷ 1000 : un lab saisi avec une
 * unité aberrante (« 38 pg/mL », bug #26) ne doit pas écraser l'axe T.
 *
 * @param {number} value  Valeur du lab
 * @param {string} unit   Unité saisie (ng/mL, ng/dL, ng/L, nmol/L…)
 * @returns {number} ng/mL
 */
export function convertTToNgMl(value, unit) {
  const u = String(unit).toLowerCase().replace(/\s/g, '');
  if (u.includes('dl')) return value / 100.0;
  if (u.includes('nmol')) return value * 0.2884;
  if (u.includes('pg')) return value / 1000.0;
  if (u.includes('µg') || u.includes('μg')) return value / 1000.0;
  if (u.includes('ng/l') || u.endsWith('/l')) return value / 1000.0;
  return value; // ng/mL (ou unité inconnue : valeur brute)
}

/**
 * L'ester « actif » à l'instant t = celui de la DERNIÈRE dose E2 ≤ t
 * (une prise de sang / un point de courbe reflète l'injection qui précède).
 * Utilisé pour choisir le k de la courbe T (k par période d'ester, v1.2.3).
 *
 * @param {object[]} treatments
 * @param {object[]} doseLogs
 * @param {number} tMs
 * @returns {string|null} clé ester, ou null (avant la 1ʳᵉ dose / ester NONE)
 */
export function activeEsterAt(treatments, doseLogs, tMs) {
  const estrogenIds = new Set(treatments.filter((tr) => tr.type === 'ESTRADIOL').map((tr) => tr.id));
  let last = null;
  for (const d of doseLogs) {
    if (d.timestamp <= tMs && estrogenIds.has(d.treatmentId)) {
      if (last === null || d.timestamp > last.timestamp) last = d;
    }
  }
  if (!last) return null;
  const tr = treatments.find((t) => t.id === last.treatmentId);
  const ester = doseEster(tr, last);
  return ester !== Esters.NONE ? ester : null;
}

/**
 * Niveau combiné (E2 + T) à un instant donné — la notification d'alerte et la
 * carte « niveau actuel » de l'accueil partagent CE calcul (même calibration :
 * `scalePerEster` + `tKPerEster`, cf §9.bis Android).
 *
 * @param {object[]} treatments
 * @param {object[]} doseLogs
 * @param {number} tMs
 * @param {TConfig} tConfig
 * @param {Object<string,number>|null} [tKPerEster=null]  k T par ester (période active)
 * @param {Object<string,number>|null} [scalePerEster=null]
 * @returns {LevelPoint}
 */
export function levelAt(treatments, doseLogs, tMs, tConfig, tKPerEster = null, scalePerEster = null) {
  const e2 = e2At(treatments, doseLogs, tMs, null, scalePerEster);
  const active = activeEsterAt(treatments, doseLogs, tMs);
  const k = (active && tKPerEster && tKPerEster[active]) ?? tConfig.k;
  return new LevelPoint(tMs, e2, testosteroneAt(e2, new TConfig(tConfig.base, tConfig.floor, k)));
}

/** Niveau actuel (now par défaut) — enveloppe de levelAt. */
export function currentLevel(treatments, doseLogs, tConfig, tKPerEster = null, nowMs = Date.now(), scalePerEster = null) {
  return levelAt(treatments, doseLogs, nowMs, tConfig, tKPerEster, scalePerEster);
}

/**
 * Courbe complète entre startMs et endMs (grille régulière stepMs, défaut 1 h).
 *  - départ clampé à la 1ʳᵉ dose (rien avant) ;
 *  - k de la T suivant l'ester ACTIF à chaque point (curseur sur les doses triées) ;
 *  - modèle via `modelOverride`, calibration via `scalePerEster` + `tKPerEster`.
 *
 * @param {object[]} treatments
 * @param {object[]} doseLogs
 * @param {number} startMs
 * @param {number} endMs
 * @param {number} [stepMs=HOUR_MS]
 * @param {TConfig} tConfig
 * @param {object} [opts={}]  { modelOverride, scalePerEster, tKPerEster }
 * @returns {LevelPoint[]}
 */
export function computeCurve(treatments, doseLogs, startMs, endMs, stepMs, tConfig, opts = {}) {
  if (stepMs === undefined) stepMs = HOUR_MS;
  const { modelOverride = null, scalePerEster = null, tKPerEster = null } = opts;
  if (treatments.length === 0 || doseLogs.length === 0 || endMs <= startMs) return [];

  const relevantTreatments = treatments.filter(
    (tr) => tr.type === 'ESTRADIOL' && doseLogs.some((d) => d.treatmentId === tr.id),
  );
  if (relevantTreatments.length === 0) return [];

  const earliestDose = doseLogs.reduce((m, d) => Math.min(m, d.timestamp), Infinity);
  const searchStart = Math.max(startMs, earliestDose);

  const batemanCache = new Map();
  relevantTreatments.forEach((tr) => batemanCache.set(tr.id, batemanParams(tr)));

  // Doses triées pour suivre l'« ester actif » le long de la grille (k de T)
  const sortedEstrogenDoses = doseLogs
    .filter((d) => relevantTreatments.some((tr) => tr.id === d.treatmentId))
    .slice()
    .sort((a, b) => a.timestamp - b.timestamp);
  let doseCursor = 0;
  let activeEster = null;

  const points = [];
  let t = searchStart;
  while (t <= endMs) {
    let e2 = 0.0;
    // Avance le curseur : la dose à t définit l'ester actif
    while (doseCursor < sortedEstrogenDoses.length && sortedEstrogenDoses[doseCursor].timestamp <= t) {
      const d = sortedEstrogenDoses[doseCursor];
      const tr = relevantTreatments.find((x) => x.id === d.treatmentId);
      activeEster = doseEster(tr, d);
      doseCursor++;
    }
    const k = (activeEster && tKPerEster && tKPerEster[activeEster]) ?? tConfig.k;
    for (const treatment of relevantTreatments) {
      const cutoff = cutoffHours(treatment);
      const p = batemanCache.get(treatment.id);
      for (const dose of doseLogs) {
        if (dose.treatmentId !== treatment.id || dose.timestamp > t) continue;
        const dtH = (t - dose.timestamp) / 3600000.0;
        if (dtH > cutoff) continue;
        const c = concentrationOfDose(treatment, dose, t, p, modelOverride);
        if (c > 0.0) {
          const scale = (scalePerEster && scalePerEster[doseEster(treatment, dose)]) ?? treatment.scaleFactor;
          e2 += c * scale;
        }
      }
    }
    points.push(new LevelPoint(t, e2, testosteroneAt(e2, new TConfig(tConfig.base, tConfig.floor, k))));
    t += stepMs;
  }
  return points;
}

// ---------------------------------------------------------------------------
// Prévision
// ---------------------------------------------------------------------------

/**
 * Génère les doses PRÉVISIONNELLES d'un traitement, à partir de sa Posologie
 * (`forecastIntervalDays`, en jours) et de la dernière dose réellement
 * enregistrée.
 *
 * Règles :
 *  - intervalle null ou ≤ 0 → aucune prévision ;
 *  - la première dose projetée suit EXACTEMENT l'intervalle après la dernière
 *    dose réelle (le rythme reste sous le contrôle de l'utilisatrice) ;
 *  - dose = standard du traitement, ester = override de la dernière injection ;
 *  - ⚠️ les créneaux déjà PASSÉS sont sautés au rythme configuré (#35, cas
 *    d'un OUBLI : simuler un créneau passé = faux pic dans l'historique) ;
 *    un RETARD décale naturellement toute la prévision (comportement voulu).
 *
 * Les doses retournées ne sont JAMAIS persistées : elles alimentent
 * computeCurve et les marqueurs du graphique quand le chip est actif.
 *
 * @param {object} treatment
 * @param {object[]} allDoseLogs
 * @param {number} toMs    Fin de la fenêtre de projection
 * @param {number} [nowMs=Date.now()]
 * @returns {object[]} DoseLog[] (objets purs, jamais en base)
 */
export function generateForecastDoses(treatment, allDoseLogs, toMs, nowMs = Date.now()) {
  const intervalDays = treatment.forecastIntervalDays;
  if (intervalDays === null || intervalDays === undefined || intervalDays <= 0.0) return [];
  const intervalMs = Math.trunc(intervalDays * 24.0 * HOUR_MS);
  if (intervalMs <= 0) return [];

  let last = null;
  for (const d of allDoseLogs) {
    if (d.treatmentId === treatment.id && d.timestamp <= nowMs) {
      if (last === null || d.timestamp > last.timestamp) last = d;
    }
  }
  if (!last) return [];

  const forecast = [];
  let t = last.timestamp + intervalMs;
  // ⚠️ Oubli d'une injection : le premier créneau théorique tombe dans le
  // PASSÉ → on avance au premier créneau STRICTEMENT FUTUR, au rythme
  // configuré (bug #35).
  while (t <= nowMs) t += intervalMs;
  while (t <= toMs) {
    forecast.push({
      id: 0, // jamais persisté
      treatmentId: treatment.id,
      timestamp: t,
      doseAmount: treatment.doseAmount,
      notes: null,
      esterType: doseEster(treatment, last),
    });
    t += intervalMs;
  }
  return forecast;
}