#!/usr/bin/env python3
"""Regenerate the two source-based MHD figures and their download data."""
from pathlib import Path
import csv
import html
import json
import xml.etree.ElementTree as ET

ROOT = Path(__file__).resolve().parent
NAVY = "#142f4b"
BLUE = "#2e79b7"
PINK = "#d95498"
GRID = "#d9e3ed"
INK = "#14253a"
MUTED = "#51677e"

NASA2011 = [
    {"mach": 4.75, "mdot_kg_s": 10.0, "eta": 0.50, "T0_K": 1180.0, "figure5_approx_MW": 69.0, "figure5_read_uncertainty_MW": 2.0},
    {"mach": 5.25, "mdot_kg_s": 10.0, "eta": 0.50, "T0_K": 1380.0, "figure5_approx_MW": 88.0, "figure5_read_uncertainty_MW": 2.0},
]
CP_KJ_KG_K = (1.0, 1.2)  # exposed sensitivity interval; exact NPSS Cp not reported at these points
NASA1980 = [
    {"area_ratio": 2.56, "gross_electric_kW": 87.5, "thermal_input_MW": 7.0},
    {"area_ratio": 4.00, "gross_electric_kW": 153.0, "thermal_input_MW": 7.0},
]


def esc(s):
    return html.escape(str(s), quote=True)


def svg_base(width, height, title, desc, content):
    return f'''<svg xmlns="http://www.w3.org/2000/svg" viewBox="0 0 {width} {height}" role="img" aria-labelledby="title desc">
<title id="title">{esc(title)}</title><desc id="desc">{esc(desc)}</desc>
<rect width="{width}" height="{height}" rx="26" fill="#f7faff"/>
{content}</svg>\n'''


def text(x, y, value, size=18, weight=400, fill=INK, anchor="start"):
    return f'<text x="{x}" y="{y}" font-family="Arial, sans-serif" font-size="{size}" font-weight="{weight}" fill="{fill}" text-anchor="{anchor}">{esc(value)}</text>'


def nasa2011_svg(mobile=False):
    w, h = (760, 840) if mobile else (1200, 650)
    left, right = (148, w - 38) if mobile else (250, 1110)
    top, bottom = (235, 610) if mobile else (180, 500)
    x = lambda v: left + (right - left) * v / 100
    bits = [text(44, 57, "A published MHD power balance to check", 27 if mobile else 30, 700)]
    bits.append(text(44, 96, "NASA/TM-2011-217136 · ambient Mach points described in the report", 16 if mobile else 18, 400, MUTED))
    bits.append(text(44, 142, "Electrical power (MW)", 18, 700))
    for tick in range(0, 101, 20):
        xx = x(tick)
        bits.append(f'<line x1="{xx:.1f}" y1="{top}" x2="{xx:.1f}" y2="{bottom}" stroke="{GRID}" stroke-width="2"/>')
        bits.append(text(xx, bottom + 32, str(tick), 16, 400, MUTED, "middle"))
    points = NASA2011
    ys = [top + 94, top + 235] if mobile else [top + 85, top + 205]
    for row, (case, yy) in enumerate(zip(points, ys)):
        low = case["eta"] * case["mdot_kg_s"] * CP_KJ_KG_K[0] * case["T0_K"] / 1000
        high = case["eta"] * case["mdot_kg_s"] * CP_KJ_KG_K[1] * case["T0_K"] / 1000
        bits.append(text(left - 25, yy + 8, f'Mach {case["mach"]:.2f}', 17, 700, INK, "end"))
        bits.append(f'<line x1="{x(low):.1f}" y1="{yy}" x2="{x(high):.1f}" y2="{yy}" stroke="{BLUE}" stroke-width="16" stroke-linecap="round"/>')
        bits.append(f'<line x1="{x(low):.1f}" y1="{yy-12}" x2="{x(low):.1f}" y2="{yy+12}" stroke="{NAVY}" stroke-width="3"/><line x1="{x(high):.1f}" y1="{yy-12}" x2="{x(high):.1f}" y2="{yy+12}" stroke="{NAVY}" stroke-width="3"/>')
        figx = x(case["figure5_approx_MW"])
        bits.append(f'<circle cx="{figx:.1f}" cy="{yy}" r="12" fill="{PINK}" stroke="#ffffff" stroke-width="4"/>')
        bits.append(text(x(high) + 13, yy - 17, f'Equation range {low:.1f}–{high:.1f}', 15, 700, BLUE))
        read = case["figure5_approx_MW"]
        read_u = case["figure5_read_uncertainty_MW"]
        bits.append(text(figx, yy + 34, f'Figure 5 ≈ {read:.0f} ± {read_u:.0f} MW', 15, 700, PINK, "middle"))
    legend_y = 675 if mobile else 565
    legend_x = 44 if mobile else left
    bits.append(f'<line x1="{legend_x}" y1="{legend_y}" x2="{legend_x+36}" y2="{legend_y}" stroke="{BLUE}" stroke-width="12" stroke-linecap="round"/>')
    bits.append(text(legend_x + 50, legend_y + 6, 'Equation with Cp sensitivity (1.0–1.2 kJ/kg K)', 15, 400, MUTED))
    marker_x = 44 if mobile else left+390
    marker_y = legend_y+35 if mobile else legend_y
    bits.append(f'<circle cx="{marker_x}" cy="{marker_y}" r="9" fill="{PINK}"/>')
    bits.append(text(marker_x+18, marker_y + 6, 'Approximate reading of published plot', 15, 400, MUTED))
    note_y = 790 if mobile else 610
    bits.append(text(44, note_y, 'Plot readings are approximate; Cp interval is an explicit sensitivity assumption.', 13, 400, MUTED))
    bits.append(f'<line x1="{left}" y1="{bottom}" x2="{right}" y2="{bottom}" stroke="{NAVY}" stroke-width="2"/>')
    return svg_base(w, h, "NASA 2011 MHD power consistency check", "For Mach 4.75 and 5.25, blue intervals show 5.90 to 7.08 MW and 6.90 to 8.28 MW calculated from the report's stated mass flow, extraction ratio and temperature using an explicit specific-heat sensitivity interval. Pink circles show approximate 70 and 88 MW readings from the report's Figure 5. This is a reproducibility question, not a corrected result.", "\n".join(bits))


def nasa1980_svg(mobile=False):
    w, h = (760, 700) if mobile else (1200, 630)
    left, right = (165, w - 38) if mobile else (260, 1110)
    top, bottom = (205, 500) if mobile else (165, 430)
    maxv = 175
    x = lambda v: left + (right - left) * v / maxv
    bits = [text(42, 58, "A measured generator output — with its boundary", 27 if mobile else 30, 700)]
    bits.append(text(42, 97, "NASA hot H₂–O₂ combustion-flow experiment · 10 atm · 0.5 kg/s · 5 T", 16 if mobile else 18, 400, MUTED))
    bits.append(text(42, 142, "Gross electrical output at generator load (kW)", 17, 700))
    for tick in range(0, 176, 25):
        xx = x(tick)
        bits.append(f'<line x1="{xx:.1f}" y1="{top}" x2="{xx:.1f}" y2="{bottom}" stroke="{GRID}" stroke-width="2"/>')
        bits.append(text(xx, bottom + 30, str(tick), 14, 400, MUTED, "middle"))
    ys = [top + 72, top + 174] if mobile else [top + 70, top + 175]
    for case, yy in zip(NASA1980, ys):
        value = case["gross_electric_kW"]
        bits.append(text(left - 20, yy + 6, f'Area ratio {case["area_ratio"]:g}:1', 16, 700, INK, "end"))
        bits.append(f'<rect x="{left}" y="{yy-18}" width="{(right-left)*value/maxv:.1f}" height="36" rx="8" fill="{BLUE}"/>')
        bits.append(text(x(value) + 12, yy + 6, f'{value:g} kW', 16, 700, BLUE))
    bits.append(f'<line x1="{left}" y1="{bottom}" x2="{right}" y2="{bottom}" stroke="{NAVY}" stroke-width="2"/>')
    note_y = 575 if mobile else 520
    bits.append(text(42, note_y, 'The report states 7 MW thermal input: gross electric / thermal = 1.25% and 2.19%.', 15, 700, INK))
    bits.append(text(42, note_y + 32, 'Not a net efficiency: magnet, cooling, seeding and other plant loads are not closed.', 14, 400, MUTED))
    return svg_base(w, h, "NASA measured gross MHD generator output", "The 1980 NASA experiment measured 87.5 and 153 kilowatts gross electrical output in two ducts at 7 megawatts thermal input. Those are generator component measurements, not net plant efficiency or thrust returned by an accelerator.", "\n".join(bits))


def main():
    for mobile in (False, True):
        suffix = "-mobile" if mobile else ""
        (ROOT / f"nasa-2011-power-check{suffix}.svg").write_text(nasa2011_svg(mobile), encoding="utf-8")
        (ROOT / f"nasa-1980-generator-output{suffix}.svg").write_text(nasa1980_svg(mobile), encoding="utf-8")
    with (ROOT / "nasa-2011-power-check.csv").open("w", newline="", encoding="utf-8") as stream:
        writer = csv.writer(stream)
        writer.writerow(["ambient_mach", "mass_flow_kg_s_from_text_approx", "extraction_ratio_from_text_approx", "stagnation_temperature_K_from_text", "Cp_sensitivity_kJ_kg_K_low", "Cp_sensitivity_kJ_kg_K_high", "equation_power_MW_low", "equation_power_MW_high", "Figure_5_approx_reading_MW", "Figure_5_read_uncertainty_MW"])
        for c in NASA2011:
            vals = [c["eta"]*c["mdot_kg_s"]*cp*c["T0_K"]/1000 for cp in CP_KJ_KG_K]
            writer.writerow([c["mach"], c["mdot_kg_s"], c["eta"], c["T0_K"], *CP_KJ_KG_K, *vals, c["figure5_approx_MW"], c["figure5_read_uncertainty_MW"]])
    with (ROOT / "nasa-1980-generator-output.csv").open("w", newline="", encoding="utf-8") as stream:
        writer = csv.writer(stream)
        writer.writerow(["duct_area_ratio", "gross_electrical_output_kW", "thermal_input_MW", "gross_output_over_thermal_input_percent"])
        for c in NASA1980:
            writer.writerow([c["area_ratio"], c["gross_electric_kW"], c["thermal_input_MW"], 100*c["gross_electric_kW"]/(c["thermal_input_MW"]*1000)])
    methods_path = ROOT / "mhd-evidence-methods.json"
    methods = json.loads(methods_path.read_text(encoding="utf-8")) if methods_path.exists() else {}
    methods["NASA_2011"] = {
        "source": "https://ntrs.nasa.gov/api/citations/20110015567/downloads/20110015567.pdf",
        "equation": "eta_N(g) = P_elec / (mdot * Cp * T0,entrance) = (T0,12 - T0,14) / T0,12, printed p. 5 (PDF p. 9), Eq. (2)",
        "text_inputs": "printed p. 9 (PDF p. 13)",
        "plot_reading": "Figure 5, printed p. 8 (PDF p. 12), visually digitized from the embedded 749x426 px graph; approximate read uncertainty ±2 MW; not tabulated raw data",
        "Cp_sensitivity_kJ_kg_K": list(CP_KJ_KG_K),
        "status": "reproducibility question; cause unknown; no source values corrected",
        "provenance_correction": "An earlier project note gave the wrong equation page; Eq. (2) is printed p. 5 / PDF p. 9.",
    }
    methods_path.write_text(json.dumps(methods, indent=2) + "\n", encoding="utf-8")
    for name in ("nasa-2011-power-check.svg", "nasa-2011-power-check-mobile.svg", "nasa-1980-generator-output.svg", "nasa-1980-generator-output-mobile.svg"):
        ET.parse(ROOT / name)
    assert [round(c["eta"]*c["mdot_kg_s"]*1.0*c["T0_K"]/1000, 2) for c in NASA2011] == [5.90, 6.90]
    assert round(100*153/(7*1000), 2) == 2.19
    print("Generated and validated 4 responsive SVGs, 2 data tables, and source methods metadata.")


if __name__ == "__main__":
    main()