#!/usr/bin/env python3
"""由 DIAMOND 双向结果构建 6 物种 1:1 直系同源表（RBH 交集）。"""
import json, sys
from pathlib import Path

ANCHOR = "GCF_024868665.1"
OTHERS = ["GCF_015846415.1", "GCA_008369825.1", "GCF_947034865.1",
          "GCF_027580225.2", "GCF_052040795.1"]
MIN_COV, MIN_PID = 50.0, 40.0


def load(path):
    best = {}
    with open(path) as fh:
        for line in fh:
            q, s, pid, ev, bs, qc, sc = line.rstrip("\n").split("\t")
            if float(qc) < MIN_COV or float(sc) < MIN_COV or float(pid) < MIN_PID:
                continue
            if q not in best or float(bs) > best[q][1]:
                best[q] = (s, float(bs))
    return {q: v[0] for q, v in best.items()}


def main(S, out):
    S = Path(S)
    rbh = {}
    stats = {}
    for acc in OTHERS:
        fwd = load(S / "orth" / f"{ANCHOR}_vs_{acc}.tsv")
        rev = load(S / "orth" / f"{acc}_vs_{ANCHOR}.tsv")
        pairs = {q: s for q, s in fwd.items() if rev.get(s) == q}
        rbh[acc] = pairs
        stats[acc] = {"forward_hits": len(fwd), "reverse_hits": len(rev), "rbh_pairs": len(pairs)}
    common = set(rbh[OTHERS[0]])
    for acc in OTHERS[1:]:
        common &= set(rbh[acc])
    sym = json.loads((S / "prep" / f"{ANCHOR}.meta.json").read_text())["protein_to_gene"]
    rows = []
    for q in sorted(common):
        row = {"gene_symbol": sym.get(q, q), ANCHOR: q}
        for acc in OTHERS:
            row[acc] = rbh[acc][q]
        rows.append(row)
    Path(out).write_text(json.dumps(
        {"anchor": ANCHOR, "species": [ANCHOR] + OTHERS, "min_query_coverage_pct": MIN_COV,
         "min_percent_identity": MIN_PID, "per_species_rbh": stats,
         "n_1to1_all_six": len(rows), "orthologs": rows}, indent=1))
    print(json.dumps(stats, indent=1))
    print("1:1 across all six species:", len(rows))


if __name__ == "__main__":
    main(*sys.argv[1:])