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

# -*- coding: utf-8 -*-
"""A翼 功率曲线机群相对筛查 (windscada M3-A, 2026-08-24).
纪律 (skill power-curve-ntf + SOP): 机舱风=self_ref ⇒ 绝对达成率 INSUFFICIENT 不产出; 只做同机型机群相对偏差(筛查非定谳);
前置 = L1 工况门 (只取'正常发电'态: 剥限电命令面+停机); 双风速计一致性作 A 类门 (风速计故障台不进曲线横比).
自校验锚: 19# 零位 -0.87° (专项确诊, 未回正) ⇒ 部分负荷区效率损失, 相对偏差应为负; 31# 已归零 ⇒ 应近中位.
边界: 风速=机舱 SecAnemo 单计 (PriAnemo 死); 机群相对口径, 绝对达成率不产出."""
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

BINS = np.arange(3.0, 15.5, 0.5)   # 机舱风 bin (m/s); 额定以上桨控吸收, 相对比意义降 → 截 15
MIN_N_BIN = 24                      # 每 bin ≥4h
ANEMO_DIV_GATE = 0.08               # 双风速计中位相对差 >8% → 风速计嫌疑, 不进横比

def per_turbine(cfg=None, since='2025-07-01', until='2026-01-01'):
    cfg = cfg or farm()
    rows, flags = [], {}
    for t in cfg['turbines']:
        try:
            d = load_10min(t, cfg, groups=['A.功率', 'A.风况', 'A.转速'])
        except Exception:
            flags[t] = '数据缺失'; continue
        d = d[(d.ts >= since) & (d.ts < until)]
        d['state'] = classify(d, cfg['rated_kw'])
        n_all = len(d)
        d = d[d['state'] == '正常发电']
        ws, p = d['tur_wtc_SecAnemo_mean'], d['grd_wtc_ActPower_mean']   # 风速源=SecAnemo (PriAnemo 全场恒0.01 死, 数据仲裁 2026-08-24); 单风速计口径
        if ws.nunique() <= 2:
            flags[t] = '风速计恒值 (A类, 不进横比)'; continue
        if len(d) < 500:
            flags[t] = flags.get(t, '') + f' 干净行不足({len(d)})'
            continue
        m12 = d[d['grd_wtc_ActPower_mean'].between(1000, 2000)]
        ws12 = float(m12['tur_wtc_SecAnemo_mean'].median()) if len(m12) >= 200 else np.nan
        b = pd.cut(ws, BINS)
        g = d.groupby(b, observed=True).agg(p=('grd_wtc_ActPower_mean', 'median'), n=('grd_wtc_ActPower_mean', 'size'))
        g = g[g['n'] >= MIN_N_BIN]
        for iv, r in g.iterrows():
            rows.append(dict(turbine=t, bin=float(iv.mid), p=float(r['p']), n=int(r['n']), ws12=ws12))
    pc = pd.DataFrame(rows)
    fleet = pc.groupby('bin')['p'].median().rename('fleet_p')
    pc = pc.join(fleet, on='bin')
    pc['dev'] = pc['p'] / pc['fleet_p'] - 1
    out = pc.groupby('turbine').apply(lambda x: float(np.average(x['dev'], weights=x['n'] * x['fleet_p'])), include_groups=False).rename('dev_w').reset_index()
    ws12 = pc.groupby('turbine')['ws12'].first()
    out = out.join(ws12.rename('ws12'), on='turbine')
    out['ws_bias'] = out['ws12'] - out['ws12'].median()   # 同功率段(1-2MW)风速反演: 偏高→曲线假低 (03# +0.58m/s 解释其-14.7% 实逮)
    out['rank'] = out['dev_w'].rank(ascending=True).astype(int)
    def verdict(r):
        if r['dev_w'] < -0.03:
            if r['ws_bias'] == r['ws_bias'] and r['ws_bias'] > 0.3:
                return '风速计偏高候选(A类)'
            return '欠发候选(筛查级)'
        if r['dev_w'] > 0.05 and r['ws_bias'] == r['ws_bias'] and r['ws_bias'] < -0.3:
            return '风速计偏低候选(A类)'
        return '—'
    out['判别'] = out.apply(verdict, axis=1)
    for t, f in flags.items():
        out.loc[out.turbine == t, '判别'] = f'A类: {f}'
    return out.sort_values('dev_w'), flags, pc

def store(cfg=None, **kw):
    cfg = cfg or farm()
    out, flags, pc = per_turbine(cfg, **kw)
    st = pathlib.Path(cfg['store']); st.mkdir(parents=True, exist_ok=True)
    out.to_parquet(st / 'powercurve_dev.parquet'); pc.to_parquet(st / 'powercurve_bins.parquet')
    return out, flags


def attribution(cfg=None, since='2025-07-01', until='2026-01-01'):
    """欠发归因路由 (筛查级): 对曲线非常态台做机制交叉.
    轴: ①桨角挂高 (1-2MW 档 PitcPosA 中位 ×fleet, 挂高=欠发经典机制) ②跨域线索 (温度/变桨/偏航登记簿由调用方并读).
    A类(风速计偏置)直接路由校风, 不入欠发。窗口与判别窗一致 (2025-H2)."""
    from src.windscada.config import farm as _farm
    from app_ETL.app_ETL_guanlan.api import load_10min
    cfg = cfg or _farm()
    rows = []
    for t in cfg['turbines']:
        try:
            d = load_10min(t, cfg, groups=['A.功率', 'B.变桨'])
            d = d[(d.ts >= since) & (d.ts < until)]
            band = d[(d['grd_wtc_ActPower_mean'] >= 1000) & (d['grd_wtc_ActPower_mean'] <= 2000)]
            rows.append(dict(turbine=t, pitch_mid=float(band['tur_wtc_PitcPosA_mean'].median()), n=len(band)))
        except Exception as e:
            rows.append(dict(turbine=t, err=str(e)[:50]))
        print(t, flush=True)
    df = pd.DataFrame(rows)
    med = df['pitch_mid'].median()
    df['pitch_dev'] = df['pitch_mid'] - med
    df.to_parquet(pathlib.Path(cfg['store']) / 'underperf_pitchmid.parquet')
    print('→ underperf_pitchmid.parquet', df.shape)
    return df