Fix the IBL Cbin extractor so that it works with 4 split shanks as well

src/spikeinterface/extractors/cbin_ibl.py

@@ -1,4 +1,5 @@
 from pathlib import Path
+import re
 import warnings
 import numpy as np
 
@@ -9,6 +10,89 @@
 from spikeinterface.core.core_tools import define_function_from_class
 
 
+def _parse_spikeglx_meta_table(meta_value):
+    """Return parenthesized entries from a SpikeGLX meta value.
+
+    Accepts either the raw string form ``"(a)(b)(c)"`` or an already-parsed list.
+    """
+    if isinstance(meta_value, str):
+        return re.findall(r"\(([^()]*)\)", meta_value)
+    return list(meta_value)
+
+
+def _parse_saved_channel_subset(subset_text):
+    if subset_text is None or subset_text == "all":
+        return None
+    channels = []
+    for entry in subset_text.split(","):
+        if ":" in entry:
+            start, stop = entry.split(":")
+            channels.extend(np.arange(int(start), int(stop) + 1))
+        else:
+            channels.append(int(entry))
+    return np.asarray(channels, dtype="int64")
+
+
+def _get_saved_channel_indices(meta):
+    subset_text = meta.get("snsSaveChanSubset_orig")
+    if subset_text is None:
+        subset_text = meta.get("snsSaveChanSubset")
+    return _parse_saved_channel_subset(subset_text)
+
+
+def _read_cbin_probe(meta_file, meta):
+    """Build a Probe for a SpikeGLX meta, honoring IBL shank-split subsets.
+    This is a specific to IBL cbin files, that sometimes do not have `~` in front of the .meta file fields.
+    """
+    if "snsSaveChanSubset_orig" not in meta:
+        return probeinterface.read_spikeglx(meta_file)
+
+    # probeinterface removed the single-shot `_read_imro_string` helper used
+    # previously (in 0.3.2). The equivalent flow is: build the full geometric
+    # probe, parse the IMRO table to a dict, resolve which contact IDs are
+    # active, and slice the full probe down to those contacts.
+    from probeinterface.neuropixels_tools import (
+        build_neuropixels_probe,
+        _parse_imro_string,
+        _get_imro_active_contact_ids,
+        _annotate_probe_with_adc_sampling_info,
+    )
+
+    imDatPrb_pn = meta.get("imDatPrb_pn")
+    # neo's read_meta_file leaves imroTbl as a raw string when the meta field
+    # has no leading '~' (the IBL cbin case noted above). _parse_imro_string
+    # requires the string form, so rebuild it defensively if neo split it.
+    imro_table_string = meta["imroTbl"]
+    if isinstance(imro_table_string, list):
+        imro_table_string = "(" + ")(".join(imro_table_string) + ")"
+
+    full_probe = build_neuropixels_probe(probe_part_number=imDatPrb_pn)
+    imro_per_channel = _parse_imro_string(imro_table_string)
+    active_contact_ids = _get_imro_active_contact_ids(imro_per_channel)
+    contact_id_to_index = {cid: i for i, cid in enumerate(full_probe.contact_ids)}
+    selected = np.array([contact_id_to_index[cid] for cid in active_contact_ids])
+    probe = full_probe.get_slice(selected)
+
+    # Attach ADC sampling annotations (num_channels_per_adc and per-contact
+    # adc_sample_order) needed by get_neuropixels_sample_shifts_from_probe for
+    # phase_shift preprocessing. Must run while contacts are still in readout
+    # channel order, i.e. before the snsSaveChanSubset_orig slice below.
+    adc_sampling_table = probe.annotations.get("adc_sampling_table")
+    _annotate_probe_with_adc_sampling_info(probe, adc_sampling_table)
+
+    saved = _get_saved_channel_indices(meta)
+    if saved is not None:
+        saved = saved[saved < probe.get_contact_count()]
+        if saved.size != probe.get_contact_count():
+            probe = probe.get_slice(saved)
+    probe.annotate(serial_number=meta.get("imDatPrb_sn", meta.get("imProbeSN", None)))
+    probe.annotate(part_number=meta.get("imDatPrb_pn", None))
+    probe.annotate(port=meta.get("imDatPrb_port", None))
+    probe.annotate(slot=meta.get("imDatPrb_slot", None))
+    probe.set_device_channel_indices(np.arange(probe.get_contact_count()))
+    return probe
+
+
 class CompressedBinaryIblExtractor(BaseRecording):
     """Load IBL data as an extractor object.
 
@@ -99,7 +183,7 @@ def __init__(self, folder_path=None, load_sync_channel=False, stream_name="ap",
         self.set_channel_offsets(offsets)
 
         if not load_sync_channel:
-            probe = probeinterface.read_spikeglx(meta_file)
+            probe = _read_cbin_probe(meta_file, meta)
 
             if probe.shank_ids is not None:
                 self.set_probe(probe, in_place=True, group_mode="by_shank")
@@ -204,9 +288,31 @@ def extract_stream_info(meta_file, meta):
     info["stream_kind"] = stream_kind
     info["stream_name"] = stream_name
     info["units"] = units
-    info["channel_names"] = [txt.split(";")[0] for txt in meta["snsChanMap"]]
+
+    chan_map_entries = _parse_spikeglx_meta_table(meta["snsChanMap"])
+    # First entry is the header tuple like "(384,0,1)"; drop it when present.
+    if chan_map_entries and "," in chan_map_entries[0] and ";" not in chan_map_entries[0]:
+        chan_map_entries = chan_map_entries[1:]
+    full_channel_names = [entry.split(";")[0] for entry in chan_map_entries]
+
+    saved_indices = _get_saved_channel_indices(meta)
+    if saved_indices is not None and len(full_channel_names) >= saved_indices.max() + 1:
+        channel_names = [full_channel_names[i] for i in saved_indices]
+    else:
+        channel_names = full_channel_names
+
+    channel_offsets = np.zeros(num_chan)
+
+    if len(channel_names) != num_chan or channel_gains.shape[0] != num_chan:
+        raise ValueError(
+            "Parsed channel metadata does not match nSavedChans: "
+            f"channel_names={len(channel_names)}, channel_gains={channel_gains.shape[0]}, "
+            f"num_chan={num_chan}"
+        )
+
+    info["channel_names"] = channel_names
     info["channel_gains"] = channel_gains
-    info["channel_offsets"] = np.zeros(info["num_chan"])
+    info["channel_offsets"] = channel_offsets
     info["has_sync_trace"] = has_sync_trace
 
     return info

src/spikeinterface/extractors/tests/test_cbin_ibl_extractors.py

@@ -1,8 +1,17 @@
 import pytest
 import unittest
 
+import numpy as np
+
 from spikeinterface.extractors.extractor_classes import CompressedBinaryIblExtractor, read_cbin_ibl
 
+from spikeinterface.extractors.cbin_ibl import (
+    _get_saved_channel_indices,
+    _parse_saved_channel_subset,
+    _parse_spikeglx_meta_table,
+    _read_cbin_probe,
+    extract_stream_info,
+)
 
 from spikeinterface.extractors.tests.common_tests import RecordingCommonTestSuite, SortingCommonTestSuite
 
@@ -13,6 +22,201 @@ class CompressedBinaryIblExtractorTest(RecordingCommonTestSuite, unittest.TestCa
     entities = []
 
 
+# --------------------------------------------------------------------------------
+# NP2.4 shank-split meta handling
+#
+# When IBL splits an NP2.4 recording into four single-shank recordings, the
+# generated .meta file has two peculiarities:
+#
+#   * none of the list-valued fields carry the leading "~" that SpikeGLX writes,
+#     so neo.rawio.spikeglxrawio.read_meta_file leaves them as raw strings;
+#   * the true acquisition channels live in "snsSaveChanSubset_orig" while
+#     "snsSaveChanSubset" only holds the local, post-split renumbering.
+#
+# probeinterface.read_spikeglx uses "snsSaveChanSubset" and therefore returns a
+# silently wrong probe for these files (sync counted as a contact, contacts
+# spread over the wrong shanks). Correcting that is what _read_cbin_probe does.
+# These tests need only a .meta file, so they run without any downloaded data.
+# --------------------------------------------------------------------------------
+
+# 16 acquired channels in blocks of four, alternating shank 0 / shank 1, so the
+# saved subset is non-contiguous exactly like a real 4-shank split. Channels
+# 4-7 and 12-15 are the shank-1 channels, i.e. electrodes 0-7 on shank 1.
+N_ACQUIRED_CHANNELS = 16
+SYNC_CHANNEL_INDEX = N_ACQUIRED_CHANNELS
+SAVED_CHANNEL_SUBSET_ORIG = "4:7,12:15,16"
+EXPECTED_CHANNEL_NAMES = ["AP4", "AP5", "AP6", "AP7", "AP12", "AP13", "AP14", "AP15", "SY0"]
+
+
+def _build_shank_split_meta_text(include_subset_orig=True):
+    """Build a minimal NP2.4 shank-split .meta text, deliberately without any "~".
+
+    Parameters
+    ----------
+    include_subset_orig : bool, default: True
+        If False, omit "snsSaveChanSubset_orig" so the meta looks like a normal
+        SpikeGLX file and _read_cbin_probe should defer to probeinterface.
+
+    Returns
+    -------
+    str
+        The full contents of the .meta file.
+    """
+    imro_entries, geom_entries, chan_map_entries = [], [], []
+    for channel in range(N_ACQUIRED_CHANNELS):
+        block, position = divmod(channel, 4)
+        shank = block % 2
+        electrode = (block // 2) * 4 + position
+        # type-24 (NP2.4) IMRO entry: (channel shank bank refid electrode)
+        imro_entries.append(f"({channel} {shank} 0 1 {electrode})")
+        # NP2.4 geometry: two columns at 27/59 um from the shank edge, 15 um row pitch
+        x_um = 27 if electrode % 2 == 0 else 59
+        geom_entries.append(f"({shank}:{x_um}:{15 * (electrode // 2)}:1)")
+        chan_map_entries.append(f"(AP{channel};{channel}:{channel})")
+    chan_map_entries.append(f"(SY0;{SYNC_CHANNEL_INDEX}:{SYNC_CHANNEL_INDEX})")
+
+    num_saved_channels = len(EXPECTED_CHANNEL_NAMES)
+    fields = {
+        "typeThis": "imec",
+        "acqApLfSy": f"{N_ACQUIRED_CHANNELS},0,1",
+        "snsApLfSy": f"{num_saved_channels - 1},0,1",
+        "nSavedChans": str(num_saved_channels),
+        # the local, post-split numbering that read_spikeglx wrongly relies on
+        "snsSaveChanSubset": f"0:{num_saved_channels - 1}",
+        "imSampRate": "30000",
+        "fileSizeBytes": str(num_saved_channels * 2 * 3000),
+        "imAiRangeMax": "0.5",
+        "imMaxInt": "8192",
+        "imChan0apGain": "80",
+        "imDatPrb_type": "24",
+        "imDatPrb_pn": "NP2010",
+        "imDatPrb_sn": "20472319942",
+        "imDatPrb_port": "1",
+        "imDatPrb_slot": "3",
+        "imroTbl": f"(24,{N_ACQUIRED_CHANNELS})" + "".join(imro_entries),
+        "snsChanMap": f"({N_ACQUIRED_CHANNELS},0,1)" + "".join(chan_map_entries),
+        "snsGeomMap": "(NP2010,4,250,70)" + "".join(geom_entries),
+    }
+    if include_subset_orig:
+        fields["snsSaveChanSubset_orig"] = SAVED_CHANNEL_SUBSET_ORIG
+
+    # NOTE: written without the leading "~" that SpikeGLX puts on list-valued
+    # fields. That omission is the whole point of this fixture.
+    return "".join(f"{key}={value}\n" for key, value in fields.items())
+
+
+def _read_meta(meta_file):
+    """Parse a .meta file the same way CompressedBinaryIblExtractor does."""
+    from neo.rawio.spikeglxrawio import read_meta_file
+
+    return read_meta_file(str(meta_file))
+
+
+@pytest.fixture
+def shank_split_meta_file(tmp_path):
+    """Path to a minimal NP2.4 shank-split meta file."""
+    meta_file = tmp_path / "np24_shank_split.ap.meta"
+    meta_file.write_text(_build_shank_split_meta_text())
+    return meta_file
+
+
+def test_parse_spikeglx_meta_table_accepts_string_and_list():
+    """The parser must cope with both the raw "~"-less string and a parsed list."""
+    assert _parse_spikeglx_meta_table("(a)(b)(c)") == ["a", "b", "c"]
+    assert _parse_spikeglx_meta_table(["a", "b", "c"]) == ["a", "b", "c"]
+
+
+def test_parse_saved_channel_subset():
+    """Ranges, singletons and the "all" sentinel."""
+    assert _parse_saved_channel_subset("all") is None
+    assert _parse_saved_channel_subset(None) is None
+    np.testing.assert_array_equal(_parse_saved_channel_subset("0:3"), [0, 1, 2, 3])
+    np.testing.assert_array_equal(_parse_saved_channel_subset("4:7,12:15,16"), [4, 5, 6, 7, 12, 13, 14, 15, 16])
+
+
+def test_saved_channel_indices_prefer_orig(shank_split_meta_file):
+    """snsSaveChanSubset_orig must take precedence over the local snsSaveChanSubset."""
+    meta = _read_meta(shank_split_meta_file)
+    np.testing.assert_array_equal(_get_saved_channel_indices(meta), [4, 5, 6, 7, 12, 13, 14, 15, 16])
+
+
+def test_meta_without_tilde_stays_a_string(shank_split_meta_file):
+    """Guards the premise of the fix: neo cannot list-parse these fields."""
+    meta = _read_meta(shank_split_meta_file)
+    assert isinstance(meta["imroTbl"], str)
+    assert isinstance(meta["snsChanMap"], str)
+
+
+def test_read_cbin_probe_selects_the_correct_shank(shank_split_meta_file):
+    """The probe must hold only the 8 shank-1 contacts, at the right positions."""
+    meta = _read_meta(shank_split_meta_file)
+    probe = _read_cbin_probe(str(shank_split_meta_file), meta)
+
+    assert probe.get_contact_count() == 8
+    assert set(probe.shank_ids) == {"1"}
+    # shank pitch 250 um, two columns 32 um apart, four rows at 15 um pitch
+    np.testing.assert_array_equal(np.unique(probe.contact_positions[:, 0]), [250.0, 282.0])
+    np.testing.assert_array_equal(np.unique(probe.contact_positions[:, 1]), [0.0, 15.0, 30.0, 45.0])
+    np.testing.assert_array_equal(probe.device_channel_indices, np.arange(8))
+    assert probe.annotations["serial_number"] == "20472319942"
+    assert probe.annotations["part_number"] == "NP2010"
+
+
+def test_read_cbin_probe_differs_from_plain_read_spikeglx(shank_split_meta_file):
+    """Regression guard: plain read_spikeglx returns a wrong probe for this meta.
+
+    It applies the local "snsSaveChanSubset" to the geometry table, so the sync
+    channel becomes a contact and the contacts land on two different shanks.
+    If probeinterface ever learns about "snsSaveChanSubset_orig" this test will
+    fail, which is the signal that _read_cbin_probe can be removed.
+    """
+    import probeinterface
+
+    meta = _read_meta(shank_split_meta_file)
+    wrong_probe = probeinterface.read_spikeglx(str(shank_split_meta_file))
+    correct_probe = _read_cbin_probe(str(shank_split_meta_file), meta)
+
+    assert wrong_probe.get_contact_count() == 9  # 8 contacts + sync wrongly counted as a contact
+    assert set(wrong_probe.shank_ids) == {"0", "1"}  # contacts wrongly spread over two shanks
+    assert correct_probe.get_contact_count() != wrong_probe.get_contact_count()
+
+
+def test_read_cbin_probe_falls_back_when_no_subset_orig(tmp_path):
+    """A normal SpikeGLX meta must be passed through to probeinterface untouched."""
+    import probeinterface
+
+    meta_file = tmp_path / "plain.ap.meta"
+    meta_file.write_text(_build_shank_split_meta_text(include_subset_orig=False))
+    meta = _read_meta(meta_file)
+
+    probe = _read_cbin_probe(str(meta_file), meta)
+    reference_probe = probeinterface.read_spikeglx(str(meta_file))
+    np.testing.assert_array_equal(probe.contact_positions, reference_probe.contact_positions)
+
+
+def test_extract_stream_info_uses_the_original_subset(shank_split_meta_file):
+    """Channel names must follow snsSaveChanSubset_orig, not the local numbering."""
+    meta = _read_meta(shank_split_meta_file)
+    info = extract_stream_info(str(shank_split_meta_file), meta)
+
+    assert info["num_chan"] == len(EXPECTED_CHANNEL_NAMES)
+    assert info["channel_names"] == EXPECTED_CHANNEL_NAMES
+    assert info["channel_gains"].shape[0] == info["num_chan"]
+    assert info["channel_offsets"].shape[0] == info["num_chan"]
+    assert info["has_sync_trace"]
+
+
+def test_extract_stream_info_raises_on_inconsistent_meta(tmp_path):
+    """The nSavedChans consistency guard must fire on a corrupted meta."""
+    meta_text = _build_shank_split_meta_text().replace("nSavedChans=9", "nSavedChans=42")
+    meta_file = tmp_path / "bad.ap.meta"
+    meta_file.write_text(meta_text)
+    meta = _read_meta(meta_file)
+
+    with pytest.raises(ValueError, match="does not match nSavedChans"):
+        extract_stream_info(str(meta_file), meta)
+
+
 # ~ def test_read_cbin_ibl():
 # ~ base_folder = '/media/samuel/dataspikesorting/DataSpikeSorting/olivier_destripe/'
 # ~ data_folder = base_folder + '4c04120d-523a-4795-ba8f-49dbb8d9f63a'