Parse DIO data and add to NWB file

src/spikegadgets_to_nwb/convert.py

@@ -2,6 +2,7 @@
 
 import pandas as pd
 
+from spikegadgets_to_nwb.convert_dios import add_dios
 from spikegadgets_to_nwb.convert_ephys import add_raw_ephys
 from spikegadgets_to_nwb.convert_position import add_position
 from spikegadgets_to_nwb.convert_rec_header import (
@@ -14,7 +15,6 @@
     add_acquisition_devices,
     add_associated_files,
     add_cameras,
-    add_dios,
     add_electrode_groups,
     add_subject,
     add_tasks,

src/spikegadgets_to_nwb/convert_analog.py

@@ -38,10 +38,17 @@ def add_analog_data(nwbfile: NWBFile, rec_file_path: list[str], **kwargs) -> Non
             analog_channel_ids.append(channel.attrib["id"])
 
     # make the data chunk iterator
+    # TODO use the stream name instead of the stream index to be more robust
     rec_dci = RecFileDataChunkIterator(
-        rec_file_path, nwb_hw_channel_order=analog_channel_ids, stream_index=0
+        rec_file_path,
+        nwb_hw_channel_order=analog_channel_ids,
+        stream_index=2,
+        is_analog=True,
     )
 
+    # add headstage channel IDs to the list of analog channel IDs
+    analog_channel_ids.extend(rec_dci.neo_io[0].multiplexed_channel_xml.keys())
+
     # (16384, 32) chunks of dtype int16 (2 bytes) is 1 MB, which is recommended
     # by studies by the NWB team.
     # could also add compression here. zstd/blosc-zstd are recommended by the NWB team, but

src/spikegadgets_to_nwb/convert_dios.py

@@ -0,0 +1,83 @@
+import numpy as np
+from pynwb import NWBFile, TimeSeries
+from pynwb.behavior import BehavioralEvents
+
+from .spike_gadgets_raw_io import SpikeGadgetsRawIO
+
+
+def _get_channel_name_map(metadata: dict) -> dict[str, str]:
+    """Parses behavioral events metadata from the yaml file
+
+    Parameters
+    ----------
+    metadata : dict
+        metadata from the yaml generator
+
+    Returns
+    -------
+    channel_name_map : dict
+        Parsed behavioral events metadata mapping hardware event name to human-readable name
+    """
+    dio_metadata = metadata["behavioral_events"]
+    channel_name_map = {}
+    for dio_event in dio_metadata:
+        channel_name_map[dio_event["description"]] = dio_event["name"]
+    return channel_name_map
+
+
+def add_dios(nwbfile: NWBFile, recfile: list[str], metadata: dict) -> None:
+    """Adds DIO event information and data to nwb file
+
+    Parameters
+    ----------
+    nwbfile : NWBFile
+        nwb file being assembled
+    recfile : list[str]
+        list of paths to rec files
+    metadata : dict
+        metadata from the yaml generator
+    """
+
+    # TODO remove redundancy with convert_ephys.py
+    neo_io = [
+        SpikeGadgetsRawIO(filename=file) for file in recfile
+    ]  # get all streams for all files
+    [neo_io.parse_header() for neo_io in neo_io]
+
+    # Make a processing module for behavior and add to the nwbfile
+    if not "behavior" in nwbfile.processing:
+        nwbfile.create_processing_module(
+            name="behavior", description="Contains all behavior-related data"
+        )
+
+    # Make BehavioralEvents object to hold DIO data
+    beh_events = BehavioralEvents(name="behavioral_events")
+
+    # Map hardware event name (encoded in `description` in metadata YAML)
+    # to a human-readable name (encoded in `name`)
+    channel_name_map = _get_channel_name_map(metadata)
+
+    # Loop through the channels from the metadata YAML and add a TimeSeries for each one
+    stream_name = "ECU_digital"
+    for channel_name in channel_name_map:
+        # merge streams from multiple files
+        all_timestamps = np.array([], dtype=np.float64)
+        all_state_changes = np.array([], dtype=np.uint8)
+        for io in neo_io:
+            timestamps, state_changes = io.get_digitalsignal(
+                stream_name, "ECU_" + channel_name
+            )
+            all_timestamps = np.concatenate((all_timestamps, timestamps))
+            all_state_changes = np.concatenate((all_state_changes, state_changes))
+
+        ts = TimeSeries(
+            name=channel_name_map[channel_name],
+            description=channel_name,
+            data=all_state_changes,
+            unit="-1",  # TODO change to "N/A",
+            timestamps=all_timestamps,  # TODO adjust timestamps
+        )
+        beh_events.add_timeseries(ts)
+
+    # Add the BehavioralEvents object to the file
+    nwbfile.processing["behavior"].add(beh_events)

src/spikegadgets_to_nwb/convert_ephys.py

@@ -21,20 +21,26 @@ def __init__(
         rec_file_path: list[str],
         nwb_hw_channel_order=[],
         conversion: float = 1.0,
-        stream_index: int = 1,
+        stream_index: int = 3,  # TODO use the stream name instead of the index
+        is_analog: bool = False,
         **kwargs,
     ):
         self.conversion = conversion
+        self.is_analog = is_analog
         self.neo_io = [
             SpikeGadgetsRawIO(filename=file) for file in rec_file_path
         ]  # get all streams for all files
         [neo_io.parse_header() for neo_io in self.neo_io]
         # TODO see what else spikeinterface does and whether it is necessary
 
-        # for now, make sure that there is only one block, one segment, and two streams
+        # for now, make sure that there is only one block, one segment, and four streams:
+        # Controller_DIO_digital
+        # ECU_digital
+        # ECU_analog
+        # trodes
         assert all([neo_io.block_count() == 1 for neo_io in self.neo_io])
         assert all([neo_io.segment_count(0) == 1 for neo_io in self.neo_io])
-        assert all([neo_io.signal_streams_count() == 2 for neo_io in self.neo_io])
+        assert all([neo_io.signal_streams_count() == 4 for neo_io in self.neo_io])
 
         self.block_index = 0
         self.seg_index = 0
@@ -64,6 +70,9 @@ def __init__(
         self.n_channel = self.neo_io[0].signal_channels_count(
             stream_index=self.stream_index
         )
+        self.n_multiplexed_channel = 0
+        if self.is_analog:
+            self.n_multiplexed_channel += len(self.neo_io[0].multiplexed_channel_xml)
 
         # order that the hw channels are in within the nwb table
         if len(nwb_hw_channel_order) == 0:  # TODO: raise error instead?
@@ -93,11 +102,18 @@ def _get_data(self, selection: Tuple[slice]) -> np.ndarray:
         # DCI will want channels 0 to X first to put into the array in that order
         # those are stored in the file as channel IDs
         # make into list form passed to neo_io
-        channel_ids = [str(x) for x in self.nwb_hw_channel_order[selection[1]]]
+        selection_list = list(selection)
+        if self.is_analog:
+            selection_list[1] = slice(
+                selection[1].start,
+                min(selection[1].stop, self.n_channel),
+                selection[1].step,
+            )
+        channel_ids = [str(x) for x in self.nwb_hw_channel_order[selection_list[1]]]
         # what global index each file starts at
         file_start_ind = np.append(np.zeros(1), np.cumsum(self.n_time))
         # the time indexes we want
-        time_index = np.arange(self._get_maxshape()[0])[selection[0]]
+        time_index = np.arange(self._get_maxshape()[0])[selection_list[0]]
         data = []
         i = time_index[0]
         while i < min(time_index[-1], self._get_maxshape()[0]):
@@ -134,12 +150,13 @@ def _get_data(self, selection: Tuple[slice]) -> np.ndarray:
                 time_index[-1] - i,  # if finished in this stream
             )
         data = (np.array(data) * self.conversion).astype("int16")
+
         return data
 
     def _get_maxshape(self) -> Tuple[int, int]:
         return (
             np.sum(self.n_time),
-            self.n_channel,
+            self.n_channel + self.n_multiplexed_channel,
         )  # TODO: Is this right for maxshape @rly
 
     def _get_dtype(self) -> np.dtype:

src/spikegadgets_to_nwb/convert_yaml.py

@@ -365,43 +365,6 @@ def add_tasks(nwbfile: NWBFile, metadata: dict) -> None:
         nwbfile.processing["tasks"].add(task)
 
 
-def add_dios(nwbfile: NWBFile, metadata: dict) -> None:
-    """Adds DIO event information and data to nwb file
-
-    Parameters
-    ----------
-    nwbfile : NWBFile
-        nwb file being assembled
-    metadata : dict
-        metadata from the yaml generator
-    """
-    # TODO: pass the dio data and include in this
-    # Make a processing module for behavior and add to the nwbfile
-    if not "behavior" in nwbfile.processing:
-        nwbfile.create_processing_module(
-            name="behavior", description="Contains all behavior-related data"
-        )
-    # Make Behavioral events object to hold DIO data
-    events = BehavioralEvents(name="behavioral_events")
-    # Loop through and add timeseries for each one
-    dio_metadata = metadata["behavioral_events"]
-    for dio_event in dio_metadata:
-        events.add_timeseries(
-            TimeSeries(
-                name=dio_event["name"],
-                description=dio_event["description"],
-                data=np.array(
-                    []
-                ),  # TODO: from rec file // self.data[dio_event['description']],
-                unit="N/A",
-                timestamps=np.array([]),
-                # TODO: data, timestamps,
-            )
-        )
-    # add it to your file
-    nwbfile.processing["behavior"].add(events)
-
-
 def add_associated_files(nwbfile: NWBFile, metadata: dict) -> None:
     """Adds associated files processing module. Reads in file referenced in metadata and stores in processing