"""F4-conf-03: harmonise five effect sizes reported under different definitions.
Reads only data/; writes out/results.json deterministically. No randomness."""
import json, math, os, re

ROOT = os.path.dirname(os.path.dirname(os.path.abspath(__file__)))
D = os.path.join(ROOT, "data")
SRC = os.path.join(D, "sources")


def read_lines(path):
    with open(path, encoding="utf-8") as f:
        return f.read().split("\n")


def find_line(lines, src, pattern):
    """Return (line_no, text) of the unique line matching regex `pattern`; print it."""
    hits = [(i + 1, l) for i, l in enumerate(lines) if re.search(pattern, l)]
    if len(hits) != 1:
        raise RuntimeError(f"{src}: pattern {pattern!r} matched {len(hits)} lines")
    n, l = hits[0]
    print(f"[{src} L{n}] {l}")
    return n, l


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


def num(s):
    """Parse a printed number, stripping stars, %, pp, $, commas."""
    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(D, "effects.json"), encoding="utf-8"))
ids = [e["id"] for e in effects]
assert ids == ["E1", "E2", "E3", "E4", "E5"]

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

harm = {}
baseline = {}

# ---------------- E1: DCR Table 6, ATT, 5 km (log points -> relative %)
s = "arxiv-2609.27732v1"
L = read_lines(os.path.join(SRC, s + ".txt"))
_, hdr = find_line(L, s, r"^\| 1 km \| 3 km \| 5 km$")
assert cells(hdr)[3] == "5 km"
n_att, att = find_line(L, s, r"^ATT \| -0\.399")
b1_str = cells(att)[3]
b1 = num(b1_str)
n_pc, pc = find_line(L, s, r"^% change \| -32\.9%")
printed_pc = cells(pc)[3]
v1 = 100 * (math.exp(b1) - 1)
baseline["E1"] = b1
harm["E1"] = {
    "value": round(v1, 4),
    "unit": "percent change relative to the counterfactual hospitalisation rate",
    "transformation": f"100*(exp(b)-1) with b = {b1_str.rstrip('*')} log points -> 100*(exp(-0.431)-1) = {v1:.4f}% (printed '% change' row: {printed_pc})",
    "harmonisable": True,
    "missing_quantity": None,
    "locator": f"{s} | Table 6 | row 'ATT' (extract line {n_att}) | column '5 km' = {b1_str}; cross-check row '% change' (line {n_pc}) column '5 km' = {printed_pc}",
}

# ---------------- E4: DCR Table 7, fatal, 5 km (log points -> absolute per 100,000)
n_fat, fat = find_line(L, s, r"^Drug-related \(fatal\) \|")
b4_str = cells(fat)[3]
b4 = num(b4_str)
baseline["E4"] = b4
rel4 = 100 * (math.exp(b4) - 1)
# Check whether any baseline level of the fatal rate is printed: search for a mean/baseline row
level_rows = [l for l in L if re.search(r"(?i)(mean|baseline|pre-treatment).*(rate|fatal)|(rate|fatal).*(mean)", l)
              and re.search(r"\d", l) and "|" in l]
print(f"[{s}] table rows printing a baseline/mean level of the fatal rate: {len(level_rows)}")
harm["E4"] = {
    "value": None,
    "unit": "hospitalisations per 100,000 inhabitants per year",
    "transformation": f"absolute change = R0*(exp(b)-1) with b = {b4_str.rstrip('*')}; R0 not printed, so not computable. Only the relative change is reachable: 100*(exp(-0.354)-1) = {rel4:.4f}%",
    "harmonisable": False,
    "missing_quantity": "baseline (counterfactual / pre-opening 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",
    "locator": f"{s} | Table 7 | row 'Drug-related (fatal)' (extract line {n_fat}) | column '5 km' = {b4_str}",
}

# ---------------- E2: RESEA Table 2, Reemployment (pp -> % of control mean)
s = "arxiv-2606.14621v1"
L = read_lines(os.path.join(SRC, s + ".txt"))
n_h2, h2 = find_line(L, s, r"^Outcome \| Treatment \| Control \|")
hc = cells(h2)
n_re, re_ = find_line(L, s, r"^Reemployment \| 75")
rc = cells(re_)
ctrl_str = rc[hc.index("Control")]
ate_cell = rc[hc.index("Average treatment effect estimate")]
ate_str = ate_cell.split(";")[0]
p = num(ctrl_str)
d = num(ate_str)
v2 = 100 * d / p
baseline["E2"] = d
harm["E2"] = {
    "value": round(v2, 4),
    "unit": "percent change relative to the control-group mean reemployment rate",
    "transformation": f"100*Delta/p with Delta = {ate_str} and p = {ctrl_str} -> 100*1.47/75.0 = {v2:.4f}%",
    "harmonisable": True,
    "missing_quantity": None,
    "locator": f"{s} | Table 2 | row 'Reemployment' (extract line {n_re}) | column 'Average treatment effect estimate' = {ate_str}; column 'Control' = {ctrl_str}",
}

# ---------------- E3: MW Table 2 elasticity x printed MW change
s = "arxiv-2505.04555v3"
L = read_lines(os.path.join(SRC, s + ".txt"))
n_el, el = find_line(L, s, r"^Employment elasticity w\.r\.t MW \|")
eta_str = cells(el)[1]
n_g, gl = find_line(L, s, r"^% MW changes \|")
g_str = cells(gl)[1]
find_line(L, s, r"about 4% minimum wage increase")
eta, g = num(eta_str), num(g_str)
v3 = 100 * eta * g
baseline["E3"] = eta
harm["E3"] = {
    "value": round(v3, 4),
    "unit": "percent change in employment",
    "transformation": f"100*eta*g with eta = {eta_str} and g = {g_str} (proportional MW change as printed) -> 100*(-0.388)*0.047 = {v3:.4f}% (first-order approximation)",
    "harmonisable": True,
    "missing_quantity": None,
    "locator": f"{s} | Table 2 | row 'Employment elasticity w.r.t MW' (extract line {n_el}) | column 'Estimate (SE)' = {eta_str}; row '% MW changes' (line {n_g}) | column 'Estimate (SE)' = {g_str}",
}

# ---------------- E5: Benin Table 3 col (4) LPM coefficient -> pp
s = "arxiv-2510.03658v1"
L = read_lines(os.path.join(SRC, s + ".txt"))
_, h5 = find_line(L, s, r"^\| \(1\) \| \(2\) \| \(3\) \| \(4\)$")
col = cells(h5).index("(4)")
n_t, tl = find_line(L, s, r"^Treatment \| -0\.0911")
b5_str = cells(tl)[col]
n_m, ml = find_line(L, s, r"^Mean Dep\. Var\. \|")
n_c, cl = find_line(L, s, r"^% Change \|")
b5 = num(b5_str)
v5 = 100 * b5
baseline["E5"] = b5
harm["E5"] = {
    "value": round(v5, 4),
    "unit": "percentage points",
    "transformation": f"100*beta with beta = {b5_str} (0-1 probability scale) -> {v5:.2f} pp (exact rescaling)",
    "harmonisable": True,
    "missing_quantity": None,
    "locator": f"{s} | Table 3 | row 'Treatment' (extract line {n_t}) | column '(4)' = {b5_str}",
}

# ---------------- internal cross-checks used to document data issues
chk_e1 = round(v1, 1) == num(printed_pc)
mean5 = num(cells(ml)[col]); pc5 = num(cells(cl)[col])
chk_e5 = b5 / mean5
print(f"check E1: computed {v1:.2f}% vs printed {printed_pc}: match={chk_e1}")
print(f"check E4 relative: {rel4:.2f}% (abstract says 35%, Sec 6.1 says ~30%)")
print(f"check E5 relative: {b5}/{mean5} = {chk_e5:.4f} vs printed % Change {pc5}")

data_issues = [
    {"issue": "arxiv-2609.27732v1: abstract states fatal drug-related hospitalisations 'fall by 35%%', but Table 7 5 km fatal ATT -0.354 implies 100*(exp(-0.354)-1) = %.1f%%, consistent with Section 6.1's 'approximately 30%%'." % rel4,
     "handling": "Internal inconsistency of the source; E4 is based on the Table 7 cell only; the abstract figure is not used."},
    {"issue": "arxiv-2609.27732v1: no baseline/counterfactual level of the fatal hospitalisation rate (per 100,000 per year) is printed in the extract.",
     "handling": "E4 declared not harmonisable (value null) with missing_quantity named; no external value substituted."},
    {"issue": "arxiv-2609.27732v1: Table A19 CS-DR column prints SE (0.218) for the same -0.431 estimate that Table 6 (5 km) prints with SE (0.210); A19 note says 640 never-treated regions while the tables print 608.",
     "handling": "Recorded only; SEs are not used (no uncertainty propagation) and E1 uses Table 6 as specified."},
    {"issue": "arxiv-2609.27732v1: Table A19 '% change' row: -0.362 -> 100*(exp(b)-1) = -30.37% printed as -30.3%; -0.149 -> -13.84% printed as -13.9% (small rounding discrepancies, possibly from unrounded coefficients).",
     "handling": "Recorded only; these cells are not inputs to any harmonised item. For E1 the computed value matches the printed -35.0%."},
    {"issue": "arxiv-2505.04555v3: Table 2 note describes 'about 4% minimum wage increase' while the table row '% MW changes' prints 0.047 (4.7%).",
     "handling": "Target definition requires the size printed in the table; g = 0.047 used. With 4% instead, the implied change would be smaller in magnitude (not reported as a result)."},
    {"issue": "arxiv-2510.03658v1: Section 4.1 text quotes a 'table 2' estimate (-0.0855, mean 0.455) not included in the extract, mislabels relative changes as 'percentage change', and interprets the relative change as lives saved per 100 deaths; Table 3 col (4) printed %% Change -0.3685 vs -0.0962/0.261 = %.4f (rounding). Footnote markers render as '^9^9 9'." % chk_e5,
     "handling": "Only Table 3 row 'Treatment', column (4) is used for E5; the target (percentage points) needs no mean, so these issues do not affect E5."},
    {"issue": "arxiv-2606.14621v1: the 'Control' column (75.0%) is a raw group mean, while the ATE (1.47 pp) adjusts for randomisation week; abstract rounds the ATE to 1.5 pp; 'p ¡ 0.001' rendering defect in text.",
     "handling": "Control-group mean used as the denominator exactly as the target definition specifies; printed table value 1.47 pp used, not the rounded abstract value."},
]

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