guanlan-src/app_algorithmModel/app_algorithmModel_guanlan/subsys/temp_nbm.py

# -*- coding: utf-8 -*-
"""B翼 温度NBM (windscada M6a, 2026-08-24). 判据出处: skill temp-nbm + 如东教训七坑设防.
方法: 发电态×功率档(500kW) 同工况 fleet 相对偏差 (ref=2025-H2 基线代 / cur=末90天), n加权.
七坑: ①恒值哨兵(契约层已排GFilB2/B3) ②传感器不可信 dev≥15K→A类非热 ③返服run-in(cur运行量不足→不可判)
④配对差(主轴承双温/发电机NDE-DE) ⑤热汇列(环温/舱温自身偏→该台ΔT类解释存疑, 单列注记)
⑥冷浸盲区(发电态过滤天然规避) ⑦相对判据配绝对量(dev以K计, 报警需绝对dev≥4K非只z)."""
import numpy as np, pandas as pd, pathlib
from src.windscada.config import farm
from app_ETL.app_ETL_guanlan.api import load_10min, contracted_cols
from ..perf.curtail import classify

PBINS = np.arange(0, 4500, 500)
SENSOR_IMPLAUSIBLE_K = 15.0
# 季节解耦闸: 环温通道自身窗间位移超此值 ⇒ 两窗季节错配, common_shift 拒绝出结论 (见其 docstring)
SEASON_MATCH_MAX_K = 1.0
PAIRS = {'主轴承双温差': ('tmp_wtc_MBearHTm_mean', 'tmp_wtc_MBearGTm_mean'),
         '发电机轴承NDE-DE(G-R)': ('tmp_wtc_GenBeGTm_mean', 'tmp_wtc_GenBeRTm_mean')}   # 官方IO表: G=NDE, R=DE (2026-08-26 定案)
HEATSINK = ['tmp_wtc_AmbieTmp_mean', 'tmp_wtc_NacelTmp_mean']

# 跨线联动跟踪名单 (2026-08-26 与振动线约定, 机器位非口头):
#   触发 = 最新月 dev≥2.5K 或 位次进前3 → temp_overview 置顶显示 [跨线通知], 见即发信箱给振动线
WATCHLIST = {
    ('WTG17', 'tmp_wtc_GenBeRTm_mean'): '17# 发电机DE(R端): 振动iso_rms_vel排1/35台, 温度2026-02起爬升(07月+2.3K第7)',
    ('WTG23', 'tmp_wtc_GenBeGTm_mean'): '23# 发电机NDE(G端): 振动单端1.5-1.9×fleet, 温度全程正常(对照台)',
    # 29# 机制=润滑断供·泵本体 (振动线v5定案, 2026-08-27更正: 非新报警非轴承损伤).
    # 温度/振动尖峰/柱塞=同一机制的果侧两观测+因侧一观测, 禁写"三源印证"(振动线④纪律采纳).
    # 柱塞两段式签名(mblub_monthly, 全纵向口径): 25-04~10单次时长由自身25-01~03的132~261s塌到34~57s→温度峰(25-05/08第1);
    # 26-07窗(至07-07)**单月突变**: 26-06单次68s完全正常→07窗1784s=自身18月区间(34~261s)的7~52倍,
    # 总打脂10.3h/天 vs 自身正常~1.4h/天, counts 6天窗应有288实测145=50% = 泵打不满挣扎态;
    # 突变后数日振动尖峰起(07-11首发在台风窗内, 窗外净尖峰07-24起). [口径勘误2026-08-27: 原"6×中位/70×
    # 对照台"为横比, 柱塞单次时长跨台天然异质(p90/p10=16~100×)横比无效, 双线一致弃用, 判读只认自身纵向]
    # 并读纪律: 温度dev与振动尖峰频次同步恶化才升级, 单独一个上升可能是工况.
    # 修泵后验收正本 = 振动线 findings.json TCM-16 (三分支决策树带量化门: A尖峰消14天/20条/max<8扣台风,
    # B柱塞回带尖峰不消才内窥, C柱塞不回带重排泵/管路; 温度不降不否证A). 本处只留指针不复制, 防载体漂移.
    # ★温度滞后条款(TCM-16③): 泵修好后温度dev允许滞后1~2月缓降 — watch通知若在修后窗内触发, 先对照TCM-16再升级.
    ('WTG29', 'tmp_wtc_MBearGTm_mean'): '29# 后主轴承(齿轮箱侧): 润滑断供机制(v5), 温度=果侧观测, 与振动尖峰并读',
    ('WTG29', 'tmp_wtc_MBearHTm_mean'): '29# 前主轴承(轮毂侧): 同上, 双端并盯',
}

def watchlist_check(cfg=None):
    """联动名单当前状态: 最新月 dev 与位次; 触发返回 [跨线通知] 行."""
    cfg = cfg or farm()
    import pandas as _pd
    f = pathlib.Path(cfg['store']) / 'watch_channels_monthly.parquet'
    if not f.exists():
        return ['[联动跟踪] watch_channels_monthly 仓缺失 — 跟踪失效, 须重建 scripts/windscada_watch_channels_build.py (守护失败必响亮)']
    df = _pd.read_parquet(f)
    out = []
    for (t, ch), why in WATCHLIST.items():
        piv = df.pivot_table(index='month', columns='turbine', values=ch)
        m = piv.index.max()
        row = piv.loc[m]
        dev = float(row[t] - row.median())
        rank = int(row.rank(ascending=False)[t])
        n = int(row.notna().sum())
        hit = dev >= 2.5 or rank <= 3
        tag = '🔔[跨线通知-发信箱给振动线]' if hit else '[联动跟踪中]'
        out.append(f"{tag} {t} {ch.replace('tmp_wtc_','').replace('_mean','')} {m}: dev{dev:+.1f}K 位次{rank}/{n} (门=dev≥2.5K∨前3) | {why}")
    return out

def build_store(cfg=None, span=None, write=True):
    """温度 NBM 比档件。`span=(起, 止)` 给了就把 `cur` 档换成**所选时间窗**（含两端，用户令 2026-09-21）；
    `ref`/`ref_sm` 两个对照窗保持不变；span 模式走 `slim10min` 窄仓且**不写盘**（不覆盖正式产物）。
    """
    cfg = cfg or farm()
    tmps = [x for x in contracted_cols(cfg, groups=['B.温度NBM']) if x.endswith('_mean')]
    rows, failed = [], []
    for t in cfg['turbines']:
        try:
            if span:
                from app_ETL.app_ETL_guanlan.api import slim
                # ★一次读、内存里切三窗（2026-09-21 性能实测）: 三窗各读一次 = 3×38 次 parquet 读，
                #   每台 261k 行 ⇒ 37 s/窗。改成"按列裁剪读一次 + 内存掩码切片" ⇒ 列少一半、读一次，
                #   口径不变（掩码仍按各窗边界，逐行同源）。
                _cols_t = ['ts', 'grd_wtc_ActPower_mean', 'tur_wtc_PowerRef_endvalue',
                           'tur_wtc_SecAnemo_mean', 'tur_wtc_GenRpm_mean',
                           *[c for c in tmps if c not in ('grd_wtc_ActPower_mean',)]]
                d = slim.load(t, cfg, span=None, columns=_cols_t)
            else:
                d = load_10min(t, cfg, groups=['A.功率', 'A.风况', 'A.转速', 'B.温度NBM'])
        except Exception as e:
            # 守护失败必响亮: 原为裸 `except: continue`, 数据根整个搬走时会**静默**产出空仓
            # (2026-09-07 实测: 源盘路径失效 → 38 台全失败 → 0 行落盘, 无任何报错, 只能手工 debug)。
            failed.append((t, f'{type(e).__name__}: {e}'))
            continue
        d['state'] = classify(d, cfg['rated_kw'])
        d = d[(d['state'].isin(['正常发电', '限电·绑定', '限电命令·未绑定'])) & (d['grd_wtc_ActPower_mean'] > 200)]
        cur_lo = str(d.ts.max() - pd.Timedelta(days=90))[:10]
        # ★ref_sm (2026-09-07): 与 cur **同历月、早一年**的季节匹配窗, 供 common_shift 做季节解耦。
        #   为什么不沿用互斥 np.where: ref_sm 与 ref 有约一周重叠, 互斥赋值会**静默削掉** ref_sm 一截
        #   (先匹配者胜)。改为每窗独立取子集 → 允许重叠, 且 **ref/cur 的内容逐行不变** (加窗零副作用)。
        #   span 模式: ref_sm 跟着**所选窗**早一年（季节匹配的是用户选的窗，不是"数据末端 90 天"）。
        if span:
            sm_lo = str((pd.Timestamp(span[0]) - pd.DateOffset(years=1)).date())
            sm_hi = str((pd.Timestamp(span[1]) - pd.DateOffset(years=1)).date())
        else:
            sm_lo = str(pd.Timestamp(cur_lo) - pd.DateOffset(years=1))[:10]
            sm_hi = str(d.ts.max() - pd.DateOffset(years=1))[:10]
        d['pbin'] = pd.cut(d['grd_wtc_ActPower_mean'], PBINS).map(lambda iv: iv.mid if pd.notna(iv) else np.nan)
        d = d.dropna(subset=['pbin'])
        # 掩码必须在 dropna **之后**建: 先建会让布尔序列索引与 d 错位, pandas 靠索引静默重对齐
        # (UserWarning: Boolean Series key will be reindexed) —— 那是会出错的地方, 不留。
        # 窗边界仍取自 state 过滤后、pbin 过滤前的 tmax (与加窗前口径一致, 保 cur 逐行不变)。
        # ★2026-09-21 用户令「判级按所选时间窗」: span 给了就把 `cur` 档换成所选窗（含两端），
        #   `ref`/`ref_sm` 仍固定（季节解耦的对照窗必须钉住，否则失去对照意义）。
        if span:
            _a, _b = pd.Timestamp(span[0]), pd.Timestamp(span[1]) + pd.Timedelta(days=1) - pd.Timedelta(seconds=1)
            cur_m = (d.ts >= _a) & (d.ts <= _b)
        else:
            cur_m = d.ts >= cur_lo
        wins = {'ref': (d.ts >= '2025-07-01') & (d.ts < '2026-01-01'),
                'cur': cur_m,
                'ref_sm': (d.ts >= sm_lo) & (d.ts <= sm_hi)}
        for wname, m in wins.items():
            dw = d[m]
            if dw.empty:                      # 早期场可能没有早一年的数据 → 该窗缺, common_shift 会 raise 而非静默
                continue
            if span:
                # ★按窗重算走**批量**分组（2026-09-21 性能实测: 逐通道 27 次 groupby/台 ⇒ 37 s/窗;
                #   一次 groupby 全通道 ⇒ ~2 s。口径逐格相同: 同 median、同 size 门 ≥18，
                #   只是行序不同 —— 默认路径保持原样不动，正式产物 byte 级不变。
                cols = [c for c in tmps if c in dw.columns]
                g = dw.groupby('pbin', observed=True)[cols].agg(['median', 'size'])
                for ch in cols:
                    sub = g[ch]
                    sub = sub[sub['size'] >= 18]
                    for pb, r in sub.iterrows():
                        rows.append(dict(turbine=t, channel=ch, win=wname, pbin=float(pb),
                                         med=float(r['median']), n=int(r['size'])))
                continue
            for ch in tmps:
                if ch not in dw.columns:
                    continue
                g = dw.groupby('pbin', observed=True)[ch].agg(['median', 'size']).reset_index()
                g = g[g['size'] >= 18]
                for _, r in g.iterrows():
                    rows.append(dict(turbine=t, channel=ch, win=wname, pbin=float(r['pbin']),
                                     med=float(r['median']), n=int(r['size'])))
        if write:
            print(t, flush=True)
    if failed:
        print(f'⚠ {len(failed)}/{len(cfg["turbines"])} 台加载失败: '
              + '; '.join(f'{t} {e}' for t, e in failed[:3])
              + (' …' if len(failed) > 3 else ''), flush=True)
    if not rows:
        raise RuntimeError(
            f'build_store 一行都没产出 ({len(failed)}/{len(cfg["turbines"])} 台加载失败) — '
            f'源目录 {cfg.get("src_10min")} 是否存在? 不静默落空仓覆盖既有 store。'
            + (f' 首个错误: {failed[0][1]}' if failed else ''))
    st = pathlib.Path(cfg['store']); st.mkdir(parents=True, exist_ok=True)
    out = pd.DataFrame(rows)
    if write:
        out.to_parquet(st / 'temp_bins.parquet')
        print('→ temp_bins.parquet', len(rows), f'(窗: {sorted({r["win"] for r in rows})})')
    return out

def registry(cfg=None, win='cur', span=None):
    """温度 NBM 判级。`span=(起, 止)` 给了就**按所选时间窗**重算 `cur` 档（用户令 2026-09-21）：
    `ref`（2025-07~2025-12 干净窗）与 `ref_sm`（同历月早一年）仍是**固定基线**——
    它们的作用是"季节解耦的对照窗"，跟着用户窗漂移就失去对照意义；页面据此如实标注。
    """
    cfg = cfg or farm()
    if span:
        tb = build_store(cfg, span=span, write=False)
    else:
        tb = pd.read_parquet(pathlib.Path(cfg['store']) / 'temp_bins.parquet')
    tb = tb[tb.win == win]
    # ★2026-09-21 实逮（服务器按窗重算时）: 窗内样本退化（如"近30日"=2026-09 只有 1 个节拍 ⇒ 比档件为空）时,
    #   `groupby(...).apply(...)` 返回的是**空 DataFrame** 而不是 Series ⇒ `.rename('dev_K')` 抛
    #   `TypeError: Index(...) must be called with a collection of some kind, 'dev_K' was passed`，
    #   把整面判级打断。退化窗应当**如实"不可判"**（样本不足），不是崩。
    if not len(tb):
        return pd.DataFrame(columns=['turbine', 'channel', 'dev_K', '判'])
    fleet = tb.groupby(['channel', 'pbin'])['med'].median().rename('fleet_med')
    tb = tb.join(fleet, on=['channel', 'pbin'])
    tb['dev'] = tb['med'] - tb['fleet_med']
    dev = tb.groupby(['turbine', 'channel']).apply(lambda x: float(np.average(x['dev'], weights=x['n'])), include_groups=False).rename('dev_K').reset_index()
    mad = dev.groupby('channel')['dev_K'].apply(lambda s: float(1.4826 * np.median(np.abs(s - s.median())))).rename('mad')
    dev = dev.join(mad, on='channel')
    dev['z'] = dev['dev_K'] / dev['mad'].clip(lower=0.3)
    nrows = tb.groupby('turbine')['n'].sum().rename('n_rows')
    dev = dev.join(nrows, on='turbine')
    def verdict(r):
        if r['n_rows'] < 500:
            return '不可判(样本不足/返服嫌疑)'
        if r['dev_K'] >= SENSOR_IMPLAUSIBLE_K:
            return '传感器候选(A类,非热)'
        if r['dev_K'] >= 4.0 and r['z'] >= 3.0:
            return '温度候选(定向查)'
        if r['dev_K'] >= 2.5 and r['z'] >= 2.0:
            return '记基线观察'
        return '—'
    dev['判'] = dev.apply(verdict, axis=1)
    dev['热汇注'] = dev['channel'].isin(HEATSINK)
    dev = dev.sort_values('dev_K', ascending=False)
    # ★批次盲声明 (实验 K): dev_K/z 是 fleet 相对量, 对全场同步平移逐值不变 ⇒ 本表**看不见批次退化**。
    #   判'—'只等于"该台不偏离机群", 不等于"该台正常"。共模轴须看 common_shift()。
    dev.attrs['batch_blind'] = {
        'blind_to': '全场同步平移 (批次退化 / 场级参考量变化)',
        'why': "dev_K = med − fleet_med; 全体同移则 dev_K 逐值不变 (实测 max|Δz|=0.00e+00)",
        'not_a_threshold_issue': True,
        'read_as': "判'—' ≠ 该台正常, 只等于'不偏离机群'; 共模轴见 temp_nbm.common_shift()",
    }
    return dev

def common_shift(cfg=None, win_cur='cur', win_ref=None):
    """共模图 — registry 缺的那另一半 (实验 K 接线 2026-09-07; 判据族见 sop.discriminators.fleet_two_component_split)。

    registry() 算 `fleet_med` 只为把它**减掉**, 自身从不被看。但 `d − fleet_med` 对全场同步平移
    **逐值恒等不变** (如东 38 台×4 通道实测 max|Δz|=0.00e+00) ⇒ **批次退化在 registry 上必然全绿**,
    且这不是阈值问题, 调阈值无效。本函数补上被丢掉的那一半: 同工况分档下 fleet_med 自身 cur vs ref 的位移。

    数据不用新采 —— `temp_bins.parquet` 里 ref(2025-07~2026-01) 与 cur(末90天) 两窗都在, registry 只读了 cur。

    ⚠⚠ **季节混杂 — 本函数在如东当前 store 上给不出结论 (2026-09-07 首跑实测, 非理论顾虑)**:
       build_store 的 ref(2025-07~2026-01) 与 cur(末90天) **不同历月**, pbin 只匹配**载荷**没匹配**环温**。
       实测 `amb_common_K = −5.5K` ⇒ 两窗差大半个季节。此时"raw 减环温"的一阶修正 (隐含各通道对环温增益≈1)
       把**全部 19 个通道同向抬成候选 = 19/19 全报 = 零分辨力** —— 真实增益逐通道不同 (热汇类≈1, 强自发热部件<1)。
       故本函数设**季节解耦闸**: |amb_common_K| ≥ SEASON_MATCH_MAX_K 时全表判 INSUFFICIENT, **拒绝出结论**,
       不发全红表让人自己悟。要真答案走 `sop.discriminators.fleet_two_component_split`
       (逐台 OLS 扣 P/环温 + 季节匹配 H0 + 时序 SPRT) 或 `nbm_residual(common_*_window=...)`。
       ⇒ 本表当前的用途 = **让批次盲这件事在 UI 上可见**, 外加 `apparent_gain` 诊断列; **不是**批次退化检出器。
       修法 (待裁, 需改 build_store 并重建 store): 加一个与 cur **同历月、早一年**的窗 (如 ref_sm), 两窗季节对齐。
    ⚠ 共模图**不能单独定因**: 全场同移既可能是批次退化, 也可能是场级参考量变化 (环温族偏/微气象/控制版本
       全场升级/限电制度切换) ⇒ 单靠本表最高 **准定论·预警**, 升定论须外部绝对锚。

    返回 DataFrame: channel / common_K_raw / common_K_net / amb_common_K / vs_mad(net 相对台间散布 MAD 的倍数)
                    / n_pbin / 判 / 热汇注。
    """
    cfg = cfg or farm()
    tb = pd.read_parquet(pathlib.Path(cfg['store']) / 'temp_bins.parquet')
    if win_ref is None:      # 优先季节匹配窗; 老 store 无 ref_sm 时退回 ref, 由季节解耦闸兜底判 INSUFFICIENT
        win_ref = 'ref_sm' if 'ref_sm' in set(tb['win'].unique()) else 'ref'
    if not {win_cur, win_ref} <= set(tb['win'].unique()):
        raise KeyError(f'temp_bins 缺窗 {win_cur}/{win_ref} (实有 {sorted(tb["win"].unique())}); '
                       '共模图是硬前置, 不静默跳过')
    # 每 (win, channel, pbin) 的 fleet 中位 = registry 里被减掉的那个量
    fm = tb.groupby(['win', 'channel', 'pbin']).agg(fleet_med=('med', 'median'), n=('n', 'sum')).reset_index()
    cur = fm[fm.win == win_cur].set_index(['channel', 'pbin'])
    ref = fm[fm.win == win_ref].set_index(['channel', 'pbin'])
    j = cur.join(ref, how='inner', lsuffix='_c', rsuffix='_r')
    j['shift'] = j['fleet_med_c'] - j['fleet_med_r']
    j['w'] = j[['n_c', 'n_r']].min(axis=1)
    raw = j.groupby('channel').apply(
        lambda x: float(np.average(x['shift'], weights=x['w'])), include_groups=False).rename('common_K_raw')
    npb = j.groupby('channel').size().rename('n_pbin')
    out = pd.concat([raw, npb], axis=1).reset_index()

    # 一阶季节修正: 减环温通道自身的共模位移 (假定增益≈1 — 粗, 见 docstring)
    amb_row = out[out.channel == 'tmp_wtc_AmbieTmp_mean']['common_K_raw']
    amb_k = float(amb_row.iloc[0]) if len(amb_row) else np.nan
    out['amb_common_K'] = amb_k
    out['common_K_net'] = out['common_K_raw'] - (amb_k if np.isfinite(amb_k) else 0.0)

    # 相对台间散布: 全场整体移了多少个"台间 MAD" — >1 即整体位移已超过 registry 能看见的散布
    dev = registry(cfg, win=win_cur)
    mad = dev.groupby('channel')['dev_K'].apply(
        lambda x: float(1.4826 * np.median(np.abs(x - x.median())))).rename('mad')
    out = out.join(mad, on='channel')
    out['vs_mad'] = out['common_K_net'] / out['mad'].clip(lower=0.3)
    # 表观增益 = 该通道随环温走了几成 (诊断量; 自发热强的部件 <1, 热汇类 ≈1)
    # ⚠ 只在环温**确有可观位移**时才算表观增益 —— 季节对齐后 amb_k→0, 除以近零分母会得出 7~8 的假增益
    #   (2026-09-07 装 ref_sm 后实测: amb_k 从 −5.54K 降到 −0.31K, gain 立刻虚涨到 7.5~8.4)。
    #   本列只是"季节错配时"的诊断量, 季节对齐时它无意义 → NaN, 不给假数。
    out['apparent_gain'] = (out['common_K_raw'] / amb_k
                            if (np.isfinite(amb_k) and abs(amb_k) >= SEASON_MATCH_MAX_K) else np.nan)

    # ★季节解耦闸 (2026-09-07 首跑实测逼出): 环温自身窗间位移 = 两窗季节错配的**实测**代理。
    #   如东实跑 amb_common_K = −5.5K ⇒ ref(7~12月) 与 cur(末90天) 差了大半个季节。此时
    #   "raw 减环温" 的一阶修正 (假定增益≈1) 会把**全部 19 通道**同向抬成候选 = 19/19 全报 = 零分辨力。
    #   故此处**不许出结论**: 全表判 INSUFFICIENT 并指向真判据, 而不是发一张全红的表让人自己悟。
    #   (失效必须落在显式状态 — 同 memory parse-failure-read-as-positive-result 的反向: 别让"全阳"冒充检出。)
    season_decoupled = bool(np.isfinite(amb_k) and abs(amb_k) < SEASON_MATCH_MAX_K)
    out.attrs['season_decoupled'] = season_decoupled
    out.attrs['amb_common_K'] = amb_k
    out.attrs['win_ref'] = win_ref      # 用了哪个参考窗 — 别让读者猜 ref 还是 ref_sm

    def verdict(r):
        if not np.isfinite(amb_k):
            return '不可判(环温通道缺, 季节无法修)'
        if not season_decoupled:
            return f'INSUFFICIENT(季节未解耦: 环温窗间移{amb_k:+.1f}K)'
        if r['channel'] == 'tmp_wtc_AmbieTmp_mean':
            return '— (季节基准本身)'
        if r['n_pbin'] < 3:
            return '不可判(共同档位<3)'
        a = abs(r['common_K_net'])
        if a >= SENSOR_IMPLAUSIBLE_K:
            return '共模巨移(口径/数据嫌疑, 非热)'
        if a >= 2.5 and abs(r['vs_mad']) >= 1.0:
            return '共模漂移候选(批次/场级, 须外部锚)'
        if a >= 1.5:
            return '记基线观察'
        return '—'
    out['判'] = out.apply(verdict, axis=1)
    out['热汇注'] = out['channel'].isin(HEATSINK)
    return out.sort_values('common_K_net', key=lambda s: s.abs(), ascending=False).reset_index(drop=True)


def pairs_table(cfg=None):
    cfg = cfg or farm()
    tb = pd.read_parquet(pathlib.Path(cfg['store']) / 'temp_bins.parquet')
    tb = tb[tb.win == 'cur']
    out = []
    for name, (a, b) in PAIRS.items():
        pa = tb[tb.channel == a].groupby('turbine')['med'].median()
        pb = tb[tb.channel == b].groupby('turbine')['med'].median()
        dd = (pa - pb).dropna().rename('pair_diff').reset_index()
        dd['pair'] = name
        med = dd['pair_diff'].median()
        dd['pair_dev'] = dd['pair_diff'] - med
        out.append(dd)
    return pd.concat(out, ignore_index=True)