"""Plot six saved Case2 results; no statistical analysis is recomputed.
Usage:
  python case2-clock-context-forest-20261012.py --csv clock-context-enrichment.csv --report clock-context-report.md --prefix case2-clock-context-forest-rerun
Dependencies: numpy, pandas, matplotlib, Pillow.
Existing output files are never overwritten. Inputs are verified by SHA-256.
"""
from pathlib import Path
import argparse, hashlib, io, re
import numpy as np
import pandas as pd
import matplotlib as mpl
mpl.use("Agg")
import matplotlib.pyplot as plt
from matplotlib.ticker import FixedLocator, FuncFormatter, NullLocator
from PIL import Image
import xml.etree.ElementTree as ET

CSV_SHA256 = "15ad30248785aff157e86c5eb35c360802e2a9260e0a4193beba943a99dfe7c1"
REPORT_SHA256 = "9428e1e32dc8e40d4d5c3fe31f8b26a8b7f77e0414eb42c8203b3fa9d933923e"
FEATURES = ["Island", "Shore", "Promoter_TSS200_or_TSS1500"]
BACKGROUNDS = ["450K_all_CpGs", "27K_shared_CpGs"]
FEATURE_LABELS = ["CpG island", "CpG shore", "Promoter\n(TSS200 or TSS1500)"]
REPORT_LABELS = ["CpG island", "CpG shore", "Promoter: TSS200 or TSS1500"]

def make_figure(csv_path, report_path, prefix):
    csv_path, report_path, prefix = map(Path, (csv_path, report_path, prefix))
    outputs = [Path(str(prefix)+ext) for ext in [".png", ".svg", "-source-check.csv"]]
    if any(p.exists() for p in outputs):
        raise FileExistsError("Refusing to overwrite an existing output; choose another prefix.")
    csv_bytes, report_bytes = csv_path.read_bytes(), report_path.read_bytes()
    assert hashlib.sha256(csv_bytes).hexdigest() == CSV_SHA256, "CSV snapshot mismatch"
    assert hashlib.sha256(report_bytes).hexdigest() == REPORT_SHA256, "Report snapshot mismatch"
    report_text = report_bytes.decode("utf-8")
    raw = pd.read_csv(io.BytesIO(csv_bytes), dtype=str, keep_default_na=False)
    raw["source_csv_line"] = np.arange(len(raw))+2
    selected = raw[raw.feature.isin(FEATURES) & raw.background.isin(BACKGROUNDS)].copy()
    assert len(selected)==6 and not selected.duplicated(["feature","background"]).any()
    expected = {(f,b) for f in FEATURES for b in BACKGROUNDS}
    assert set(zip(selected.feature,selected.background)) == expected
    assert "conditional maximum-likelihood odds ratios" in report_text
    assert "exact conditional 95% confidence intervals" in report_text
    assert "proxy, not the original 21,369-site training universe" in report_text
    assert "TSS200 **or** TSS1500" in report_text

    # The CSV is the only numeric plotting input. The report is used only for checks.
    indexed = selected.set_index(["feature","background"])
    records = []
    for fi, f in enumerate(FEATURES):
        lines=[l for l in report_text.splitlines() if l.startswith("| "+REPORT_LABELS[fi]+" |")]
        assert len(lines)==1
        cells=[s.strip() for s in lines[0].strip("|").split("|")]
        for bi,b in enumerate(BACKGROUNDS):
            row=indexed.loc[(f,b)]
            numbers=np.array([float(row[k]) for k in ["conditional_odds_ratio","ci95_low","ci95_high"]])
            estimate, lower, upper=numbers
            assert np.isfinite(numbers).all() and (numbers>0).all(), "Cannot faithfully show interval on log axis"
            assert lower <= estimate <= upper
            report_cell=cells[3 if bi==0 else 5]
            reported=re.fullmatch(r"([0-9.]+) \(([0-9.]+)[–-]([0-9.]+)\)",report_cell)
            assert reported and [f"{v:.2f}" for v in numbers] == list(reported.groups())
            assert row.clock_n == "353"
            records.append({
                "feature":f, "feature_label":FEATURE_LABELS[fi].replace("\n"," "),
                "background":b,
                "background_label":"Full 450K background" if bi==0 else "27K-shared proxy background",
                "source_csv_line":row.source_csv_line,
                "conditional_odds_ratio":row.conditional_odds_ratio,
                "ci95_low":row.ci95_low, "ci95_high":row.ci95_high,
                "clock_n":row.clock_n, "nonclock_n":row.nonclock_n,
                "report_or_ci_text":report_cell, "report_matches_csv_at_display_precision":True,
                "positive_finite_log_interval":True,
                "csv_sha256":CSV_SHA256, "report_sha256":REPORT_SHA256})
    checked=pd.DataFrame(records)
    # Explicit portable styling matches the Notebook figure-style size ladder.
    mpl.rcParams.update({
        "font.family":"DejaVu Sans", "font.size":11, "axes.labelsize":11,
        "axes.titlesize":11, "xtick.labelsize":9, "ytick.labelsize":11,
        "legend.fontsize":10, "axes.linewidth":0.8, "svg.fonttype":"none",
        "svg.hashsalt":"case2-clock-context", "figure.facecolor":"white",
        "axes.facecolor":"white", "text.color":"#222222", "axes.labelcolor":"#222222",
        "xtick.color":"#222222", "ytick.color":"#222222",
        "savefig.facecolor":"white", "axes.grid":False, "savefig.bbox":None})
    fig=plt.figure(figsize=(7.2,4.5), dpi=300)
    ax=fig.add_axes([0.32,0.32,0.64,0.42])
    fig.text(0.03,0.97,"Horvath clock: enrichment by comparison background",
             ha="left",va="top",fontsize=11)
    colors=["#0072B2","#D55E00"]
    markers=["o","s"]
    handles=[]
    checks=[]
    for bi,b in enumerate(BACKGROUNDS):
        group=checked[checked.background.eq(b)].set_index("feature").loc[FEATURES]
        center=group.conditional_odds_ratio.astype(float).to_numpy()
        lo=group.ci95_low.astype(float).to_numpy()
        hi=group.ci95_high.astype(float).to_numpy()
        y=np.array([2.,1.,0.])+(0.14 if bi==0 else -0.14)
        nc=int(group.nonclock_n.iloc[0])
        assert group.nonclock_n.nunique()==1
        label=("Full 450K background" if bi==0 else "27K-shared proxy background")+f" ({nc:,} non-clock CpGs)"
        h=ax.errorbar(center,y,xerr=np.vstack([center-lo,hi-center]),fmt=markers[bi],
                      markersize=5.5,color=colors[bi],ecolor=colors[bi],
                      elinewidth=1.4,capsize=3,capthick=1.2,
                      markeredgewidth=0.8,linestyle="none",label=label,zorder=3)
        handles.append(h)
        h.lines[0].set_gid("estimates-"+b)
        for collection in h.lines[2]: collection.set_gid("intervals-"+b)
        actual_centers=np.asarray(h.lines[0].get_xdata(),dtype=float)
        segments=h.lines[2][0].get_segments()
        for i,f in enumerate(FEATURES):
            endpoints=np.asarray(segments[i])[:,0]
            assert actual_centers[i] == center[i]
            np.testing.assert_allclose(endpoints,[lo[i],hi[i]],rtol=1e-14,atol=0)
            checks.append((f,b,True))
    ax.set_xscale("log")
    bounds=checked[["ci95_low","ci95_high"]].astype(float).to_numpy()
    low,high=bounds.min(),bounds.max()
    pad=np.exp((np.log(high)-np.log(low))*0.08)
    ax.set_xlim(low/pad,high*pad)
    ax.set_ylim(-0.42,2.42)
    ticks=[v for v in [0.5,0.75,1,1.5,2,3,4,5,6] if low/pad<=v<=high*pad]
    ax.xaxis.set_major_locator(FixedLocator(ticks))
    ax.xaxis.set_major_formatter(FuncFormatter(lambda v,p:f"{v:g}"))
    ax.xaxis.set_minor_locator(NullLocator())
    ax.set_yticks([2,1,0],FEATURE_LABELS)
    ax.tick_params(axis="y",length=0,pad=12)
    ax.tick_params(axis="x",length=4,pad=5)
    for side in ["top","right","left"]:ax.spines[side].set_visible(False)
    ax.spines["bottom"].set_color("#777777")
    ax.axvline(1,color="#777777",linestyle=(0,(4,3)),linewidth=1,zorder=1,gid="reference-OR-1")
    ax.set_xlabel("Conditional odds ratio (log scale)",labelpad=10)
    ax.grid(axis="x",color="#E8E8E8",linewidth=0.6,zorder=0)
    fig.legend(handles=handles,loc="upper left",bbox_to_anchor=(0.025,0.895),
               frameon=False,borderaxespad=0,handlelength=2,handletextpad=0.6,labelspacing=0.65)
    assert checked.clock_n.nunique()==1
    clock_count=int(checked.clock_n.iloc[0])
    fig.text(0.03,0.135,
             f"{clock_count:,} clock CpGs; points: conditional odds ratios; whiskers: matching 95% confidence intervals.",
             ha="left",va="top",fontsize=9)
    fig.text(0.03,0.085,
             "27K-shared is a proxy, not the exact training-eligible universe.",
             ha="left",va="top",fontsize=9)
    fig.text(0.03,0.035,
             "OR > 1 indicates enrichment; dashed line: OR = 1. Intervals are not multiplicity-adjusted.",
             ha="left",va="top",fontsize=9)
    fig.canvas.draw()
    renderer=fig.canvas.get_renderer()
    texts=[(t,t.get_window_extent(renderer)) for t in fig.findobj(mpl.text.Text)
           if t.get_visible() and t.get_text().strip()]
    overlaps=[(a.get_text(),b.get_text()) for i,(a,ba) in enumerate(texts)
              for b,bb in texts[i+1:] if ba.overlaps(bb)]
    outside=[t.get_text() for t,box in texts if not (
        fig.bbox.x0<=box.x0 and box.x1<=fig.bbox.x1 and
        fig.bbox.y0<=box.y0 and box.y1<=fig.bbox.y1)]
    assert not overlaps, f"Overlapping text: {overlaps}"
    assert not outside, f"Clipped text: {outside}"
    assert ax.get_xscale()=="log"
    assert all(ax.get_xlim()[0]<float(r.ci95_low)<float(r.ci95_high)<ax.get_xlim()[1]
               for r in checked.itertuples())
    checked["plotted_estimate_and_endpoints_match_csv"]=True
    checked["labels_match_report"]=True
    checked["all_labels_inside_canvas"]=True
    checked["no_text_overlaps"]=True
    # Exclusive-create mode protects every existing output.
    for path,fmt in zip(outputs[:2],["png","svg"]):
        with path.open("xb") as out:
            fig.savefig(out,format=fmt,dpi=300)
    with Image.open(outputs[0]) as im:
        assert im.size == (2160,1350)
        assert abs(im.info["dpi"][0]-300)<0.1
    # Inspect only relevant SVG label text and the six interval/estimate group IDs.
    svg=ET.parse(outputs[1]).getroot()
    svg_text=" ".join("".join(e.itertext()) for e in svg.iter()
                      if e.tag.endswith("}text"))
    for term in ["CpG island","CpG shore","Promoter","TSS200 or TSS1500",
                 "27K-shared proxy background","not the exact training-eligible universe",
                 "Conditional odds ratio (log scale)"]:
        assert term in svg_text, term
    ids={e.attrib.get("id") for e in svg.iter()}
    assert "reference-OR-1" in ids
    assert all("intervals-"+b in ids and "estimates-"+b in ids for b in BACKGROUNDS)
    checked["png_and_svg_structure_validated"]=True
    with outputs[2].open("x",encoding="utf-8",newline="") as out:
        checked.to_csv(out,index=False)
    saved=pd.read_csv(outputs[2],dtype=str,keep_default_na=False)
    assert len(saved)==6
    for r in saved.itertuples():
        original=indexed.loc[(r.feature,r.background)]
        for field in ["conditional_odds_ratio","ci95_low","ci95_high"]:
            assert getattr(r,field)==original[field], (r.feature,r.background,field)
    assert hashlib.sha256(csv_path.read_bytes()).hexdigest()==CSV_SHA256
    assert hashlib.sha256(report_path.read_bytes()).hexdigest()==REPORT_SHA256
    return fig,checked,{"rows":6,"all_requested_intervals_representable":True,
                       "labels_and_values_checked":True,"inputs_unchanged":True,
                       "png_pixels":[2160,1350],"png_dpi":300,
                       "matplotlib_version":mpl.__version__,"outputs":[str(p) for p in outputs]}

if __name__=="__main__":
    ap=argparse.ArgumentParser(description=__doc__)
    ap.add_argument("--csv",type=Path,default=Path("clock-context-enrichment.csv"))
    ap.add_argument("--report",type=Path,default=Path("clock-context-report.md"))
    ap.add_argument("--prefix",type=Path,default=Path("case2-clock-context-forest-rerun"))
    args=ap.parse_args()
    figure,table,validation=make_figure(args.csv,args.report,args.prefix)
    print(validation)
    plt.close(figure)
