Repository: SpikeInterface/spikeinterface
Pull request #1008 (2022-10-14 10:53:18)

Fixes to plot_timeseries
- Fix bug of ordering channels **before** getting traces with channel ids
- Added `return_scaled` option in `plot_timeseries`
- Simplified `order_channels_by_depth` by first oreding x and then y (instead of ordering by distance to an arbitrary origin)
- Fixed channel selection in ipywidgets when channels are ordered
- Flipped y when order by depth both for line and map mode: now deepest channel is at the bottom
- Set default cmap for map mode to more  intuitive `RdBu_r` (from `RdBu`): blue is negative, red is positive
- For sortingview backend, do not display timeseries in jupyter (since it doesn't display)


Patch:
diff --git a/spikeinterface/core/__init__.py b/spikeinterface/core/__init__.py
index b56db747a..459045b60 100644
--- a/spikeinterface/core/__init__.py
+++ b/spikeinterface/core/__init__.py
@@ -58,7 +58,7 @@ from .core_tools import write_binary_recording, write_to_h5_dataset_format, writ
     write_python
 from .job_tools import ensure_n_jobs, ensure_chunk_size, ChunkRecordingExecutor
 from .recording_tools import (get_random_data_chunks, get_channel_distances, get_closest_channels, 
-                              get_noise_levels, get_chunk_with_margin)
+                              get_noise_levels, get_chunk_with_margin, order_channels_by_depth)
 from .waveform_tools import extract_waveforms_to_buffers
 from .snippets_tools import snippets_from_sorting
 
diff --git a/spikeinterface/core/recording_tools.py b/spikeinterface/core/recording_tools.py
index 48cf999c0..28e450492 100644
--- a/spikeinterface/core/recording_tools.py
+++ b/spikeinterface/core/recording_tools.py
@@ -200,3 +200,28 @@ def get_chunk_with_margin(rec_segment, start_frame, end_frame,
             right_margin = margin
 
     return traces_chunk, left_margin, right_margin
+
+
+def order_channels_by_depth(recording, channel_ids=None):
+    """
+    Order channels by depth, by first ordering the x-axis, and then the y-axis.
+
+    Parameters
+    ----------
+    recording : BaseRecording
+        The input recording
+    channel_ids : list/array or None
+        If given, a subset of channels to order locations for
+
+    Returns
+    -------
+    order : np.array
+        Array with sorted indices
+    """
+    locations = recording.get_channel_locations()
+    channel_inds = recording.ids_to_indices(channel_ids)
+    locations = locations[channel_inds, :]
+
+    order = np.lexsort((locations[:, 0], locations[:, 1]))
+
+    return order
diff --git a/spikeinterface/widgets/ipywidgets/timeseries.py b/spikeinterface/widgets/ipywidgets/timeseries.py
index 6d29aba6b..002aed2e5 100644
--- a/spikeinterface/widgets/ipywidgets/timeseries.py
+++ b/spikeinterface/widgets/ipywidgets/timeseries.py
@@ -3,12 +3,12 @@ import numpy as np
 import matplotlib.pyplot as plt
 import ipywidgets.widgets as widgets
 
+from ...core import order_channels_by_depth
 
 from .base_ipywidgets import IpywidgetsPlotter
 from .utils import make_timeseries_controller, make_channel_controller, make_scale_controller
 
 from ..timeseries import TimeseriesWidget, _get_trace_list
-from ..utils import order_channels_by_depth
 from ..matplotlib.timeseries import TimeseriesPlotter as MplTimeseriesPlotter
 
 from IPython.display import display
@@ -32,7 +32,7 @@ class TimeseriesPlotter(IpywidgetsPlotter):
         with plt.ioff():
             output = widgets.Output()
             with output:
-                fig, ax = plt.subplots(figsize=(ratios[1] * width_cm * cm, height_cm * cm))
+                fig, ax = plt.subplots(figsize=(0.9 * ratios[1] * width_cm * cm, height_cm * cm))
                 plt.show()
 
         t_start = 0.
@@ -89,6 +89,7 @@ class PlotUpdater:
         self.controller = controller
 
         self.recordings = data_plot['recordings']
+        self.return_scaled = data_plot['return_scaled']
         self.next_data_plot = data_plot.copy()
         self.list_traces = None
 
@@ -153,13 +154,16 @@ class PlotUpdater:
         data_plot = self.next_data_plot
         
         if retrieve_traces:
-            channel_ids = self.recordings[list(self.recordings.keys())[0]].channel_ids[channel_indices]
+            all_channel_ids = self.recordings[list(self.recordings.keys())[0]].channel_ids
+            if self.data_plot["order"] is not None:
+                all_channel_ids = all_channel_ids[self.data_plot["order"]]
+            channel_ids = all_channel_ids[channel_indices]
             if self.data_plot['order_channel_by_depth']:
                 order = order_channels_by_depth(self.rec0, channel_ids)
             else:
                 order = None
             times, list_traces, frame_range, channel_ids = _get_trace_list(self.recordings, channel_ids, time_range,
-                                                                           segment_index, order)
+                                                                           segment_index, order, self.return_scaled)
             self.list_traces = list_traces
         else:
             times = data_plot['times']
diff --git a/spikeinterface/widgets/matplotlib/timeseries.py b/spikeinterface/widgets/matplotlib/timeseries.py
index 78a9dbaa1..27cf9aaa9 100644
--- a/spikeinterface/widgets/matplotlib/timeseries.py
+++ b/spikeinterface/widgets/matplotlib/timeseries.py
@@ -27,14 +27,14 @@ class TimeseriesPlotter(MplPlotter):
 
             for layer_key, traces in zip(dp.layer_keys, dp.list_traces):
                 for i, chan_id in enumerate(dp.channel_ids):
-                    offset = dp.vspacing * (n - 1 - i)
+                    offset = dp.vspacing * i
                     color = dp.colors[layer_key][chan_id]
                     ax.plot(dp.times, offset + traces[:, i], color=color)
                 ax.get_lines()[-1].set_label(layer_key)
 
             if dp.show_channel_ids:
                 ax.set_yticks(np.arange(n) * dp.vspacing)
-                channel_labels = np.array([str(chan_id) for chan_id in dp.channel_ids])[::-1]
+                channel_labels = np.array([str(chan_id) for chan_id in dp.channel_ids])
                 ax.set_yticklabels(channel_labels)
             ax.set_xlim(*dp.time_range)
             ax.set_ylim(-dp.vspacing, dp.vspacing * n)
@@ -49,7 +49,7 @@ class TimeseriesPlotter(MplPlotter):
             clim = list(dp.clims.values())[0]
             extent = (dp.time_range[0], dp.time_range[1], min_y, max_y)
             im = ax.imshow(dp.list_traces[0].T, interpolation='nearest',
-                           origin='upper', aspect='auto', extent=extent, cmap=dp.cmap)
+                           origin='lower', aspect='auto', extent=extent, cmap=dp.cmap)
 
             im.set_clim(*clim)
 
@@ -58,7 +58,7 @@ class TimeseriesPlotter(MplPlotter):
 
             if dp.show_channel_ids:
                 ax.set_yticks(np.linspace(min_y, max_y, n) + (max_y - min_y) / n * 0.5)
-                channel_labels = np.array([str(chan_id) for chan_id in dp.channel_ids])[::-1]
+                channel_labels = np.array([str(chan_id) for chan_id in dp.channel_ids])
                 ax.set_yticklabels(channel_labels)
 
 TimeseriesPlotter.register(TimeseriesWidget)
diff --git a/spikeinterface/widgets/sortingview/timeseries.py b/spikeinterface/widgets/sortingview/timeseries.py
index 68e2a31a8..95a83b562 100644
--- a/spikeinterface/widgets/sortingview/timeseries.py
+++ b/spikeinterface/widgets/sortingview/timeseries.py
@@ -45,6 +45,8 @@ class TimeseriesPlotter(SortingviewPlotter):
 
         self.set_view(view_ts)
 
+        # timeseries currently doesn't display on the jupyter backend
+        backend_kwargs["display"] = False
         self.handle_display_and_url(view_ts, **backend_kwargs)
         return view_ts
 
diff --git a/spikeinterface/widgets/timeseries.py b/spikeinterface/widgets/timeseries.py
index 14f0070ab..93ca48578 100644
--- a/spikeinterface/widgets/timeseries.py
+++ b/spikeinterface/widgets/timeseries.py
@@ -1,8 +1,8 @@
 import numpy as np
 
-from ..core import BaseRecording
+from ..core import BaseRecording, order_channels_by_depth
 from .base import BaseWidget
-from .utils import get_some_colors, order_channels_by_depth
+from .utils import get_some_colors
 
 
 class TimeseriesWidget(BaseWidget):
@@ -28,6 +28,8 @@ class TimeseriesWidget(BaseWidget):
             * 'line' : classical for low channel count
             * 'map' : for high channel count use color heat map
             * 'auto' : auto switch depending the channel count <32ch
+    return_scaled: bool
+        If True and recording.has_scaled(), it plots the scaled traces. Default False
     cmap: str default 'RdBu'
         matplotlib colormap used in mode 'map'
     show_channel_ids: bool
@@ -55,7 +57,7 @@ class TimeseriesWidget(BaseWidget):
     
 
     def __init__(self, recording, segment_index=None, channel_ids=None, order_channel_by_depth=False,
-                 time_range=None, mode='auto', cmap='RdBu', show_channel_ids=False,
+                 time_range=None, mode='auto', return_scaled=False, cmap='RdBu_r', show_channel_ids=False,
                  color_groups=False, color=None, clim=None, tile_size=512, seconds_per_row=0.2, 
                  with_colorbar=True, add_legend=True, backend=None, **backend_kwargs):
         if isinstance(recording, BaseRecording):
@@ -107,7 +109,7 @@ class TimeseriesWidget(BaseWidget):
         cmap = cmap
         
         times, list_traces, frame_range, channel_ids = _get_trace_list(recordings, channel_ids, time_range, 
-                                                                       segment_index, order)
+                                                                       segment_index, order, return_scaled)
         
         # stat for auto scaling done on the first layer
         traces0 = list_traces[0]
@@ -185,39 +187,45 @@ class TimeseriesWidget(BaseWidget):
             order_channel_by_depth=order_channel_by_depth,
             order=order,
             tile_size=tile_size,
-            num_timepoints_per_row=int(seconds_per_row * fs)
+            num_timepoints_per_row=int(seconds_per_row * fs),
+            return_scaled=return_scaled
         )
 
         BaseWidget.__init__(self, plot_data, backend=backend, **backend_kwargs)
 
 
-def _get_trace_list(recordings, channel_ids, time_range, segment_index, order=None):
+def _get_trace_list(recordings, channel_ids, time_range, segment_index, order=None,
+                    return_scaled=False):
     # function also used in ipywidgets plotter
     k0 = list(recordings.keys())[0]
     rec0 = recordings[k0]
 
     fs = rec0.get_sampling_frequency()
     
+    if return_scaled:
+        assert all(rec.has_scaled() for rec in recordings.values()), \
+            ("Some recording layers do not have scaled traces. Use `return_scaled=False`")
     frame_range = (time_range * fs).astype('int64')
     a_max = rec0.get_num_frames(segment_index=segment_index)
     frame_range = np.clip(frame_range, 0, a_max)
     time_range = frame_range / fs
     times = np.arange(frame_range[0], frame_range[1]) / fs
 
-    if order is not None:
-        channel_ids = np.array(channel_ids)[order]
-
     list_traces = []
     for rec_name, rec in recordings.items():
         traces = rec.get_traces(
             segment_index=segment_index,
             channel_ids=channel_ids,
             start_frame=frame_range[0],
-            end_frame=frame_range[1]
+            end_frame=frame_range[1],
+            return_scaled=return_scaled
         )
 
         if order is not None:
             traces = traces[:, order]
         list_traces.append(traces)
 
+    if order is not None:
+        channel_ids = np.array(channel_ids)[order]
+
     return times, list_traces, frame_range, channel_ids
\ No newline at end of file
diff --git a/spikeinterface/widgets/utils.py b/spikeinterface/widgets/utils.py
index dc23e444a..b000447a3 100644
--- a/spikeinterface/widgets/utils.py
+++ b/spikeinterface/widgets/utils.py
@@ -80,7 +80,7 @@ def array_to_image(data,
     Useful for visualizing data before/after preprocessing
 
     Params
-    =======
+    ------
     data : np.array
         2D numpy array
     colormap : str 
@@ -95,7 +95,7 @@ def array_to_image(data,
         Ratio of row spacing to overall channel height
 
     Returns
-    ========
+    -------
     output_image : 3D numpy array
 
     """
@@ -122,7 +122,7 @@ def array_to_image(data,
     scaled_data[scaled_data < clim[0]] = clim[0]
     scaled_data[scaled_data > clim[1]] = clim[1]
     scaled_data = (scaled_data - clim[0]) / np.ptp(clim)
-    a = (cmap(scaled_data.T)[:, :, :3]*255).astype(np.uint8)  # colorize and convert to uint8
+    a = np.flip((cmap(scaled_data.T)[:, :, :3]*255).astype(np.uint8), axis=0)  # colorize and convert to uint8
 
     num_rows = int(np.ceil(num_timepoints / num_timepoints_per_row))
 
@@ -175,16 +175,3 @@ def array_to_image(data,
             output_image = np.frombuffer(image.tobytes(), dtype=np.uint8).reshape(output_image.shape)
 
     return output_image
-
-
-def order_channels_by_depth(recording, channel_ids):
-    import scipy.spatial
-    locations = recording.get_channel_locations()
-    channel_inds = recording.ids_to_indices(channel_ids)
-    locations = locations[channel_inds, :]
-    origin = np.array([np.max(locations[:, 0]), np.min(locations[:, 1])])[None, :]
-    dist = scipy.spatial.distance.cdist(locations, origin, metric='euclidean')
-    dist = dist[:, 0]
-    order = np.argsort(dist)
-
-    return order


Test patch:
diff --git a/spikeinterface/core/tests/test_recording_tools.py b/spikeinterface/core/tests/test_recording_tools.py
index 8a2af551f..296fc2dfb 100644
--- a/spikeinterface/core/tests/test_recording_tools.py
+++ b/spikeinterface/core/tests/test_recording_tools.py
@@ -3,7 +3,8 @@ import numpy as np
 from spikeinterface.core.testing_tools import generate_recording
 
 from spikeinterface.core.recording_tools import (get_random_data_chunks, get_chunk_with_margin,
-                                                 get_closest_channels, get_channel_distances, get_noise_levels)
+                                                 get_closest_channels, get_channel_distances, get_noise_levels,
+                                                 order_channels_by_depth)
 
 
 def test_get_random_data_chunks():
@@ -82,9 +83,27 @@ def test_get_chunk_with_margin():
     # fig, ax = plt.subplots()
     # ax.plot(traces_windowed[:, 0], color='r', ls='--')
     # plt.show()
-    
+
+
+def test_order_channels_by_depth():
+    rec = generate_recording(num_channels=10, sampling_frequency=1000., durations=[10.], 
+                             set_probe=False)
+    locations = np.zeros((10, 2))
+    locations[:, 1] = np.arange(10) // 2 * 50
+    locations[:, 0] = np.arange(10) % 2* 30
+
+    locations_copy = locations.copy()
+    locations_copy = locations_copy[::-1]
+    rec.set_channel_locations(locations_copy)
+
+    order = order_channels_by_depth(rec)
+
+    assert np.array_equal(locations, locations_copy[order])
+
+
 
 if __name__ == '__main__':
     # test_get_random_data_chunks()
-    test_get_closest_channels()
+    # test_get_closest_channels()
     # test_get_noise_levels()
+    test_order_channels_by_depth()