guanlan-src/app_ontology/app_ontology_guanlan/sop/analysis_kit.py

# -*- coding: utf-8 -*-
"""analysis_kit — 分析层原子动作库 (2026-07-26)

**定位**: discriminators.py = 高层判别器 (nbm_residual/fleet_sector_devz…, 复用良好 25-36处);
本模块 = 其下的**原子动作层** —— 每个分析脚本都要写、每次都可能写歪的那些小函数。

**为什么建**: 如东一场 14 个 per 场脚本只 1 个 import 库; 4 个并行 agent 各自重写 MAD-z
(一个脚本内实现 19 次)、各自定义发电态 → 口径微差 → 独立审逮到"方向对、绝对值飘"。
**每个函数都是踩坑后定下的纪律**, 收进库 = 把纪律从文档里的话变成代码里的默认行为。

**已有工具索引 (别重复造! 这些在 discriminators.py 里)**:
  near_equal_columns      — 重复列/别名检测 (逐值近等+线性别名双路; ★禁只用 corr, 如东油压列误弃戒)
  reference_channel_liveness — 参考通道 per台存活校验 (淄川舱外温冻结戒)
  fleet_sector_devz       — 控扇区机群相对功率 devz
  nbm_residual            — 温度 NBM 残差
  availability_calibers   — A_WIL/A_PBA 多口径可用率 (★横比须传独立风列, 阜山 A11 戒)
  curtail_vs_fault_sync   — 限电 vs 故障同步性三证
  mast_usability / ntf_correction / cp_operating_peak / matched_load_overheat …

用法: `from src.sop.analysis_kit import mad_z, gen_mask, drop_partial_periods, ...`
"""
from __future__ import annotations

import numpy as np
import pandas as pd

MADK = 1.4826

__all__ = ['mad_z', 'near_miss', 'expected_fp', 'gen_mask', 'drop_partial_periods',
           'fleet_relative_dev', 'liveness', 'caliber_stamp', 'sensitivity_sweep',
           'theil_sen', 'episode_segments', 'changepoint_events', 'trajectory_shape',
           'repair_split_eval']


# ---------------------------------------------------------------- 统计原子
def mad_z(s, *, center=None, scale_floor=None, resolution=None):
    """稳健 z = (x − median) / (1.4826·MAD)。

    ★铁律 (阜山 A12 戒): **必须对全精度值算**, 禁 mad_z(round(x, 2)) ——
    门边台 z 偏移可达 0.07, 直接翻转 gate 成员 (A12 TSR −2.52 → −2.45 伪影)。
    round 只用于显示, 不用于判定。MAD=0 时返回全 0 (不炸)。

    ★★铁律二 (如东 34# 戒, 2026-08-07): **台间高度一致时 z 会虚高** ——
    分母 MADK·MAD 趋零, 微小绝对差被放大成巨大 z。实证: 高速轴窗内 stddev
    台间 MAD 仅 0.03, 34# 与 fleet 绝对差 **0.066℃** 却得 z=**−8.2**(本轮最强信号),
    复核后是假阳 —— 该测点**量化步长 1℃**, 0.066℃ 无物理意义。
      `resolution=`: 给通道的量化步长/噪声底, **绝对差 < resolution 的一律置 0**
                     (比设常数地板更物理: 闸绑定测量分辨率而非拍脑袋)
      `scale_floor=`: 直接给 MAD 下限 (SOP §7 缺口⑤ 原设想; 温控紧列 MAD~0.04℃ 型)
    ⚠适用边界 (如东润滑响应戒): 本闸**只约束单点/单窗比较**, 不约束大样本聚合量 ——
      180 次事件的中位差 0.59K 虽 < 1℃ 步长, 但大数定律下有效, 那种场景勿套本闸。
    """
    s = pd.Series(s, dtype='float64')
    med = float(s.median()) if center is None else float(center)
    mad = float((s - med).abs().median())
    scale = MADK * mad
    if scale_floor is not None:
        scale = max(scale, float(scale_floor))
    if scale <= 0:
        return s * 0.0
    z = (s - med) / scale
    if resolution is not None:
        z = z.where((s - med).abs() >= float(resolution), 0.0)
    return z


def near_miss(z, *, gate=2.5, band=(2.0, 2.5)):
    """门边披露: 返回 {'gate': 破门集, 'near_miss': 刀刃带集}。

    ★铁律 (精确集断言錯型, 已复发 3 次): 任何 `==[]` / `全台<阈` / `恰 N 台` 的
    空集/全称/精确计数断言, 若集合来自硬 cutoff, **交付前必带 near_miss**, 且
    绝对断言须配 sensitivity_sweep 多口径核 —— 换口径即翻的断言不可印给业主。
    """
    z = pd.Series(z, dtype='float64')
    a = z.abs()
    return {'gate': sorted(z.index[a > gate].tolist()),
            'near_miss': sorted(z.index[(a > band[0]) & (a <= band[1])].tolist()),
            'gate_thresh': gate, 'near_miss_band': list(band)}


def expected_fp(n_compare, alpha=0.046):
    """RULE-4: 筛查必附期望假阳数 = N_比较 × α。

    比较族按 **metric × 台 总数**计 (逐台筛查的多重空间是台数, 非 metric 数)。
    α 默认 0.046 ≈ |z|>2 双尾。返回 {'n_compare', 'alpha', 'expected_fp'}。
    """
    return {'n_compare': int(n_compare), 'alpha': alpha,
            'expected_fp': round(float(n_compare) * alpha, 2)}


def sensitivity_sweep(fn, values, *, label='caliber'):
    """口径敏感性核: 对一组口径值跑同一判定, 返回 {值: 结果} —— 供门边绝对断言使用。

    fn(v) 应返回可比对象 (如 flag 集合 / z 值)。用法:
        sweep = sensitivity_sweep(lambda c: set(redflag(cutin=c)), [2.5, 3.0, 3.5, 4.0])
    结果集随口径变 → 断言必须写成"口径敏感 (值域 …)", 禁写单点绝对结论。
    """
    out = {}
    for v in values:
        try:
            r = fn(v)
            out[str(v)] = sorted(r) if isinstance(r, (set, list, tuple)) else r
        except Exception as e:                      # 口径不适用如实报, 不静默吞
            out[str(v)] = f'ERROR: {type(e).__name__}: {e}'
    stable = len({str(x) for x in out.values()}) == 1
    return {'sweep_by_' + label: out, 'stable_across_calibers': stable,
            'note': '结果随口径变 → 断言须带口径限定, 禁单点绝对结论' if not stable else '口径稳健'}


# ---------------------------------------------------------------- 数据面原子
def gen_mask(df, *, power_col, rated_kw, min_frac=0.0125, setpoint_col=None,
             setpoint_min=None, extra=None):
    """发电态掩码 (统一口径, 防各脚本各写一版)。

    P > min_frac·rated (默认 1.25% ≈ 阜山 50/4000 量级) [∧ setpoint ≥ setpoint_min]。
    setpoint_col 给出时同时剥限电/降档样本 (阜山 pset≥2000 型)。
    返回 (mask, caliber_dict) —— caliber 直接塞进产物 json (见 caliber_stamp)。

    ⚠ **`setpoint ≥ 常数` 不等于"剥限电" —— 用前必先验该 setpoint 列是不是调度封顶**
      (2026-07-26 独立审逮, 原 docstring 的"如东 PowerRef≥3990 型"举例**是错的, 已删**):
      如东 `wtc_PowerRef_endvalue` **跟随可用功率**而非调度封顶 —— 实测 corr(PowerRef, 风速)=0.484,
      且逐风速 bin 的 PowerRef 中位紧贴该 bin 的 ActPower 中位 (5-6m/s: 535 vs 320;
      9-10m/s: 2234 vs 1943; >12m/s: 4000 vs 3986)。故 `PowerRef<3990` ≈ **"运行在额定以下"**,
      发电态占比 66.8% 只是海上场的低于额定占空比, **不是限电率**; 真限电只在**额定以上**可辨
      (>12m/s 段该占比降到 33.1%)。
      **正确判据**: setpoint 显著低于**该风速下的可达功率**, 或 `setpoint<额定 ∧ ws>额定风速`。
      误用后果 (本轮实证): 气动闸只剩 33.2% 样本, 风速中位 9.12→10.86 m/s,
      Region-2 (桨距<2°) 样本 −68% —— 恰好闸掉了 Cp 工作峰所在段。
    """
    m = df[power_col] > min_frac * rated_kw
    cal = {'power_col': power_col, 'rated_kw': rated_kw, 'min_frac': min_frac,
           'gen_thresh_kw': round(min_frac * rated_kw, 2)}
    if setpoint_col is not None and setpoint_min is not None:
        m = m & (df[setpoint_col] >= setpoint_min)
        cal.update({'setpoint_col': setpoint_col, 'setpoint_min': setpoint_min,
                    'note': '同时剥限电/降档样本'})
    if extra is not None:
        m = m & extra
        cal['extra_filter'] = str(extra.name) if hasattr(extra, 'name') else 'custom'
    cal['n_rows_kept'] = int(m.sum())
    return m, cal


def drop_partial_periods(df, *, time_col, freq='M', min_frac=0.8, expected_rows=None,
                         per_tid=None):
    """剔不完整周期 (★阜山残月戒: 2024-01 仅 13h 进月度回归 → E_mast σ 被搅 5 倍,
    A2 假下行 −14.8→+0.25, A15 假近门 −7.02→−0.59)。

    按 freq 分组, 保留行数 ≥ min_frac × 该周期期望行数的周期。
    expected_rows=None 时用各周期行数中位数作期望 (自适应, 免手填采样率)。
    per_tid=<台号列>: 多台长表**必传** —— 否则 fleet 合计行数会掩盖单台残月
    (且台数变化的周期会被误判完整)。返回 (df_filtered, caliber_dict)。
    """
    if per_tid is not None:                          # 多台长表: 按 (周期) 计每台平均行数
        n_tid = df[per_tid].nunique()
        p = df[time_col].dt.to_period(freq)
        cnt = (p.value_counts() / max(n_tid, 1)).round(0)
    else:
        p = df[time_col].dt.to_period(freq)
        cnt = p.value_counts()
    exp = float(cnt.median()) if expected_rows is None else float(expected_rows)
    keep = cnt[cnt >= min_frac * exp].index
    cal = {'freq': freq, 'min_frac': min_frac, 'expected_rows_per_period': round(exp, 1),
           'periods_kept': len(keep), 'periods_dropped': sorted(str(x) for x in cnt.index if x not in set(keep))}
    return df[p.isin(keep)], cal


def fleet_relative_dev(df, *, tid_col, time_col, value_col, freq='M', ref='median',
                       min_tids_per_period=3):
    """机群相对偏差 (扣同期共模) —— 趋势/横比的标准前处理。

    per (周期, 台) 聚合 → 减该周期 fleet 参考 (median/mean) → 相对偏差 %。
    ★为什么必须扣共模: 阜山单塔山地 ratio 季节摆 0.130 会冒充漂移 (13/18 台假 flag)。
    返回 (wide_df[台×周期 的相对偏差%], caliber_dict)。
    """
    p = df[time_col].dt.to_period(freq)
    piv = df.assign(_p=p).groupby([tid_col, '_p'])[value_col].median().unstack()
    # ★稀台周期守卫 (2026-07-26 轴1+5 重构逮): 某周期只剩 1-2 台时, 扣共模会把该台
    # 自身钉成 0 (单台) 或造 ±对称伪点 (两台) → 趋势斜率被伪造点污染 (如东 38B PV
    # 斜率 0.72 → 剔稀台月后 4.81)。周期有效台数 < min_tids_per_period 的列整列置 NaN。
    n_per_period = piv.notna().sum(axis=0)
    thin = n_per_period[n_per_period < min_tids_per_period].index.tolist()
    fleet = piv.median(axis=0) if ref == 'median' else piv.mean(axis=0)
    rel = piv.sub(fleet, axis=1).div(fleet, axis=1) * 100.0
    if thin:
        rel[thin] = np.nan
    return rel, {'freq': freq, 'ref': ref, 'value_col': value_col,
                 'min_tids_per_period': min_tids_per_period,
                 'thin_periods_nulled': [str(x) for x in thin],
                 'note': '扣同期 fleet 共模后的相对偏差(%); 趋势斜率须在此基础上算; 稀台周期已置NaN防伪造点'}


def liveness(df, *, tid_col, cols=None, kind='continuous', work_thresh=None,
             stuck_min_run=None, time_col=None):
    """per台×per列 存活扫描 (★阜山 A12 戒: 全场聚合 nonzero% 掩盖单台死列)。

    ★★ kind 必须按列型选 (如东两次误杀戒, SOP §2.2 CT-2 断言③):
      - 'continuous' 连续量 (温度/风速/功率): nonzero% + nuniq + std
      - 'event'      事件/计数/状态位列: **低 nuniq + 大量零值是正常形态非死列** →
                     只报取值枚举与非零行数, 不下死列判定 (须调用方与功率/状态行为交叉)
    work_thresh: 工作区间阈 (如转速 >10rpm) → 判活阈与工作阈**双报** (只报前者会低估死度)。
    stuck_min_run + time_col: 卡滞检测 (最长恒值 run; ★如东 03E 卡值达数月,
    nuniq 判据抓不到, 只有 run-length 逮得到)。

    ★★ 卡滞结果必须回喂给下游趋势/漂移计算 (2026-07-26 kit 重构对拍二次逮):
      ① 卡值未必是 0 —— 如东实测 0.010 (03E 连 6 月/10F 连 5 月中位恒 0.010),
         `x > 0` 型过滤完全无效, 须按实测卡值设阈 (查 stuck_value/月度中位);
      ② **只滤卡值行不够** —— 卡滞占前 5-7 个月时, 剩余月份的"从卡滞恢复"本身
         造成假斜率 (03E 滤行后仍 28.7%/yr 居首, 压过真候选 06F 的 4.05%);
         → **卡滞台须整体排除出漂移/趋势评估** (其漂移不可评估), 并在产物里
         显式列出被排除台及其假斜率, 禁静默丢弃。
    """
    cols = cols or [c for c in df.columns
                    if c not in (tid_col, time_col) and pd.api.types.is_numeric_dtype(df[c])]
    rows = []
    for tid, g in df.groupby(tid_col):
        for c in cols:
            s = g[c]
            r = {'tid': tid, 'col': c, 'kind': kind, 'n': len(s),
                 'notna_pct': round(100 * s.notna().mean(), 2), 'nuniq': int(s.nunique())}
            if kind == 'event':
                vc = s.dropna().value_counts().head(8)
                r.update({'values': {float(k): int(v) for k, v in vc.items()},
                          'nonzero_rows': int((s.fillna(0) != 0).sum()),
                          'verdict': 'EVENT_COL_需行为交叉判活 (禁用连续量口径判死)'})
            else:
                nz = float((s.fillna(0) != 0).mean())
                r.update({'nonzero_pct': round(100 * nz, 2), 'std': float(s.std() or 0)})
                if work_thresh is not None:
                    r['work_thresh'] = work_thresh
                    r['above_work_pct'] = round(100 * float((s > work_thresh).mean()), 2)
                r['verdict'] = ('DEAD' if (r['nuniq'] <= 1 or r['std'] == 0) else
                                'SPARSE' if nz < 0.01 else 'ALIVE')
            if stuck_min_run and time_col is not None and len(s) > 1:
                # ★NaN-run 与恒值-run 必须分扫 (2026-07-26 三路重构独立报同一 bug):
                # 混计会把断档误标 stuck (如东 24C 8.2天缺失段假阳), 且掩盖真恒值段。
                v = s.to_numpy()
                isna = pd.isna(v)
                mx_val, mx_val_v, mx_nan, i = 0, None, 0, 0
                while i < len(v):
                    j = i
                    if isna[i]:
                        while j + 1 < len(v) and isna[j + 1]:
                            j += 1
                        if j - i + 1 > mx_nan:
                            mx_nan = j - i + 1
                    else:
                        while j + 1 < len(v) and (not isna[j + 1]) and v[j + 1] == v[i]:
                            j += 1
                        if j - i + 1 > mx_val:
                            mx_val, mx_val_v = j - i + 1, v[i]
                    i = j + 1
                r['max_const_run'] = int(mx_val)          # 真恒值 run (不含 NaN)
                r['max_nan_run'] = int(mx_nan)            # 断档 run (单列报, 不当 stuck)
                if mx_val >= stuck_min_run:
                    r['stuck_flag'] = True
                    r['stuck_value'] = float(mx_val_v) if mx_val_v is not None else None
                if mx_nan >= stuck_min_run:
                    r['gap_flag'] = True                  # 断档≠卡滞: 处置不同(补数 vs 换传感器)
            rows.append(r)
    return pd.DataFrame(rows)


# ---------------------------------------------------------------- 交付原子
def caliber_stamp(**calibers):
    """口径戳: 把各步骤 caliber_dict 汇总成产物 json 的 `_caliber` 节。

    ★铁律 (独立审 R7 实测): 多 agent 数字"方向全一致但绝对值屡有口径差" ——
    无口径定义则跨场复用时"方向对、数字飘"。每个判据的口径 (列/阈/窗/台集/参考系)
    必随产物落盘。
    """
    return {'_caliber': calibers,
            '_caliber_note': '口径定义随产物落盘 (SOP 元层纪律); 复现须按此口径'}

# ------------------------------------------------- 轨迹/趋势原子 (2026-08-07 如东回灌)
def theil_sen(y, x=None):
    """Theil-Sen 稳健斜率 (公开原子; 解 SOP §7 缺口②)。

    原 `discriminators._theil_sen(y)` 私有且签名只收 y (隐含等间距索引), 外部脚本
    无法传自定义 x → 本轮如东扫老化趋势时直接踩到 TypeError。此处公开并支持 x。

    ⚠**铁律 (SOP §7 缺口④原文): Theil 对阶跃会摊成假渐进** —— 断崖型失效若用
    Theil 拟合, 会报出一个"缓慢下降"的斜率, 掩盖真实的突变时刻。**用前必先过
    `trajectory_shape()` 判形态**: 只有 shape='drift' 才可用斜率描述。
    实证 (如东 10#): 月度 V_cyc 4.02→3.06 用 Theil 读作"缓降", 实为 2026-06-29
    **单周** 4.73→0.76 的断崖。
    """
    y = np.asarray(pd.Series(y, dtype='float64').dropna())
    n = len(y)
    if n < 3:
        return 0.0
    x = np.arange(n, dtype='float64') if x is None else np.asarray(x, dtype='float64')[-n:]
    i, j = np.triu_indices(n, k=1)
    dx = x[j] - x[i]
    ok = dx != 0
    if not ok.any():
        return 0.0
    return float(np.median((y[j][ok] - y[i][ok]) / dx[ok]))


def episode_segments(flag, *, min_len=1):
    """连续段分割 (解 SOP §7 缺口⑥): 布尔序列 → [(start_idx, end_idx, length), ...]。

    ★为什么需要 (如东 27# 戒, 2026-08-07): 只看"超限周占比"会把
    **一次长发作**与**多次短发作**混为一谈 —— 二者机制完全不同
    (前者=持续劣化, 后者=间歇性)。27# 41 周中 15 周超限(37%), 分割后是
    **5 个发作段**(长度 1/5/7/1/1) → 判为反复发作型, 而非单次事件。
    """
    f = pd.Series(flag).fillna(False).astype(bool).reset_index(drop=True)
    out, start = [], None
    for i, v in enumerate(f):
        if v and start is None:
            start = i
        elif not v and start is not None:
            if i - start >= min_len:
                out.append((start, i - 1, i - start))
            start = None
    if start is not None and len(f) - start >= min_len:
        out.append((start, len(f) - 1, len(f) - start))
    return out


def changepoint_events(s, *, k_mad=6.0, min_abs=None, rel_to=None):
    """轨迹突变点 (解 SOP §7 缺口③): 返回 [(idx, before, after, delta), ...]。

    单步差分超过 k_mad 倍**自身差分稳健尺度**即记一次突变。

    ★★铁律一 (如东 10# 戒): **粒度决定能否看见突变** —— 同一台 V_cyc,
    月粒度读作"缓慢下滑"(4.02→3.06), 周粒度才看出是 **2026-06-29 单周
    4.73→0.76(−84%)** 的断崖。**断崖型失效必须用周粒度复核, 月度是钝的**。
    ★★铁律二 (如东 27# 戒): 返回的 idx 是**差分序列的位置**, 对应
    "变化发生后的那一点"。本轮曾把它当成峰值周去做交叉核对, 查到的却是
    相邻的正常周, 得出"孤立通道"的**假阴性**。**定位事件请用 idxmax(|dev|),
    不要用差分索引**。
      `min_abs=`: 绝对幅度门 (配合 mad_z 的 resolution 思路, 防高一致序列虚报)
      `rel_to=`: 给基准值时, min_abs 按其比例解释
    """
    s = pd.Series(s, dtype='float64').dropna()
    if len(s) < 4:
        return []
    d = s.diff().dropna()
    mad = float((d - d.median()).abs().median()) * MADK
    if mad <= 0:
        # ★阶跃序列的 diff 大部分为 0 ⇒ MAD(diff)=0, 直接返回会漏掉唯一的那次突变
        # (本测试逼出的真 bug, 2026-08-07)。退回非零 diff 的稳健尺度。
        # MAD(diff)=0 ⇒ "绝大多数时刻无变化" ⇒ 任何非零变化都是突变。
        # (曾误用"非零 diff 的中位"做尺度 → 唯一那次突变自己就是中位, 永远检不出)
        out0 = []
        for pos0, (idx0, dv0) in enumerate(d.items()):
            if dv0 == 0:
                continue
            if min_abs is not None and abs(dv0) < (min_abs * (abs(rel_to) if rel_to else abs(float(s.median()))) if rel_to else min_abs):
                continue
            out0.append((idx0, float(s.iloc[pos0]), float(s.iloc[pos0 + 1]), float(dv0)))
        return out0
    base = float(abs(rel_to)) if rel_to else float(abs(s.median()))
    out = []
    for pos, (idx, dv) in enumerate(d.items()):
        if abs(dv - d.median()) <= k_mad * mad:
            continue
        if min_abs is not None and abs(dv) < (min_abs * base if rel_to else min_abs):
            continue
        out.append((idx, float(s.iloc[pos]), float(s.iloc[pos + 1]), float(dv)))
    return out


def trajectory_shape(s, *, gate=None, tail=6, head=6):
    """轨迹形状分型 (解 SOP §7 缺口④): 阶跃/反复发作/漂移/平台/正常。

    ★为什么形态比数值重要 (如东 2026-08-07 全轮教训): **"现在多差"不是紧迫性依据,
    "还在不在变"才是**。同为蓄能器失效: 已跌到底且斜率归零的台**损失已发生、
    再快也追不回**(可排期); 仍在下降的台**每早一周处置就少损失一周**(该抢)。
    本轮按"现值高低"排的派工序列因此是错的, 改按形态重排。

    返回 dict: shape ∈ {'cliff_then_plateau','recurrent','drifting','flat_offset','normal'}
    + n_episodes / last_mean / delta / slope。
      shape='drifting' 时斜率才有意义 (见 theil_sen 的 Theil-对阶跃警告)。
    ★派工排序 override (如东 U6 外审双席 2026-08-08): 形态排序 (还在变>反复>已到底) 之上
      再加一层 — **已触保护阈随时跳机的台排最前** (绝对水平接近/超过保护整定值时,
      形态再"平"也压不过"随时会停机"; 排序=override(近保护阈) → 形态 → 现值)。
    """
    s = pd.Series(s, dtype='float64').dropna()
    if len(s) < 8:
        return {'shape': 'insufficient', 'n': len(s)}
    cps = changepoint_events(s, k_mad=6.0)
    last, first = float(s.iloc[-tail:].mean()), float(s.iloc[:head].mean())
    slope = theil_sen(s.values)
    g = gate if gate is not None else float(s.median() + 3 * MADK * (s - s.median()).abs().median())
    eps = episode_segments(s > g)
    tail_slope = theil_sen(s.iloc[-tail:].values)
    # ★优先级按 派工紧迫性 排 (如东 2026-08-07 对拍逼出的修正):
    #   还在变的最急 → 反复发作次之 → 已到底可排期 → 稳态偏置 → 正常。
    #   曾把 recurrent 放最前, 致尾斜率 −1.94 (仍在急降) 的 7# 被判 recurrent, 掩盖紧迫性。
    scale = MADK * float((s - s.median()).abs().median())
    if scale <= 0:
        scale = abs(float(s.median())) * 0.02 or 1e-9
    if abs(tail_slope) > 0.15 * scale:
        shape = 'drifting'
    elif len(eps) >= 2:
        shape = 'recurrent'
    elif cps:
        shape = 'cliff_then_plateau'
    elif abs(last - first) > scale:
        shape = 'flat_offset'
    else:
        shape = 'normal'
    return {'shape': shape, 'n_episodes': len(eps), 'n_changepoints': len(cps),
            'last_mean': last, 'first_mean': first, 'delta': last - first,
            'slope': slope, 'tail_slope': tail_slope}


def repair_split_eval(s, *, event, exclude=None, min_side=3):
    """换件/检修前后分段对比原子 (维修后评估 known-answer 回测; 如东 U4 五台换件回灌 2026-08-08)。

    s: 周期聚合 Series (月度中位等), index 为 'YYYY-MM' 字符串或可与 event 比较的时间标签。
    event: 换件时点 (同 index 类型)。event 当期**恒剔** (换件月混两态); exclude 再剔污染期
    (如充脂月污染油样/泵秒数)。

    ★铁律一 (更换是最好的实验): 换件在窗内 = 天然 known-answer 机制检验 —
      换后回归基准 = 根因在被换件上; 换后不回归/有残余 = 根因另有所在
      (如东 14# 换泵即愈=泵侧根因 vs 26# 清废脂不愈=排油侧根因; 19# 换泵后残余+20%=排油背压另存)。
    ★铁律二 (memory degradation-timewindow): 换件在窗前 (数据窗全在 event 后) → verdict='all_post',
      **全窗都是"换后状态", 禁把全窗统计当"换前基线"**; 判"现在坏不坏"必现窗。
      对偶: 数据窗全在 event 前 → 'all_pre' (换后效果本窗不可评)。
    ★铁律三 (不设 magic 阈): improved/worsened 的判定阈由调用方按通道物理定 (RULE-3),
      本原子只给 pre/post 中位 + delta + ratio; 方向语义 (低=好 or 高=好) 调用方持有。

    返回 {pre_median, post_median, delta, ratio, n_pre, n_post, verdict, excluded}
    verdict ∈ {'pre_post_available', 'all_post', 'all_pre', 'insufficient'}。
    """
    s = pd.Series(s).dropna()
    excl = {str(event)} | {str(x) for x in (exclude or [])}
    s = s[~s.index.astype(str).isin(excl)]
    pre = s[s.index < event]
    post = s[s.index > event]
    if len(s) == 0 or (len(pre) < min_side and len(post) < min_side):
        return {'pre_median': None, 'post_median': None, 'delta': None, 'ratio': None,
                'n_pre': int(len(pre)), 'n_post': int(len(post)),
                'verdict': 'insufficient', 'excluded': sorted(excl)}
    if len(pre) < min_side:
        return {'pre_median': None, 'post_median': round(float(post.median()), 3), 'delta': None,
                'ratio': None, 'n_pre': int(len(pre)), 'n_post': int(len(post)),
                'verdict': 'all_post', 'excluded': sorted(excl),
                'note': '换件在窗前: 全窗=换后状态, 禁当"换前基线"'}
    if len(post) < min_side:
        return {'pre_median': round(float(pre.median()), 3), 'post_median': None, 'delta': None,
                'ratio': None, 'n_pre': int(len(pre)), 'n_post': int(len(post)),
                'verdict': 'all_pre', 'excluded': sorted(excl),
                'note': '换件在窗末/窗外: 换后效果本窗不可评'}
    pm, qm = float(pre.median()), float(post.median())
    return {'pre_median': round(pm, 3), 'post_median': round(qm, 3),
            'delta': round(qm - pm, 3), 'ratio': round(qm / pm, 3) if pm else None,
            'n_pre': int(len(pre)), 'n_post': int(len(post)),
            'verdict': 'pre_post_available', 'excluded': sorted(excl)}