# -*- coding: utf-8 -*-
"""
fig_dgcup2020a.py —— 2020 电工杯 A 题范文配图（24 张 SVG）
============================================================
数据全部来自 gen_dgcup2020a.gen_dgcup2020a()，与范文正文、附录同源，
保证"图 / 文 / 附录 / 源码"四方一致。

运行： python3 fig_dgcup2020a.py
输出： ../assets/problems/dgcup2020a/figs/dgcup2020a-{1,2,3}-fig{1..8}.svg
"""
import os
import math
import gen_dgcup2020a as G
from _svg import (_fig, _save, _txt, _bar, _grouped_bar, _line, _scatter,
                  _pie, _heatmap, _network, _flow,
                  C_ACC, C_RED, C_GREEN, C_AMB, C_PUR, C_CYAN, C_TEAL, C_PINK,
                  C_MUT, PALETTE)

HERE = os.path.dirname(os.path.abspath(__file__))
OUT = os.path.normpath(os.path.join(HERE, "..", "assets", "problems", "dgcup2020a", "figs"))
D = G.gen_dgcup2020a()
COND = G.COND
HKEY = (3, 5, 7, 9, 11, 13)


def _ds(seq, m):
    """把序列均匀降采样到不超过 m 个点。"""
    n = len(seq)
    if n <= m:
        return list(seq)
    step = max(1, n // m)
    return [seq[i] for i in range(0, n, step)]


def _wf(cd, npts=61, tmax=0.02):
    """由权威频谱解析重建标称化电流波形（p.u. 基波幅值）。"""
    q = D["pq"][cd]
    pts = []
    for k in range(npts):
        t = tmax * k / (npts - 1)
        v = math.sin(2 * math.pi * G.F0 * t)
        for h in HKEY:
            v += q["hr"][h] / 100.0 * math.sin(2 * math.pi * G.F0 * h * t + 0.45 * h)
        pts.append((t * 1000.0, v))
    return pts


# ============================ 范文一：电能质量 ============================
def p1f1():
    return _flow("图1 高铁牵引供电系统运行数据分析总体建模路线", [
        ("数据预处理", "20kHz 录波分窗、去趋势"),
        ("电能质量评估", "序分量/不平衡度/THD"),
        ("负荷辨识预测", "脉冲提取+分类+回归"),
        ("等值建模校验", "戴维南+RLC+误差"),
    ])


def p1f2():
    return _line("图2 三种运行工况馈线电流标称化波形（0~20 ms，基波幅值为 1 p.u.）",
                 [("空载", C_GREEN, _wf("空载")),
                  ("牵引", C_RED, _wf("牵引")),
                  ("制动", C_ACC, _wf("制动"))],
                 xlabel="时间 t / ms（标称化电流 p.u.）", ymin=-1.6, ymax=1.6)


def p1f3():
    vals = [[D["pq"][c]["I1"], D["pq"][c]["I2"], D["pq"][c]["I0"]] for c in COND]
    return _grouped_bar("图3 三种工况 220 kV 侧电流对称分量（A）",
                        COND, ["正序 I1", "负序 I2", "零序 I0"], vals,
                        [C_ACC, C_RED, C_GREEN], fmt="%.3f")


def p1f4():
    vals = [[D["pq"][c]["eps_I"], D["pq"][c]["eps_U"] * 10.0] for c in COND]
    return _grouped_bar("图4 电流不平衡度与电压不平衡度（电压值已放大 10 倍便于同轴显示，%）",
                        COND, ["电流不平衡度", "电压不平衡度×10"], vals,
                        [C_PUR, C_AMB], fmt="%.2f")


def p1f5():
    M = [[D["pq"][c]["hr"][h] for h in HKEY] for c in COND]
    return _heatmap("图5 三种工况馈线电流各次谐波含有率热力图（%）",
                    COND, ["%d 次" % h for h in HKEY], M, fmt="%.2f")


def p1f6():
    vals = [D["pq"][c]["THD_TR"] for c in COND]
    return _bar("图6 三种工况馈线电流总谐波畸变率 THD_I（%，国标限值参考 8%）",
                COND, vals, [C_GREEN, C_ACC, C_RED], fmt="%.2f")


def p1f7():
    R = D["rec"]
    pts = [(t, R["P"][t]) for t in range(0, len(R["P"]), 4)]
    return _line("图7 全程 200 s 牵引侧有功功率曲线（含牵引—惰行—制动全过程）",
                 [("有功功率 P", C_ACC, pts)],
                 xlabel="时间 t / s（80 s 为牵引窗口，200 s 为制动窗口）",
                 ymin=-3.5, ymax=10.0,
                 xticks=[(0, "0"), (40, "40"), (80, "80"), (120, "120"),
                         (160, "160"), (200, "200")])


def p1f8():
    vals = [[abs(D["pq"][c]["P"]), D["pq"][c]["Q"], D["pq"][c]["S"]] for c in COND]
    return _grouped_bar("图8 三种工况牵引侧有功、无功与视在容量（MW / Mvar / MVA）",
                        COND, ["|有功 P|", "无功 Q", "视在 S"], vals,
                        [C_RED, C_AMB, C_TEAL], fmt="%.3f")


# ============================ 范文二：节能与预测 ============================
def p2f1():
    hm = [(h, D["hour_mean"][h]) for h in range(24)]
    hp = [(h, D["hour_pk"][h]) for h in range(24)]
    return _line("图1 全天 86400 点牵引负荷的逐小时平均功率与小时峰值（MW）",
                 [("逐小时平均", C_ACC, hm), ("小时峰值", C_RED, hp)],
                 xlabel="时刻 h / 时", ymin=0.0, ymax=12.0,
                 xticks=[(0, "0"), (4, "4"), (8, "8"), (12, "12"),
                         (16, "16"), (20, "20"), (23, "23")])


def p2f2():
    E = D["en"]
    used = E["E_bk"] * E["eta0"]
    lost = E["E_lost"]
    net = E["E_tr"] - used
    return _pie("图2 全天电量构成（kWh）：净取用电量 / 邻车吸收再生 / 弃用再生",
                ["净取用电量", "邻车吸收再生", "弃用再生"], [net, used, lost],
                [C_ACC, C_GREEN, C_RED])


def p2f3():
    A, B = D["en"]["A"], D["en"]["B"]
    vals = [[A["share"] * 100.0, A["eta"] * 100.0, A["cut"]],
            [B["share"] * 100.0, B["eta"] * 100.0, B["cut"]]]
    return _grouped_bar("图3 两种再生制动能量回收方案的份额、效率与节电率对比（%）",
                        ["方案A 车载储能", "方案B 地面储能"],
                        ["可利用份额", "综合效率", "节电率"], vals,
                        [C_CYAN, C_TEAL, C_GREEN], fmt="%.2f")


def p2f4():
    A, B = D["en"]["A"], D["en"]["B"]
    return _bar("图4 两方案投资与年节约电费对比（万元）",
                ["A 投资", "A 年节约", "B 投资", "B 年节约"],
                [A["inv"], A["save_y"], B["inv"], B["save_y"]],
                [C_MUT, C_GREEN, C_PUR, C_ACC], fmt="%.1f")


def p2f5():
    A, B = D["en"]["A"], D["en"]["B"]
    pa, pb = [], []
    p = 0.45
    while p <= 0.851:
        pa.append((p, A["inv"] / (A["rec"] * p / 10000.0 * 365.0)))
        pb.append((p, B["inv"] / (B["rec"] * p / 10000.0 * 365.0)))
        p += 0.05
    return _line("图5 电价灵敏度：电价 0.45~0.85 元/kWh 对两方案静态回收期的影响（年）",
                 [("方案A 车载储能", C_RED, pa), ("方案B 地面储能", C_ACC, pb)],
                 xlabel="电价 / (元·kWh⁻¹)", ymin=3.0, ymax=10.0,
                 xticks=[(0.45, "0.45"), (0.55, "0.55"), (0.63, "0.63"),
                         (0.75, "0.75"), (0.85, "0.85")])


def p2f6():
    M = D["mdl"]
    pts = _ds(list(zip(M["y"], M["yh"])), 60)
    ref = [(3.0, 3.0), (10.5, 10.5)]
    return _scatter("图6 牵引峰值功率回归模型拟合效果（R²=%.4f，RMSE=%.4f MW）"
                    % (M["R2"], M["RMSE"]),
                    [("实测—拟合散点", C_ACC, pts), ("理想对角线", C_RED, ref)],
                    xlabel="实测峰值 / MW（纵轴为模型拟合峰值 / MW）",
                    xmin=3.0, xmax=10.5, ymin=3.0, ymax=10.5)


def p2f7():
    lab = {(0, 0): "8编组下行", (0, 1): "8编组上行",
           (1, 0): "16编组下行", (1, 1): "16编组上行"}
    grp = {}
    segs = D["segs"]
    trs = D["trs"]
    for sg in segs:
        best = min(trs, key=lambda t: abs(t["t"] - sg["t"]))
        k = (best["g16"], best["up"])
        grp.setdefault(k, []).append((sg["pk"], sg["ratio"]))
    ser = []
    for i, k in enumerate(sorted(grp)):
        ser.append((lab[k], PALETTE[i], _ds(grp[k], 24)))
    return _scatter("图7 牵引脉冲特征空间（峰值—制动深度比）与四类聚簇分离效果",
                    ser, xlabel="牵引峰值 / MW（纵轴为制动峰值/牵引峰值）",
                    xmin=3.0, xmax=10.8, ymin=0.15, ymax=0.60)


def p2f8():
    P = D["prd"]
    a = [(h, P["hh_a"][h]) for h in range(6, 24)]
    p = [(h, P["hh_p"][h]) for h in range(6, 24)]
    return _line("图8 逐小时牵引电量实测与预测对比（kWh，MAPE=%.2f%%）" % P["h_mape"],
                 [("实测电量", C_ACC, a), ("模型预测", C_RED, p)],
                 xlabel="时刻 h / 时", ymin=0.0,
                 xticks=[(6, "6"), (9, "9"), (12, "12"), (15, "15"),
                         (18, "18"), (21, "21"), (23, "23")])


# ============================ 范文三：等值建模 ============================
def p3f1():
    return _flow("图1 牵引供电系统等值建模四步流程", [
        ("工况分窗", "空载/牵引/制动"),
        ("外特性辨识", "U–I 最小二乘"),
        ("参数等效", "R–L–C 分布参数"),
        ("误差校验", "谐波截断+压降"),
    ])


def p3f2():
    Q = D["eq"]
    pts = _ds(list(zip(Q["ii"], Q["uu"])), 60)
    ln = [(0.0, Q["E0"] * 1000.0), (350.0, Q["E0"] * 1000.0 - Q["Zeq"] * 350.0)]
    return _scatter("图2 牵引母线电压—电流外特性与戴维南等值拟合（R²=%.4f）" % Q["R2"],
                    [("逐秒实测点", C_ACC, pts), ("拟合直线端点", C_RED, ln)],
                    xlabel="馈线电流 I / A（纵轴为母线电压 U / V）",
                    xmin=0.0, xmax=350.0, ymin=26950.0, ymax=27560.0)


def p3f3():
    Q = D["eq"]
    pts = [(h, z) for h, z in Q["Zh"] if 5 <= h <= 25]
    return _line("图3 牵引网谐波阻抗频率特性（串联谐振点 %.1f Hz ≈ %.2f 次）"
                 % (Q["f_res"], Q["h_res_exact"]),
                 [("等值阻抗模 |Z(h)|", C_PUR, pts)],
                 xlabel="谐波次数 h（纵轴为阻抗模 / Ω）", ymin=0.0, ymax=270.0,
                 xticks=[(5, "5"), (9, "9"), (12, "12"), (15, "15"),
                         (19, "19"), (25, "25")])


def p3f4():
    nodes = [("s", "电源", C_MUT), ("t", "变压器", C_ACC), ("f", "馈线", C_CYAN),
             ("c", "接触网", C_TEAL), ("v", "动车组", C_RED), ("r", "钢轨", C_AMB),
             ("b", "储能", C_GREEN)]
    edges = [(0, 1, 3.0), (1, 2, 2.6), (2, 3, 2.4), (3, 4, 2.8),
             (4, 5, 2.2), (5, 1, 1.8), (1, 6, 1.4), (4, 6, 1.0)]
    return _network("图4 高铁牵引供电系统等值网络拓扑（节点间连线粗细表示功率通道权重）",
                    nodes, edges)


def p3f5():
    vals = [D["eq"]["chk"][c]["err"] for c in COND]
    return _bar("图5 三种工况等值模型电流有效值相对误差（%，工程可接受阈值 5%）",
                COND, vals, [C_RED, C_GREEN, C_ACC], fmt="%.2f", ymax=6.0)


def p3f6():
    ch = D["eq"]["chk"]
    vals = [[ch[c]["err_h"], ch[c]["err_f"], ch[c]["err"]] for c in COND]
    return _grouped_bar("图6 等值模型误差来源分解：谐波截断误差、工频压降误差与合成误差（%）",
                        COND, ["谐波截断", "工频压降", "合成误差"], vals,
                        [C_PUR, C_AMB, C_RED], fmt="%.2f")


def p3f7():
    Q = D["eq"]
    L = Q["L"] / 1000.0
    pts = []
    c = 1.0
    while c <= 3.61:
        f = 1.0 / (2.0 * math.pi * math.sqrt(L * (c + Q["C_line"]) * 1e-6))
        pts.append((c, f / G.F0))
        c += 0.25
    return _line("图7 补偿电容容量灵敏度：补偿电容 1.0~3.5 μF 对串联谐振次数的影响",
                 [("谐振次数 h_res", C_TEAL, pts)],
                 xlabel="补偿电容 C_comp / μF（纵轴为谐振次数）",
                 ymin=6.0, ymax=18.0,
                 xticks=[(1.0, "1.0"), (1.85, "1.85"), (2.5, "2.5"), (3.5, "3.5")])


def p3f8():
    M, Q = D["mdl"], D["eq"]
    labs = ["外特性R²", "回归R²", "编组辨识", "方向辨识", "六类准确", "提取成功"]
    vals = [Q["R2"] * 100.0, M["R2"] * 100.0, M["acc_g"], M["acc_d"],
            M["acc4"], D["extract_rate"]]
    return _bar("图8 全题模型可信度指标汇总（%，R² 已折算为百分数）",
                labs, vals, [C_ACC, C_CYAN, C_TEAL, C_GREEN, C_PUR, C_AMB],
                fmt="%.2f", ymax=110.0)


JOBS = [
    (1, 1, p1f1), (1, 2, p1f2), (1, 3, p1f3), (1, 4, p1f4),
    (1, 5, p1f5), (1, 6, p1f6), (1, 7, p1f7), (1, 8, p1f8),
    (2, 1, p2f1), (2, 2, p2f2), (2, 3, p2f3), (2, 4, p2f4),
    (2, 5, p2f5), (2, 6, p2f6), (2, 7, p2f7), (2, 8, p2f8),
    (3, 1, p3f1), (3, 2, p3f2), (3, 3, p3f3), (3, 4, p3f4),
    (3, 5, p3f5), (3, 6, p3f6), (3, 7, p3f7), (3, 8, p3f8),
]


def main():
    if not os.path.isdir(OUT):
        os.makedirs(OUT)
    ok = 0
    bad = 0
    for paper, idx, fn in JOBS:
        name = "dgcup2020a-%d-fig%d.svg" % (paper, idx)
        try:
            _save(os.path.join(OUT, name), fn())
            ok += 1
        except Exception as e:
            bad += 1
            print("  失败 %s : %s" % (name, e))
    print("输出目录 %s" % OUT)
    print("生成 SVG %d 张，失败 %d 张" % (ok, bad))


if __name__ == "__main__":
    main()