guanlan-src/scripts/rudong_fusion_run.py

#!/usr/bin/env python3
# -*- coding: utf-8 -*-
r"""六层链 · `fusion` 步（**逆向工程实现**，2026-09-17 用户令 #2）。

## 这个文件是怎么来的

包里的编排壳 `src/windcms/pipeline.py::analyze` 会调 `scripts/rudong_fusion_run.py`，而该脚本没随包
（六层链四步脚本都缺）。用户令 #2 问"能不能根据旧包产物逆向工程推导出来"——本文件是**可对拍的那部分**
的答案：`fleet_scalar_z.parquet` 的标准答案在盘上（随包件），于是口径可以**反推 + 逐值对拍**。

## 反推出来的口径（三条都逐值验过）

    val        = 窗内 median(scalar_value)，按 (turbine, sensor_name, meas_name, condition_key) 分组
    n          = 该组记录数
    fleet_med  = 跨机组的 median(val)（同一 sensor/meas/bin 下）
    z          = (val - fleet_med) / (1.4826 × MAD(val))      ← 稳健 z 分数（MAD 抗离群）

验证（`--verify`，样件 outputs/<场>/m5_cms_tcm/fleet_scalar_z.parquet，8,887 行）：
    样例 WTG01 / Gear_HS_generator_side / CrestFactor / WPS-ActivePower 0-1600,
      样件: val=5.66453  fleet_med=5.91444  z=-0.210156  n=7
      本器: val=5.66453  fleet_med=5.91444  z=-0.210153  n=7
    尺度候选实测: std=1.44859(z=-0.1725) · MAD=0.80209(z=-0.3116) · IQR/1.349=1.46899(z=-0.1701)
                  **1.4826×MAD=1.18918(z=-0.2102) ← 命中**

## 用法

    python scripts/rudong_fusion_run.py                    # 算并落盘 (默认窗口 = 配置里的首窗)
    python scripts/rudong_fusion_run.py --verify           # 与盘上样件逐值对拍 (不改任何件)
    python scripts/rudong_fusion_run.py --out <路径>       # 写到别处 (对拍/试验用)

★ 本器**只实现已被标准答案验证过的那一件**。`fusion_38.csv` 与 `model_run_l6.parquet` 从未随过包,
没有标准答案 ⇒ 不在本器里猜着写（本项目硬规矩: 允许响亮降级, 不许造数; 见
`docs/振动六层链_接口规格与缺口_v0.1.md`）。

退出码: 0 成功/对拍通过 · 5 对拍不一致 · 2 找不到输入
"""
from __future__ import annotations

import argparse
import pathlib
import sys

ROOT = pathlib.Path(__file__).resolve().parents[1]
sys.path.insert(0, str(ROOT))
from src import paths as P                                              # noqa: E402


def _window_index(win: str) -> pathlib.Path:
    """窗索引在哪: w0127 = 首窗特例 (索引在 m5 根), 其余在 windows/<w>/index.parquet。"""
    m5 = P.m5()
    if win == 'w0127' and (m5 / 'tcm_index.parquet').is_file():
        return m5 / 'tcm_index.parquet'
    return m5 / 'windows' / win / 'index.parquet'


def compute(win: str = 'w0127', min_n: int = 3):
    """→ 与样件同构的 DataFrame（sensor, meas, bin, turbine, val, fleet_med, z, n）。"""
    import numpy as np
    import pandas as pd
    ix = _window_index(win)
    if not ix.is_file():
        raise SystemExit(f'[X] 窗索引不存在: {P.rel(ix)} —— 先跑 scripts/rudong_tcm_index.py')
    d = pd.read_parquet(ix, columns=['turbine', 'sensor_name', 'meas_name', 'condition_key', 'scalar_value', 'y_unit'])
    d['val'] = pd.to_numeric(d['scalar_value'], errors='coerce')
    d = d.dropna(subset=['val', 'turbine', 'sensor_name', 'meas_name', 'condition_key'])
    # ★ 量纲闸 (反推自样件, 逐值验过): 只留振动量纲的标量 —— 1(无量纲指标) / m/s / m/s²。
    #   排除 %(Disk Usage / Memory Usage)、RPM(rms_rawRPM_DC) 与全部波形/谱名(FFT_/Time_/Env_/Cep_)。
    #   加上这一条之后, 样件 8,887 行的 val / fleet_med / z **逐值 100% 一致**。
    d = d[d['y_unit'].astype(str).isin({'1', 'm/s', 'm/s²'})]
    g = (d.groupby(['sensor_name', 'meas_name', 'condition_key', 'turbine'])['val']
         .agg(val='median', n='size').reset_index())
    g = g[g['n'] >= min_n]      # ★ 样件口径: 每组至少 3 条记录 (反推自 n 分布 min=3)
    # 跨机组的稳健基线: 中位数 + 1.4826×MAD (与样件逐值一致)
    gm = g.groupby(['sensor_name', 'meas_name', 'condition_key'])['val'].median().rename('fleet_med')
    g = g.merge(gm, on=['sensor_name', 'meas_name', 'condition_key'], how='left')
    mad = (g.assign(dev=(g['val'] - g['fleet_med']).abs())
           .groupby(['sensor_name', 'meas_name', 'condition_key'])['dev'].median()
           .mul(1.4826).rename('scale'))
    g = g.merge(mad, on=['sensor_name', 'meas_name', 'condition_key'], how='left')
    g['z'] = (g['val'] - g['fleet_med']) / g['scale']
    out = g.rename(columns={'sensor_name': 'sensor', 'meas_name': 'meas', 'condition_key': 'bin'})
    return out[['sensor', 'meas', 'bin', 'turbine', 'val', 'fleet_med', 'z', 'n']].sort_values(
        ['sensor', 'meas', 'bin', 'turbine']).reset_index(drop=True)


def verify(win: str = 'w0127', sample: pathlib.Path | None = None) -> int:
    """与盘上样件逐值对拍（严格: 同键同值, 容差 1e-4）。"""
    import pandas as pd
    sp = pathlib.Path(sample) if sample else (P.m5() / 'fleet_scalar_z.parquet')
    if not sp.is_file():
        print(f'[X] 没有样件可对拍: {P.rel(sp)}')
        return 2
    got = compute(win)
    want = pd.read_parquet(sp)
    key = ['sensor', 'meas', 'bin', 'turbine']
    m = want.merge(got, on=key, how='outer', suffixes=('_样件', '_本器'), indicator=True)
    both = m[m['_merge'] == 'both']
    only_w = int((m['_merge'] == 'left_only').sum())
    only_g = int((m['_merge'] == 'right_only').sum())

    def near(col, tol=1e-4):
        a, b = pd.to_numeric(both[f'{col}_样件'], errors='coerce'), pd.to_numeric(both[f'{col}_本器'], errors='coerce')
        ok = (a - b).abs() <= tol
        return int(ok.sum()), int(len(ok) - ok.sum())

    print(f'== fleet_scalar_z 逐值对拍 · 窗 {win} ==')
    print(f'   样件 {len(want)} 行 · 本器 {len(got)} 行 · 同键 {len(both)} 行'
          f' · 仅样件 {only_w} · 仅本器 {only_g}')
    bad_any = 0
    for col in ('val', 'fleet_med', 'z', 'n'):
        if f'{col}_样件' not in both.columns:
            print(f'   {col}: (样件无此列)')
            continue
        ok, bad = near(col, tol=1e-4 if col != 'n' else 0.5)
        bad_any += bad
        print(f'   {col:10s} 一致 {ok:5d} / {len(both):5d}   不一致 {bad}')
    ok_all = (len(both) == len(want) == len(got) and bad_any == 0)
    print(f'   结论: {"逐值完全一致（口径已复现）" if ok_all else "有差异 —— 见上, 不要拿本器产物替换样件"}')
    return 0 if ok_all else 5


def main() -> int:
    ap = argparse.ArgumentParser(description='六层链 fusion 步（逆向工程实现: fleet_scalar_z）')
    ap.add_argument('--window', default='w0127', help='用哪个窗算（默认首窗 w0127，与样件同源）')
    ap.add_argument('--verify', action='store_true', help='与盘上样件逐值对拍，不写盘')
    ap.add_argument('--sample', default=None, help='--verify 用哪个样件（默认 outputs/<场>/m5_cms_tcm/fleet_scalar_z.parquet）')
    ap.add_argument('--out', default=None, help='输出路径（默认写回 m5/fleet_scalar_z.parquet）')
    ap.add_argument('--min-n', type=int, default=3, help='每组最少记录数（样件口径 = 3）')
    a = ap.parse_args()
    if a.verify:
        return verify(a.window, pathlib.Path(a.sample) if a.sample else None)
    df = compute(a.window, a.min_n)
    out = pathlib.Path(a.out) if a.out else (P.m5() / 'fleet_scalar_z.parquet')
    out.parent.mkdir(parents=True, exist_ok=True)
    df.to_parquet(out, index=False)
    print(f'已写 {P.rel(out)}: {len(df)} 行 × {len(df.columns)} 列（窗 {a.window}）')
    # 自登记 (产物来源自登记: 谁算的谁登记)
    # ★2026-09-18 修: 原先按 `record(store, rel, builder=…, by=…)` 逐件传参, 而
    #   src.derived_manifest.record 的签名是 `record(store_root, files: dict, by: str)` ——
    #   于是**登记从来没成功过**(被 except 吞成一行 [i] 提示): 产物在盘上、台账里没有来路。
    #   与 vib_raw_build.py 的调用形态对齐(那处传的是 dict)。
    try:
        from src import derived_manifest as DM
        DM.record(P.out_root(),
                  {out.relative_to(P.out_root()).as_posix():
                   'scripts/rudong_fusion_run.py (val=median(scalar); z=(val-fleet_med)/(1.4826*MAD); '
                   '逆向工程口径, 与随包样件逐值对拍 val 列 100%)'},
                  by='rudong_fusion_run')
        print('   已自登记 → _derived_manifest.json')
    except Exception as e:
        print(f'   [i] 自登记跳过: {type(e).__name__}: {e}')
    return 0


if __name__ == '__main__':
    for _s in (sys.stdout, sys.stderr):
        try:
            _s.reconfigure(errors='replace')
        except Exception:
            pass
    sys.exit(main())