# -*- coding: utf-8 -*-
"""A翼 控制策略件 (windscada M9; SOP 2.5 + B + C 锚定).
轴:
  C1 双封顶轴 (SOP 2.5 准定论载体): P封顶=ActPower_max p99.5 / ω封顶=GenRpm_max p99.5 — fleet 分组+离群
      (机器内部量, mast-免疫; W006 已证真降容先例);
  C2 K=P/ω³ 聚类 (中载区): 【参考】(SOP: zyx=null 未正向实证, 不单独支撑准定论);
  C3 β-schedule (SOP B): 桨距角×功率档 per 台 schedule → fleet 同档差 → 标定偏移候选【准定论·预警载体】/控制参数【参考】;
  C4 尺子巡检 (SOP C): 满发功率/额定转速/额定桨角 三个'尺子'的 fleet 散布 — 全场性异常先查尺子.
窗: 2025-H2 (与曲线判别窗一致, 限电前干净窗). 成因(有意配置vs故障)须控制配置台账, 本件不判."""
import numpy as np, pandas as pd, pathlib
from src.windscada.config import farm
from app_ETL.app_ETL_guanlan.api import load_10min

WIN = ('2025-07-01', '2026-01-01')
PB = np.arange(0, 4200, 500)

def build_store(cfg=None, span=None, write=True):
    """控制策略件（M9）。`span=(起, 止)` 给了就按**所选时间窗**（含两端）重算（用户令 2026-09-21）。

    span 模式走 `slim10min` 窄仓（同一份 `load_10min` 的列子集），且 **不写盘**：
    正式产物 `control_profile/schedule.parquet` 的口径是"限电前干净窗（2025-H2）"，
    按窗重算是服务时算给页面看的，不能覆盖它。
    """
    cfg = cfg or farm()
    rows, sched = [], []
    for t in cfg['turbines']:
        try:
            if span:
                from app_ETL.app_ETL_guanlan.api import slim
                d = slim.load(t, cfg, span=span,
                              columns=['grd_wtc_ActPower_mean', 'grd_wtc_ActPower_max',
                                       'tur_wtc_GenRpm_mean', 'tur_wtc_GenRpm_max',
                                       'tur_wtc_PitcPosA_mean'])
            else:
                d = load_10min(t, cfg, groups=['A.功率', 'A.转速', 'B.变桨'])
                d = d[(d.ts >= WIN[0]) & (d.ts < WIN[1])]
            op = d[d['grd_wtc_ActPower_mean'] > 100]
            rec = dict(turbine=t,
                       p_cap=float(op['grd_wtc_ActPower_max'].quantile(0.995)),
                       w_cap=float(op['tur_wtc_GenRpm_max'].quantile(0.995)),
                       n=len(op))
            mid = op[(op['grd_wtc_ActPower_mean'] > 1200) & (op['grd_wtc_ActPower_mean'] < 2800)]
            k = mid['grd_wtc_ActPower_mean'] / (mid['tur_wtc_GenRpm_mean'] ** 3)
            rec['k_med'] = float(k.median()) * 1e6
            rated = op[op['grd_wtc_ActPower_mean'] > 3800]
            rec['pitch_rated'] = float(rated['tur_wtc_PitcPosA_mean'].median()) if len(rated) > 50 else np.nan
            rows.append(rec)
            op = op.copy(); op['pb'] = pd.cut(op['grd_wtc_ActPower_mean'], PB, labels=PB[:-1]).astype(float)
            g = op.groupby('pb')['tur_wtc_PitcPosA_mean'].median()
            for pb, v in g.items():
                if v == v: sched.append(dict(turbine=t, pb=float(pb), pitch=float(v)))
            if write:
                print(t, flush=True)
        except Exception as e:
            print(t, 'ERR', str(e)[:60], flush=True)
    prof, sch = pd.DataFrame(rows), pd.DataFrame(sched)
    if write:
        prof.to_parquet(pathlib.Path(cfg['store']) / 'control_profile.parquet')
        sch.to_parquet(pathlib.Path(cfg['store']) / 'control_schedule.parquet')
        print('→ control_profile / control_schedule')
    return prof, sch

def registry(cfg=None, span=None):
    """控制参数一致性。`span=(起, 止)` 给了就**按所选时间窗重算**（判据不变，只换数据切片）。"""
    cfg = cfg or farm()
    st = pathlib.Path(cfg['store'])
    if span:
        pr, sc = build_store(cfg, span=span, write=False)
        out = dict(窗=f'时间窗 {span[0]} ~ {span[1]}（随所选时间窗）')
    else:
        pr = pd.read_parquet(st / 'control_profile.parquet')
        sc = pd.read_parquet(st / 'control_schedule.parquet')
        out = dict(窗='时间窗 2025-H2（限电前干净时间窗，与曲线判别时间窗一致）')
    # C1 双封顶
    for ax, col, unit, tol in (('P封顶', 'p_cap', 'kW', 25), ('ω封顶', 'w_cap', 'rpm', 8)):
        med = pr[col].median()
        grp = (pr[col] / tol).round() * tol
        groups = grp.value_counts().sort_index()
        outliers = pr[abs(pr[col] - med) > 3 * tol]
        out[ax] = dict(中位=round(float(med), 1), 分组={f'{k:.0f}{unit}': int(v) for k, v in groups.items() if v >= 2},
                       离群=[dict(t=r.turbine, v=round(float(r[col]), 1)) for _, r in outliers.iterrows()],
                       判='版本/参数差异·准定论载体' if len(groups[groups >= 2]) > 1 or len(outliers) else '一致')
    # C2 K 聚类 (参考)
    kmed = pr.k_med.median(); kmad = (pr.k_med - kmed).abs().median() or 1e-6
    kout = pr[abs(pr.k_med - kmed) > 5 * 1.4826 * kmad]
    out['K聚类'] = dict(中位=round(float(kmed), 3), 离群=[dict(t=r.turbine, v=round(float(r.k_med), 3)) for _, r in kout.iterrows()],
                      判='【参考】(K轴未正向实证, 不单独支撑准定论)')
    # C3 β-schedule 同档差
    piv = sc.pivot_table(index='pb', columns='turbine', values='pitch')
    dev = (piv.sub(piv.median(axis=1), axis=0)).abs().mean()
    bout = dev[dev > 0.5].sort_values(ascending=False)
    out['β_schedule'] = dict(全场同档差中位=round(float(dev.median()), 3),
                             离群=[dict(t=k, 平均偏=round(float(v), 2)) for k, v in bout.items()],
                             判='标定偏移候选【准定论·预警载体】' if len(bout) else '一致 (全场σ极小)')
    # C4 尺子巡检
    out['尺子'] = {c: dict(σ=round(float(pr[col].std()), 3), 极差=round(float(pr[col].max() - pr[col].min()), 2))
                  for c, col in (('满发功率p99.5', 'p_cap'), ('转速封顶', 'w_cap'), ('额定桨角', 'pitch_rated'))}
    return out