guanlan-src/app_algorithmModel/app_algorithmModel_guanlan/perf/availability.py

# -*- coding: utf-8 -*-
"""A翼 可用率与损失账 (windscada M6b). 口径纪律:
①时间可用率=非(故障类停机)占比, '发电态占比≠可利用率'(低发可能低风) — 只出状态账+时间可用率, WIL口径需独立风不产出;
②损失=Σ签号(潜在−实发) 分态归集, 禁逐行clip(整流噪声实逮16.8%假损失); 配零本底自检(正常态签号和≈0);
③潜在功率=本台2025-H2基线曲线(bin插值, 机舱风口径); ④月度损失≤额定×时长 物理上限核; ⑤停机拆 调度令/其他."""
import numpy as np, pandas as pd, pathlib
from src.windscada.config import farm
from app_ETL.app_ETL_guanlan.api import load_10min
from .curtail import classify, dispatch_intervals, refine_stop

def _curve(cfg):
    pc = pd.read_parquet(pathlib.Path(cfg['store']) / 'powercurve_bins.parquet')
    return {t: g.set_index('bin')['p'].sort_index() for t, g in pc.groupby('turbine')}

def build(cfg=None):
    cfg = cfg or farm()
    curves = _curve(cfg)
    iv = dispatch_intervals(cfg)
    rows = []
    for t in cfg['turbines']:
        try:
            d = load_10min(t, cfg, groups=['A.功率', 'A.风况', 'A.转速'])
        except Exception:
            continue
        d['state'] = classify(d, cfg['rated_kw'])
        d = refine_stop(d, t, iv)
        cv = curves.get(t)
        if cv is None or len(cv) < 8:
            continue
        ws = d['tur_wtc_SecAnemo_mean'].clip(lower=0)
        pot = np.interp(ws, cv.index.values, cv.values, left=0, right=float(cv.values[-8:].mean()))
        pot = np.where(ws < 3.0, 0.0, pot)   # 切入下潜在=0
        act = d['grd_wtc_ActPower_mean'].fillna(0).values
        d['loss_kwh'] = (pot - act) / 6.0    # 10min → kWh, 签号不clip
        d['month'] = d.ts.dt.strftime('%Y-%m')
        g = d.groupby(['month', 'state']).agg(rows_=('state', 'size'), loss=('loss_kwh', 'sum'), act=('grd_wtc_ActPower_mean', lambda x: float(x.fillna(0).sum() / 6))).reset_index()
        g['turbine'] = t
        rows.append(g)
        print(t, flush=True)
    out = pd.concat(rows, ignore_index=True)
    st = pathlib.Path(cfg['store']); out.to_parquet(st / 'loss_monthly.parquet')
    print('→ loss_monthly.parquet', len(out))

def summary(cfg=None, month_from='2026-01'):
    """口径 (2026-08-26 外审后钉紧): 时间可用率 = 1 − 停机/(总时长−调度令−低风待机);
    低风待机(SecAnemo<3m/s 切入以下)已剔出分母 — 与 IEC 61400-26 '风况在限内'方向一致
    (完整 IEC 口径还需高风切出/环境排除, 本件未做, 故仍标注非合同级)。"""
    cfg = cfg or farm()
    d = pd.read_parquet(pathlib.Path(cfg['store']) / 'loss_monthly.parquet')
    d = d[d.month >= month_from]
    # 零本底自检: 正常发电态 签号损失/实发 ≈ 0
    norm = d[d.state == '正常发电']
    base = float(norm.loss.sum() / max(norm.act.sum(), 1))
    hours = d.groupby('state')['rows_'].sum() / 6
    loss = d[d.state != '正常发电'].groupby('state')['loss'].sum() / 1000   # MWh
    stop_fault = d[d.state == '停机']
    # 实现对齐口径声明 (2026-08-26 自逮: docstring 钉的分母剔项此前未实现, 流通的 95.5% 系全分母产物):
    #   分母 = 总 − 调度令 − 低风待机
    denom = d.rows_.sum() - d[d.state.isin(['停机(调度令)', '低风待机'])].rows_.sum()
    avail_t = 1 - (stop_fault.rows_.sum() / max(denom, 1))
    # 物理上限核
    cap_ok = all(d.groupby(['turbine', 'month'])['loss'].sum() / 1000 <= cfg['rated_kw'] * 24 * 31 / 1000 + 1)
    per_t = d[d.state.isin(['停机', '限电·绑定'])].groupby(['turbine', 'state'])['loss'].sum().unstack(fill_value=0) / 1000
    return dict(零本底=round(base, 4), 时长h={k: round(float(v), 0) for k, v in hours.items()},
                损失MWh={k: round(float(v), 1) for k, v in loss.items()}, 时间可用率_不含调度令=round(float(avail_t), 4),
                物理上限核=cap_ok, 台级损失top=per_t.sum(axis=1).nlargest(5).round(1).to_dict())