control.py 5.9 KB

1234567891011121314151617181920212223242526272829303132333435363738394041424344454647484950515253545556575859606162636465666768697071727374757677787980818283848586878889909192939495969798
  1. # -*- coding: utf-8 -*-
  2. """A翼 控制策略件 (windscada M9; SOP 2.5 + B + C 锚定).
  3. 轴:
  4. C1 双封顶轴 (SOP 2.5 准定论载体): P封顶=ActPower_max p99.5 / ω封顶=GenRpm_max p99.5 — fleet 分组+离群
  5. (机器内部量, mast-免疫; W006 已证真降容先例);
  6. C2 K=P/ω³ 聚类 (中载区): 【参考】(SOP: zyx=null 未正向实证, 不单独支撑准定论);
  7. C3 β-schedule (SOP B): 桨距角×功率档 per 台 schedule → fleet 同档差 → 标定偏移候选【准定论·预警载体】/控制参数【参考】;
  8. C4 尺子巡检 (SOP C): 满发功率/额定转速/额定桨角 三个'尺子'的 fleet 散布 — 全场性异常先查尺子.
  9. 窗: 2025-H2 (与曲线判别窗一致, 限电前干净窗). 成因(有意配置vs故障)须控制配置台账, 本件不判."""
  10. import numpy as np, pandas as pd, pathlib
  11. from src.windscada.config import farm
  12. from app_ETL.app_ETL_guanlan.api import load_10min
  13. WIN = ('2025-07-01', '2026-01-01')
  14. PB = np.arange(0, 4200, 500)
  15. def build_store(cfg=None, span=None, write=True):
  16. """控制策略件(M9)。`span=(起, 止)` 给了就按**所选时间窗**(含两端)重算(用户令 2026-09-21)。
  17. span 模式走 `slim10min` 窄仓(同一份 `load_10min` 的列子集),且 **不写盘**:
  18. 正式产物 `control_profile/schedule.parquet` 的口径是"限电前干净窗(2025-H2)",
  19. 按窗重算是服务时算给页面看的,不能覆盖它。
  20. """
  21. cfg = cfg or farm()
  22. rows, sched = [], []
  23. for t in cfg['turbines']:
  24. try:
  25. if span:
  26. from app_ETL.app_ETL_guanlan.api import slim
  27. d = slim.load(t, cfg, span=span,
  28. columns=['grd_wtc_ActPower_mean', 'grd_wtc_ActPower_max',
  29. 'tur_wtc_GenRpm_mean', 'tur_wtc_GenRpm_max',
  30. 'tur_wtc_PitcPosA_mean'])
  31. else:
  32. d = load_10min(t, cfg, groups=['A.功率', 'A.转速', 'B.变桨'])
  33. d = d[(d.ts >= WIN[0]) & (d.ts < WIN[1])]
  34. op = d[d['grd_wtc_ActPower_mean'] > 100]
  35. rec = dict(turbine=t,
  36. p_cap=float(op['grd_wtc_ActPower_max'].quantile(0.995)),
  37. w_cap=float(op['tur_wtc_GenRpm_max'].quantile(0.995)),
  38. n=len(op))
  39. mid = op[(op['grd_wtc_ActPower_mean'] > 1200) & (op['grd_wtc_ActPower_mean'] < 2800)]
  40. k = mid['grd_wtc_ActPower_mean'] / (mid['tur_wtc_GenRpm_mean'] ** 3)
  41. rec['k_med'] = float(k.median()) * 1e6
  42. rated = op[op['grd_wtc_ActPower_mean'] > 3800]
  43. rec['pitch_rated'] = float(rated['tur_wtc_PitcPosA_mean'].median()) if len(rated) > 50 else np.nan
  44. rows.append(rec)
  45. op = op.copy(); op['pb'] = pd.cut(op['grd_wtc_ActPower_mean'], PB, labels=PB[:-1]).astype(float)
  46. g = op.groupby('pb')['tur_wtc_PitcPosA_mean'].median()
  47. for pb, v in g.items():
  48. if v == v: sched.append(dict(turbine=t, pb=float(pb), pitch=float(v)))
  49. if write:
  50. print(t, flush=True)
  51. except Exception as e:
  52. print(t, 'ERR', str(e)[:60], flush=True)
  53. prof, sch = pd.DataFrame(rows), pd.DataFrame(sched)
  54. if write:
  55. prof.to_parquet(pathlib.Path(cfg['store']) / 'control_profile.parquet')
  56. sch.to_parquet(pathlib.Path(cfg['store']) / 'control_schedule.parquet')
  57. print('→ control_profile / control_schedule')
  58. return prof, sch
  59. def registry(cfg=None, span=None):
  60. """控制参数一致性。`span=(起, 止)` 给了就**按所选时间窗重算**(判据不变,只换数据切片)。"""
  61. cfg = cfg or farm()
  62. st = pathlib.Path(cfg['store'])
  63. if span:
  64. pr, sc = build_store(cfg, span=span, write=False)
  65. out = dict(窗=f'时间窗 {span[0]} ~ {span[1]}(随所选时间窗)')
  66. else:
  67. pr = pd.read_parquet(st / 'control_profile.parquet')
  68. sc = pd.read_parquet(st / 'control_schedule.parquet')
  69. out = dict(窗='时间窗 2025-H2(限电前干净时间窗,与曲线判别时间窗一致)')
  70. # C1 双封顶
  71. for ax, col, unit, tol in (('P封顶', 'p_cap', 'kW', 25), ('ω封顶', 'w_cap', 'rpm', 8)):
  72. med = pr[col].median()
  73. grp = (pr[col] / tol).round() * tol
  74. groups = grp.value_counts().sort_index()
  75. outliers = pr[abs(pr[col] - med) > 3 * tol]
  76. out[ax] = dict(中位=round(float(med), 1), 分组={f'{k:.0f}{unit}': int(v) for k, v in groups.items() if v >= 2},
  77. 离群=[dict(t=r.turbine, v=round(float(r[col]), 1)) for _, r in outliers.iterrows()],
  78. 判='版本/参数差异·准定论载体' if len(groups[groups >= 2]) > 1 or len(outliers) else '一致')
  79. # C2 K 聚类 (参考)
  80. kmed = pr.k_med.median(); kmad = (pr.k_med - kmed).abs().median() or 1e-6
  81. kout = pr[abs(pr.k_med - kmed) > 5 * 1.4826 * kmad]
  82. out['K聚类'] = dict(中位=round(float(kmed), 3), 离群=[dict(t=r.turbine, v=round(float(r.k_med), 3)) for _, r in kout.iterrows()],
  83. 判='【参考】(K轴未正向实证, 不单独支撑准定论)')
  84. # C3 β-schedule 同档差
  85. piv = sc.pivot_table(index='pb', columns='turbine', values='pitch')
  86. dev = (piv.sub(piv.median(axis=1), axis=0)).abs().mean()
  87. bout = dev[dev > 0.5].sort_values(ascending=False)
  88. out['β_schedule'] = dict(全场同档差中位=round(float(dev.median()), 3),
  89. 离群=[dict(t=k, 平均偏=round(float(v), 2)) for k, v in bout.items()],
  90. 判='标定偏移候选【准定论·预警载体】' if len(bout) else '一致 (全场σ极小)')
  91. # C4 尺子巡检
  92. out['尺子'] = {c: dict(σ=round(float(pr[col].std()), 3), 极差=round(float(pr[col].max() - pr[col].min()), 2))
  93. for c, col in (('满发功率p99.5', 'p_cap'), ('转速封顶', 'w_cap'), ('额定桨角', 'pitch_rated'))}
  94. return out