P1: file parsers, sample sheet, discovery, CLI inspect

Readers for IDAT v3, BPM (locus entries v4-v8), CSV manifest, EGT
(data blocks 5-9) and GTC v3-v5 (plus a layout-preserving GTC writer);
sample-sheet validation (IN-001..IN-014), IDAT discovery and pairing
(IN-020..IN-022), message catalogue and exit codes (CLI-014/016),
JSON Lines logging (CLI-007/008), atomic writes (CLI-011), run manifest,
and `simurg inspect`.

Validation (docs/validation/P1.md): GTC, BPM and EGT readers agree
exactly with Illumina BeadArrayFiles on GSA v4 / GDA v1 files; GTCs
round-trip byte-identically; BPM and CSV manifests are consistent.

Co-Authored-By: Claude Opus 5.5 <[email protected]>
Claude-Session: https://claude.ai/code/session_01814GKVzZVL3FuF3DeiZXhv

tests/oracle/test_p1_real_files.py

@@ -0,0 +1,182 @@
+"""Phase P1 acceptance on real files (plan §18 P1): parsers agree field-by-field with Illumina's BeadArrayFiles
+(BSD-2) and GTCs round-trip byte-identically. Needs ``scripts/fetch_data.sh p0``; GTC oracles are produced with
+gtc2vcf ``idat2gtc`` (see docs/validation/P1.md).
+"""
+
+from __future__ import annotations
+
+import numpy as np
+import pytest
+from tests.conftest import (
+    BIOC_IDATS,
+    GDA_EGT,
+    GDA_SHEET,
+    GSA_BPM,
+    GSA_CSV,
+    GSA_EGT,
+    GSA_IDATS,
+    ORACLE_GTCS,
+    need,
+)
+
+from simurg.core.errors import FormatError
+from simurg.formats.bpm import read_bpm
+from simurg.formats.csv_manifest import compare_manifests, read_csv_manifest
+from simurg.formats.discovery import discover_idats
+from simurg.formats.egt import read_egt
+from simurg.formats.gtc import encode_gtc, read_gtc
+from simurg.formats.idat import check_pair, read_idat
+from simurg.formats.samplesheet import parse_samplesheet
+
+pytestmark = pytest.mark.data
+
+
+@pytest.fixture(scope="module")
+def gsa_bpm():  # type: ignore[no-untyped-def]
+    need(GSA_BPM)
+    return read_bpm(GSA_BPM)
+
+
+def test_gsa_idats_all_parse_and_pair() -> None:
+    need(GSA_IDATS)
+    d = discover_idats(GSA_IDATS)
+    assert len(d.pairs) == 36 and all(p.complete for p in d.pairs.values())
+    ref = None
+    for pair in d.pairs.values():
+        g, r = read_idat(pair.grn), read_idat(pair.red)  # type: ignore[arg-type]
+        check_pair(g, r)
+        assert g.barcode == pair.barcode and g.position == pair.position
+        assert g.scan_time_iso  # every demo IDAT has a scan entry
+        ids = np.sort(g.address_ids)
+        if ref is None:
+            ref = ids
+        assert np.array_equal(ids, ref)
+
+
+def test_bioc_idats() -> None:
+    need(BIOC_IDATS)
+    ok = {
+        "4019585376_B_Red.idat": 381079,
+        "4019585376_B_Red.idat.gz": 381079,
+        "5723646052_R02C02_Grn.idat": 622399,
+        "5723646052_R02C02_Grn.idat.gz": 622399,
+    }
+    for name, n in ok.items():
+        assert read_idat(BIOC_IDATS / name).n_addresses == n
+    for name in ("4343238080_A_Grn.idat", "6016741005_A_Grn.idat"):
+        with pytest.raises(FormatError) as e:
+            read_idat(BIOC_IDATS / name)
+        assert e.value.code == "E_IDAT_ENCRYPTED"
+
+
+def test_idat_addresses_cover_manifest(gsa_bpm) -> None:  # type: ignore[no-untyped-def]  # precursor of IN-027
+    need(GSA_IDATS)
+    pair = next(iter(discover_idats(GSA_IDATS).pairs.values()))
+    ids = read_idat(pair.grn).address_ids  # type: ignore[arg-type]
+    wanted = np.concatenate([gsa_bpm.address_a, gsa_bpm.address_b[gsa_bpm.address_b != 0]])
+    assert np.isin(wanted, ids).mean() == 1.0
+
+
+def test_gtc_matches_beadarrayfiles_and_roundtrips(beadarrayfiles) -> None:  # type: ignore[no-untyped-def]
+    need(ORACLE_GTCS)
+    files = sorted(ORACLE_GTCS.glob("*.gtc"))
+    assert files
+    GenotypeCalls = beadarrayfiles.GenotypeCalls
+    for f in files:
+        g = read_gtc(f)
+        assert encode_gtc(g) == f.read_bytes(), f"round-trip differs: {f}"
+        b = GenotypeCalls(str(f))
+        assert np.array_equal(g.raw_x, b.get_raw_x_intensities())
+        assert np.array_equal(g.raw_y, b.get_raw_y_intensities())
+        assert g.genotypes.tolist() == list(b.get_genotypes())
+        assert np.array_equal(g.genotype_scores, b.get_genotype_scores(), equal_nan=True)
+        assert np.array_equal(g.b_allele_freqs, b.get_ballele_freqs(), equal_nan=True)
+        assert np.array_equal(g.logr_ratios, b.get_logr_ratios(), equal_nan=True)
+        assert g.base_call_strings() == [x.decode() for x in b.get_base_calls()]
+        assert (g.call_rate, g.logr_dev, g.gc10, g.gc50) == (
+            b.get_call_rate(),
+            b.get_logr_dev(),
+            b.get_gc10(),
+            b.get_gc50(),
+        )
+        assert g.gender == b.get_gender()
+        assert (g.num_calls, g.num_no_calls, g.num_intensity_only) == (
+            b.get_num_calls(),
+            b.get_num_no_calls(),
+            b.get_num_intensity_only(),
+        )
+        assert g.percentiles_x.tolist() == list(b.get_percentiles_x())
+        xf = [
+            (t.offset_x, t.offset_y, t.scale_x, t.scale_y, t.shear, t.theta) for t in b.get_normalization_transforms()
+        ]
+        mine = [
+            tuple(float(x[k]) for k in ("offset_x", "offset_y", "scale_x", "scale_y", "shear", "theta"))
+            for x in g.normalization_transforms
+        ]
+        assert mine == [tuple(float(v) for v in t) for t in xf]
+
+
+def test_gtc_norm_transform_count_matches_bpm_bins(gsa_bpm) -> None:  # type: ignore[no-untyped-def]
+    need(ORACLE_GTCS)
+    g = read_gtc(sorted(ORACLE_GTCS.glob("*.gtc"))[0])
+    assert len(g.normalization_transforms) == len(np.unique(gsa_bpm.norm_bin_ids()))
+    assert g.num_snps == gsa_bpm.num_loci
+
+
+@pytest.mark.slow
+def test_bpm_matches_beadarrayfiles(gsa_bpm, beadarrayfiles) -> None:  # type: ignore[no-untyped-def]
+    b = beadarrayfiles.BeadPoolManifest(str(GSA_BPM))
+    m = gsa_bpm
+    assert m.names == b.names
+    assert m.snps == b.snps
+    assert m.chroms == b.chroms
+    assert m.map_infos.tolist() == b.map_infos
+    assert m.address_a.tolist() == b.addresses
+    assert m.assay_types.tolist() == b.assay_types
+    assert m.normalization_lookups().tolist() == b.normalization_lookups
+    rs = beadarrayfiles.RefStrand
+    assert [rs.from_string(s) for s in m.ref_strands] == b.ref_strands
+    assert m.manifest_name == b.manifest_name and m.control_config == b.control_config
+
+
+@pytest.mark.slow
+@pytest.mark.parametrize("egt_path", [GSA_EGT, GDA_EGT], ids=["gsa", "gda"])
+def test_egt_matches_beadarrayfiles(egt_path, beadarrayfiles) -> None:  # type: ignore[no-untyped-def]
+    need(egt_path)
+    e = read_egt(egt_path)
+    with open(egt_path, "rb") as fh:
+        b = beadarrayfiles.ClusterFile.read_cluster_file(fh)
+    assert e.num_records == len(b.records) if hasattr(b, "records") else True
+    idx = e.name_index()
+    for name in e.names[:: max(1, e.num_records // 5000)]:  # ~5000 evenly spaced records
+        rec = b.get_record(name)
+        i = idx[name]
+        for k, stats in enumerate((rec.aa_cluster_stats, rec.ab_cluster_stats, rec.bb_cluster_stats)):
+            assert e.n[i, k] == stats.N
+            assert e.theta_mean[i, k] == np.float32(stats.theta_mean)
+            assert e.theta_dev[i, k] == np.float32(stats.theta_dev)
+            assert e.r_mean[i, k] == np.float32(stats.r_mean)
+            assert e.r_dev[i, k] == np.float32(stats.r_dev)
+        assert e.intensity_threshold[i] == np.float32(rec.intensity_threshold)
+        assert e.gentrain_score[i] == np.float32(rec.cluster_score.total_score)
+        assert e.cluster_sep[i] == np.float32(rec.cluster_score.cluster_separation)
+        assert e.addresses[i] == rec.address
+
+
+@pytest.mark.slow
+def test_csv_manifest_consistent_with_bpm(gsa_bpm) -> None:  # type: ignore[no-untyped-def]  # IN-042
+    need(GSA_CSV)
+    c = read_csv_manifest(GSA_CSV)
+    assert compare_manifests(gsa_bpm, c) == []
+    # Assay type inferred from CSV sequences (fallback for BPM locus entries < v6, gtc2vcf rule) agrees with the
+    # BPM's explicit value for 99.71% of Infinium I loci on GSA v4; the BPM is authoritative (known_unknowns K-003).
+    inf1 = gsa_bpm.address_b != 0
+    agree = (c.assay_types[inf1] == gsa_bpm.assay_types[inf1]).mean()
+    assert agree >= 0.995
+
+
+def test_gda_demo_sheet() -> None:
+    need(GDA_SHEET)
+    s = parse_samplesheet(GDA_SHEET)
+    assert len(s.rows) == 36 and s.has_sample_name
+    assert {r.values["Parent1"] for r in s.rows} >= {"NA12891"}