Merge branch 'master' into neuronexus

Author: zm711

neo/rawio/plexon2rawio/plexon2rawio.py

@@ -53,7 +53,7 @@ class Plexon2RawIO(BaseRawIO):
     pl2_dll_file_path: str | Path | None, default: None
         The path to the necessary dll for loading pl2 files
         If None will find correct dll for architecture and if it does not exist will download it
-    reading_attempts: int, default: 15
+    reading_attempts: int, default: 25
         Number of attempts to read the file before raising an error
         This opening process is somewhat unreliable and might fail occasionally. Adjust this higher
         if you encounter problems in opening the file.
@@ -92,7 +92,7 @@ class Plexon2RawIO(BaseRawIO):
     extensions = ["pl2"]
     rawmode = "one-file"
 
-    def __init__(self, filename, pl2_dll_file_path=None, reading_attempts=15):
+    def __init__(self, filename, pl2_dll_file_path=None, reading_attempts=25):
 
         # signals, event and spiking data will be cached
         # cached signal data can be cleared using `clear_analogsignal_cache()()`
@@ -196,6 +196,7 @@ def _parse_header(self):
             "FP": "FPl-Low Pass Filtered",
             "SP": "SPKC-High Pass Filtered",
             "AI": "AI-Auxiliary Input",
+            "AIF": "AIF-Auxiliary Input Filtered",
         }
 
         unique_stream_ids = np.unique(signal_channels["stream_id"])
@@ -209,17 +210,17 @@ def _parse_header(self):
 
         signal_streams = np.array(signal_streams, dtype=_signal_stream_dtype)
 
-        self.stream_id_samples = {}
-        self.stream_index_to_stream_id = {}
+        self._stream_id_samples = {}
+        self._stream_index_to_stream_id = {}
         for stream_index, stream_id in enumerate(signal_streams["id"]):
             # Keep a mapping from stream_index to stream_id
-            self.stream_index_to_stream_id[stream_index] = stream_id
+            self._stream_index_to_stream_id[stream_index] = stream_id
 
             # We extract the number of samples for each stream
             mask = signal_channels["stream_id"] == stream_id
             signal_num_samples = np.unique(channel_num_samples[mask])
             assert signal_num_samples.size == 1, "All channels in a stream must have the same number of samples"
-            self.stream_id_samples[stream_id] = signal_num_samples[0]
+            self._stream_id_samples[stream_id] = signal_num_samples[0]
 
         # pre-loading spike channel_data for later usage
         self._spike_channel_cache = {}
@@ -231,7 +232,14 @@ def _parse_header(self):
             if not (schannel_info.m_ChannelEnabled and schannel_info.m_ChannelRecordingEnabled):
                 continue
 
-            for channel_unit_id in range(schannel_info.m_NumberOfUnits):
+            # In a PL2 spike channel header, the field "m_NumberOfUnits" denotes the number
+            # of units to which spikes detected on that channel have been assigned. It does
+            # not account for unsorted spikes, i.e., spikes that have not been assigned to 
+            # a unit. It is Plexon's convention to assign unsorted spikes to channel_unit_id=0, 
+            # and sorted spikes to channel_unit_id's 1, 2, 3...etc. Therefore, for a given value of 
+            # m_NumberOfUnits, there are m_NumberOfUnits+1 channel_unit_ids to consider - 1 
+            # unsorted channel_unit_id (0) + the m_NumberOfUnits sorted channel_unit_ids. 
+            for channel_unit_id in range(schannel_info.m_NumberOfUnits+1):
                 unit_name = f"{schannel_info.m_Name.decode()}.{channel_unit_id}"
                 unit_id = f"unit{schannel_info.m_Channel}.{channel_unit_id}"
                 wf_units = schannel_info.m_Units
@@ -386,8 +394,8 @@ def _segment_t_stop(self, block_index, seg_index):
         return float(end_time / self.pl2reader.pl2_file_info.m_TimestampFrequency)
 
     def _get_signal_size(self, block_index, seg_index, stream_index):
-        stream_id = self.stream_index_to_stream_id[stream_index]
-        num_samples = int(self.stream_id_samples[stream_id])
+        stream_id = self._stream_index_to_stream_id[stream_index]
+        num_samples = int(self._stream_id_samples[stream_id])
         return num_samples
 
     def _get_signal_t_start(self, block_index, seg_index, stream_index):

neo/rawio/plexonrawio.py

@@ -43,6 +43,7 @@
     _event_channel_dtype,
 )
 
+from neo.core.baseneo import NeoReadWriteError
 
 class PlexonRawIO(BaseRawIO):
     extensions = ["plx"]
@@ -230,9 +231,19 @@ def _parse_header(self):
                 self._data_blocks[bl_type][chan_id] = data_block
 
         # signals channels
-        sig_channels = []
-        all_sig_length = []
         source_id = []
+        
+        # Scanning sources and populating signal channels at the same time. Sources have to have
+        # same sampling rate and number of samples to belong to one stream.
+        signal_channels = []
+        channel_num_samples = []
+        
+        # We will build the stream ids based on the channel prefixes
+        # The channel prefixes are the first characters of the channel names which have the following format:
+        # WB{number}, FPX{number}, SPKCX{number}, AI{number}, etc
+        # We will extract the prefix and use it as stream id
+        regex_prefix_pattern = r"^\D+"  # Match any non-digit character at the beginning of the string
+        
         if self.progress_bar:
             chan_loop = trange(nb_sig_chan, desc="Parsing signal channels", leave=True)
         else:
@@ -245,7 +256,7 @@ def _parse_header(self):
             if length == 0:
                 continue  # channel not added
             source_id.append(h["SrcId"])
-            all_sig_length.append(length)
+            channel_num_samples.append(length)
             sampling_rate = float(h["ADFreq"])
             sig_dtype = "int16"
             units = ""  # I don't know units
@@ -258,61 +269,60 @@ def _parse_header(self):
                     0.5 * (2 ** global_header["BitsPerSpikeSample"]) * h["Gain"] * h["PreampGain"]
                 )
             offset = 0.0
-            stream_id = "0"  # This is overwritten later
-            sig_channels.append((name, str(chan_id), sampling_rate, sig_dtype, units, gain, offset, stream_id))
+            channel_prefix = re.match(regex_prefix_pattern, name).group(0)
+            stream_id = channel_prefix
+            
+            signal_channels.append((name, str(chan_id), sampling_rate, sig_dtype, units, gain, offset, stream_id))
 
-        sig_channels = np.array(sig_channels, dtype=_signal_channel_dtype)
+        signal_channels = np.array(signal_channels, dtype=_signal_channel_dtype)
 
-        if sig_channels.size == 0:
+        if signal_channels.size == 0:
             signal_streams = np.array([], dtype=_signal_stream_dtype)
 
         else:
             # Detect streams
-            all_sig_length = np.asarray(all_sig_length)
-
-            # names are WB{number}, FPX{number}, SPKCX{number}, AI{number}, etc
-            pattern = r"^\D+"  # Match any non-digit character at the beginning of the string
-            channels_prefixes = np.asarray([re.match(pattern, name).group(0) for name in sig_channels["name"]])
-            buffer_stream_groups = set(zip(channels_prefixes, sig_channels["sampling_rate"], all_sig_length))
-
-            # There are explanations of the streams based on channel names
-            # provided by a Plexon Engineer, see here:
+            channel_num_samples = np.asarray(channel_num_samples)
+            # We are using channel prefixes as stream_ids
+            # The meaning of the channel prefixes was provided by a Plexon Engineer, see here:
             # https://github.com/NeuralEnsemble/python-neo/pull/1495#issuecomment-2184256894
-            channel_prefix_to_stream_name = {
+            stream_id_to_stream_name = {
                 "WB": "WB-Wideband",
-                "FP": "FPl-Low Pass Filtered ",
+                "FP": "FPl-Low Pass Filtered",
                 "SP": "SPKC-High Pass Filtered",
                 "AI": "AI-Auxiliary Input",
+                "AIF": "AIF-Auxiliary Input Filtered",
             }
 
-            # Using a mapping to ensure consistent order of stream_index
-            channel_prefix_to_stream_id = {
-                "WB": "0",
-                "FP": "1",
-                "SP": "2",
-                "AI": "3",
-            }
-
+            unique_stream_ids = np.unique(signal_channels["stream_id"])
             signal_streams = []
-            self._signal_length = {}
-            self._sig_sampling_rate = {}
-
-            for stream_index, (channel_prefix, sr, length) in enumerate(buffer_stream_groups):
-                # The users of plexon can modify the prefix of the channel names (e.g. `my_prefix` instead of `WB`). This is not common but in that case
-                # We assign the channel_prefix both as stream_name and stream_id
-                stream_name = channel_prefix_to_stream_name.get(channel_prefix, channel_prefix)
-                stream_id = channel_prefix_to_stream_id.get(channel_prefix, channel_prefix)
-
-                mask = (sig_channels["sampling_rate"] == sr) & (all_sig_length == length)
-                sig_channels["stream_id"][mask] = stream_id
-
-                self._sig_sampling_rate[stream_index] = sr
-                self._signal_length[stream_index] = length
-
+            for stream_id in unique_stream_ids:
+                # We are using the channel prefixes as ids
+                # The users of plexon can modify the prefix of the channel names (e.g. `my_prefix` instead of `WB`).
+                # In that case we use the channel prefix both as stream id and name
+                stream_name = stream_id_to_stream_name.get(stream_id, stream_id)
                 signal_streams.append((stream_name, stream_id))
 
             signal_streams = np.array(signal_streams, dtype=_signal_stream_dtype)
 
+            self._stream_id_samples = {}
+            self._stream_id_sampling_frequency = {}
+            self._stream_index_to_stream_id = {}
+            for stream_index, stream_id in enumerate(signal_streams["id"]):
+                # Keep a mapping from stream_index to stream_id
+                self._stream_index_to_stream_id[stream_index] = stream_id
+
+                mask = signal_channels["stream_id"] == stream_id
+                
+                signal_num_samples = np.unique(channel_num_samples[mask])
+                if signal_num_samples.size > 1:
+                    raise NeoReadWriteError(f"Channels in stream {stream_id} don't have the same number of samples")
+                self._stream_id_samples[stream_id] = signal_num_samples[0]
+
+                signal_sampling_frequency = np.unique(signal_channels[mask]["sampling_rate"])
+                if signal_sampling_frequency.size > 1:
+                    raise NeoReadWriteError(f"Channels in stream {stream_id} don't have the same sampling frequency")
+                self._stream_id_sampling_frequency[stream_id] = signal_sampling_frequency[0]
+                
         self._global_ssampling_rate = global_header["ADFrequency"]
 
         # Determine number of units per channels
@@ -374,7 +384,7 @@ def _parse_header(self):
             "nb_block": 1,
             "nb_segment": [1],
             "signal_streams": signal_streams,
-            "signal_channels": sig_channels,
+            "signal_channels": signal_channels,
             "spike_channels": spike_channels,
             "event_channels": event_channels,
         }
@@ -392,28 +402,31 @@ def _segment_t_start(self, block_index, seg_index):
 
     def _segment_t_stop(self, block_index, seg_index):
         t_stop = float(self._last_timestamps) / self._global_ssampling_rate
-        if hasattr(self, "_signal_length"):
-            for stream_index in self._signal_length.keys():
-                t_stop_sig = self._signal_length[stream_index] / self._sig_sampling_rate[stream_index]
+        if hasattr(self, "__stream_id_samples"):
+            for stream_id in self._stream_id_samples.keys():
+                t_stop_sig = self._stream_id_samples[stream_id] / self._stream_id_sampling_frequency[stream_id]
                 t_stop = max(t_stop, t_stop_sig)
         return t_stop
 
     def _get_signal_size(self, block_index, seg_index, stream_index):
-        return self._signal_length[stream_index]
+        stream_id = self._stream_index_to_stream_id[stream_index]
+        return self._stream_id_samples[stream_id]
 
     def _get_signal_t_start(self, block_index, seg_index, stream_index):
         return 0.0
 
     def _get_analogsignal_chunk(self, block_index, seg_index, i_start, i_stop, stream_index, channel_indexes):
+        signal_channels = self.header["signal_channels"]
+        signal_streams = self.header["signal_streams"]
+        stream_id = signal_streams[stream_index]["id"]
+
         if i_start is None:
             i_start = 0
         if i_stop is None:
-            i_stop = self._signal_length[stream_index]
+            i_stop = self._stream_id_samples[stream_id]
+
 
-        signal_channels = self.header["signal_channels"]
-        signal_streams = self.header["signal_streams"]
 
-        stream_id = signal_streams[stream_index]["id"]
         mask = signal_channels["stream_id"] == stream_id
         signal_channels = signal_channels[mask]
         if channel_indexes is not None:

neo/test/rawiotest/common_rawio_test.py

@@ -59,14 +59,15 @@ class BaseTestRawIO:
 
     local_test_dir = get_local_testing_data_folder()
 
-    def setUp(self):
+    @classmethod
+    def setUpClass(cls):
         """
-        Set up the test fixture.  This is run for every test
+        Set up the test fixture. This is run once before all tests.
         """
-        self.shortname = self.rawioclass.__name__.lower().replace("rawio", "")
+        cls.shortname = cls.rawioclass.__name__.lower().replace("rawio", "")
 
-        if HAVE_DATALAD and self.use_network:
-            for remote_path in self.entities_to_download:
+        if HAVE_DATALAD and cls.use_network:
+            for remote_path in cls.entities_to_download:
                 download_dataset(repo=repo_for_test, remote_path=remote_path)
         else:
             raise unittest.SkipTest("Requires datalad download of data from the web")

neo/test/rawiotest/test_maxwellrawio.py

@@ -1,5 +1,4 @@
 import unittest
-import os
 
 from neo.rawio.maxwellrawio import MaxwellRawIO, auto_install_maxwell_hdf5_compression_plugin
 from neo.test.rawiotest.common_rawio_test import BaseTestRawIO
@@ -9,6 +8,7 @@ class TestMaxwellRawIO(
     BaseTestRawIO,
     unittest.TestCase,
 ):
+
     rawioclass = MaxwellRawIO
     entities_to_download = ["maxwell"]
     entities_to_test = files_to_test = [
@@ -18,7 +18,6 @@ class TestMaxwellRawIO(
 
     def setUp(self):
         auto_install_maxwell_hdf5_compression_plugin(force_download=False)
-        BaseTestRawIO.setUp(self)
 
 
 if __name__ == "__main__":

neo/test/rawiotest/test_plexon2rawio.py

@@ -9,6 +9,7 @@
 
 from neo.test.rawiotest.common_rawio_test import BaseTestRawIO
 
+from numpy.testing import assert_equal
 
 try:
     from neo.rawio.plexon2rawio.pypl2 import pypl2lib
@@ -25,8 +26,45 @@ class TestPlexon2RawIO(
 ):
     rawioclass = Plexon2RawIO
     entities_to_download = ["plexon"]
-    entities_to_test = ["plexon/4chDemoPL2.pl2"]
+    entities_to_test = ["plexon/4chDemoPL2.pl2",
+                        "plexon/NC16FPSPKEVT_1m.pl2"
+                        ]
+    
 
+    def test_check_enabled_flags(self):
+        """
+        This test loads a 1-minute PL2 file with 16 channels' each
+        of field potential (FP), and spike (SPK) data. The channels 
+        cycle through 4 possible combinations of m_ChannelEnabled 
+        and m_ChannelRecordingEnabled - (True, True), (True, False), 
+        (False, True), and (False, False). With 16 channels for each 
+        source, each combination of flags occurs 4 times. Only the 
+        first combination (True, True) causes data to be recorded to
+        disk. Therefore, we expect the following channels to be loaded by
+        Neo: FP01, FP05, FP09, FP13, SPK01, SPK05, SPK09, and SPK13.
+        
+        Note: the file contains event (EVT) data as well. Although event
+        channel headers do contain m_ChannelEnabled and m_ChannelRecording-
+        Enabled flags, the UI for recording PL2 files does not expose any
+        controls by which these flags can be changed from (True, True). 
+        Therefore, no test for event channels is necessary here.
+        """
+
+        # Load data from NC16FPSPKEVT_1m.pl2, a 1-minute PL2 recording containing
+        # 16-channels' each of field potential (FP), spike (SPK), and event (EVT)
+        # data.
+        reader = Plexon2RawIO(filename=self.get_local_path("plexon/NC16FPSPKEVT_1m.pl2"))
+        reader.parse_header()
+
+        # Check that the names of the loaded signal channels match what we expect
+        signal_channel_names = reader.header["signal_channels"]["name"].tolist()
+        expected_signal_channel_names = ["FP01","FP05","FP09","FP13"]
+        assert_equal(signal_channel_names, expected_signal_channel_names)
+
+        # Check that the names of the loaded spike channels match what we expect
+        spike_channel_names = reader.header["spike_channels"]["name"].tolist()
+        expected_spike_channel_names = ["SPK01.0", "SPK05.0", "SPK09.0", "SPK13.0"]
+        assert_equal(spike_channel_names, expected_spike_channel_names)
 
 if __name__ == "__main__":
     unittest.main()