Merge branch 'master' into add_brw_4.x_sparse

Author: zm711

neo/rawio/blackrockrawio.py

@@ -276,8 +276,8 @@ def _parse_header(self):
             self.internal_unit_ids = []  # pair of chan['packet_id'], spikes['unit_class_nb']
             for i in range(len(self.__nev_ext_header[b"NEUEVWAV"])):
 
-                # electrode_id values are stored at uint16 which can overflow when 
-                # multiplying by 1000 below. We convert to a regular python int which 
+                # electrode_id values are stored at uint16 which can overflow when
+                # multiplying by 1000 below. We convert to a regular python int which
                 # won't overflow
                 channel_id = int(self.__nev_ext_header[b"NEUEVWAV"]["electrode_id"][i])
 

neo/rawio/intanrawio.py

@@ -263,7 +263,7 @@ def _parse_header(self):
                         t_start = copy.copy(file_mmap[0][0])
                     else:  # empty file
                         t_start = 0
-                    stream_prop = (info["sampling_rate"], n_packets, t_start)
+                    stream_prop = (float(info["sampling_rate"]), int(n_packets), float(t_start))
                     if stream_prop not in stream_props:
                         stream_props[stream_prop] = {"stream_id": len(stream_props), "filenames": [filename]}
                     else:

neo/rawio/neuralynxrawio/neuralynxrawio.py

@@ -227,7 +227,9 @@ def _parse_header(self):
             if num_ephy_channels % num_chan_per_chip == 0:
                 all_hw_chans = [int(schan.attrib["hwChan"]) for trode in sconf for schan in trode]
                 missing_hw_chans = set(range(num_ephy_channels)) - set(all_hw_chans)
-                channel_ids = self._produce_ephys_channel_ids(num_ephy_channels_xml, num_chan_per_chip, missing_hw_chans)
+                channel_ids = self._produce_ephys_channel_ids(
+                    num_ephy_channels_xml, num_chan_per_chip, missing_hw_chans
+                )
                 raw_channel_ids = False
             else:
                 raw_channel_ids = True

neo/rawio/spikegadgetsrawio.py

@@ -376,7 +376,6 @@ def get_segment_tuple(info):
     for info in info_list:
         info["seg_index"] = segment_tuple_to_segment_index[get_segment_tuple(info)]
 
-
     for info in info_list:
         # device_kind is imec, nidq
         if info.get("device_kind") == "imec":

neo/rawio/spikeglxrawio.py

@@ -165,7 +165,7 @@ def random_segment():
         rec_datetime=random_datetime(),
         **random_annotations(4),
     )
-    
+
     n_sigs = random.randint(0, 5)
     for i in range(n_sigs):
         seg.analogsignals.append(random_signal())

neo/test/generate_datasets.py

@@ -23,6 +23,9 @@ class TestIntanIO(
         "intan/intan_fpc_rhs_test_240329_091637/info.rhs",  # Format one-file-per-channel
         "intan/intan_fps_rhs_test_240329_091536/info.rhs",  # Format one-file-per-signal
         "intan/rhd_fpc_multistim_240514_082044/info.rhd",  # Multiple digital channels one-file-per-channel rhd
+        "intan/rhs_stim_data_single_file_format/intanTestFile.rhs", # header-attached rhs data with stimulus current 
+        "intan/test_fcs_dc_250327_154333/info.rhs", # this is an example of only having dc amp rather than amp files
+        #"intan/test_fpc_stim_250327_151617/info.rhs", # wrong files Heberto will fix
     ]
 
 

neo/test/iotest/test_intanio.py

@@ -21,6 +21,10 @@ class TestIntanRawIO(
         "intan/intan_fpc_rhs_test_240329_091637/info.rhs",  # Format one-file-per-channel
         "intan/intan_fps_rhs_test_240329_091536/info.rhs",  # Format one-file-per-signal
         "intan/rhd_fpc_multistim_240514_082044/info.rhd",  # Multiple digital channels one-file-per-channel rhd
+        "intan/rhs_stim_data_single_file_format/intanTestFile.rhs", # header-attached rhs data with stimulus current
+        "intan/test_fcs_dc_250327_154333/info.rhs", # this is an example of only having dc amp rather than amp files
+        "intan/test_fpc_stim_250327_151617/info.rhs", # wrong files names Heberto will fix naimgin in the future
+
     ]
 
     def test_annotations(self):
@@ -67,7 +71,7 @@ def test_annotations(self):
         )
         np.testing.assert_array_equal(signal_array_annotations["board_stream_num"][:10], [0, 0, 0, 0, 0, 0, 0, 0, 0, 0])
 
-    def test_correct_reading_one_file_per_channel(self):
+    def test_correct_reading_one_file_per_channel_amplifiers(self):
         "Issue: https://github.com/NeuralEnsemble/python-neo/issues/1599"
         # Test reading of one-file-per-channel format file. The channels should match the raw files
         file_path = Path(self.get_local_path("intan/intan_fpc_test_231117_052630/info.rhd"))
@@ -81,11 +85,134 @@ def test_correct_reading_one_file_per_channel(self):
         amplifier_file_paths = [path for path in folder_path.iterdir() if "amp" in path.name]
         channel_names = [path.name[4:-4] for path in amplifier_file_paths]
 
+        amplifier_stream_index = 0
         for channel_name, amplifier_file_path in zip(channel_names, amplifier_file_paths):
             data_raw = np.fromfile(amplifier_file_path, dtype=np.int16).squeeze()
-            data_from_neo = intan_reader.get_analogsignal_chunk(channel_ids=[channel_name], stream_index=0).squeeze()
+            data_from_neo = intan_reader.get_analogsignal_chunk(channel_ids=[channel_name], stream_index=amplifier_stream_index).squeeze()
             np.testing.assert_allclose(data_raw, data_from_neo)
 
+    def test_correct_reading_one_file_per_channel_rhs_stim(self):
+        "Zach request for testing that the channel order is correct"
+        # Test reading of one-file-per-channel format file. The channels should match the raw files
+        file_path = Path(self.get_local_path("intan/test_fpc_stim_250327_151617/info.rhs"))
+        intan_reader = IntanRawIO(filename=file_path)
+        intan_reader.parse_header()
+
+        # This should be the folder where the files of all the channels are stored
+        folder_path = file_path.parent
+
+        # The paths for the stim channels are stim-A-000.dat, stim-A-001.dat, stim-A-002.dat, 
+        # Whereas the ids are A-001_STIM, A-002_STIM, A-003_STIM, etc
+        stim_file_paths = [path for path in folder_path.iterdir() if "stim" in path.name]
+        channel_ids = [f"{p.stem[5:]}_STIM" for p in stim_file_paths]
+        
+        stim_stream_index  = 2
+        for channel_id, amplifier_file_path in zip(channel_ids, stim_file_paths):
+            data_raw = np.fromfile(amplifier_file_path, dtype=np.uint16)
+            decoded_data = intan_reader._decode_current_from_stim_data(data_raw, 0, data_raw.shape[0])
+            data_from_neo = intan_reader.get_analogsignal_chunk(channel_ids=[channel_id], stream_index=stim_stream_index).squeeze()
+            np.testing.assert_allclose(decoded_data, data_from_neo)
+
+
+    def test_correct_decoding_of_stimulus_current(self):
+        # See https://github.com/NeuralEnsemble/python-neo/pull/1660 for discussion
+        # See https://gin.g-node.org/NeuralEnsemble/ephy_testing_data/src/master/intan/README.md#rhs_stim_data_single_file_format 
+        # For a description of the data 
+        
+        file_path = Path(self.get_local_path("intan/rhs_stim_data_single_file_format/intanTestFile.rhs"))
+        intan_reader = IntanRawIO(filename=file_path)
+        intan_reader.parse_header()
+        
+        signal_streams = intan_reader.header['signal_streams']
+        stream_ids = signal_streams['id'].tolist()
+        stream_index = stream_ids.index('11')
+        sampling_rate = intan_reader.get_signal_sampling_rate(stream_index=stream_index)
+        sig_chunk = intan_reader.get_analogsignal_chunk(stream_index=stream_index, channel_ids=["D-016_STIM"])
+        final_stim = intan_reader.rescale_signal_raw_to_float(sig_chunk, stream_index=stream_index, channel_ids=["D-016_STIM"])
+
+        # This contains only the first pulse and I got this by visual inspection
+        data_to_test = final_stim[200:250]
+        
+        positive_pulse_size = np.max(data_to_test).item()
+        negative_pulse_size = np.min(data_to_test).item()
+
+        expected_value = 60 * 1e-6# 60 microamperes
+
+        # Assert is close float
+        assert np.isclose(positive_pulse_size, expected_value)
+        assert np.isclose(negative_pulse_size, -expected_value)
+
+        # Check that negative pulse is leading
+        argmin = np.argmin(data_to_test)
+        argmax = np.argmax(data_to_test)
+        assert argmin < argmax
+        
+        # Check that the negative pulse is 200 us long
+        negative_pulse_frames = np.where(data_to_test > 0)[0]
+        number_of_negative_frames = negative_pulse_frames.size
+        duration_of_negative_pulse = number_of_negative_frames / sampling_rate
+
+        expected_duration = 200 * 1e-6  # 400 microseconds / 2 
+        assert np.isclose(duration_of_negative_pulse, expected_duration)
+        
+        # Check that the positive pulse is 200 us long
+        positive_pulse_frames = np.where(data_to_test > 0)[0]
+        number_of_positive_frames = positive_pulse_frames.size
+        duration_of_positive_pulse = number_of_positive_frames / sampling_rate
+        expected_duration = 200 * 1e-6  # 400 microseconds / 2
+
+        assert np.isclose(duration_of_positive_pulse, expected_duration)
+
+
+    def test_correct_decoding_of_stimulus_current(self):
+        # See https://github.com/NeuralEnsemble/python-neo/pull/1660 for discussion
+        # See https://gin.g-node.org/NeuralEnsemble/ephy_testing_data/src/master/intan/README.md#rhs_stim_data_single_file_format 
+        # For a description of the data 
+        
+        file_path = Path(self.get_local_path("intan/rhs_stim_data_single_file_format/intanTestFile.rhs"))
+        intan_reader = IntanRawIO(filename=file_path)
+        intan_reader.parse_header()
+        
+        signal_streams = intan_reader.header['signal_streams']
+        stream_ids = signal_streams['id'].tolist()
+        stream_index = stream_ids.index('11')
+        sampling_rate = intan_reader.get_signal_sampling_rate(stream_index=stream_index)
+        sig_chunk = intan_reader.get_analogsignal_chunk(stream_index=stream_index, channel_ids=["D-016_STIM"])
+        final_stim = intan_reader.rescale_signal_raw_to_float(sig_chunk, stream_index=stream_index, channel_ids=["D-016_STIM"])
+
+        # This contains only the first pulse and I got this by visual inspection
+        data_to_test = final_stim[200:250]
+        
+        positive_pulse_size = np.max(data_to_test).item()
+        negative_pulse_size = np.min(data_to_test).item()
+
+        expected_value = 60 * 1e-6# 60 microamperes
+
+        # Assert is close float
+        assert np.isclose(positive_pulse_size, expected_value)
+        assert np.isclose(negative_pulse_size, -expected_value)
+
+        # Check that negative pulse is leading
+        argmin = np.argmin(data_to_test)
+        argmax = np.argmax(data_to_test)
+        assert argmin < argmax
+        
+        # Check that the negative pulse is 200 us long
+        negative_pulse_frames = np.where(data_to_test > 0)[0]
+        number_of_negative_frames = negative_pulse_frames.size
+        duration_of_negative_pulse = number_of_negative_frames / sampling_rate
+
+        expected_duration = 200 * 1e-6  # 400 microseconds / 2 
+        assert np.isclose(duration_of_negative_pulse, expected_duration)
+        
+        # Check that the positive pulse is 200 us long
+        positive_pulse_frames = np.where(data_to_test > 0)[0]
+        number_of_positive_frames = positive_pulse_frames.size
+        duration_of_positive_pulse = number_of_positive_frames / sampling_rate
+        expected_duration = 200 * 1e-6  # 400 microseconds / 2
+
+        assert np.isclose(duration_of_positive_pulse, expected_duration)
+
 
 if __name__ == "__main__":
     unittest.main()

neo/test/rawiotest/test_intanrawio.py

@@ -21,11 +21,12 @@ class TestSpikeGadgetsRawIO(
     def test_parse_header_missing_channels(self):
 
         file_path = Path(self.get_local_path("spikegadgets/SL18_D19_S01_F01_BOX_SLP_20230503_112642_stubbed.rec"))
-        reader = SpikeGadgetsRawIO(filename = file_path)
+        reader = SpikeGadgetsRawIO(filename=file_path)
         reader.parse_header()
 
         assert_array_equal(
-            reader.header['signal_channels']['id'],
+            reader.header["signal_channels"]["id"],
+            # fmt: off
             [
                 'ECU_Ain1', 'ECU_Ain2', 'ECU_Ain3', 'ECU_Ain4', 'ECU_Ain5', 'ECU_Ain6',
                 'ECU_Ain7', 'ECU_Ain8', 'ECU_Aout1', 'ECU_Aout2', 'ECU_Aout3', 'ECU_Aout4', '0',
@@ -46,6 +47,7 @@ def test_parse_header_missing_channels(self):
                 '218', '250', '27', '59', '91', '123', '155', '187', '219', '251', '28', '60', '92',
                 '156', '188', '220', '252', '29', '61', '93', '125', '157', '189', '221', '253', '30',
                 '62', '94', '158', '190', '222', '254', '31', '63', '95', '127', '159', '191', '223',
-                '255',
+                '255'
             ]
+            # fmt: on
         )