Materialize electrode-table columns before indexing in the NWB recording extractor

src/spikeinterface/extractors/nwbextractors.py

@@ -750,7 +750,9 @@ def _fetch_recording_segment_info_backend(self, file, cache, load_time_vector, s
 
         # If channel names are present, use them as channel_ids instead of the electrode ids
         if "channel_name" in electrode_table_columns:
-            channel_names = electrodes_table["channel_name"]
+            # Materialize the column to numpy first: the electrodes region can reference rows in
+            # any order, but h5py only fancy-indexes with strictly increasing indices.
+            channel_names = electrodes_table["channel_name"][:]
             channel_ids = channel_names[electrodes_indices]
             # Decode if bytes with utf-8
             channel_ids = [x.decode("utf-8") if isinstance(x, bytes) else x for x in channel_ids]
@@ -774,15 +776,17 @@ def _fetch_locations_and_groups(self, electrodes_table, electrodes_indices):
             if "rel_y" in electrodes_table:
                 ndim = 3 if "rel_z" in electrodes_table else 2
                 locations = np.zeros((self.get_num_channels(), ndim), dtype=float)
-                locations[:, 0] = electrodes_table["rel_x"][electrodes_indices]
-                locations[:, 1] = electrodes_table["rel_y"][electrodes_indices]
+                # Materialize each column to numpy first: the electrodes region can reference rows in
+                # any order, but h5py only fancy-indexes with strictly increasing indices.
+                locations[:, 0] = electrodes_table["rel_x"][:][electrodes_indices]
+                locations[:, 1] = electrodes_table["rel_y"][:][electrodes_indices]
                 if "rel_z" in electrodes_table:
-                    locations[:, 2] = electrodes_table["rel_z"][electrodes_indices]
+                    locations[:, 2] = electrodes_table["rel_z"][:][electrodes_indices]
 
         # Channel groups
         groups = None
         if "group_name" in electrodes_table:
-            groups = electrodes_table["group_name"][electrodes_indices][:]
+            groups = electrodes_table["group_name"][:][electrodes_indices]
         if groups is not None:
             groups = np.array([x.decode("utf-8") if isinstance(x, bytes) else x for x in groups])
         return locations, groups
@@ -809,7 +813,9 @@ def _fetch_other_properties(self, electrodes_table, electrodes_indices, columns)
                 continue
             else:
                 column_name = rename_properties.get(column, column)
-                properties[column_name] = electrodes_table[column][electrodes_indices]
+                # Materialize the column to numpy first: the electrodes region can reference rows in
+                # any order, but h5py only fancy-indexes with strictly increasing indices.
+                properties[column_name] = electrodes_table[column][:][electrodes_indices]
 
         return properties
 
@@ -833,7 +839,8 @@ def _fetch_main_properties_pynwb(self):
         # Channel offsets
         offset = self.electrical_series.offset if hasattr(self.electrical_series, "offset") else 0
         if offset == 0 and "offset" in electrodes_table:
-            offset = electrodes_table["offset"].data[electrodes_indices]
+            # See note in _fetch_locations_and_groups: materialize to numpy before reordering.
+            offset = electrodes_table["offset"].data[:][electrodes_indices]
         offsets = offset * 1e6
 
         locations, groups = self._fetch_locations_and_groups(electrodes_table, electrodes_indices)
@@ -865,7 +872,8 @@ def _fetch_main_properties_backend(self):
         # Channel offsets
         offset = data_attributes["offset"] if "offset" in data_attributes else 0
         if offset == 0 and "offset" in electrodes_table:
-            offset = electrodes_table["offset"][electrodes_indices]
+            # See note in _fetch_locations_and_groups: materialize to numpy before reordering.
+            offset = electrodes_table["offset"][:][electrodes_indices]
         offsets = offset * 1e6
 
         # Channel locations and groups

src/spikeinterface/extractors/tests/test_nwbextractors.py

@@ -254,6 +254,63 @@ def test_nwb_extractor_offset_from_electrodes_table(generate_nwbfile, use_pynwb)
     assert np.array_equal(extracted_offsets_uV, expected_offsets_uV)
 
 
+@pytest.mark.parametrize("use_pynwb", [True, False])
+def test_nwb_extractor_electrodes_region_out_of_order(tmp_path, use_pynwb):
+    """An ElectricalSeries may reference its electrodes in any order (e.g. channels reordered by
+    depth during processing). h5py rejects fancy indexing with non-increasing indices, so reading
+    such a file used to raise "Indexing elements must be in increasing order" (GH-4619)."""
+    from pynwb import NWBHDF5IO
+    from pynwb.ecephys import ElectricalSeries
+    from pynwb.testing.mock.file import mock_NWBFile
+    from pynwb.testing.mock.device import mock_Device
+    from pynwb.testing.mock.ecephys import mock_ElectrodeGroup
+
+    nwbfile = mock_NWBFile()
+    device = mock_Device(name="probe")
+    nwbfile.add_device(device)
+    nwbfile.add_electrode_column(name="channel_name", description="channel name")
+    nwbfile.add_electrode_column(name="rel_x", description="rel_x")
+    nwbfile.add_electrode_column(name="rel_y", description="rel_y")
+    nwbfile.add_electrode_column(name="property", description="A property")
+    nwbfile.add_electrode_column(name="offset", description="offset")
+    electrode_group = mock_ElectrodeGroup(device=device)
+    nwbfile.add_electrode_group(electrode_group)
+
+    num_electrodes = 5
+    for index in range(num_electrodes):
+        nwbfile.add_electrode(
+            id=index,
+            group=electrode_group,
+            location="brain",
+            channel_name=f"ch{index}",
+            rel_x=float(index),
+            rel_y=float(index),
+            property=f"prop{index}",
+            offset=float(index),
+        )
+
+    # The region is deliberately not in increasing order.
+    region = [4, 2, 0, 3, 1]
+    electrode_region = nwbfile.create_electrode_table_region(region=region, description="electrodes")
+    data = np.random.default_rng(0).random(size=(100, num_electrodes))
+    electrical_series = ElectricalSeries(name="ElectricalSeries", data=data, electrodes=electrode_region, rate=30_000.0)
+    nwbfile.add_acquisition(electrical_series)
+
+    nwbfile_path = tmp_path / "out_of_order.nwb"
+    with NWBHDF5IO(str(nwbfile_path), mode="w") as io:
+        io.write(nwbfile)
+
+    recording = NwbRecordingExtractor(
+        nwbfile_path, electrical_series_path="acquisition/ElectricalSeries", use_pynwb=use_pynwb
+    )
+
+    # Everything pulled from the electrodes table must follow the region order, not the table order.
+    assert np.array_equal(recording.channel_ids, np.array([f"ch{i}" for i in region]))
+    assert np.array_equal(recording.get_channel_locations(), np.array([[float(i), float(i)] for i in region]))
+    assert np.array_equal(recording.get_property("property"), np.array([f"prop{i}" for i in region]))
+    assert np.array_equal(recording.get_channel_offsets(), np.array([float(i) for i in region]) * 1e6)
+
+
 @pytest.mark.parametrize("use_pynwb", [True, False])
 def test_nwb_extractor_offset_from_series(generate_nwbfile, use_pynwb):
     """Test that the offset is retrieved from the ElectricalSeries if it is present."""