Merge pull request #949 from samuelgarcia/rawio_improvement

Rawio improvement

Author: JuliaSprenger

doc/source/rawio.rst

@@ -77,11 +77,11 @@ Then browse the internal header and display information::
     nb_block: 1
     nb_segment:  [1]
     signal_channels: [V1]
-    unit_channels: [Wspk1u, Wspk2u, Wspk4u, Wspk5u ... Wspk29u Wspk30u Wspk31u Wspk32u]
+    spike_channels: [Wspk1u, Wspk2u, Wspk4u, Wspk5u ... Wspk29u Wspk30u Wspk31u Wspk32u]
     event_channels: []
 
 You get the number of blocks and segments per block. You have information
-about channels: **signal_channels**, **unit_channels**, **event_channels**.
+about channels: **signal_channels**, **spike_channels**, **event_channels**.
 
 All this information is internally available in the *header* dict::
 
@@ -91,7 +91,7 @@ All this information is internally available in the *header* dict::
     event_channels []
     nb_segment [1]
     nb_block 1
-    unit_channels [('Wspk1u', 'ch1#0', '',  0.00146484,  0., 0,  30000.)
+    spike_channels [('Wspk1u', 'ch1#0', '',  0.00146484,  0., 0,  30000.)
     ('Wspk2u', 'ch2#0', '',  0.00146484,  0., 0,  30000.)
     ...
 
@@ -141,7 +141,7 @@ Inspect units channel. Each channel gives a SpikeTrain for each Segment.
 Note that for many formats a physical channel can have several units after spike
 sorting. So the nb_unit could be more than physical channel or signal channels.
 
-    >>> nb_unit = reader.unit_channels_count()
+    >>> nb_unit = reader.spike_channels_count()
     >>> print('nb_unit', nb_unit)
     nb_unit 30
     >>> for unit_index in range(nb_unit):

examples/read_files_neo_rawio.py

@@ -31,8 +31,8 @@
 print(float_sigs.shape, float_sigs.dtype)
 print(sampling_rate, t_start, units)
 
-# Count unit and spike per units
-nb_unit = reader.unit_channels_count()
+# Count units and spikes per unit
+nb_unit = reader.spike_channels_count()
 print('nb_unit', nb_unit)
 for unit_index in range(nb_unit):
     nb_spike = reader.spike_count(block_index=0, seg_index=0, unit_index=unit_index)

neo/io/basefromrawio.py

@@ -134,31 +134,24 @@ def read_block(self, block_index=0, lazy=False,
 
         bl = Block(**bl_annotations)
 
-        # Group for AnalogSignals
+        # Group for AnalogSignals coming from signal_streams
         if create_group_across_segment['AnalogSignal']:
-            all_channels = self.header['signal_channels']
-            channel_indexes_list = self.get_group_signal_channel_indexes()
-            sig_groups = []
-            for channel_index in channel_indexes_list:
-                for i, (ind_within, ind_abs) in self._make_signal_channel_subgroups(
-                        channel_index, signal_group_mode=signal_group_mode).items():
-                        group = Group(name='AnalogSignal group {}'.format(i))
-                        # @andrew @ julia @michael : do we annotate group across segment with this arrays ?
-                        group.annotate(ch_names=all_channels[ind_abs]['name'].astype('U'))  # ??
-                        group.annotate(channel_ids=all_channels[ind_abs]['id'])  # ??
-                        bl.groups.append(group)
-                        sig_groups.append(group)
+            signal_streams = self.header['signal_streams']
+            sub_streams = self.get_sub_signal_streams(signal_group_mode)
+            sub_stream_groups = []
+            for sub_stream in sub_streams:
+                stream_index, inner_stream_channels, name = sub_stream
+                group = Group(name=name, stream_id=signal_streams[stream_index]['id'])
+                bl.groups.append(group)
+                sub_stream_groups.append(group)
 
         if create_group_across_segment['SpikeTrain']:
-            unit_channels = self.header['unit_channels']
+            spike_channels = self.header['spike_channels']
             st_groups = []
-            for c in range(unit_channels.size):
+            for c in range(spike_channels.size):
                 group = Group(name='SpikeTrain group {}'.format(c))
-                group.annotate(unit_name=unit_channels[c]['name'])
-                group.annotate(unit_id=unit_channels[c]['id'])
-                unit_annotations = self.raw_annotations['unit_channels'][c]
-                unit_annotations = check_annotations(unit_annotations)
-                group.annotate(**unit_annotations)
+                group.annotate(unit_name=spike_channels[c]['name'])
+                group.annotate(unit_id=spike_channels[c]['id'])
                 bl.groups.append(group)
                 st_groups.append(group)
 
@@ -183,7 +176,7 @@ def read_block(self, block_index=0, lazy=False,
         for seg in bl.segments:
             if create_group_across_segment['AnalogSignal']:
                 for c, anasig in enumerate(seg.analogsignals):
-                    sig_groups[c].add(anasig)
+                    sub_stream_groups[c].add(anasig)
 
             if create_group_across_segment['SpikeTrain']:
                 for c, sptr in enumerate(seg.spiketrains):
@@ -231,38 +224,35 @@ def read_segment(self, block_index=0, seg_index=0, lazy=False,
             signal_group_mode = self._prefered_signal_group_mode
 
         # annotations
-        seg_annotations = dict(self.raw_annotations['blocks'][block_index]['segments'][seg_index])
-        for k in ('signals', 'units', 'events'):
+        seg_annotations = self.raw_annotations['blocks'][block_index]['segments'][seg_index].copy()
+        for k in ('signals', 'spikes', 'events'):
             seg_annotations.pop(k)
         seg_annotations = check_annotations(seg_annotations)
 
         seg = Segment(index=seg_index, **seg_annotations)
 
         # AnalogSignal
-        signal_channels = self.header['signal_channels']
-        if signal_channels.size > 0:
-            channel_indexes_list = self.get_group_signal_channel_indexes()
-            for channel_indexes in channel_indexes_list:
-                for i, (ind_within, ind_abs) in self._make_signal_channel_subgroups(
-                        channel_indexes,
-                        signal_group_mode=signal_group_mode).items():
-                    # make a proxy...
-                    anasig = AnalogSignalProxy(rawio=self, global_channel_indexes=ind_abs,
-                                    block_index=block_index, seg_index=seg_index)
-
-                    if not lazy:
-                        # ... and get the real AnalogSIgnal if not lazy
-                        anasig = anasig.load(time_slice=time_slice, strict_slicing=strict_slicing)
-                        # TODO magnitude_mode='rescaled'/'raw'
-
-                    anasig.segment = seg
-                    seg.analogsignals.append(anasig)
+        signal_streams = self.header['signal_streams']
+        sub_streams = self.get_sub_signal_streams(signal_group_mode)
+        for sub_stream in sub_streams:
+            stream_index, inner_stream_channels, name = sub_stream
+            anasig = AnalogSignalProxy(rawio=self, stream_index=stream_index,
+                            inner_stream_channels=inner_stream_channels,
+                            block_index=block_index, seg_index=seg_index)
+            anasig.name = name
+
+            if not lazy:
+                # ... and get the real AnalogSignal if not lazy
+                anasig = anasig.load(time_slice=time_slice, strict_slicing=strict_slicing)
+
+            anasig.segment = seg
+            seg.analogsignals.append(anasig)
 
         # SpikeTrain and waveforms (optional)
-        unit_channels = self.header['unit_channels']
-        for unit_index in range(len(unit_channels)):
+        spike_channels = self.header['spike_channels']
+        for spike_channel_index in range(len(spike_channels)):
             # make a proxy...
-            sptr = SpikeTrainProxy(rawio=self, unit_index=unit_index,
+            sptr = SpikeTrainProxy(rawio=self, spike_channel_index=spike_channel_index,
                                                 block_index=block_index, seg_index=seg_index)
 
             if not lazy:
@@ -286,7 +276,7 @@ def read_segment(self, block_index=0, seg_index=0, lazy=False,
                 seg.events.append(e)
             elif event_channels['type'][chan_ind] == b'epoch':
                 e = EpochProxy(rawio=self, event_channel_index=chan_ind,
-                                        block_index=block_index, seg_index=seg_index)
+                               block_index=block_index, seg_index=seg_index)
                 if not lazy:
                     e = e.load(time_slice=time_slice, strict_slicing=strict_slicing)
                 e.segment = seg
@@ -295,37 +285,50 @@ def read_segment(self, block_index=0, seg_index=0, lazy=False,
         seg.create_many_to_one_relationship()
         return seg
 
-    def _make_signal_channel_subgroups(self, channel_indexes,
-                                       signal_group_mode='group-by-same-units'):
+    def get_sub_signal_streams(self, signal_group_mode='group-by-same-units'):
         """
-        For some RawIO channel are already splitted in groups.
-        But in any cases, channel need to be splitted again in sub groups
-        because they do not have the same units.
-
-        They can also be splitted one by one to match previous behavior for
-        some IOs in older version of neo (<=0.5).
+        When signal streams don't have homogeneous SI units across channels,
+        they have to be split in sub streams to construct AnalogSignal objects with unique units.
 
-        This method aggregate signal channels with same units or split them all.
+        For backward compatibility (neo version <= 0.5) sub-streams can also be 
+        used to generate one AnalogSignal per channel.
         """
-        all_channels = self.header['signal_channels']
-        if channel_indexes is None:
-            channel_indexes = np.arange(all_channels.size, dtype=int)
-        channels = all_channels[channel_indexes]
-
-        groups = collections.OrderedDict()
-        if signal_group_mode == 'group-by-same-units':
-            all_units = np.unique(channels['units'])
-
-            for i, unit in enumerate(all_units):
-                ind_within, = np.nonzero(channels['units'] == unit)
-                ind_abs = channel_indexes[ind_within]
-                groups[i] = (ind_within, ind_abs)
-
-        elif signal_group_mode == 'split-all':
-            for i, chan_index in enumerate(channel_indexes):
-                ind_within = [i]
-                ind_abs = channel_indexes[ind_within]
-                groups[i] = (ind_within, ind_abs)
-        else:
-            raise (NotImplementedError)
-        return groups
+        signal_streams = self.header['signal_streams']
+        signal_channels = self.header['signal_channels']
+
+        sub_streams = []
+        for stream_index in range(len(signal_streams)):
+            stream_id = signal_streams[stream_index]['id']
+            stream_name = signal_streams[stream_index]['name']
+            mask = signal_channels['stream_id'] == stream_id
+            channels = signal_channels[mask]
+            if signal_group_mode == 'group-by-same-units':
+                # this does not keep the original order
+                _, idx = np.unique(channels['units'], return_index=True)
+                all_units = channels['units'][np.sort(idx)]
+
+                if len(all_units) == 1:
+                    # no substream
+                    #  None iwill be transform as slice later
+                    inner_stream_channels = None
+                    name = stream_name
+                    sub_stream = (stream_index, inner_stream_channels, name)
+                    sub_streams.append(sub_stream)
+                else:
+                    for units in all_units:
+                        inner_stream_channels, = np.nonzero(channels['units'] == units)
+                        chan_names = channels[inner_stream_channels]['name']
+                        name = 'Channels: (' + ' '.join(chan_names) + ')'
+                        sub_stream = (stream_index, inner_stream_channels, name)
+                        sub_streams.append(sub_stream)
+            elif signal_group_mode == 'split-all':
+                # mimic all neo <= 0.5 behavior
+                for i, channel in enumerate(channels):
+                    inner_stream_channels = [i]
+                    name = channels[i]['name']
+                    sub_stream = (stream_index, inner_stream_channels, name)
+                    sub_streams.append(sub_stream)
+            else:
+                raise (NotImplementedError)
+
+        return sub_streams