"""Rebuild the NIST TN 2006 study figures from its public companion data."""
from __future__ import annotations

import csv
import hashlib
import html
import json
import math
from pathlib import Path

HERE = Path(__file__).resolve().parent
RESULTS = HERE / "TN2006-ILS-results.csv"
ATTRIBUTES = HERE / "TN2006-ILS-results-attributes.csv"
README = HERE / "nist-tn2006-dataset-readme.txt"
STRESS_STRAIN = HERE / "TN2006-stress-strain.csv"
STRESS_STRAIN_ATTRIBUTES = HERE / "TN2006-stress-strain-attributes.csv"
RATIOS = HERE / "nist-tn2006-table5-ratios.csv"
EXPECTED_HASHES = {
    RESULTS: "ea3923cec9af659e3d691ac7490ef1bc393477445dd51975a5b5b5524131fdee",
    ATTRIBUTES: "37f9a313a4ec17af9fc4d2414910a74791f009994eaad57d6d898fd63ae35cd7",
    README: "97a66b5a195535b9787c927a36eac84eb1e874fa687942a55e3c4b8c8478e21a",
    STRESS_STRAIN: "2ec86b9d1257bd6528ebcdd1095169decc48da6c23cec9718478ac4702db7da4",
    STRESS_STRAIN_ATTRIBUTES: "ae69b8a22ff407036dd4a742db56271830cbbed8a0a1e289e7737c0b039c6cf8",
}
EXPECTED_RATIOS = {
    "Tensile strength": 2.7,
    "0.2% offset yield strength": 4.9,
    "Elongation": 3.5,
    "Elastic modulus": 0.5,
}
COLORS = ["#5b3c91", "#007b92", "#c03a88", "#b05c16", "#32735b"]


def sha256(path: Path) -> str:
    h = hashlib.sha256()
    with path.open("rb") as f:
        for block in iter(lambda: f.read(1024 * 1024), b""):
            h.update(block)
    return h.hexdigest()


def check_sources() -> dict[str, str]:
    actual = {}
    for path, expected in EXPECTED_HASHES.items():
        actual[path.name] = sha256(path)
        if actual[path.name] != expected:
            raise ValueError(f"Source checksum mismatch for {path.name}")
    return actual


def load_data() -> tuple[list[dict[str, str]], list[dict[str, str]]]:
    with RESULTS.open(newline="", encoding="utf-8-sig") as f:
        all_rows = list(csv.DictReader(f))
    plotted = [r for r in all_rows if int(r["Lab"]) in range(1, 6)]
    if len(plotted) != 34 or any(not r["Sy"] for r in plotted):
        raise ValueError("Expected 34 reported yield-strength results from Labs 1–5")
    with ATTRIBUTES.open(newline="", encoding="utf-8-sig") as f:
        attributes = list(csv.DictReader(f))
    sy = next((r for r in attributes if r["Variable"] == "Sy"), None)
    if not sy or sy["Units"] != "MPa" or "0.2% offset yield strength" not in sy["Description"]:
        raise ValueError("The source attribute definition for Sy changed")
    counts = {lab: sum(int(r["Lab"]) == lab for r in plotted) for lab in range(1, 6)}
    if counts != {1: 8, 2: 6, 3: 6, 4: 7, 5: 7}:
        raise ValueError(f"Unexpected specimen counts: {counts}")
    with RATIOS.open(newline="", encoding="utf-8-sig") as f:
        ratios = list(csv.DictReader(f))
    ratio_map = {r["property"]: float(r["between_to_within_variability_ratio"]) for r in ratios}
    if ratio_map != EXPECTED_RATIOS:
        raise ValueError("Table 5 transcriptions differ from the audited source values")
    return plotted, ratios


def svg_header(width: int, height: int, title: str, desc: str) -> list[str]:
    return [
        f'<svg xmlns="http://www.w3.org/2000/svg" viewBox="0 0 {width} {height}" role="img" aria-labelledby="title desc">',
        f'<title id="title">{html.escape(title)}</title>',
        f'<desc id="desc">{html.escape(desc)}</desc>',
        '<rect width="100%" height="100%" fill="#fff"/>',
        '<style>text{font-family:Arial,sans-serif;fill:#192c42}.grid{stroke:#dbe1e8;stroke-width:1}.axis{stroke:#536477;stroke-width:1.5}.minor{font-size:14px;fill:#536477}.label{font-size:16px;font-weight:600}.heading{font-size:24px;font-weight:700}.sub{font-size:15px;fill:#536477}.value{font-size:14px;font-weight:700}</style>',
    ]


def write_lab_chart(rows: list[dict[str, str]], phone: bool = False) -> None:
    vals = [(int(r["Lab"]), r["Specimen"], float(r["Sy"])) for r in rows]
    out = []
    if not phone:
        w, h = 920, 520
        left, right, top, bottom = 90, 42, 110, 92
        ymin, ymax = 790.0, 960.0
        x0, x1 = left, w - right
        y0, y1 = top, h - bottom
        sy_y = lambda v: y1 - (v - ymin) / (ymax - ymin) * (y1 - y0)
        out += svg_header(w, h, "Reported yield strength by laboratory", "Thirty-four NIST test-level results for five laboratories in the NIST TN 2006 same-process comparison. Each dot is one reported 0.2 percent offset yield strength in megapascals.")
        out += [f'<text class="heading" x="{left}" y="40">Reported yield strength by laboratory</text>', f'<text class="sub" x="{left}" y="67">Five laboratories in the NIST comparison · one dot per reported test · MPa</text>']
        for tick in range(800, 951, 25):
            yy = sy_y(tick)
            out += [f'<line class="grid" x1="{x0}" y1="{yy:.2f}" x2="{x1}" y2="{yy:.2f}"/>', f'<text class="minor" x="{left-14}" y="{yy+5:.2f}" text-anchor="end">{tick}</text>']
        out += [f'<line class="axis" x1="{x0}" y1="{y0}" x2="{x0}" y2="{y1}"/>', f'<line class="axis" x1="{x0}" y1="{y1}" x2="{x1}" y2="{y1}"/>', f'<text class="label" x="24" y="{(y0+y1)/2}" text-anchor="middle" transform="rotate(-90 24 {(y0+y1)/2})">0.2% offset yield strength (MPa)</text>']
        centers = [x0 + (i + 0.5) * (x1 - x0) / 5 for i in range(5)]
        for lab, cx in enumerate(centers, 1):
            group = [r for r in vals if r[0] == lab]
            for j, (_, specimen, val) in enumerate(group):
                jitter = (j - (len(group) - 1) / 2) * 8
                yy = sy_y(val)
                out += [f'<circle cx="{cx+jitter:.1f}" cy="{yy:.2f}" r="6.5" fill="{COLORS[lab-1]}" fill-opacity=".82" stroke="#fff" stroke-width="1.5"><title>Lab {lab}, specimen {html.escape(specimen)}: {val:.5f} MPa</title></circle>']
            out += [f'<text class="label" x="{cx:.1f}" y="{y1+28}" text-anchor="middle">Lab {lab}</text>', f'<text class="minor" x="{cx:.1f}" y="{y1+52}" text-anchor="middle">n = {len(group)}</text>']
        out += [f'<text class="minor" x="{x0}" y="{h-14}">Each dot is a reported test result; horizontal spacing is only to separate overlapping marks.</text>', '</svg>']
        target = HERE / "nist-tn2006-yield-by-lab.svg"
    else:
        w, h = 420, 610
        left, right, top, bottom = 98, 25, 104, 74
        xmin, xmax = 790.0, 960.0
        x0, x1 = left, w - right
        ytop, ybottom = top, h - bottom
        xx = lambda v: x0 + (v - xmin) / (xmax - xmin) * (x1 - x0)
        out += svg_header(w, h, "Reported yield strength by laboratory", "Thirty-four NIST test-level results for five laboratories in the NIST TN 2006 same-process comparison. Each dot is one reported 0.2 percent offset yield strength in megapascals.")
        out += ['<text class="heading" x="20" y="39">Reported yield strength</text>', '<text class="sub" x="20" y="66">Five laboratories · one dot per test · MPa</text>']
        for tick in range(800, 951, 25):
            xt = xx(tick)
            out += [f'<line class="grid" x1="{xt:.2f}" y1="{ytop-10}" x2="{xt:.2f}" y2="{ybottom}"/>', f'<text class="minor" x="{xt:.2f}" y="{ytop-20}" text-anchor="middle">{tick}</text>']
        out += [f'<line class="axis" x1="{x0}" y1="{ytop-10}" x2="{x0}" y2="{ybottom}"/>']
        for lab in range(1, 6):
            group = [r for r in vals if r[0] == lab]
            cy = ytop + (lab - 0.5) * (ybottom - ytop) / 5
            out += [f'<text class="label" x="{left-12}" y="{cy+5:.1f}" text-anchor="end">Lab {lab}</text>']
            for j, (_, specimen, val) in enumerate(group):
                jitter = (j - (len(group) - 1) / 2) * 7
                cx = xx(val)
                out += [f'<circle cx="{cx:.2f}" cy="{cy+jitter:.1f}" r="6.5" fill="{COLORS[lab-1]}" fill-opacity=".82" stroke="#fff" stroke-width="1.5"><title>Lab {lab}, specimen {html.escape(specimen)}: {val:.5f} MPa</title></circle>']
            out += [f'<text class="minor" x="{x1}" y="{cy+5:.1f}">n={len(group)}</text>']
        out += [f'<line class="axis" x1="{x0}" y1="{ybottom}" x2="{x1}" y2="{ybottom}"/>', f'<text class="label" x="{(x0+x1)/2}" y="{h-24}" text-anchor="middle">0.2% offset yield strength (MPa)</text>', '</svg>']
        target = HERE / "nist-tn2006-yield-by-lab-phone.svg"
    target.write_text("\n".join(out), encoding="utf-8")


def write_ratio_chart(ratios: list[dict[str, str]], phone: bool = False) -> None:
    values = [(r["property"], float(r["between_to_within_variability_ratio"])) for r in ratios]
    if not phone:
        w, h = 900, 430
        left, right, top, bottom = 250, 90, 95, 82
        xmax = 5.5
        x0, x1 = left, w - right
        ytop, ybottom = top, h - bottom
        xx = lambda v: x0 + v / xmax * (x1 - x0)
        out = svg_header(w, h, "Between-laboratory variability relative to within-laboratory variability", "Four ratios transcribed from Table 5 of NIST TN 2006. A ratio of one is the reference where between-laboratory and median within-laboratory variability are equal. The values are 2.7 for tensile strength, 4.9 for yield strength, 3.5 for elongation and 0.5 for elastic modulus.")
        out += [f'<text class="heading" x="{left}" y="39">How much did laboratories differ?</text>', f'<text class="sub" x="{left}" y="66">Between-laboratory variability ÷ median within-laboratory variability</text>']
        for tick in range(0, 6):
            xt = xx(tick)
            out += [f'<line class="grid" x1="{xt:.2f}" y1="{ytop-10}" x2="{xt:.2f}" y2="{ybottom}"/>', f'<text class="minor" x="{xt:.2f}" y="{ybottom+27}" text-anchor="middle">{tick}×</text>']
        ref = xx(1)
        out += [f'<line x1="{ref:.2f}" y1="{ytop-13}" x2="{ref:.2f}" y2="{ybottom}" stroke="#192c42" stroke-width="2" stroke-dasharray="6 5"/>']
        for i, (name, val) in enumerate(values):
            yy = ytop + (i + 0.5) * (ybottom - ytop) / 4
            out += [f'<text class="label" x="{left-16}" y="{yy+6:.1f}" text-anchor="end">{html.escape(name)}</text>', f'<rect x="{x0}" y="{yy-17:.1f}" width="{xx(val)-x0:.2f}" height="34" rx="4" fill="{COLORS[i]}"/>', f'<text class="value" x="{xx(val)+10:.2f}" y="{yy+5:.1f}">{val:.1f}×</text>']
        out += [f'<text class="minor" x="{(x0+x1)/2}" y="{h-18}" text-anchor="middle">Ratio of variability estimates · 1× means equal</text>', '</svg>']
        target = HERE / "nist-tn2006-between-within.svg"
    else:
        w, h = 420, 510
        left, right, top, bottom = 165, 48, 128, 75
        xmax = 5.5
        x0, x1 = left, w - right
        ytop, ybottom = top, h - bottom
        xx = lambda v: x0 + v / xmax * (x1 - x0)
        out = svg_header(w, h, "Between-laboratory variability relative to within-laboratory variability", "Four ratios transcribed from Table 5 of NIST TN 2006. A ratio of one is the reference where between-laboratory and median within-laboratory variability are equal. The values are 2.7 for tensile strength, 4.9 for yield strength, 3.5 for elongation and 0.5 for elastic modulus.")
        out += ['<text class="heading" x="18" y="37">How much did labs differ?</text>', '<text class="sub" x="18" y="65">Between-lab ÷ typical within-lab variability</text>']
        for tick in range(0, 6):
            xt = xx(tick)
            out += [f'<line class="grid" x1="{xt:.2f}" y1="{ytop-13}" x2="{xt:.2f}" y2="{ybottom}"/>', f'<text class="minor" x="{xt:.2f}" y="{ybottom+25}" text-anchor="middle">{tick}×</text>']
        ref = xx(1)
        out += [f'<line x1="{ref:.2f}" y1="{ytop-13}" x2="{ref:.2f}" y2="{ybottom}" stroke="#192c42" stroke-width="2" stroke-dasharray="6 5"/>']
        for i, (name, val) in enumerate(values):
            yy = ytop + (i + 0.5) * (ybottom - ytop) / 4
            short = ["Tensile", "Yield strength", "Elongation", "Elastic modulus"][i]
            out += [f'<text class="minor" x="{left-8}" y="{yy+5:.1f}" text-anchor="end">{short}</text>', f'<rect x="{x0}" y="{yy-18:.1f}" width="{xx(val)-x0:.2f}" height="36" rx="4" fill="{COLORS[i]}"/>', f'<text class="value" x="{xx(val)+7:.2f}" y="{yy+5:.1f}">{val:.1f}×</text>']
        out += [f'<text class="minor" x="{(x0+x1)/2}" y="{h-20}" text-anchor="middle">1× means equal variability</text>', '</svg>']
        target = HERE / "nist-tn2006-between-within-phone.svg"
    target.write_text("\n".join(out), encoding="utf-8")


def main() -> None:
    hashes = check_sources()
    rows, ratios = load_data()
    write_lab_chart(rows)
    write_lab_chart(rows, phone=True)
    write_ratio_chart(ratios)
    write_ratio_chart(ratios, phone=True)
    summary = {}
    for lab in range(1, 6):
        values = [float(r["Sy"]) for r in rows if int(r["Lab"]) == lab]
        summary[str(lab)] = {
            "n": len(values),
            "mean_MPa": sum(values) / len(values),
            "min_MPa": min(values),
            "max_MPa": max(values),
        }
    result = {
        "study": "NIST TN 2006, interlaboratory mechanical-property study of cobalt-chromium alloy",
        "classification": "published study data; educational plot generated by DMGB Systems",
        "primary_report": "https://doi.org/10.6028/NIST.TN.2006",
        "dataset_doi": "https://doi.org/10.18434/T4/1502566",
        "plotted_property": "Sy: 0.2% offset yield strength",
        "plotted_units": "MPa",
        "subset": "Laboratories 1 through 5, the five-laboratory comparison in NIST TN 2006 Table 5",
        "published_rows_plotted": len(rows),
        "rows_by_lab": {k: v["n"] for k, v in summary.items()},
        "computed_summary_MPa": summary,
        "transcribed_table5_ratios": EXPECTED_RATIOS,
        "source_sha256": hashes,
        "checks": {
            "nist_source_hashes_match_publication_record": True,
            "Sy_units_and_definition_match_attribute_file": True,
            "all_34_reported_values_plotted": True,
            "counts_by_lab_match_audited_data": True,
            "table5_ratio_transcription_matches_audited_pdf_text": True,
        },
        "limitations": [
            "The five-laboratory model comparison is a subset of an eight-laboratory study.",
            "These results concern the alloy, process subset, specimens and methods in this study, not additive manufacturing generally.",
            "Between/within ratios are model-based study estimates, not ratios of material strength values.",
            "The plotted ILS CSV contains published test-level property results, not the original force-sensor time series; the separate stress–strain source CSV is provided as a download.",
            "This plot does not qualify a component or rank laboratories.",
        ],
    }
    (HERE / "nist-tn2006-study.json").write_text(json.dumps(result, indent=2) + "\n", encoding="utf-8")
    print(json.dumps({"rows": len(rows), "labs": {k: v["n"] for k, v in summary.items()}, "ratios": EXPECTED_RATIOS, "source_sha256": hashes}, indent=2))


if __name__ == "__main__":
    main()