"""F4-conf-03: harmonise five printed effect sizes onto target definitions.
Reads only data/; writes out/results.json deterministically. No randomness.
Run from workspace root: ~/wp29/venv/bin/python out/analysis.py
"""
import json, math, os, re

DATA = "data"
SRC = os.path.join(DATA, "sources")


def lines_of(fname):
    with open(os.path.join(SRC, fname), encoding="utf-8") as f:
        return f.read().split("\n")


def find_line(fname, startswith, after=None):
    """Return (1-based line number, text) of first line starting with `startswith`
    (optionally after a line containing `after`)."""
    L = lines_of(fname)
    start = 0
    if after is not None:
        start = next(i for i, l in enumerate(L) if after in l)
    for i in range(start, len(L)):
        if L[i].startswith(startswith):
            print(f"[{fname}:{i+1}] {L[i]}")
            return i + 1, L[i]
    raise ValueError(f"{startswith!r} not found in {fname}")


def cells(line):
    return [c.strip() for c in line.split("|")]


def num(s):
    """Parse a printed number, dropping stars, %, 'pp' and thousands separators."""
    m = re.search(r"-?\d[\d,]*\.?\d*", s.replace("−", "-"))
    return float(m.group(0).replace(",", ""))


def wc_l(path):
    with open(path, "rb") as f:
        return f.read().count(b"\n")


effects = json.load(open(os.path.join(DATA, "effects.json"), encoding="utf-8"))
assert [e["id"] for e in effects] == ["E1", "E2", "E3", "E4", "E5"]
H, baseline = {}, {}

# ---------- E1: Armendariz Table 6, ATT, 5 km (log points -> % change) ----------
f1 = "arxiv-2609.27732v1.txt"
_, hdr = find_line(f1, "| 1 km | 3 km | 5 km", after="Table 6:")
ln, row = find_line(f1, "ATT |", after="Table 6:")
c = cells(row); col = cells(hdr).index("5 km")
b1_txt = c[col]; b1 = num(b1_txt)
lnp, rowp = find_line(f1, "% change |", after="Table 6:")
pc_printed = cells(rowp)[col]
v1 = 100.0 * (math.exp(b1) - 1.0)
H["E1"] = dict(value=round(v1, 4), unit="percent (relative change vs counterfactual rate)",
               transformation=f"100*(exp(b)-1) with b = {b1_txt} log points -> 100*(exp({b1})-1) = {v1:.4f}%",
               harmonisable=True, missing_quantity=None,
               locator=f"arxiv-2609.27732v1 | Table 6 | row 'ATT' (line {ln}) | column '5 km': {b1_txt}; "
                       f"cross-check row '% change' (line {lnp}) column '5 km': {pc_printed}")
baseline["E1"] = {"printed_value": b1_txt, "printed_unit": "log points", "author_printed_conversion": pc_printed}

# ---------- E2: Li et al. Table 2, Reemployment ATE (pp -> % of control mean) ----------
f2 = "arxiv-2606.14621v1.txt"
_, hdr2 = find_line(f2, "Outcome | Treatment | Control")
ln2, row2 = find_line(f2, "Reemployment |", after="Outcome | Treatment | Control")
h2, c2 = cells(hdr2), cells(row2)
ctrl_txt = c2[h2.index("Control")]
ate_cell = c2[h2.index("Average treatment effect estimate")]
ate_txt = ate_cell.split(";")[0]
p_ctrl, d2 = num(ctrl_txt), num(ate_txt)
v2 = 100.0 * d2 / p_ctrl
H["E2"] = dict(value=round(v2, 4), unit="percent (relative change vs control-group mean)",
               transformation=f"100*Delta/p with Delta = {ate_txt}, p = {ctrl_txt} -> 100*{d2}/{p_ctrl} = {v2:.4f}%",
               harmonisable=True, missing_quantity=None,
               locator=f"arxiv-2606.14621v1 | Table 2 | row 'Reemployment' (line {ln2}) | column 'Average treatment effect estimate': {ate_txt}; "
                       f"column 'Control': {ctrl_txt}")
baseline["E2"] = {"printed_value": ate_txt, "printed_unit": "percentage points", "author_printed_conversion": None}

# ---------- E3: Kanayama et al. Table 2, elasticity x printed MW change ----------
f3 = "arxiv-2505.04555v3.txt"
ln3, row3 = find_line(f3, "Employment elasticity w.r.t MW |")
ln3g, row3g = find_line(f3, "% MW changes |")
find_line(f3, "Note: The table reports the effects of a minimum wage increase")
eta_txt, g_txt = cells(row3)[1], cells(row3g)[1]
eta, g = num(eta_txt), num(g_txt)
v3 = 100.0 * eta * g
H["E3"] = dict(value=round(v3, 4), unit="percent (relative change in employment)",
               transformation=f"100*eta*g with eta = {eta_txt}, g = {g_txt} (proportional MW change) -> 100*({eta})*({g}) = {v3:.4f}%",
               harmonisable=True, missing_quantity=None,
               locator=f"arxiv-2505.04555v3 | Table 2 | row 'Employment elasticity w.r.t MW' (line {ln3}) | column 'Estimate (SE)': {eta_txt}; "
                       f"row '% MW changes' (line {ln3g}) | column 'Estimate (SE)': {g_txt}")
baseline["E3"] = {"printed_value": eta_txt, "printed_unit": "elasticity", "author_printed_conversion": None}

# ---------- E4: Armendariz Table 7, fatal, 5 km (log points -> absolute per 100,000) ----------
_, hdr7 = find_line(f1, "Outcome | 1 km | 3 km | 5 km", after="Table 7:")
ln4, row4 = find_line(f1, "Drug-related (fatal) |", after="Table 7:")
b4_txt = cells(row4)[cells(hdr7).index("5 km")]; b4 = num(b4_txt)
rel4 = 100.0 * (math.exp(b4) - 1.0)
# Search the extract for any printed level of the fatal hospitalisation rate (mean / baseline).
txt1 = "\n".join(lines_of(f1))
level_hits = re.findall(r"(?i)(mean|baseline level|pre-treatment (mean|level)|average rate)[^\n]{0,80}per 100,000", txt1)
print(f"[E4] search for a printed level of the fatal rate per 100,000: {len(level_hits)} hit(s)")
H["E4"] = dict(value=None, unit="hospitalisations per 100,000 inhabitants per year",
               transformation=(f"Would require: Delta = Y0*(exp(b)-1) with b = {b4_txt} log points and Y0 the counterfactual "
                               f"fatal rate per 100,000 per year. Only the relative change 100*(exp({b4})-1) = {rel4:.4f}% "
                               f"is computable from the extract; Y0 is not printed."),
               harmonisable=False,
               missing_quantity=("the counterfactual (or pre-treatment / control-group mean) level of the fatal drug-related "
                                 "hospitalisation rate per 100,000 inhabitants per year for the 5 km treated regions "
                                 "(not printed in the extract; Tables 6, 7, A19 print only log-point ATTs, SEs and counts)"),
               locator=f"arxiv-2609.27732v1 | Table 7 | row 'Drug-related (fatal)' (line {ln4}) | column '5 km': {b4_txt}")
baseline["E4"] = {"printed_value": b4_txt, "printed_unit": "log points", "author_printed_conversion": None,
                  "computable_relative_change_percent": round(rel4, 4)}

# ---------- E5: Tossou Table 3, Treatment, column (4) (0-1 -> pp) ----------
f5 = "arxiv-2510.03658v1.txt"
_, hdr5 = find_line(f5, "| (1) | (2) | (3) | (4)")
ln5, row5 = find_line(f5, "Treatment |", after="Table 3:")
b5_txt = cells(row5)[cells(hdr5).index("(4)")]; b5 = num(b5_txt)
v5 = 100.0 * b5
H["E5"] = dict(value=round(v5, 4), unit="percentage points",
               transformation=f"100*beta with beta = {b5_txt} (0-1 probability scale) -> 100*({b5}) = {v5:.4f} pp",
               harmonisable=True, missing_quantity=None,
               locator=f"arxiv-2510.03658v1 | Table 3 | row 'Treatment' (line {ln5}) | column '(4)': {b5_txt}")
baseline["E5"] = {"printed_value": b5_txt, "printed_unit": "probability (0-1)", "author_printed_conversion": None}

# ---------- consistency evidence for data_issues ----------
find_line(f1, "This paper estimates the causal effect")  # abstract (fatal 'fall by 35%')
find_line(f1, "The sub-outcome results point to two channels")  # Sec 6.1 (~30%)
_, a19 = find_line(f1, "| (0.218)", after="Table A19: Robustness")
find_line(f1, "never-treated comparison regions (640 regions)")
_, a19pc = find_line(f1, "% change |", after="Table A19: Robustness")
_, a19b = find_line(f1, "ATT |", after="Table A19: Robustness")
a19_check = {h: (round(100 * (math.exp(num(b)) - 1), 1), p) for h, b, p in
             zip(["CS-DR", "Sun & Abraham", "Synthetic DiD", "Stacked DiD"], cells(a19b)[1:], cells(a19pc)[1:])}
print("[A19 check] computed vs printed % change:", a19_check)
_, t5pc = find_line(f5, "% Change |"); _, t5m = find_line(f5, "Mean Dep. Var. |")
t5_check = round(b5 / num(cells(t5m)[4]), 4)
print(f"[Tossou check] -0.0962/0.261 = {t5_check} vs printed {cells(t5pc)[4]}")
find_line(f2, "The results for the short-term outcomes")  # 'p ¡ 0.001' rendering

data_issues = [
    {"issue": f"arxiv-2609.27732v1: abstract says fatal drug-related hospitalisations 'fall by 35%', but Table 7 (5 km) prints {b4_txt} "
              f"log points, i.e. {rel4:.1f}%; Section 6.1 says 'approximately 30%'.",
     "handling": "Table 7 value used for E4; abstract figure not used."},
    {"issue": "arxiv-2609.27732v1: Table A19 prints SE (0.218) for the CS-DR 5 km ATT -0.431, whereas Table 6 prints (0.210) for the same estimate; "
              "A19 note says 640 never-treated regions while the table body and Tables 6-7 say 608.",
     "handling": "Not used: no standard errors are propagated and region counts do not enter any transformation. Reported only."},
    {"issue": f"arxiv-2609.27732v1: Table A19 '% change' for Synthetic DiD and Stacked DiD differ by 0.1 from 100*(exp(b)-1) on the printed coefficients "
              f"(computed vs printed: {a19_check}).",
     "handling": "Not used (E1 uses Table 6, where printed -35.0% matches the computed value). Attributed to rounding of coefficients."},
    {"issue": "arxiv-2609.27732v1: no level (mean, baseline or counterfactual) of the fatal drug-related hospitalisation rate per 100,000 is printed, "
              "and the time unit of the panel ('per year') is not stated explicitly.",
     "handling": "E4 declared not harmonisable (value null); missing quantity named; no external value substituted."},
    {"issue": f"arxiv-2505.04555v3: Table 2 note describes 'about 4% minimum wage increase' while the row '% MW changes' prints {g_txt} (4.7%).",
     "handling": f"Task instructs use of the size printed in the table: g = {g_txt}. With 0.04 the result would differ; noted in the report."},
    {"issue": f"arxiv-2510.03658v1: abstract/text quote an effect of 0.0855 relative to a mean of 0.455 from 'table 2' (not in the extract); "
              f"Table 3 col (4) '% Change' prints {cells(t5pc)[4]} vs {t5_check} recomputed from printed coefficient and mean (rounding); footnote 11 text is garbled.",
     "handling": "Only Table 3, row Treatment, column (4) is used for E5, as specified in effects.json."},
    {"issue": "arxiv-2606.14621v1: rendering defect 'p ¡ 0.001' in the short-term results paragraph; the control mean (75.0%) is the raw aggregate "
              "control-group mean, whereas the ATE (1.47 pp) is block (randomisation-week) adjusted.",
     "handling": "Rendering defect irrelevant to E2. Control mean used as the counterfactual level p as the target definition requires; caveat stated."},
]

n_input_rows = {"effects.json": len(effects), "README.md": wc_l(os.path.join(DATA, "README.md"))}
for fn in sorted(os.listdir(SRC)):
    n_input_rows[fn] = wc_l(os.path.join(SRC, fn))

n_h = sum(1 for v in H.values() if v["value"] is not None)
n_nh = sum(1 for v in H.values() if v["harmonisable"] is False)
results = {
    "task_id": "F4-conf-03",
    "harmonised": H,
    "n_harmonised": n_h,
    "n_not_harmonisable": n_nh,
    "feasibility": "complete",
    "feasibility_reason": None,
    "data_issues": data_issues,
    "n_input_rows": n_input_rows,
    "n_analysed": len(effects),
    "random_seed": None,
    "baseline_reported": baseline,
}
os.makedirs("out", exist_ok=True)
with open(os.path.join("out", "results.json"), "w", encoding="utf-8") as f:
    json.dump(results, f, indent=2, ensure_ascii=False)
print("\nHarmonised values:")
for k, v in H.items():
    print(f"  {k}: value={v['value']} [{v['unit']}] harmonisable={v['harmonisable']}")
print(f"n_harmonised={n_h} n_not_harmonisable={n_nh} feasibility=complete")