From 0546cfab3973e547d7e11783b53d94c8779aaffc Mon Sep 17 00:00:00 2001 From: Alessio Buccino Date: Mon, 28 Sep 2026 15:04:04 +0200 Subject: [PATCH 1/5] fix: docs --- examples/tutorials/core/plot_4_sorting_analyzer.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/examples/tutorials/core/plot_4_sorting_analyzer.py b/examples/tutorials/core/plot_4_sorting_analyzer.py index 7c89088c09..3587aec178 100644 --- a/examples/tutorials/core/plot_4_sorting_analyzer.py +++ b/examples/tutorials/core/plot_4_sorting_analyzer.py @@ -33,6 +33,8 @@ # First let's generate a simulated recording and sorting recording, sorting = generate_ground_truth_recording() +# let's make the recording JSON serializable +recording = recording.save(folder="recording_binary") print(recording) print(sorting) From e5df6213a7b4270a4762bc7662cdd12bcd2b8e4d Mon Sep 17 00:00:00 2001 From: Alessio Buccino Date: Wed, 7 Oct 2026 12:42:38 +0200 Subject: [PATCH 2/5] fix: group property propagation in channel slice --- src/spikeinterface/core/channelslice.py | 2 +- .../core/tests/test_channelslicerecording.py | 11 +++++++++++ 2 files changed, 12 insertions(+), 1 deletion(-) diff --git a/src/spikeinterface/core/channelslice.py b/src/spikeinterface/core/channelslice.py index 4630e0be24..bfb46da619 100644 --- a/src/spikeinterface/core/channelslice.py +++ b/src/spikeinterface/core/channelslice.py @@ -58,7 +58,6 @@ def __init__(self, recording, channel_ids=None, renamed_channel_ids=None): self.add_recording_segment(sub_segment) # copy annotation and properties - recording.copy_metadata(self, only_main=False, ids=self._channel_ids) self._parent = recording # change the wiring of the probe @@ -67,6 +66,7 @@ def __init__(self, recording, channel_ids=None, renamed_channel_ids=None): sliced_probegroup = parent_probegroup.get_slice(self._parent_channel_indices) sliced_probegroup.set_global_device_channel_indices(np.arange(len(self._channel_ids))) self.set_probegroup(sliced_probegroup) + recording.copy_metadata(self, only_main=False, ids=self._channel_ids) # update dump dict self._kwargs = { diff --git a/src/spikeinterface/core/tests/test_channelslicerecording.py b/src/spikeinterface/core/tests/test_channelslicerecording.py index c563888e07..da255d6953 100644 --- a/src/spikeinterface/core/tests/test_channelslicerecording.py +++ b/src/spikeinterface/core/tests/test_channelslicerecording.py @@ -1,6 +1,7 @@ import shutil from pathlib import Path +from numpy import rec import pytest import numpy as np @@ -129,5 +130,15 @@ def test_select_channels_preserves_probe_metadata(): assert probes[0].annotations.get("manufacturer") == "vendor_Y" +def test_channel_slice_retain_original_groups(): + rec = generate_recording(num_channels=8, durations=[1.0]) # has a probe attached + rec.set_property("group", [0, 0, 1, 1, 2, 2, 3, 3]) + for g, sub in rec.split_by("group").items(): + assert np.all(sub.get_property("group") == g) + + rec_sliced = rec.remove_channels(remove_channel_ids=["1", "3", "5", "7"]) + assert np.all(rec_sliced.get_property("group") == [0, 1, 2, 3]) + + if __name__ == "__main__": test_ChannelSliceRecording() From cc4999da72069e0ec0d287305a16c54becf56b14 Mon Sep 17 00:00:00 2001 From: Alessio Buccino Date: Wed, 7 Oct 2026 14:49:48 +0200 Subject: [PATCH 3/5] feat: destroy snippets --- doc/api.rst | 5 - doc/get_started/install_sorters.rst | 1 - doc/get_started/quickstart.rst | 2 +- doc/modules/core.rst | 71 +- doc/modules/preprocessing.rst | 15 - src/spikeinterface/core/__init__.py | 9 +- src/spikeinterface/core/baserecording.py | 683 ++++++++++++++++-- .../core/baserecordingsnippets.py | 637 ---------------- src/spikeinterface/core/basesnippets.py | 294 -------- src/spikeinterface/core/channelslice.py | 116 +-- src/spikeinterface/core/generate.py | 75 -- src/spikeinterface/core/npyfoldersnippets.py | 103 --- .../core/npysnippetsextractor.py | 166 ----- src/spikeinterface/core/numpyextractors.py | 136 ---- src/spikeinterface/core/snippets_tools.py | 66 -- .../core/tests/test_basesnippets.py | 198 ----- .../core/tests/test_channelslicerecording.py | 4 - .../core/tests/test_npyfoldersnippets.py | 33 - .../core/tests/test_npysnippetsextractor.py | 25 - .../extractors/extractor_classes.py | 16 - .../extractors/waveclussnippetstextractors.py | 154 ---- src/spikeinterface/generation/__init__.py | 1 - src/spikeinterface/preprocessing/__init__.py | 3 - .../preprocessing/align_snippets.py | 113 --- .../tests/test_align_snippets.py | 39 - .../external/tests/test_waveclus_snippets.py | 24 - .../sorters/external/waveclus_snippets.py | 217 ------ .../external/waveclus_snippets_master.m | 78 -- src/spikeinterface/sorters/runsorter.py | 1 - src/spikeinterface/sorters/sorterlist.py | 2 - .../sorters/tests/common_tests.py | 55 -- 31 files changed, 639 insertions(+), 2703 deletions(-) delete mode 100644 src/spikeinterface/core/baserecordingsnippets.py delete mode 100644 src/spikeinterface/core/basesnippets.py delete mode 100644 src/spikeinterface/core/npyfoldersnippets.py delete mode 100644 src/spikeinterface/core/npysnippetsextractor.py delete mode 100644 src/spikeinterface/core/snippets_tools.py delete mode 100644 src/spikeinterface/core/tests/test_basesnippets.py delete mode 100644 src/spikeinterface/core/tests/test_npyfoldersnippets.py delete mode 100644 src/spikeinterface/core/tests/test_npysnippetsextractor.py delete mode 100644 src/spikeinterface/extractors/waveclussnippetstextractors.py delete mode 100644 src/spikeinterface/preprocessing/align_snippets.py delete mode 100644 src/spikeinterface/preprocessing/tests/test_align_snippets.py delete mode 100644 src/spikeinterface/sorters/external/tests/test_waveclus_snippets.py delete mode 100644 src/spikeinterface/sorters/external/waveclus_snippets.py delete mode 100644 src/spikeinterface/sorters/external/waveclus_snippets_master.m diff --git a/doc/api.rst b/doc/api.rst index c0ace59d3b..9e8e5ef2c7 100755 --- a/doc/api.rst +++ b/doc/api.rst @@ -10,8 +10,6 @@ spikeinterface.core :members: .. autoclass:: BaseSorting :members: - .. autoclass:: BaseSnippets - :members: .. autoclass:: BaseEvent :members: .. autoclass:: SortingAnalyzer @@ -27,10 +25,8 @@ spikeinterface.core .. autoclass:: BinaryRecordingExtractor .. autoclass:: BinaryFolderRecording .. autoclass:: NumpyFolderSorting - .. autoclass:: NpyFolderSnippets .. autoclass:: NumpyRecording .. autoclass:: NumpySorting - .. autoclass:: NumpySnippets .. autoclass:: AppendSegmentRecording .. autoclass:: ConcatenateSegmentRecording .. autoclass:: SelectSegmentRecording @@ -422,7 +418,6 @@ Core .. autofunction:: generate_recording .. autofunction:: generate_sorting - .. autofunction:: generate_snippets .. autofunction:: generate_templates .. autofunction:: generate_recording_by_size .. autofunction:: generate_ground_truth_recording diff --git a/doc/get_started/install_sorters.rst b/doc/get_started/install_sorters.rst index fa3d883da9..75de4921f3 100644 --- a/doc/get_started/install_sorters.rst +++ b/doc/get_started/install_sorters.rst @@ -239,7 +239,6 @@ Waveclus ^^^^^^^^ * Matlab -* Also supports Snippets (waveform cutouts) objects (:py:class:`~spikeinterface.core.BaseSnippets`) * Url: https://github.com/csn-le/wave_clus/wiki * Authors: Fernando Chaure, Hernan Rey and Rodrigo Quian Quiroga * Installation requires Matlab:: diff --git a/doc/get_started/quickstart.rst b/doc/get_started/quickstart.rst index 3fbd113ba7..3d36bb45e5 100644 --- a/doc/get_started/quickstart.rst +++ b/doc/get_started/quickstart.rst @@ -262,7 +262,7 @@ installed .. parsed-literal:: - Available sorters ['combinato', 'hdsort', 'herdingspikes', 'ironclust', 'kilosort', 'kilosort2', 'kilosort2_5', 'kilosort3', 'kilosort4', 'klusta', 'mountainsort4', 'mountainsort5', 'pykilosort', 'simple', 'spykingcircus', 'spykingcircus2', 'tridesclous', 'tridesclous2', 'waveclus', 'waveclus_snippets', 'yass'] + Available sorters ['combinato', 'hdsort', 'herdingspikes', 'ironclust', 'kilosort', 'kilosort2', 'kilosort2_5', 'kilosort3', 'kilosort4', 'klusta', 'mountainsort4', 'mountainsort5', 'pykilosort', 'simple', 'spykingcircus', 'spykingcircus2', 'tridesclous', 'tridesclous2', 'waveclus', 'yass'] Installed sorters ['mountainsort4', 'simple', 'spykingcircus2', 'tridesclous', 'tridesclous2'] diff --git a/doc/modules/core.rst b/doc/modules/core.rst index 297f22e945..40b6b090ee 100644 --- a/doc/modules/core.rst +++ b/doc/modules/core.rst @@ -11,9 +11,7 @@ The core classes are: :py:class:`~spikeinterface.core.BaseRecording` (for raw da :py:class:`~spikeinterface.core.BaseSorting` (for spike-sorted data), and :py:class:`~spikeinterface.core.SortingAnalyzer` (for postprocessing, quality metrics, and waveform extraction). -There are additional classes to allow to retrieve events (:py:class:`~spikeinterface.core.BaseEvent`) and to -handle unsorted waveform cutouts, or *snippets*, which are recorded by some acquisition systems -(:py:class:`~spikeinterface.core.BaseSnippets`). +In addition, this module includes a base class to retrieve events (:py:class:`~spikeinterface.core.BaseEvent`). All classes support: @@ -483,52 +481,6 @@ with events from two channels: ### NOTE ### # 'segment_index' is required for multi-segment objects - -Snippets --------- - -The :py:class:`~spikeinterface.core.BaseSnippets` class serves as basis for all :code:`Snippets` -classes (currently only :py:class:`~spikeinterface.core.NumpySnippets` and -:code:`WaveClusSnippetsExtractor` are implemented). - -It represents unsorted waveform cutouts. Some acquisition systems, in fact, allow users to set a -threshold and only record the times at which a peak was detected and the waveform cut out around -the peak. - -.. note:: - - While we support this class (mainly for legacy formats), this approach is a bad practice - and is highly discouraged! Most modern spike sorters, in fact, require the raw traces to perform - template matching to recover spikes! - -Here we assume :code:`snippets` is a :py:class:`~spikeinterface.core.BaseSnippets` object -with 16 channels: - -.. code-block:: python - - channel_ids = snippets.channel_ids - num_channels = snippets.get_num_channels() - # retrieve number of snippets - num_snippets = snippets.get_num_snippets(segment_index=0) - ### NOTE ### - # 'segment_index' is required for multi-segment objects - # retrieve total number of snippets across segments - total_snippets = snippets.get_total_snippets() - - # retrieve snippet size - nbefore = snippets.nbefore # samples before peak - nsamples_per_snippet = snippets.snippet_len # total - nafter = nsamples_per_snippet - nbefore # samples after peak - - # retrieve sample/frame indices - frames = snippets.get_frames(segment_index=0) - # retrieve snippet cutouts - snippet_cutouts = snippets.get_snippets(segment_index=0) - # retrieve snippet cutouts on first 4 channels - snippet_cutouts_slice = snippets.get_snippets(channel_ids=channel_ids[:4], - segment_index=0) - - Handling probes --------------- @@ -648,15 +600,14 @@ Saving, loading, and compression -------------------------------- The Base SpikeInterface objects (:py:class:`~spikeinterface.core.BaseRecording`, -:py:class:`~spikeinterface.core.BaseSorting`, and -:py:class:`~spikeinterface.core.BaseSnippets`) hold full information about their history to maintain provenance. +:py:class:`~spikeinterface.core.BaseSorting`) hold full information about their history to maintain provenance. Each object is in fact internally represented as a dictionary (:code:`si_object.to_dict()`) which can be used to re-instantiate the object from scratch (this is true for all objects except in-memory ones, see :ref:`in_memory`). The :code:`save()` function allows to easily store SI objects to a folder on disk. :py:class:`~spikeinterface.core.BaseRecording` objects are stored in binary (.raw) or `Zarr `__ (.zarr) format and -:py:class:`~spikeinterface.core.BaseSorting` and :py:class:`~spikeinterface.core.BaseSnippets` object in numpy (.npz) +:py:class:`~spikeinterface.core.BaseSorting` object in numpy (.npy) or Zarr (.zarr) format. format. With the actual data, the :code:`save()` function also stores the provenance dictionary and all the properties and annotations associated to the object. The save function also supports parallel processing to speed up the writing process. @@ -673,10 +624,6 @@ This saving/loading features enables us to store SpikeInterface objects efficien recording_bin = recording.save(folder="recording", **job_kwargs) # save recording to folder in zarr format (.zarr is appended automatically) recording_zarr = recording.save(folder="recording", format="zarr", **job_kwargs) - # save snippets to NPZ - snippets_saved = snippets.save(folder="snippets") - # save sorting to NPZ - sorting_saved = sorting.save(folder="sorting") **NOTE:** the Zarr format by default applies data compression with :code:`Blosc.Zstandard` codec with BIT shuffling. Any other Zarr-compatible compressors and filters can be applied using the :code:`compressor` and :code:`filters` @@ -749,8 +696,7 @@ in-memory objects (for example, for testing a new method) or "manually" add some workflow. In order to do this, one can use the :code:`Numpy*` classes, :py:class:`~spikeinterface.core.NumpyRecording`, -:py:class:`~spikeinterface.core.NumpySorting`, :py:class:`~spikeinterface.core.NumpyEvent`, and -:py:class:`~spikeinterface.core.NumpySnippets`. These object behave exactly like normal SpikeInterface objects, +:py:class:`~spikeinterface.core.NumpySorting`, and :py:class:`~spikeinterface.core.NumpyEvent`. These object behave exactly like normal SpikeInterface objects, but they are not bound to a file. Also note the class :py:class:`~spikeinterface.core.SharedMemorySorting` which is very similar to @@ -814,7 +760,7 @@ Any sorting object can be transformed into a :py:class:`~spikeinterface.core.Num Manipulating objects: slicing, aggregating ------------------------------------------- -:py:class:`~spikeinterface.core.BaseRecording` (and :py:class:`~spikeinterface.core.BaseSnippets`) +:py:class:`~spikeinterface.core.BaseRecording` and :py:class:`~spikeinterface.core.BaseSorting` objects can be sliced on the time or channel/unit axis. These operations are completely lazy, as there is no data duplication. After slicing or aggregating, @@ -949,7 +895,7 @@ They are useful to make examples, tests, and small demos: .. code-block:: python - from spikeinterface.core import generate_recording, generate_sorting, generate_snippets + from spikeinterface.core import generate_recording, generate_sorting # recording with 2 segments and 4 channels recording = generate_recording(num_channels=4, sampling_frequency=30000., @@ -959,11 +905,6 @@ They are useful to make examples, tests, and small demos: sorting = generate_sorting(num_units=5, sampling_frequency=30000., durations=[10.325, 3.5], firing_rates=15, refractory_period_ms=1.5) - # snippets of 60 samples on 2 channels from 5 units - snippets = generate_snippets(nbefore=20, nafter=40, num_channels=2, - sampling_frequency=30000., durations=[10.325, 3.5], - set_probe=True, num_units=5) - There are also some more advanced functions to generate sorting objects with varioues "mistakes" (mainly for testing purposes): diff --git a/doc/modules/preprocessing.rst b/doc/modules/preprocessing.rst index 765996ad80..13bf385d58 100644 --- a/doc/modules/preprocessing.rst +++ b/doc/modules/preprocessing.rst @@ -532,21 +532,6 @@ Of course, these pipelines can be enhanced and customized using other available - -Preprocessing on Snippets -------------------------- - -Some preprocessing steps are available also for :py:class:`~spikeinterface.core.BaseSnippets` objects: - -align_snippets() -^^^^^^^^^^^^^^^^ - -This function aligns waveform snippets. - -* :py:func:`~spikeinterface.preprocessing.align_snippets()` - - - References ---------- diff --git a/src/spikeinterface/core/__init__.py b/src/spikeinterface/core/__init__.py index 72222753f7..7f30078a6d 100644 --- a/src/spikeinterface/core/__init__.py +++ b/src/spikeinterface/core/__init__.py @@ -1,8 +1,6 @@ from .baserecording import BaseRecording, BaseRecordingSegment from .basesorting import BaseSorting, BaseSortingSegment, SpikeVectorSortingSegment from .baseevent import BaseEvent, BaseEventSegment -from .basesnippets import BaseSnippets, BaseSnippetsSegment -from .baserecordingsnippets import BaseRecordingSnippets from .loading import load @@ -15,16 +13,13 @@ NumpySorting, SharedMemorySorting, NumpyEvent, - NumpySnippets, ) from .zarrextractors import read_zarr, read_zarr_array, get_default_zarr_compressor from .binaryfolder import BinaryFolderRecording, read_binary_folder from .sortingfolder import NumpyFolderSorting, NpzFolderSorting, read_numpy_sorting_folder, read_npz_folder -from .npysnippetsextractor import NpySnippetsExtractor, read_npy_snippets -from .npyfoldersnippets import NpyFolderSnippets, read_npy_snippets_folder # utility extractors (equivalent to OLD subrecording/subsorting) -from .channelslice import ChannelSliceRecording, ChannelSliceSnippets +from .channelslice import ChannelSliceRecording from .unitsselectionsorting import UnitsSelectionSorting from .frameslicerecording import FrameSliceRecording from .frameslicesorting import FrameSliceSorting @@ -38,7 +33,6 @@ generate_sorting, add_synchrony_to_sorting, create_sorting_npz, - generate_snippets, synthesize_random_firings, inject_some_duplicate_units, inject_some_split_units, @@ -121,7 +115,6 @@ ) from .waveform_tools import extract_waveforms_to_buffers, estimate_templates, estimate_templates_with_accumulator -from .snippets_tools import snippets_from_sorting # waveform extractor # Important not for compatibility!! diff --git a/src/spikeinterface/core/baserecording.py b/src/spikeinterface/core/baserecording.py index 772c026ccc..9ff382b3b3 100644 --- a/src/spikeinterface/core/baserecording.py +++ b/src/spikeinterface/core/baserecording.py @@ -3,21 +3,22 @@ from pathlib import Path import numpy as np -from probeinterface import read_probeinterface, write_probeinterface +from probeinterface import Probe, ProbeGroup, select_axes from .time_series import TimeSeriesSegment, TimeSeries -from .baserecordingsnippets import BaseRecordingSnippets +from .base import BaseExtractor from .core_tools import convert_bytes_to_str, convert_seconds_to_str from .job_tools import split_job_kwargs +from .recording_tools import _set_group_property_based_on_probegroup, check_probe_do_not_overlap -class BaseRecording(BaseRecordingSnippets, TimeSeries): +class BaseRecording(BaseExtractor, TimeSeries): """ Abstract class representing several a multichannel timeseries (or block of raw ephys traces). Internally handle list of RecordingSegment """ - _main_annotations = BaseRecordingSnippets._main_annotations + ["is_filtered"] + _main_annotations = BaseExtractor._main_annotations + ["is_filtered"] _main_properties = [ "group", "gain_to_uV", @@ -38,10 +39,11 @@ class BaseRecording(BaseRecordingSnippets, TimeSeries): ] def __init__(self, sampling_frequency: float, channel_ids: list, dtype): - BaseRecordingSnippets.__init__( - self, channel_ids=channel_ids, sampling_frequency=sampling_frequency, dtype=dtype - ) + BaseExtractor.__init__(self, channel_ids) TimeSeries.__init__(self) + self._sampling_frequency = float(sampling_frequency) + self._dtype = np.dtype(dtype) + self._probegroup = None # initialize main annotation and properties self.annotate(is_filtered=False) @@ -167,6 +169,31 @@ def __sub__(self, other): return SubtractRecordings(self, other) + @property + def channel_ids(self): + return self._main_ids + + @property + def sampling_frequency(self): + return self._sampling_frequency + + @property + def dtype(self): + return self._dtype + + def get_sampling_frequency(self): + return self._sampling_frequency + + def get_channel_ids(self): + return self._main_ids + + def get_num_channels(self): + return len(self.get_channel_ids()) + + def get_dtype(self): + return self._dtype + + # TimeSeries zone @property def segments(self) -> list["BaseRecordingSegment"]: """List of recording segments.""" @@ -305,6 +332,495 @@ def get_data(self, start_frame: int, end_frame: int, segment_index: int | None = def get_shape(self, segment_index: int | None = None) -> tuple[int, ...]: return (self.get_num_samples(segment_index=segment_index), self.get_num_channels()) + # Main properties + def is_filtered(self): + # the is_filtered is handle with annotation + return self._annotations.get("is_filtered", False) + + def set_channel_groups(self, groups, channel_ids=None): + if "probes" in self._annotations: + warn("set_channel_groups() destroys the probe description. Using set_probe() is preferable") + self._annotations.pop("probes") + self.set_property("group", groups, ids=channel_ids) + + def get_channel_groups(self, channel_ids=None): + groups = self.get_property("group", ids=channel_ids) + return groups + + def clear_channel_groups(self, channel_ids=None): + if channel_ids is None: + n = self.get_num_channels() + else: + n = len(channel_ids) + groups = np.zeros(n, dtype="int64") + self.set_property("group", groups, ids=channel_ids) + + def set_channel_gains(self, gains, channel_ids=None): + if np.isscalar(gains): + gains = [gains] * self.get_num_channels() + self.set_property("gain_to_uV", gains, ids=channel_ids) + + def get_channel_gains(self, channel_ids=None): + return self.get_property("gain_to_uV", ids=channel_ids) + + def set_channel_offsets(self, offsets, channel_ids=None): + if np.isscalar(offsets): + offsets = [offsets] * self.get_num_channels() + self.set_property("offset_to_uV", offsets, ids=channel_ids) + + def get_channel_offsets(self, channel_ids=None): + return self.get_property("offset_to_uV", ids=channel_ids) + + def get_channel_property(self, channel_id, key): + values = self.get_property(key) + v = values[self.id_to_index(channel_id)] + return v + + # Probe + def has_scaleable_traces(self) -> bool: + if self.get_property("gain_to_uV") is None or self.get_property("offset_to_uV") is None: + return False + else: + return True + + def has_probe(self) -> bool: + # probe group is saved and loaded to binary/zarr, so we don't need to check for legacy "contact_vector" property + return self._probegroup is not None + + def has_3d_probe(self) -> bool: + if self.has_probe(): + probe = self.get_probegroup().probes[0] + return probe.ndim == 3 + else: + return False + + def has_channel_location(self) -> bool: + return self.has_probe() + + def remove_probe(self): + """ + Removes probe information + """ + self._probegroup = None + + def set_probe( + self, + probe: Probe, + group_mode: Literal["auto", "by_probe", "by_shank", "by_side"] = "auto", + in_place: bool | None = None, + ) -> None: + """ + Attach a Probe object to a recording. + + Parameters + ---------- + probe: Probe + The probe to be attached to the recording + group_mode: "auto" | "by_probe" | "by_shank" | "by_side", default: "auto" + How to add the "group" property. + "auto" is the best splitting possible that can be all at once when multiple probes, multiple shanks + and two sides are present. + in_place: (deprecated) bool | None, default: None + Deprecated argument to indicate whether to modify the recording in place + or return a new recording. The function is always in place now. + Use the `recording.select_channels_with_probegroup()` method instead of `in_place=False` + to return a new recording with a channel selection to match the probe/probegroup. + + Notes + ----- + Internally, this will construct a ProbeGroup with the probe and call `set_probegroup()`. + """ + assert isinstance(probe, Probe), "The input must be a Probe object" + probegroup = ProbeGroup() + probegroup.add_probe(probe) + # TODO: remove return in 0.106.0 after removing in_place argument + return self.set_probegroup(probegroup, group_mode=group_mode, in_place=in_place) + + def set_probegroup( + self, + probegroup: ProbeGroup, + group_mode: Literal["auto", "by_probe", "by_shank", "by_side"] = "auto", + in_place: bool | None = None, + check_overlap: bool = True, + ) -> None: + """ + Attach a ProbeGroup or dict to a recording. + For this Probe.device_channel_indices is used to link contacts to recording channels. + After removing unconnected contacts, the number of connected contacts must match the + number of channels in the recording. If this is not the case, use the `recording.select_with_probegroup()` + method instead to return a new recording with a channel selection to match the probe/probegroup. + + Note: The probe order of the probegroup is not kept. Channel ids are re-ordered to match the channel_ids of the recording. + + Parameters + ---------- + probe_or_probegroup: ProbeGroup, or dict + The probe(s) to be attached to the recording + group_mode: "auto" | "by_probe" | "by_shank" | "by_side", default: "auto" + How to add the "group" property. + "auto" is the best splitting possible that can be all at once when multiple probes, multiple shanks and two sides are present. + in_place: (deprecated) bool | None, default: None + Deprecated argument to indicate whether to modify the recording in place + or return a new recording. The function is always in place now. + Use the `recording.select_channels_with_probegroup()` method instead of `in_place=False` + to return a new recording with a channel selection to match the probe/probegroup. + check_overlap: bool, default: True + If True, check that the probes in the probegroup do not overlap in space. + This should be set to False when aggregating recordings whose probes share + the same physical space (e.g. channels split by group from a single probe), + where contact positions are unique but probe bounding boxes may overlap. + """ + if in_place is not None: + warnings.warn( + "The 'in_place' argument is deprecated and will be removed in version 0.106.0. " + "The `set_probe/probegroup()` are now always in place; please remove the in_place argument.", + FutureWarning, + stacklevel=2, + ) + if not in_place: + return self.select_channels_with_probegroup(probegroup, group_mode=group_mode) + + if check_overlap and len(probegroup.probes) > 0: + check_probe_do_not_overlap(probegroup.probes) + + probegroup_sorted = self._get_probegroup_based_on_device_channel_indices(probegroup) + + if probegroup_sorted.get_contact_count() != self.get_num_channels(): + raise ValueError( + "The probe/probegroup must have the same number of connected contacts " + f"as the number of channels as the recording, but the probe has {probegroup.get_contact_count()} " + f"connected channels and the recording has {self.get_num_channels()} channels. " + "Use the `recording.select_channels_with_probegroup()` method instead to return a new recording with " + "a channel selection to match the probe/probegroup." + ) + + device_channel_indices = probegroup_sorted.get_global_device_channel_indices()["device_channel_indices"] + if not np.array_equal(device_channel_indices, np.arange(self.get_num_channels())): + raise ValueError( + "`device_channel_indices` is wrong! " + "It should contain only values [0...n-1] after ordering, " + f"but they are: {device_channel_indices}" + ) + + # probegroup_sorted.set_global_device_channel_indices(np.arange(probegroup_sorted.get_contact_count())) + self._probegroup = probegroup_sorted + + # Handle and set channel groups + _set_group_property_based_on_probegroup(self, probegroup_sorted, group_mode=group_mode) + + def select_channels_with_probe( + self, probe: Probe, group_mode: Literal["auto", "by_probe", "by_shank", "by_side"] = "auto" + ) -> "BaseRecording": + """ + Returns a new recording with channels selected based on the probe. + + Parameters + ---------- + probe: Probe + The probe to be used for channel selection + group_mode: "auto" | "by_probe" | "by_shank" | + "by_side", default: "auto" + How to add the "group" property. + "auto" is the best splitting possible that can be all at once when multiple probes, multiple shanks and two sides are present. + + Returns + ------- + sub_recording: BaseRecording + A view of the recording (ChannelSlice or clone or itself) + """ + assert isinstance(probe, Probe), "The input must be a Probe object" + probegroup = ProbeGroup() + probegroup.add_probe(probe) + return self.select_channels_with_probegroup(probegroup, group_mode=group_mode) + + def select_channels_with_probegroup( + self, probegroup: ProbeGroup, group_mode: Literal["auto", "by_probe", "by_shank", "by_side"] = "auto" + ) -> "BaseRecording": + """ + Selects channels based on the given ProbeGroup and returns a new recording with the selected channels. + + Parameters + ---------- + probegroup: ProbeGroup + The probegroup to be used for channel selection + group_mode: "auto" | "by_probe" | "by_shank" | + "by_side", default: "auto" + How to add the "group" property. + "auto" is the best splitting possible that can be all at once when multiple probes, multiple shanks + and two sides are present. + + Returns + ------- + sub_recording: BaseRecording + A view of the recording (ChannelSlice or clone or itself) + """ + probegroup_sorted = self._get_probegroup_based_on_device_channel_indices(probegroup) + if probegroup_sorted.get_contact_count() > 0: + sorted_dci = probegroup_sorted.get_global_device_channel_indices()["device_channel_indices"] + new_channel_ids = self.channel_ids[sorted_dci] + probegroup_sorted.set_global_device_channel_indices(np.arange(len(new_channel_ids))) + if np.array_equal(new_channel_ids, self.channel_ids): + sub_recording = self.clone() + else: + sub_recording = self.select_channels(new_channel_ids) + sub_recording._probegroup = probegroup_sorted + _set_group_property_based_on_probegroup(sub_recording, probegroup_sorted, group_mode=group_mode) + else: + sub_recording = self.select_channels([]) # empty recording + sub_recording._probegroup = ProbeGroup() # empty probegroup + return sub_recording + + def _get_probegroup_based_on_device_channel_indices(self, probegroup: ProbeGroup) -> ProbeGroup: + """ + Returns a new probegroup sorted based on their device_channel_indices. + This is useful to ensure that the probes are ordered correctly when attached to a recording. + Also checks that the device_channel_indices are consistent with the recording channel count and + contacts are unique across probes in the probegroup. + + Parameters + ---------- + probegroup : ProbeGroup + The probegroup to be sorted. + + Returns + ------- + ProbeGroup + The sorted probegroup. + """ + if not isinstance(probegroup, ProbeGroup): + raise ValueError("The input must be a ProbeGroup or dict") + + assert all( + probe.device_channel_indices is not None for probe in probegroup.probes + ), "Probe must have device_channel_indices" + + # Remove unconnected contacts and slice the probe group accordingly + device_channel_indices = probegroup.get_global_device_channel_indices()["device_channel_indices"] + keep_indices = np.flatnonzero(device_channel_indices >= 0) + if len(keep_indices) < len(device_channel_indices): + if len(keep_indices) == 0: + device_channel_indices = np.array([], dtype="int64") + else: + probegroup = probegroup.get_slice(keep_indices) + device_channel_indices = device_channel_indices[keep_indices] + + if len(device_channel_indices) > 0: + # Check consistency of device_channel_indices with the recording channel count + number_of_device_channel_indices = np.max(list(device_channel_indices) + [0]) + if number_of_device_channel_indices >= self.get_num_channels(): + error_msg = ( + f"The given Probe either has 'device_channel_indices' that does not match channel count \n" + f"{len(device_channel_indices)} vs {self.get_num_channels()} \n" + f"or it's max index {number_of_device_channel_indices} is the same as the number of channels {self.get_num_channels()} \n" + f"If using all channels remember that python is 0-indexed so max device_channel_index should be {self.get_num_channels() - 1} \n" + f"device_channel_indices are the following: {device_channel_indices} \n" + f"recording channels are the following: {self.get_channel_ids()} \n" + ) + raise ValueError(error_msg) + # Now slice the probe using the device channel indices to match the recording channel_ids + order = np.argsort(device_channel_indices) + probegroup = probegroup.get_slice(order) + else: + warn( + "No connected channels in the probegroup! " + "The probegroup will be attached but no channel will be selected." + ) + probegroup = ProbeGroup() # empty probegroup + + return probegroup + + def get_probe(self): + probes = self.get_probes() + assert len(probes) == 1, "There are several probe use .get_probes() or get_probegroup()" + return probes[0] + + def get_probes(self): + probegroup = self.get_probegroup() + return probegroup.probes + + def get_probegroup(self): + if self._probegroup is None: + raise ValueError("There is no Probe attached to this recording. Use set_probe(...) to attach one.") + return self._probegroup + + def create_dummy_probe_from_locations(self, locations, shape="circle", shape_params={"radius": 1}, axes="xy"): + """ + Creates a "dummy" probe based on locations. + + Parameters + ---------- + locations : np.array + Array with channel locations (num_channels, ndim) [ndim can be 2 or 3] + shape : str, default: "circle" + Electrode shapes + shape_params : dict, default: {"radius": 1} + Shape parameters + axes : str, default: "xy" + If ndim is 3, indicates the axes that define the plane of the electrodes + + Returns + ------- + probe : Probe + The created probe + """ + ndim = locations.shape[1] + probe = Probe(ndim=2) + if ndim == 3: + locations_2d = select_axes(locations, axes) + else: + locations_2d = locations + probe.set_contacts(locations_2d, shapes=shape, shape_params=shape_params) + probe.set_device_channel_indices(np.arange(self.get_num_channels())) + + if ndim == 3: + probe = probe.to_3d(axes=axes) + + return probe + + def set_dummy_probe_from_locations(self, locations, shape="circle", shape_params={"radius": 1}, axes="xy"): + """ + Sets a "dummy" probe based on locations. + + Parameters + ---------- + locations : np.array + Array with channel locations (num_channels, ndim) [ndim can be 2 or 3] + shape : str, default: "circle" + Electrode shapes + shape_params : dict, default: {"radius": 1} + Shape parameters + axes : "xy" | "yz" | "xz", default: "xy" + If ndim is 3, indicates the axes that define the plane of the electrodes + """ + probe = self.create_dummy_probe_from_locations( + np.array(locations), shape=shape, shape_params=shape_params, axes=axes + ) + self.set_probe(probe) + + def set_channel_locations(self, locations, channel_ids=None): + warnings.warn( + ( + "set_channel_locations() is deprecated and will be removed in version 0.106.0. " + "If you want to set probe information, use `set_dummy_probe_from_locations()`." + ), + FutureWarning, + stacklevel=2, + ) + self.set_dummy_probe_from_locations(locations, axes="xy") + + def get_channel_locations( + self, + channel_ids: list | np.ndarray | tuple | None = None, + axes: Literal["xy", "yz", "xz", "xyz"] = "xy", + ) -> np.ndarray: + """ + Get the physical locations of specified channels. + + Parameters + ---------- + channel_ids : array-like, optional + The IDs of the channels for which to retrieve locations. If None, retrieves locations + for all available channels. Default is None. + axes : "xy" | "yz" | "xz" | "xyz", default: "xy" + The spatial axes to return, specified as a string (e.g., "xy", "xyz"). Default is "xy". + + Returns + ------- + np.ndarray + A 2D or 3D array of shape (n_channels, n_dimensions) containing the locations of the channels. + The number of dimensions depends on the `axes` argument (e.g., 2 for "xy", 3 for "xyz"). + """ + if channel_ids is None: + channel_ids = self.get_channel_ids() + channel_indices = self.ids_to_indices(channel_ids) + if not self.has_probe(): + raise ValueError("get_channel_locations(..) needs a probe to be attached to the recording") + probegroup = self.get_probegroup() + contact_positions = probegroup.get_global_contact_positions() + return select_axes(contact_positions, axes)[channel_indices] + + def is_probe_3d(self) -> bool: + if not self.has_probe(): + raise ValueError("is_probe_3d() needs a probe to be attached to the recording") + probegroup = self.get_probegroup() + return probegroup.ndim == 3 + + def clear_channel_locations(self, channel_ids=None): + warnings.warn( + ( + "clear_channel_locations() is deprecated and will be removed in version 0.106.0. " + "If you want to remove probe information, use `remove_probe()`." + ), + FutureWarning, + stacklevel=2, + ) + self.remove_probe() + + def planarize(self, axes: str = "xy"): + """ + Returns a Recording with a 2D probe from one with a 3D probe + + Parameters + ---------- + axes : "xy" | "yz" |"xz", default: "xy" + The axes to keep + + Returns + ------- + BaseRecording + The recording with 2D positions + """ + assert self.has_3d_probe(), "The 'planarize' function needs a recording with 3d locations" + assert len(axes) == 2, "You need to specify 2 dimensions (e.g. 'xy', 'zy')" + + probe2d = self.get_probe().to_2d(axes=axes) + recording2d = self.clone() + recording2d.set_probe(probe2d) + + return recording2d + + def split_by(self, property="group", outputs="dict"): + """ + Splits object based on a certain property (e.g. "group") + + Parameters + ---------- + property : str, default: "group" + The property to use to split the object, default: "group" + outputs : "dict" | "list", default: "dict" + Whether to return a dict or a list + + Returns + ------- + dict or list + A dict or list with grouped objects based on property + + Raises + ------ + ValueError + Raised when property is not present + """ + assert outputs in ("list", "dict") + values = self.get_property(property) + if values is None: + raise ValueError(f"property {property} is not set") + + if outputs == "list": + recordings = [] + elif outputs == "dict": + recordings = {} + for value in np.unique(values).tolist(): + (inds,) = np.nonzero(values == value) + new_channel_ids = self.channel_ids[inds] + subrec = self.select_channels(new_channel_ids) + subrec.set_annotation("split_by_property", value=property) + if outputs == "list": + recordings.append(subrec) + elif outputs == "dict": + recordings[value] = subrec + return recordings + + # Save and metadata handling and propagation def save(self, format="binary", verbose: bool = False, **save_kwargs): """ Save a `BaseRecording` object to a specified format: @@ -416,7 +932,10 @@ def save(self, format="binary", verbose: bool = False, **save_kwargs): return cached def _extra_metadata_to_dict(self, dump_dict): - super()._extra_metadata_to_dict(dump_dict) + # Save probe + if self.has_probe(): + probegroup = self.get_probegroup() + dump_dict["probegroup"] = probegroup.to_dict() # Add times_kwargs if the recording has been modified in memory (e.g. by set_times / shift_times / reset_times) if self._time_info_modified: @@ -426,7 +945,11 @@ def _extra_metadata_to_dict(self, dump_dict): dump_dict["times_kwargs"].append(times_kwargs) def _extra_metadata_from_dict(self, dump_dict): - super()._extra_metadata_from_dict(dump_dict) + # Load probe and handle backward-compatibility with legacy "contact_vector"/"location" property + if "probegroup" in dump_dict: + # this is for SI>=0.105.0 + probegroup = dump_dict["probegroup"] + self._probegroup = ProbeGroup.from_dict(probegroup) if "times_kwargs" in dump_dict: # When serializing, dump timestamps information because this could have been @@ -437,6 +960,52 @@ def _extra_metadata_from_dict(self, dump_dict): self.segments[segment_index]._t_start = times_kwargs["t_start"] self.segments[segment_index]._time_vector = times_kwargs["time_vector"] + def _extra_metadata_copy(self, other): + if self._probegroup is not None: + other._probegroup = self._probegroup.copy() + + def _handle_extractor_backward_compatibility(self): + """ + This handles backward compatibility for recordings that were saved with older versions of spikeinterface. + + Options: + + 1. "contact_vector" property: This was used in versions < 0.105.0 to store probe information, when saved to + pickle + 2. "location" property: This was used in versions < 0.105.0 to store probe information, when saved to JSON + (no contact_vector saved) + 3. probe annotation: probe annotations and contours were saved as recording properties in versions < 0.105.0, + but now they are saved in the probegroup. This method will copy the annotations and the contour to the probes + in the the probegroup and remove the annotations from the recording. + """ + if self._probegroup is None: + check_for_probes_info = False + if "contact_vector" in self.get_property_keys(): + # this is for SI<0.105.0 and from pickle + contact_vector = self.get_property("contact_vector") + probegroup = ProbeGroup.from_numpy(contact_vector) + self._probegroup = probegroup + check_for_probes_info = True + elif "location" in self.get_property_keys(): + # this is for SI<0.105.0 and from JSON (no contact_vector saved) + locations = self.get_property("location") + self.set_dummy_probe_from_locations(locations) + check_for_probes_info = True + + if check_for_probes_info: + for i, probe in enumerate(self._probegroup.probes): + if "probes_info" in self._annotations: + probe_dict = self._annotations["probes_info"][i] + probe.annotations.update(probe_dict) + if f"probe_{i}_planar_contour" in self._annotations: + contour = self.get_annotation(f"probe_{i}_planar_contour") + if contour is not None: + probe.set_planar_contour(contour) + self.delete_annotation(f"probe_{i}_planar_contour") + if "probes_info" in self._annotations: + self._annotations.pop("probes_info") + + # Utility methods for channel/time manipulation def select_channels(self, channel_ids: list | np.ndarray | tuple) -> "BaseRecording": """ Returns a new recording object with a subset of channels. @@ -452,26 +1021,21 @@ def select_channels(self, channel_ids: list | np.ndarray | tuple) -> "BaseRecord return ChannelSliceRecording(self, channel_ids) - def rename_channels(self, new_channel_ids: list | np.ndarray | tuple) -> "BaseRecording": + def remove_channels(self, remove_channel_ids): """ - Returns a new recording object with renamed channel ids. + Returns a new object with removed channels. - Note that this method does not modify the current recording and instead returns a new recording object. Parameters ---------- - new_channel_ids : list or np.array or tuple - The new channel ids. They are mapped positionally to the old channel ids. - """ - from .channelslice import ChannelSliceRecording - - assert len(new_channel_ids) == self.get_num_channels(), ( - "new_channel_ids must have the same length as the " "number of channels in the recording" - ) - - return ChannelSliceRecording(self, renamed_channel_ids=new_channel_ids) + remove_channel_ids : np.array or list + The list of channels to remove - def _remove_channels(self, remove_channel_ids): + Returns + ------- + BaseRecording + The object with removed channels + """ from .channelslice import ChannelSliceRecording recording_channel_ids = self.get_channel_ids() @@ -485,6 +1049,43 @@ def _remove_channels(self, remove_channel_ids): sub_recording = ChannelSliceRecording(self, new_channel_ids) return sub_recording + def select_segments(self, segment_indices): + """ + Return a new object with the segments specified by "segment_indices". + + Parameters + ---------- + segment_indices : list of int + List of segment indices to keep in the returned recording + + Returns + ------- + BaseRecording + The onject with the selected segments + """ + from .segmentutils import SelectSegmentRecording + + return SelectSegmentRecording(self, segment_indices=segment_indices) + + def rename_channels(self, new_channel_ids: list | np.ndarray | tuple) -> "BaseRecording": + """ + Returns a new recording object with renamed channel ids. + + Note that this method does not modify the current recording and instead returns a new recording object. + + Parameters + ---------- + new_channel_ids : list or np.array or tuple + The new channel ids. They are mapped positionally to the old channel ids. + """ + from .channelslice import ChannelSliceRecording + + assert len(new_channel_ids) == self.get_num_channels(), ( + "new_channel_ids must have the same length as the " "number of channels in the recording" + ) + + return ChannelSliceRecording(self, renamed_channel_ids=new_channel_ids) + def frame_slice(self, start_frame: int | None, end_frame: int | None) -> "BaseRecording": """ Returns a new recording with sliced frames. Note that this operation is not in place. @@ -557,35 +1158,12 @@ def time_slice(self, start_time: float | None, end_time: float | None) -> "BaseR return self.frame_slice(start_frame=start_frame, end_frame=end_frame) - def _select_segments(self, segment_indices): - from .segmentutils import SelectSegmentRecording - - return SelectSegmentRecording(self, segment_indices=segment_indices) - - def get_channel_locations( - self, - channel_ids: list | np.ndarray | tuple | None = None, - axes: Literal["xy", "yz", "xz", "xyz"] = "xy", - ) -> np.ndarray: - """ - Get the physical locations of specified channels. - - Parameters - ---------- - channel_ids : array-like, optional - The IDs of the channels for which to retrieve locations. If None, retrieves locations - for all available channels. Default is None. - axes : "xy" | "yz" | "xz" | "xyz", default: "xy" - The spatial axes to return, specified as a string (e.g., "xy", "xyz"). Default is "xy". + def astype(self, dtype, round: bool | None = None): + from spikeinterface.preprocessing.astype import astype - Returns - ------- - np.ndarray - A 2D or 3D array of shape (n_channels, n_dimensions) containing the locations of the channels. - The number of dimensions depends on the `axes` argument (e.g., 2 for "xy", 3 for "xyz"). - """ - return super().get_channel_locations(channel_ids=channel_ids, axes=axes) + return astype(self, dtype=dtype, round=round) + # Binary compatibility def is_binary_compatible(self) -> bool: """ Checks if the recording is "binary" compatible. @@ -648,11 +1226,6 @@ def binary_compatible_with( # good job you pass all crucible return True - def astype(self, dtype, round: bool | None = None): - from spikeinterface.preprocessing.astype import astype - - return astype(self, dtype=dtype, round=round) - class BaseRecordingSegment(TimeSeriesSegment): """ diff --git a/src/spikeinterface/core/baserecordingsnippets.py b/src/spikeinterface/core/baserecordingsnippets.py deleted file mode 100644 index 11eb909b59..0000000000 --- a/src/spikeinterface/core/baserecordingsnippets.py +++ /dev/null @@ -1,637 +0,0 @@ -from pathlib import Path -from typing import Literal -import warnings -import numpy as np - -from probeinterface import Probe, ProbeGroup, write_probeinterface, read_probeinterface, select_axes - -from .base import BaseExtractor -from .recording_tools import _set_group_property_based_on_probegroup, check_probe_do_not_overlap - -from warnings import warn - - -class BaseRecordingSnippets(BaseExtractor): - """ - Mixin that handles all probe and channel operations - """ - - def __init__(self, sampling_frequency: float, channel_ids: list[str, int], dtype: np.dtype): - BaseExtractor.__init__(self, channel_ids) - self._sampling_frequency = float(sampling_frequency) - self._dtype = np.dtype(dtype) - self._probegroup = None - - @property - def channel_ids(self): - return self._main_ids - - @property - def sampling_frequency(self): - return self._sampling_frequency - - @property - def dtype(self): - return self._dtype - - def get_sampling_frequency(self): - return self._sampling_frequency - - def get_channel_ids(self): - return self._main_ids - - def get_num_channels(self): - return len(self.get_channel_ids()) - - def get_dtype(self): - return self._dtype - - def has_scaleable_traces(self) -> bool: - if self.get_property("gain_to_uV") is None or self.get_property("offset_to_uV") is None: - return False - else: - return True - - def has_probe(self) -> bool: - # probe group is saved and loaded to binary/zarr, so we don't need to check for legacy "contact_vector" property - return self._probegroup is not None - - def has_3d_probe(self) -> bool: - if self.has_probe(): - probe = self.get_probegroup().probes[0] - return probe.ndim == 3 - else: - return False - - def has_channel_location(self) -> bool: - return self.has_probe() - - def is_filtered(self): - # the is_filtered is handle with annotation - return self._annotations.get("is_filtered", False) - - def remove_probe(self): - """ - Removes probe information - """ - self._probegroup = None - - def set_probe( - self, - probe: Probe, - group_mode: Literal["auto", "by_probe", "by_shank", "by_side"] = "auto", - in_place: bool | None = None, - ) -> None: - """ - Attach a Probe object to a recording. - - Parameters - ---------- - probe: Probe - The probe to be attached to the recording - group_mode: "auto" | "by_probe" | "by_shank" | "by_side", default: "auto" - How to add the "group" property. - "auto" is the best splitting possible that can be all at once when multiple probes, multiple shanks - and two sides are present. - in_place: (deprecated) bool | None, default: None - Deprecated argument to indicate whether to modify the recording in place - or return a new recording. The function is always in place now. - Use the `recording.select_channels_with_probegroup()` method instead of `in_place=False` - to return a new recording with a channel selection to match the probe/probegroup. - - Notes - ----- - Internally, this will construct a ProbeGroup with the probe and call `set_probegroup()`. - """ - assert isinstance(probe, Probe), "The input must be a Probe object" - probegroup = ProbeGroup() - probegroup.add_probe(probe) - # TODO: remove return in 0.106.0 after removing in_place argument - return self.set_probegroup(probegroup, group_mode=group_mode, in_place=in_place) - - def set_probegroup( - self, - probegroup: ProbeGroup, - group_mode: Literal["auto", "by_probe", "by_shank", "by_side"] = "auto", - in_place: bool | None = None, - check_overlap: bool = True, - ) -> None: - """ - Attach a ProbeGroup or dict to a recording. - For this Probe.device_channel_indices is used to link contacts to recording channels. - After removing unconnected contacts, the number of connected contacts must match the - number of channels in the recording. If this is not the case, use the `recording.select_with_probegroup()` - method instead to return a new recording with a channel selection to match the probe/probegroup. - - Note: The probe order of the probegroup is not kept. Channel ids are re-ordered to match the channel_ids of the recording. - - Parameters - ---------- - probe_or_probegroup: ProbeGroup, or dict - The probe(s) to be attached to the recording - group_mode: "auto" | "by_probe" | "by_shank" | "by_side", default: "auto" - How to add the "group" property. - "auto" is the best splitting possible that can be all at once when multiple probes, multiple shanks and two sides are present. - in_place: (deprecated) bool | None, default: None - Deprecated argument to indicate whether to modify the recording in place - or return a new recording. The function is always in place now. - Use the `recording.select_channels_with_probegroup()` method instead of `in_place=False` - to return a new recording with a channel selection to match the probe/probegroup. - check_overlap: bool, default: True - If True, check that the probes in the probegroup do not overlap in space. - This should be set to False when aggregating recordings whose probes share - the same physical space (e.g. channels split by group from a single probe), - where contact positions are unique but probe bounding boxes may overlap. - """ - if in_place is not None: - warnings.warn( - "The 'in_place' argument is deprecated and will be removed in version 0.106.0. " - "The `set_probe/probegroup()` are now always in place; please remove the in_place argument.", - FutureWarning, - stacklevel=2, - ) - if not in_place: - return self.select_channels_with_probegroup(probegroup, group_mode=group_mode) - - if check_overlap and len(probegroup.probes) > 0: - check_probe_do_not_overlap(probegroup.probes) - - probegroup_sorted = self._get_probegroup_based_on_device_channel_indices(probegroup) - - if probegroup_sorted.get_contact_count() != self.get_num_channels(): - raise ValueError( - "The probe/probegroup must have the same number of connected contacts " - f"as the number of channels as the recording, but the probe has {probegroup.get_contact_count()} " - f"connected channels and the recording has {self.get_num_channels()} channels. " - "Use the `recording.select_channels_with_probegroup()` method instead to return a new recording with " - "a channel selection to match the probe/probegroup." - ) - - device_channel_indices = probegroup_sorted.get_global_device_channel_indices()["device_channel_indices"] - if not np.array_equal(device_channel_indices, np.arange(self.get_num_channels())): - raise ValueError( - "`device_channel_indices` is wrong! " - "It should contain only values [0...n-1] after ordering, " - f"but they are: {device_channel_indices}" - ) - - # probegroup_sorted.set_global_device_channel_indices(np.arange(probegroup_sorted.get_contact_count())) - self._probegroup = probegroup_sorted - - # Handle and set channel groups - _set_group_property_based_on_probegroup(self, probegroup_sorted, group_mode=group_mode) - - def select_channels_with_probe( - self, probe: Probe, group_mode: Literal["auto", "by_probe", "by_shank", "by_side"] = "auto" - ) -> "BaseRecordingSnippets": - """ - Returns a new recording with channels selected based on the probe. - - Parameters - ---------- - probe: Probe - The probe to be used for channel selection - group_mode: "auto" | "by_probe" | "by_shank" | - "by_side", default: "auto" - How to add the "group" property. - "auto" is the best splitting possible that can be all at once when multiple probes, multiple shanks and two sides are present. - - Returns - ------- - sub_recording: BaseRecording - A view of the recording (ChannelSlice or clone or itself) - """ - assert isinstance(probe, Probe), "The input must be a Probe object" - probegroup = ProbeGroup() - probegroup.add_probe(probe) - return self.select_channels_with_probegroup(probegroup, group_mode=group_mode) - - def select_channels_with_probegroup( - self, probegroup: ProbeGroup, group_mode: Literal["auto", "by_probe", "by_shank", "by_side"] = "auto" - ) -> "BaseRecordingSnippets": - """ - Selects channels based on the given ProbeGroup and returns a new recording with the selected channels. - - Parameters - ---------- - probegroup: ProbeGroup - The probegroup to be used for channel selection - group_mode: "auto" | "by_probe" | "by_shank" | - "by_side", default: "auto" - How to add the "group" property. - "auto" is the best splitting possible that can be all at once when multiple probes, multiple shanks - and two sides are present. - - Returns - ------- - sub_recording: BaseRecording - A view of the recording (ChannelSlice or clone or itself) - """ - probegroup_sorted = self._get_probegroup_based_on_device_channel_indices(probegroup) - if probegroup_sorted.get_contact_count() > 0: - sorted_dci = probegroup_sorted.get_global_device_channel_indices()["device_channel_indices"] - new_channel_ids = self.channel_ids[sorted_dci] - probegroup_sorted.set_global_device_channel_indices(np.arange(len(new_channel_ids))) - if np.array_equal(new_channel_ids, self.channel_ids): - sub_recording = self.clone() - else: - sub_recording = self.select_channels(new_channel_ids) - sub_recording._probegroup = probegroup_sorted - _set_group_property_based_on_probegroup(sub_recording, probegroup_sorted, group_mode=group_mode) - else: - sub_recording = self.select_channels([]) # empty recording - sub_recording._probegroup = ProbeGroup() # empty probegroup - return sub_recording - - def _get_probegroup_based_on_device_channel_indices(self, probegroup: ProbeGroup) -> ProbeGroup: - """ - Returns a new probegroup sorted based on their device_channel_indices. - This is useful to ensure that the probes are ordered correctly when attached to a recording. - Also checks that the device_channel_indices are consistent with the recording channel count and - contacts are unique across probes in the probegroup. - - Parameters - ---------- - probegroup : ProbeGroup - The probegroup to be sorted. - - Returns - ------- - ProbeGroup - The sorted probegroup. - """ - if not isinstance(probegroup, ProbeGroup): - raise ValueError("The input must be a ProbeGroup or dict") - - assert all( - probe.device_channel_indices is not None for probe in probegroup.probes - ), "Probe must have device_channel_indices" - - # Remove unconnected contacts and slice the probe group accordingly - device_channel_indices = probegroup.get_global_device_channel_indices()["device_channel_indices"] - keep_indices = np.flatnonzero(device_channel_indices >= 0) - if len(keep_indices) < len(device_channel_indices): - if len(keep_indices) == 0: - device_channel_indices = np.array([], dtype="int64") - else: - probegroup = probegroup.get_slice(keep_indices) - device_channel_indices = device_channel_indices[keep_indices] - - if len(device_channel_indices) > 0: - # Check consistency of device_channel_indices with the recording channel count - number_of_device_channel_indices = np.max(list(device_channel_indices) + [0]) - if number_of_device_channel_indices >= self.get_num_channels(): - error_msg = ( - f"The given Probe either has 'device_channel_indices' that does not match channel count \n" - f"{len(device_channel_indices)} vs {self.get_num_channels()} \n" - f"or it's max index {number_of_device_channel_indices} is the same as the number of channels {self.get_num_channels()} \n" - f"If using all channels remember that python is 0-indexed so max device_channel_index should be {self.get_num_channels() - 1} \n" - f"device_channel_indices are the following: {device_channel_indices} \n" - f"recording channels are the following: {self.get_channel_ids()} \n" - ) - raise ValueError(error_msg) - # Now slice the probe using the device channel indices to match the recording channel_ids - order = np.argsort(device_channel_indices) - probegroup = probegroup.get_slice(order) - else: - warn( - "No connected channels in the probegroup! " - "The probegroup will be attached but no channel will be selected." - ) - probegroup = ProbeGroup() # empty probegroup - - return probegroup - - def get_probe(self): - probes = self.get_probes() - assert len(probes) == 1, "There are several probe use .get_probes() or get_probegroup()" - return probes[0] - - def get_probes(self): - probegroup = self.get_probegroup() - return probegroup.probes - - def get_probegroup(self): - if self._probegroup is None: - raise ValueError("There is no Probe attached to this recording. Use set_probe(...) to attach one.") - return self._probegroup - - def _extra_metadata_copy(self, other): - if self._probegroup is not None: - other._probegroup = self._probegroup.copy() - - def _extra_metadata_from_dict(self, dump_dict): - # load probe and handle backward-compatibility with legacy "contact_vector"/"location" property - if "probegroup" in dump_dict: - # this is for SI>=0.105.0 - probegroup = dump_dict["probegroup"] - self._probegroup = ProbeGroup.from_dict(probegroup) - - def _extra_metadata_to_dict(self, dump_dict): - # save probe - if self.has_probe(): - probegroup = self.get_probegroup() - dump_dict["probegroup"] = probegroup.to_dict() - - def _handle_extractor_backward_compatibility(self): - """ - This handles backward compatibility for recordings that were saved with older versions of spikeinterface. - - Options: - - 1. "contact_vector" property: This was used in versions < 0.105.0 to store probe information, when saved to - pickle - 2. "location" property: This was used in versions < 0.105.0 to store probe information, when saved to JSON - (no contact_vector saved) - 3. probe annotation: probe annotations and contours were saved as recording properties in versions < 0.105.0, - but now they are saved in the probegroup. This method will copy the annotations and the contour to the probes - in the the probegroup and remove the annotations from the recording. - """ - if self._probegroup is None: - check_for_probes_info = False - if "contact_vector" in self.get_property_keys(): - # this is for SI<0.105.0 and from pickle - contact_vector = self.get_property("contact_vector") - probegroup = ProbeGroup.from_numpy(contact_vector) - self._probegroup = probegroup - check_for_probes_info = True - elif "location" in self.get_property_keys(): - # this is for SI<0.105.0 and from JSON (no contact_vector saved) - locations = self.get_property("location") - self.set_dummy_probe_from_locations(locations) - check_for_probes_info = True - - if check_for_probes_info: - for i, probe in enumerate(self._probegroup.probes): - if "probes_info" in self._annotations: - probe_dict = self._annotations["probes_info"][i] - probe.annotations.update(probe_dict) - if f"probe_{i}_planar_contour" in self._annotations: - contour = self.get_annotation(f"probe_{i}_planar_contour") - if contour is not None: - probe.set_planar_contour(contour) - self.delete_annotation(f"probe_{i}_planar_contour") - if "probes_info" in self._annotations: - self._annotations.pop("probes_info") - - def create_dummy_probe_from_locations(self, locations, shape="circle", shape_params={"radius": 1}, axes="xy"): - """ - Creates a "dummy" probe based on locations. - - Parameters - ---------- - locations : np.array - Array with channel locations (num_channels, ndim) [ndim can be 2 or 3] - shape : str, default: "circle" - Electrode shapes - shape_params : dict, default: {"radius": 1} - Shape parameters - axes : str, default: "xy" - If ndim is 3, indicates the axes that define the plane of the electrodes - - Returns - ------- - probe : Probe - The created probe - """ - ndim = locations.shape[1] - probe = Probe(ndim=2) - if ndim == 3: - locations_2d = select_axes(locations, axes) - else: - locations_2d = locations - probe.set_contacts(locations_2d, shapes=shape, shape_params=shape_params) - probe.set_device_channel_indices(np.arange(self.get_num_channels())) - - if ndim == 3: - probe = probe.to_3d(axes=axes) - - return probe - - def set_dummy_probe_from_locations(self, locations, shape="circle", shape_params={"radius": 1}, axes="xy"): - """ - Sets a "dummy" probe based on locations. - - Parameters - ---------- - locations : np.array - Array with channel locations (num_channels, ndim) [ndim can be 2 or 3] - shape : str, default: "circle" - Electrode shapes - shape_params : dict, default: {"radius": 1} - Shape parameters - axes : "xy" | "yz" | "xz", default: "xy" - If ndim is 3, indicates the axes that define the plane of the electrodes - """ - probe = self.create_dummy_probe_from_locations( - np.array(locations), shape=shape, shape_params=shape_params, axes=axes - ) - self.set_probe(probe) - - def set_channel_locations(self, locations, channel_ids=None): - warnings.warn( - ( - "set_channel_locations() is deprecated and will be removed in version 0.106.0. " - "If you want to set probe information, use `set_dummy_probe_from_locations()`." - ), - FutureWarning, - stacklevel=2, - ) - self.set_dummy_probe_from_locations(locations, axes="xy") - - def get_channel_locations(self, channel_ids=None, axes: str = "xy") -> np.ndarray: - if channel_ids is None: - channel_ids = self.get_channel_ids() - channel_indices = self.ids_to_indices(channel_ids) - if not self.has_probe(): - raise ValueError("get_channel_locations(..) needs a probe to be attached to the recording") - probegroup = self.get_probegroup() - contact_positions = probegroup.get_global_contact_positions() - return select_axes(contact_positions, axes)[channel_indices] - - def is_probe_3d(self) -> bool: - if not self.has_probe(): - raise ValueError("is_probe_3d() needs a probe to be attached to the recording") - probegroup = self.get_probegroup() - return probegroup.ndim == 3 - - def clear_channel_locations(self, channel_ids=None): - warnings.warn( - ( - "clear_channel_locations() is deprecated and will be removed in version 0.106.0. " - "If you want to remove probe information, use `remove_probe()`." - ), - FutureWarning, - stacklevel=2, - ) - self.remove_probe() - - def set_channel_groups(self, groups, channel_ids=None): - if "probes" in self._annotations: - warn("set_channel_groups() destroys the probe description. Using set_probe() is preferable") - self._annotations.pop("probes") - self.set_property("group", groups, ids=channel_ids) - - def get_channel_groups(self, channel_ids=None): - groups = self.get_property("group", ids=channel_ids) - return groups - - def clear_channel_groups(self, channel_ids=None): - if channel_ids is None: - n = self.get_num_channels() - else: - n = len(channel_ids) - groups = np.zeros(n, dtype="int64") - self.set_property("group", groups, ids=channel_ids) - - def set_channel_gains(self, gains, channel_ids=None): - if np.isscalar(gains): - gains = [gains] * self.get_num_channels() - self.set_property("gain_to_uV", gains, ids=channel_ids) - - def get_channel_gains(self, channel_ids=None): - return self.get_property("gain_to_uV", ids=channel_ids) - - def set_channel_offsets(self, offsets, channel_ids=None): - if np.isscalar(offsets): - offsets = [offsets] * self.get_num_channels() - self.set_property("offset_to_uV", offsets, ids=channel_ids) - - def get_channel_offsets(self, channel_ids=None): - return self.get_property("offset_to_uV", ids=channel_ids) - - def get_channel_property(self, channel_id, key): - values = self.get_property(key) - v = values[self.id_to_index(channel_id)] - return v - - def planarize(self, axes: str = "xy"): - """ - Returns a Recording with a 2D probe from one with a 3D probe - - Parameters - ---------- - axes : "xy" | "yz" |"xz", default: "xy" - The axes to keep - - Returns - ------- - BaseRecording - The recording with 2D positions - """ - assert self.has_3d_probe(), "The 'planarize' function needs a recording with 3d locations" - assert len(axes) == 2, "You need to specify 2 dimensions (e.g. 'xy', 'zy')" - - probe2d = self.get_probe().to_2d(axes=axes) - recording2d = self.clone() - recording2d.set_probe(probe2d) - - return recording2d - - def select_channels(self, channel_ids): - """ - Returns a new object with sliced channels. - - Parameters - ---------- - channel_ids : np.array or list - The list of channels to keep - - Returns - ------- - BaseRecordingSnippets - The object with sliced channels - """ - raise NotImplementedError - - def remove_channels(self, remove_channel_ids): - """ - Returns a new object with removed channels. - - - Parameters - ---------- - remove_channel_ids : np.array or list - The list of channels to remove - - Returns - ------- - BaseRecordingSnippets - The object with removed channels - """ - return self._remove_channels(remove_channel_ids) - - def frame_slice(self, start_frame, end_frame): - """ - Returns a new object with sliced frames. - - Parameters - ---------- - start_frame : int - The start frame - end_frame : int - The end frame - - Returns - ------- - BaseRecordingSnippets - The object with sliced frames - """ - raise NotImplementedError - - def select_segments(self, segment_indices): - """ - Return a new object with the segments specified by "segment_indices". - - Parameters - ---------- - segment_indices : list of int - List of segment indices to keep in the returned recording - - Returns - ------- - BaseRecordingSnippets - The onject with the selected segments - """ - return self._select_segments(segment_indices) - - def split_by(self, property="group", outputs="dict"): - """ - Splits object based on a certain property (e.g. "group") - - Parameters - ---------- - property : str, default: "group" - The property to use to split the object, default: "group" - outputs : "dict" | "list", default: "dict" - Whether to return a dict or a list - - Returns - ------- - dict or list - A dict or list with grouped objects based on property - - Raises - ------ - ValueError - Raised when property is not present - """ - assert outputs in ("list", "dict") - values = self.get_property(property) - if values is None: - raise ValueError(f"property {property} is not set") - - if outputs == "list": - recordings = [] - elif outputs == "dict": - recordings = {} - for value in np.unique(values).tolist(): - (inds,) = np.nonzero(values == value) - new_channel_ids = self.channel_ids[inds] - subrec = self.select_channels(new_channel_ids) - subrec.set_annotation("split_by_property", value=property) - if outputs == "list": - recordings.append(subrec) - elif outputs == "dict": - recordings[value] = subrec - return recordings diff --git a/src/spikeinterface/core/basesnippets.py b/src/spikeinterface/core/basesnippets.py deleted file mode 100644 index 7678f0323b..0000000000 --- a/src/spikeinterface/core/basesnippets.py +++ /dev/null @@ -1,294 +0,0 @@ -import numpy as np -from warnings import warn - -from copy import deepcopy - -from pathlib import Path - -from .base import BaseSegment -from .baserecordingsnippets import BaseRecordingSnippets - - -class BaseSnippets(BaseRecordingSnippets): - """ - Abstract class representing several multichannel snippets. - """ - - _main_properties = ["group", "gain_to_uV", "offset_to_uV"] - _main_features = [] - - def __init__(self, sampling_frequency: float, nbefore: int | None, snippet_len: int, channel_ids: list, dtype): - warn( - "`BaseSnippets` is deprecated and will be removed in version 0.106.0." - "Only continuous recordings with `BaseRecording` will be supported.", - FutureWarning, - stacklevel=2, - ) - BaseRecordingSnippets.__init__( - self, channel_ids=channel_ids, sampling_frequency=sampling_frequency, dtype=dtype - ) - self._nbefore = nbefore - self._snippet_len = snippet_len - - self._snippets_segments: list[BaseSnippetsSegment] = [] - # initialize main annotation and properties - - def __repr__(self): - clsname = self.__class__.__name__ - nchan = self.get_num_channels() - nseg = self.get_num_segments() - sf_khz = self.get_sampling_frequency() / 1000.0 - txt = f"{clsname}: {nchan} channels - {nseg} segments - {sf_khz:0.1f}kHz \n snippet_len:{self._snippet_len} before peak:{self._nbefore}" - return txt - - @property - def nbefore(self): - return self._nbefore - - @property - def snippet_len(self): - return self._snippet_len - - def get_num_segments(self): - return len(self._snippets_segments) - - def add_snippets_segment(self, snippets_segment): - # todo: check channel count and sampling frequency - self._snippets_segments.append(snippets_segment) - snippets_segment.set_parent_extractor(self) - - @property - def nafter(self): - if self._nbefore is None: - return None - return self._snippet_len - self._nbefore - - @property - def snippet_len(self): - return self._snippet_len - - def get_num_snippets(self, segment_index=None): - segment_index = self._check_segment_index(segment_index) - return self._snippets_segments[segment_index].get_num_snippets() - - def get_total_snippets(self): - s = 0 - for segment_index in range(self.get_num_segments()): - s += self.get_num_snippets(segment_index) - return s - - def is_aligned(self): - return self._nbefore is not None - - def get_num_segments(self): - return len(self._snippets_segments) - - def get_frames(self, indices=None, segment_index: int | None = None): - segment_index = self._check_segment_index(segment_index) - spts = self._snippets_segments[segment_index] - return spts.get_frames(indices) - - def get_snippets( - self, - indices=None, - segment_index: int | None = None, - channel_ids: list | None = None, - return_in_uV: bool = False, - ): - """ - Return the snippets, optionally for a subset of samples and/or channels - - Parameters - ---------- - indices : list[int], default: None - Indices of the snippets to return. If None, all snippets are returned. - segment_index : int | None, default: None - The segment index to get snippets from. If snippets is multi-segment, it is required. - channel_ids : list | None, default: None - The channel ids. If None, all channels are used. - return_in_uV : bool, default: False - If True and the snippets has scaling (gain_to_uV and offset_to_uV properties), - snippets are scaled to uV - - Returns - ------- - np.array - The snippets (num_snippets, num_samples, num_channels) - """ - segment_index = self._check_segment_index(segment_index) - spts = self._snippets_segments[segment_index] - channel_indices = self.ids_to_indices(channel_ids, prefer_slice=True) - wfs = spts.get_snippets(indices, channel_indices=channel_indices) - - if return_in_uV: - if not self.has_scaleable_traces(): - raise ValueError( - "These snippets do not support return_in_uV=True (need gain_to_uV and offset_" "to_uV properties)" - ) - else: - gains = self.get_property("gain_to_uV") - offsets = self.get_property("offset_to_uV") - gains = gains[channel_indices].astype("float32") - offsets = offsets[channel_indices].astype("float32") - wfs = wfs.astype("float32") * gains + offsets - return wfs - - def get_snippets_from_frames( - self, - segment_index: int | None = None, - start_frame: int | None = None, - end_frame: int | None = None, - channel_ids: list | None = None, - return_in_uV: bool = False, - ): - """ - Return the snippets from frames, optionally for a subset of samples and/or channels - - Parameters - ---------- - segment_index : int | None, default: None - The segment index to get snippets from. If snippets is multi-segment, it is required. - start_frame : int | None, default: None - The start frame. If None, 0 is used. - end_frame : int | None, default: None - The end frame. If None, the number of samples in the segment is used. - channel_ids : list | None, default: None - The channel ids. If None, all channels are used. - return_in_uV : bool, default: False - If True and the snippets has scaling (gain_to_uV and offset_to_uV properties), - snippets are scaled to uV - - Returns - ------- - np.array - The snippets (num_snippets, num_samples, num_channels) - """ - segment_index = self._check_segment_index(segment_index) - spts = self._snippets_segments[segment_index] - indices = spts.frames_to_indices(start_frame, end_frame) - - return self.get_snippets(indices, channel_ids=channel_ids, return_in_uV=return_in_uV) - - def select_channels(self, channel_ids: list | np.ndarray | tuple) -> "BaseSnippets": - from .channelslice import ChannelSliceSnippets - - return ChannelSliceSnippets(self, channel_ids) - - def _remove_channels(self, remove_channel_ids): - from .channelslice import ChannelSliceSnippets - - new_channel_ids = self.channel_ids[~np.isin(self.channel_ids, remove_channel_ids)] - sub_recording = ChannelSliceSnippets(self, new_channel_ids) - return sub_recording - - def _select_segments(self, segment_indices): - from .segmentutils import SelectSegmentSnippets - - return SelectSegmentSnippets(self, segment_indices=segment_indices) - - def save(self, format="npy", **save_kwargs): - """ - Save a `BaseSnippets` object to a specified format: - - * "npy" - * "memory" - """ - if format == "npy": - from spikeinterface.core.npyfoldersnippets import NpyFolderSnippets - - folder = save_kwargs["folder"] - folder = Path(folder) - - cached = NpyFolderSnippets.write_snippets(self, folder, dtype=save_kwargs.get("dtype", None)) - - elif format == "memory": - snippets_list = [] - spikesframes_list = [] - for i in range(self.get_num_segments()): - spikesframes_list.append(self.get_frames(segment_index=i)) - snippets_list.append(self.get_snippets(segment_index=i)) - - from .numpyextractors import NumpySnippets - - cached = NumpySnippets( - snippets_list=snippets_list, - spikesframes_list=spikesframes_list, - sampling_frequency=self.get_sampling_frequency(), - nbefore=self.nbefore, - channel_ids=self.channel_ids, - ) - if self.has_probe(): - probegroup = self.get_probegroup() - cached.set_probegroup(probegroup) - else: - raise ValueError(f"format {format} not supported") - - return cached - - def get_times(self): - return self.get_frames() / self.sampling_frequency - - -class BaseSnippetsSegment(BaseSegment): - """ - Abstract class representing multichannel snippets - """ - - def __init__(self): - BaseSegment.__init__(self) - - def get_snippets( - self, - indices, - channel_indices: list | None = None, - ) -> np.ndarray: - """ - Return the snippets, optionally for a subset of samples and/or channels - - Parameters - ---------- - indices : list[int] - Indices of the snippets to return - channel_indices : list | None, default: None - Indices of channels to return, or all channels if None - - Returns - ------- - snippets : np.ndarray - Array of snippets, num_snippets x num_samples x num_channels - """ - raise NotImplementedError - - def get_num_snippets(self): - """Returns the number of snippets in this segment - - Returns: - SampleIndex : Number of snippets in the segment - """ - raise NotImplementedError - - def get_frames(self, indices): - """Returns the frames of the snippets in this segment - - Returns: - SampleIndex : Number of samples in the segment - """ - raise NotImplementedError - - def frames_to_indices(self, start_frame: int | None = None, end_frame: int | None = None): - """ - Return the slice of snippets - - Parameters - ---------- - start_frame : int | None, default: None - start sample index, or zero if None - end_frame : int | None, default: None - end_sample, or number of samples if None - - Returns - ------- - snippets : slice - slice of selected snippets - """ - raise NotImplementedError diff --git a/src/spikeinterface/core/channelslice.py b/src/spikeinterface/core/channelslice.py index bfb46da619..a307381921 100644 --- a/src/spikeinterface/core/channelslice.py +++ b/src/spikeinterface/core/channelslice.py @@ -1,7 +1,6 @@ import numpy as np from .baserecording import BaseRecording, BaseRecordingSegment -from .basesnippets import BaseSnippets, BaseSnippetsSegment class ChannelSliceRecording(BaseRecording): @@ -59,6 +58,7 @@ def __init__(self, recording, channel_ids=None, renamed_channel_ids=None): # copy annotation and properties self._parent = recording + recording.copy_metadata(self, only_main=False, ids=self._channel_ids) # change the wiring of the probe if self._parent.has_probe(): @@ -66,7 +66,8 @@ def __init__(self, recording, channel_ids=None, renamed_channel_ids=None): sliced_probegroup = parent_probegroup.get_slice(self._parent_channel_indices) sliced_probegroup.set_global_device_channel_indices(np.arange(len(self._channel_ids))) self.set_probegroup(sliced_probegroup) - recording.copy_metadata(self, only_main=False, ids=self._channel_ids) + # Reset channel groups to original ones to avoid remapping + self.set_channel_groups(recording.get_channel_groups()[self._parent_channel_indices]) # update dump dict self._kwargs = { @@ -111,114 +112,3 @@ def get_traces( parent_indices = self._parent_channel_indices[channel_indices] traces = self._parent_recording_segment.get_traces(start_frame, end_frame, parent_indices) return traces - - -class ChannelSliceSnippets(BaseSnippets): - """ - Class to slice a Snippets object based on channel_ids. - - Do not use this class directly but use `snippets.select_channels(...)` - - """ - - def __init__(self, parent_snippets, channel_ids=None, renamed_channel_ids=None): - if channel_ids is None: - channel_ids = parent_snippets.get_channel_ids() - if renamed_channel_ids is None: - renamed_channel_ids = channel_ids - - self._parent_snippets = parent_snippets - self._channel_ids = np.asarray(channel_ids) - self._renamed_channel_ids = np.asarray(renamed_channel_ids) - - parents_chan_ids = self._parent_snippets.get_channel_ids() - - # some checks - assert all( - chan_id in parents_chan_ids for chan_id in self._channel_ids - ), "ChannelSliceSnippets : channel ids are not all in parents" - assert len(self._channel_ids) == len( - self._renamed_channel_ids - ), "ChannelSliceSnippets: renamed channel_ids must be the same size" - assert ( - self._channel_ids.size == np.unique(self._channel_ids).size - ), "ChannelSliceSnippets : channel_ids are not unique" - - sampling_frequency = parent_snippets.get_sampling_frequency() - - BaseSnippets.__init__( - self, - sampling_frequency=sampling_frequency, - nbefore=parent_snippets.nbefore, - snippet_len=parent_snippets.snippet_len, - channel_ids=self._renamed_channel_ids, - dtype=parent_snippets.get_dtype(), - ) - - self._parent_channel_indices = parent_snippets.ids_to_indices(self._channel_ids) - - # link recording segment - for parent_segment in self._parent_snippets._snippets_segments: - sub_segment = ChannelSliceSnippetsSegment(parent_segment, self._parent_channel_indices) - self.add_snippets_segment(sub_segment) - - # copy annotation and properties - parent_snippets.copy_metadata(self, only_main=False, ids=self._channel_ids) - - # change the wiring of the probe - if self._parent_snippets.has_probe(): - parent_probegroup = self._parent_snippets.get_probegroup() - sliced_probegroup = parent_probegroup.get_slice(self._parent_channel_indices) - sliced_probegroup.set_global_device_channel_indices(np.arange(len(self._channel_ids))) - self.set_probegroup(sliced_probegroup) - - # update dump dict - self._kwargs = { - "parent_snippets": parent_snippets, - "channel_ids": channel_ids, - "renamed_channel_ids": renamed_channel_ids, - } - - -class ChannelSliceSnippetsSegment(BaseSnippetsSegment): - """ - Class to return a channel-sliced segment snippets. - """ - - def __init__(self, parent_snippets_segment, parent_channel_indices): - BaseSnippetsSegment.__init__(self) - self._parent_snippets_segment = parent_snippets_segment - self._parent_channel_indices = parent_channel_indices - - def get_num_snippets(self) -> int: - return self._parent_snippets_segment.get_num_snippets() - - def frames_to_indices(self, start_frame: int | None = None, end_frame: int | None = None): - return self._parent_snippets_segment.frames_to_indices(start_frame, end_frame) - - def get_frames(self, indices=None): - return self._parent_snippets_segment.get_frames(indices) - - def get_snippets( - self, - indices: list[int], - channel_indices: list | None = None, - ) -> np.ndarray: - """ - Return the snippets, optionally for a subset of samples and/or channels - - Parameters - ---------- - indices: list[int] - Indices of the snippets to return - channel_indices: list | None, default: None - Indices of channels to return, or all channels if None - - Returns - ------- - snippets: np.ndarray - Array of snippets, num_snippets x num_samples x num_channels - """ - parent_indices = self._parent_channel_indices[channel_indices] - snippets = self._parent_snippets_segment.get_snippets(indices, parent_indices) - return snippets diff --git a/src/spikeinterface/core/generate.py b/src/spikeinterface/core/generate.py index 28ae361fc0..5e8c352e44 100644 --- a/src/spikeinterface/core/generate.py +++ b/src/spikeinterface/core/generate.py @@ -11,7 +11,6 @@ from probeinterface import Probe, generate_linear_probe, generate_multi_columns_probe from spikeinterface.core import BaseRecording, BaseRecordingSegment, BaseSorting -from .snippets_tools import snippets_from_sorting from .core_tools import define_function_from_class, ms_to_samples, _ensure_seed @@ -602,81 +601,7 @@ def create_sorting_npz(num_seg, file_path): np.savez(file_path, **d) -def generate_snippets( - nbefore=20, - nafter=44, - num_channels=2, - wf_folder=None, - sampling_frequency=30000.0, - durations=[10.325, 3.5], #  in s for 2 segments - set_probe=True, - ndim=2, - num_units=5, - empty_units=None, - **job_kwargs, -): - """ - Generates a synthetic Snippets object. - - Parameters - ---------- - nbefore : int, default: 20 - Number of samples before the peak. - nafter : int, default: 44 - Number of samples after the peak. - num_channels : int, default: 2 - Number of channels. - wf_folder : str | Path | None, default: None - Optional folder to save the waveform snippets. If None, snippets are in memory. - sampling_frequency : float, default: 30000.0 - The sampling frequency of the snippets in Hz. - ndim : int, default: 2 - The number of dimensions of the probe. - num_units : int, default: 5 - The number of units. - empty_units : list | None, default: None - A list of units that will have no spikes. - durations : List[float], default: [10.325, 3.5] - The duration in seconds of each segment in the recording. - The number of segments is determined by the length of this list. - set_probe : bool, default: True - If true, attaches probe to the returned snippets object - **job_kwargs : dict, default: None - Job keyword arguments for `snippets_from_sorting` - - Returns - ------- - snippets : NumpySnippets - The snippets object. - sorting : NumpySorting - The associated sorting object. - """ - recording = generate_recording( - durations=durations, - num_channels=num_channels, - sampling_frequency=sampling_frequency, - ndim=ndim, - set_probe=set_probe, - ) - - sorting = generate_sorting( - num_units=num_units, sampling_frequency=sampling_frequency, durations=durations, empty_units=empty_units - ) - - snippets = snippets_from_sorting( - recording=recording, sorting=sorting, nbefore=nbefore, nafter=nafter, wf_folder=wf_folder, **job_kwargs - ) - - if set_probe: - probe = recording.get_probe() - snippets.set_probe(probe) - - return snippets, sorting - - ## spiketrain zone ## - - def _ensure_firing_rates(firing_rates, num_units, seed): if isinstance(firing_rates, tuple): rng = np.random.default_rng(seed=seed) diff --git a/src/spikeinterface/core/npyfoldersnippets.py b/src/spikeinterface/core/npyfoldersnippets.py deleted file mode 100644 index 31dee04078..0000000000 --- a/src/spikeinterface/core/npyfoldersnippets.py +++ /dev/null @@ -1,103 +0,0 @@ -from pathlib import Path -import json - -from copy import deepcopy - -from probeinterface import read_probeinterface, write_probeinterface - -from .npysnippetsextractor import NpySnippetsExtractor -from .core_tools import ( - define_function_from_class, - make_paths_absolute, - load_properties_from_folder, - save_properties_to_folder, -) - - -class NpyFolderSnippets(NpySnippetsExtractor): - """ - NpyFolderSnippets is an internal format used in spikeinterface. - It is a NpySnippetsExtractor + metadata contained in a folder. - - It is created with the function: `snippets.save(format="npy", folder="/myfolder")` - - Parameters - ---------- - folder_path : str or Path - The path to the folder - - Returns - ------- - snippets : NpyFolderSnippets - The snippets - """ - - mode = "folder" - name = "npyfolder" - - def __init__(self, folder_path): - folder_path = Path(folder_path) - - with open(folder_path / "npy.json", "r") as f: - d = json.load(f) - - if not d["class"].endswith(".NpySnippetsExtractor"): - raise ValueError("This folder is not a binary spikeinterface folder") - - assert d["relative_paths"] - - d = make_paths_absolute(d, folder_path) - - NpySnippetsExtractor.__init__(self, **d["kwargs"]) - - probe_file = folder_path / "probegroup.json" - if probe_file.is_file(): - self._probegroup = read_probeinterface(probe_file) - - load_properties_from_folder(folder_path / "properties", self) - - self._kwargs = dict(folder_path=str(Path(folder_path).absolute())) - self._bin_kwargs = d["kwargs"] - - @staticmethod - def write_snippets(snippets, folder, dtype=None): - - folder = Path(folder) - - if dtype is None: - dtype = snippets.dtype - - file_paths = [folder / f"traces_cached_seg{i}.npy" for i in range(snippets.get_num_segments())] - - if dtype is None: - dtype = snippets.get_dtype() - - # This is weird but for backward compatibility - # maybe this can be removed - NpySnippetsExtractor.write_snippets(snippets=snippets, file_paths=file_paths, dtype=dtype) - cached = NpySnippetsExtractor( - file_paths=file_paths, - sampling_frequency=snippets.get_sampling_frequency(), - channel_ids=snippets.get_channel_ids(), - nbefore=snippets.nbefore, - gain_to_uV=snippets.get_channel_gains(), - offset_to_uV=snippets.get_channel_offsets(), - ) - cached.dump(folder / "npy.json", relative_to=folder) - - save_properties_to_folder(folder / "properties", snippets) - - if snippets.has_probe(): - probegroup = snippets.get_probegroup() - write_probeinterface(folder / "probegroup.json", probegroup) - - cached = NpyFolderSnippets(folder_path=folder) - # important backward compatibility : annoations are handled (sadly) only is this file - # so we need to set then here (sad hack) - cached._annotations = deepcopy({k: snippets._annotations[k] for k in snippets._annotations.keys()}) - cached.dump(folder / "si_folder.json", relative_to=folder) - - return cached - - -read_npy_snippets_folder = define_function_from_class(source_class=NpyFolderSnippets, name="read_npy_snippets_folder") diff --git a/src/spikeinterface/core/npysnippetsextractor.py b/src/spikeinterface/core/npysnippetsextractor.py deleted file mode 100644 index 88c2ea9de3..0000000000 --- a/src/spikeinterface/core/npysnippetsextractor.py +++ /dev/null @@ -1,166 +0,0 @@ -from pathlib import Path - -import numpy as np - -from spikeinterface.core.basesnippets import BaseSnippets, BaseSnippetsSegment -from .core_tools import define_function_from_class - - -class NpySnippetsExtractor(BaseSnippets): - """ - Dead simple and super light format based on the NPY numpy format. - - It is in fact an archive of several .npy format. - All spike are store in two columns maner index+labels - """ - - mode = "file" - name = "npy" - - def __init__( - self, file_paths, sampling_frequency, channel_ids=None, nbefore=None, gain_to_uV=None, offset_to_uV=None - ): - if not isinstance(file_paths, list): - file_paths = [file_paths] - - num_segments = len(file_paths) - data = np.load(file_paths[0], mmap_mode="r") - - if channel_ids is None: - channel_ids = np.arange(data["snippet"].shape[2]) - - BaseSnippets.__init__( - self, - sampling_frequency, - nbefore=nbefore, - snippet_len=data["snippet"].shape[1], - dtype=data["snippet"].dtype, - channel_ids=channel_ids, - ) - - for i in range(num_segments): - snp_segment = NpySnippetsSegment(file_paths[i]) - self.add_snippets_segment(snp_segment) - - if gain_to_uV is not None: - self.set_channel_gains(gain_to_uV) - - if offset_to_uV is not None: - self.set_channel_offsets(offset_to_uV) - - self._kwargs = { - "file_paths": [str(Path(f).absolute()) for f in file_paths], - "sampling_frequency": sampling_frequency, - "channel_ids": channel_ids, - "nbefore": nbefore, - "gain_to_uV": gain_to_uV, - "offset_to_uV": offset_to_uV, - } - - @staticmethod - def write_snippets(snippets, file_paths, dtype=None): - """ - Save snippet extractor in binary .npy format. - - Parameters - ---------- - snippets: SnippetsExtractor - The snippets extractor object to be saved in .npy format - file_paths: str - The paths to the files. - dtype: None, str or dtype - Typecode or data-type to which the snippets will be cast. - {} - """ - if not isinstance(file_paths, list): - file_paths = [file_paths] - if dtype is None: - dtype = snippets.dtype - assert len(file_paths) == snippets.get_num_segments() - snippets_t = np.dtype( - [("frame", np.int64), ("snippet", dtype, (snippets.snippet_len, snippets.get_num_channels()))] - ) - - for i in range(snippets.get_num_segments()): - n = snippets.get_num_snippets(i) - arr = np.empty(n, dtype=snippets_t, order="F") - arr["frame"] = snippets.get_frames(segment_index=i) - arr["snippet"] = snippets.get_snippets(segment_index=i).astype(dtype, copy=False) - file_paths[i].parent.mkdir(parents=True, exist_ok=True) - np.save(file_paths[i], arr) - - -class NpySnippetsSegment(BaseSnippetsSegment): - def __init__(self, file): - BaseSnippetsSegment.__init__(self) - - npy = np.load(file, mmap_mode="r") - self._snippets = npy["snippet"] - self._spikestimes = npy["frame"] - - def get_snippets( - self, - indices: list[int], - channel_indices: list | None = None, - ) -> np.ndarray: - """ - Return the snippets, optionally for a subset of samples and/or channels - - Parameters - ---------- - indices: list[int] - Indices of the snippets to return, or all if None - channel_indices: list | None, default: None - Indices of channels to return, or all channels if None - - Returns - ------- - snippets: np.ndarray - Array of snippets, num_snippets x num_samples x num_channels - """ - if indices is None: - return self._snippets[:, :, channel_indices] - return self._snippets[indices, :, channel_indices] - - def get_num_snippets(self): - return self._spikestimes.shape[0] - - def frames_to_indices(self, start_frame: int | None = None, end_frame: int | None = None): - """ - Return the slice of snippets - - Parameters - ---------- - start_frame: int | None, default: None - start sample index, or zero if None - end_frame: int | None, default: None - end_sample, or number of samples if None - - Returns - ------- - snippets: slice - slice of selected snippets - """ - # must be implemented in subclass - if start_frame is None: - init = 0 - else: - init = np.searchsorted(self._spikestimes, start_frame, side="left") - if end_frame is None: - endi = self._spikestimes.shape[0] - else: - endi = np.searchsorted(self._spikestimes, end_frame, side="left") - return slice(init, endi, 1) - - def get_frames(self, indices=None): - """Returns the frames of the snippets in this segment - - Returns: - SampleIndex: Number of samples in the segment - """ - if indices is None: - return self._spikestimes - raise self._spikestimes[indices] - - -read_npy_snippets = define_function_from_class(source_class=NpySnippetsExtractor, name="read_npy_snippets") diff --git a/src/spikeinterface/core/numpyextractors.py b/src/spikeinterface/core/numpyextractors.py index 7f5f5321b9..b62cfc584e 100644 --- a/src/spikeinterface/core/numpyextractors.py +++ b/src/spikeinterface/core/numpyextractors.py @@ -7,8 +7,6 @@ SpikeVectorSortingSegment, BaseEvent, BaseEventSegment, - BaseSnippets, - BaseSnippetsSegment, ) from .base import minimum_spike_dtype from .core_tools import make_shared_array @@ -624,137 +622,3 @@ def get_events(self, channel_id, start_time, end_time): if end_time is not None: events = events[times <= end_time] return events - - -class NumpySnippets(BaseSnippets): - """ - In memory recording. - Contrary to previous version this class does not handle npy files. - - Parameters - ---------- - snippets_list : list of array or array (if mono segment) - The snippets to instantiate a mono or multisegment basesnippet - spikesframes_list : list of array or array (if mono segment) - Frame of each snippet - sampling_frequency : float - The sampling frequency in Hz - - channel_ids : list - An optional list of channel_ids. If None, linear channels are assumed - """ - - def __init__(self, snippets_list, spikesframes_list, sampling_frequency, nbefore=None, channel_ids=None): - if isinstance(snippets_list, list): - assert all(isinstance(e, np.ndarray) for e in snippets_list), "must give a list of numpy array" - else: - assert isinstance(snippets_list, np.ndarray), "must give a list of numpy array" - snippets_list = [snippets_list] - if isinstance(spikesframes_list, list): - assert all(isinstance(e, np.ndarray) for e in spikesframes_list), "must give a list of numpy array" - else: - assert isinstance(spikesframes_list, np.ndarray), "must give a list of numpy array" - spikesframes_list = [spikesframes_list] - - dtype = snippets_list[0].dtype - assert all(dtype == ts.dtype for ts in snippets_list) - - if channel_ids is None: - channel_ids = np.arange(snippets_list[0].shape[2]) - else: - channel_ids = np.asarray(channel_ids) - assert channel_ids.size == snippets_list[0].shape[2] - BaseSnippets.__init__( - self, - sampling_frequency, - nbefore=nbefore, - snippet_len=snippets_list[0].shape[1], - channel_ids=channel_ids, - dtype=dtype, - ) - - self._serializability["memory"] = False - self._serializability["json"] = False - self._serializability["pickle"] = False - - for snippets, spikesframes in zip(snippets_list, spikesframes_list): - snp_segment = NumpySnippetsSegment(snippets, spikesframes) - self.add_snippets_segment(snp_segment) - - self._kwargs = { - "snippets_list": snippets_list, - "spikesframes_list": spikesframes_list, - "nbefore": nbefore, - "sampling_frequency": sampling_frequency, - "channel_ids": channel_ids, - } - - -class NumpySnippetsSegment(BaseSnippetsSegment): - def __init__(self, snippets, spikesframes): - BaseSnippetsSegment.__init__(self) - self._snippets = snippets - self._spikestimes = spikesframes - - def get_snippets( - self, - indices, - channel_indices: list | None = None, - ) -> np.ndarray: - """ - Return the snippets, optionally for a subset of samples and/or channels - - Parameters - ---------- - indices : list[int] - Indices of the snippets to return - channel_indices : list | None, default: None - Indices of channels to return, or all channels if None - - Returns - ------- - snippets : np.ndarray - Array of snippets, num_snippets x num_samples x num_channels - """ - if indices is None: - return self._snippets[:, :, channel_indices] - return self._snippets[indices, :, channel_indices] - - def get_num_snippets(self): - return self._spikestimes.shape[0] - - def frames_to_indices(self, start_frame: int | None = None, end_frame: int | None = None): - """ - Return the slice of snippets - - Parameters - ---------- - start_frame : int | None, default: None - start sample index, or zero if None - end_frame : int | None, default: None - end_sample, or number of samples if None - Returns - ------- - snippets : slice - slice of selected snippets - """ - # must be implemented in subclass - if start_frame is None: - init = 0 - else: - init = np.searchsorted(self._spikestimes, start_frame, side="left") - if end_frame is None: - endi = self._spikestimes.shape[0] - else: - endi = np.searchsorted(self._spikestimes, end_frame, side="left") - return slice(init, endi, 1) - - def get_frames(self, indices=None): - """Returns the frames of the snippets in this segment - - Returns: - SampleIndex : Number of samples in the segment - """ - if indices is None: - return self._spikestimes - raise self._spikestimes[indices] diff --git a/src/spikeinterface/core/snippets_tools.py b/src/spikeinterface/core/snippets_tools.py deleted file mode 100644 index 499e4f91f3..0000000000 --- a/src/spikeinterface/core/snippets_tools.py +++ /dev/null @@ -1,66 +0,0 @@ -from .job_tools import fix_job_kwargs -from .waveform_tools import extract_waveforms_to_buffers -from .numpyextractors import NumpySnippets - - -def snippets_from_sorting(recording, sorting, nbefore=20, nafter=44, wf_folder=None, **job_kwargs): - """ - Extract snippets from recording and sorting instances - - Parameters - ---------- - recording: BaseRecording - The recording to get snippets from - sorting: BaseSorting - The sorting to get the frames from - nbefore: int - N samples before spike - nafter: int - N samples after spike - wf_folder: None, str or path - Folder to save npy files, if None shared_memory will be used to extract waveforms - Returns - ------- - snippets: NumpySnippets - Snippets extractor created - """ - job_kwargs = fix_job_kwargs(job_kwargs) - spikes = sorting.to_spike_vector(concatenated=False) - - peaks2 = sorting.to_spike_vector() - peaks2["unit_index"] = 0 - - if wf_folder is None: - mode = "shared_memory" - folder = None - else: - mode = "memmap" - folder = wf_folder - - wfs_arrays = extract_waveforms_to_buffers( - recording, - peaks2, - [0], - nbefore, - nafter, - mode=mode, - return_in_uV=False, - folder=folder, - dtype=recording.get_dtype(), - sparsity_mask=None, - copy=True, - **job_kwargs, - ) - wfs = [] - for i in range(recording.get_num_segments()): - wfs.append(wfs_arrays[0][peaks2["segment_index"] == i, :, :]) # extract class zero - - nse = NumpySnippets( - snippets_list=wfs, - spikesframes_list=[s["sample_index"] for s in spikes], - sampling_frequency=recording.get_sampling_frequency(), - nbefore=nbefore, - channel_ids=recording.get_channel_ids(), - ) - nse.set_property("gain_to_uV", recording.get_property("gain_to_uV")) - return nse diff --git a/src/spikeinterface/core/tests/test_basesnippets.py b/src/spikeinterface/core/tests/test_basesnippets.py deleted file mode 100644 index 585672a1ea..0000000000 --- a/src/spikeinterface/core/tests/test_basesnippets.py +++ /dev/null @@ -1,198 +0,0 @@ -""" -test for BaseSnippets are done with NumpySnippets. -but check only for BaseRecording general methods. -""" - -from pathlib import Path -import pytest -import numpy as np -from numpy.testing import assert_raises - -from probeinterface import Probe -from spikeinterface.core import generate_snippets, load, NumpySnippets -from spikeinterface.core.npysnippetsextractor import NpySnippetsExtractor -from spikeinterface.core.base import BaseExtractor - - -def test_BaseSnippets(create_cache_folder): - cache_folder = create_cache_folder - duration = [4, 3] - num_channels = 3 - nbefore = 20 - nafter = 44 - - nse, sorting = generate_snippets( - durations=duration, num_channels=num_channels, nbefore=nbefore, nafter=nafter, wf_folder=None - ) - num_seg = len(duration) - file_paths = [cache_folder / f"test_base_snippets_{i}.npy" for i in range(num_seg)] - - NpySnippetsExtractor.write_snippets(nse, file_paths) - snippets = NpySnippetsExtractor( - file_paths, - sampling_frequency=nse.get_sampling_frequency(), - channel_ids=nse.channel_ids, - nbefore=nse.nbefore, - gain_to_uV=nse.get_property("gain_to_uV"), - offset_to_uV=nse.get_property("offset_to_uV"), - ) - - assert snippets.get_num_segments() == len(duration) - assert snippets.get_num_channels() == num_channels - - assert np.all(snippets.ids_to_indices(["0", "1", "2"]) == [0, 1, 2]) - assert np.all(snippets.ids_to_indices(["0", "1", "2"], prefer_slice=True) == slice(0, 3, None)) - - # annotations / properties - snippets.annotate(gre="ta") - assert snippets.get_annotation("gre") == "ta" - - snippets.set_channel_groups([0, 0, 1]) - snippets.set_property("quality", [1.0, 3.3, np.nan]) - values = snippets.get_property("quality") - assert np.all( - values[:2] - == [ - 1.0, - 3.3, - ] - ) - - # missing property - snippets.set_property("string_property", ["ciao", "bello"], ids=["0", "1"]) - values = snippets.get_property("string_property") - assert values[2] == "" - - # setting an different type raises an error - assert_raises( - Exception, - snippets.set_property, - key="string_property_nan", - values=["hola", "chabon"], - ids=["0", "1"], - missing_value=np.nan, - ) - - # int properties without missing values raise an error - assert_raises(Exception, snippets.set_property, key="int_property", values=[5, 6], ids=[1, 2]) - - snippets.set_property("int_property", [5, 6], ids=["1", "2"], missing_value=200) - values = snippets.get_property("int_property") - assert values.dtype.kind == "i" - - times0 = snippets.get_frames(segment_index=0) - - seg0_times = sorting.to_spike_vector(concatenated=False)[0]["sample_index"] - - assert np.array_equal(seg0_times, times0) - - # dump/load dict - d = snippets.to_dict() - snippets2 = BaseExtractor.from_dict(d) - snippets3 = load(d) - - # dump/load json - snippets.dump_to_json(cache_folder / "test_BaseSnippets.json") - snippets2 = load(cache_folder / "test_BaseSnippets.json") - snippets3 = load(cache_folder / "test_BaseSnippets.json") - - # dump/load pickle - snippets.dump_to_pickle(cache_folder / "test_BaseSnippets.pkl") - snippets2 = load(cache_folder / "test_BaseSnippets.pkl") - snippets3 = load(cache_folder / "test_BaseSnippets.pkl") - - # dump/load dict - relative - d = snippets.to_dict(relative_to=cache_folder, recursive=True) - snippets2 = BaseExtractor.from_dict(d, base_folder=cache_folder) - snippets3 = load(d, base_folder=cache_folder) - - # dump/load json - snippets.dump_to_json(cache_folder / "test_BaseSnippets_rel.json", relative_to=cache_folder) - snippets2 = load(cache_folder / "test_BaseSnippets_rel.json", base_folder=cache_folder) - snippets3 = load(cache_folder / "test_BaseSnippets_rel.json", base_folder=cache_folder) - - # cache to npy - folder = cache_folder / "simple_snippets" - print(folder) - snippets.save(format="npy", folder=folder) - snippets2 = load(folder) - assert "quality" in snippets2.get_property_keys() - values = snippets2.get_property("quality") - assert values[0] == 1.0 - assert values[1] == 3.3 - assert np.isnan(values[2]) - - groups = snippets2.get_channel_groups() - assert np.array_equal(groups, [0, 0, 1]) - - # but also possible - snippets3 = load(cache_folder / "simple_snippets") - - # cache to memory - snippets4 = snippets3.save(format="memory") - - waveform4 = snippets4.get_snippets(segment_index=0) - waveform = snippets.get_snippets(segment_index=0) - assert np.array_equal(waveform4, waveform) - - # set/get Probe only 2 channels - probe = Probe(ndim=2) - positions = [[0.0, 0.0], [0.0, 15.0], [0, 30.0]] - probe.set_contacts(positions=positions, shapes="circle", shape_params={"radius": 5}) - probe.set_device_channel_indices([2, -1, 0]) - probe.create_auto_shape() - - snippets_p = snippets.select_channels_with_probe(probe, group_mode="auto") - positions2 = snippets_p.get_channel_locations() - assert np.array_equal(positions2, [[0, 30.0], [0.0, 0.0]]) - - probe2 = snippets_p.get_probe() - positions3 = probe2.contact_positions - assert np.array_equal(positions2, positions3) - - assert np.array_equal(probe2.device_channel_indices, [0, 1]) - - # test save with probe - folder = cache_folder / "simple_snippets3" - snippets2 = snippets_p.save(folder=folder) - snippets2 = load(folder) - probe2 = snippets2.get_probe() - assert np.array_equal(probe2.contact_positions, [[0, 30.0], [0.0, 0.0]]) - positions2 = snippets_p.get_channel_locations() - assert np.array_equal(positions2, [[0, 30.0], [0.0, 0.0]]) - wavefroms2 = snippets2.get_snippets(segment_index=0) - assert np.array_equal(wavefroms2, snippets_p.get_snippets(segment_index=0)) - - # test return_scale - sampling_frequency = 30000 - - spikesframes_list = np.arange(100) - snippets_list = np.empty((100, 64, 5), dtype="int16") - - nse_int16 = NumpySnippets( - snippets_list=snippets_list, spikesframes_list=spikesframes_list, sampling_frequency=30000 - ) - - assert nse_int16.get_dtype() == "int16" - - waveforms_int16 = nse_int16.get_snippets() - assert waveforms_int16.dtype == "int16" - - # return_in_uV raises error when no gain_to_uV/offset_to_uV properties - with pytest.raises(ValueError): - waveforms_float32 = nse_int16.get_snippets(return_in_uV=True) - - # Set properties and test return_in_uV parameter - nse_int16.set_property("gain_to_uV", [0.195] * 5) - nse_int16.set_property("offset_to_uV", [0.0] * 5) - - waveforms_float32 = nse_int16.get_snippets(return_in_uV=True) - assert waveforms_float32.dtype == "float32" - - -if __name__ == "__main__": - import tempfile - - tmp_path = Path(tempfile.mkdtemp()) - - test_BaseSnippets(tmp_path) diff --git a/src/spikeinterface/core/tests/test_channelslicerecording.py b/src/spikeinterface/core/tests/test_channelslicerecording.py index da255d6953..0f8247e85e 100644 --- a/src/spikeinterface/core/tests/test_channelslicerecording.py +++ b/src/spikeinterface/core/tests/test_channelslicerecording.py @@ -1,7 +1,3 @@ -import shutil -from pathlib import Path - -from numpy import rec import pytest import numpy as np diff --git a/src/spikeinterface/core/tests/test_npyfoldersnippets.py b/src/spikeinterface/core/tests/test_npyfoldersnippets.py deleted file mode 100644 index ebf56e3985..0000000000 --- a/src/spikeinterface/core/tests/test_npyfoldersnippets.py +++ /dev/null @@ -1,33 +0,0 @@ -import pytest - -from pathlib import Path -import shutil - -from spikeinterface.core import load -from spikeinterface.core import generate_snippets - - -@pytest.fixture(scope="module") -def cache_folder_creation(tmp_path_factory): - cache_folder = tmp_path_factory.mktemp("cache_folder") - return cache_folder - - -def test_NpyFolderSnippets(cache_folder_creation): - - cache_folder = cache_folder_creation - snippets, _ = generate_snippets(num_channels=10, durations=[2.0, 1.0]) - folder = cache_folder / "npy_folder_1" - - if folder.is_dir(): - shutil.rmtree(folder) - - saved_snippets = snippets.save(folder=folder) - print(snippets) - - loaded_snippets = load(folder) - print(loaded_snippets) - - -if __name__ == "__main__": - test_NpyFolderSnippets() diff --git a/src/spikeinterface/core/tests/test_npysnippetsextractor.py b/src/spikeinterface/core/tests/test_npysnippetsextractor.py deleted file mode 100644 index c3fbfcd885..0000000000 --- a/src/spikeinterface/core/tests/test_npysnippetsextractor.py +++ /dev/null @@ -1,25 +0,0 @@ -import pytest -from pathlib import Path - -from spikeinterface.core import NpySnippetsExtractor -from spikeinterface.core import generate_snippets - - -def test_NpySnippetsExtractor(create_cache_folder): - cache_folder = create_cache_folder - segment_durations = [2, 5] - sampling_frequency = 30000 - file_path = [cache_folder / f"test_NpySnippetsExtractor_{i}.npy" for i in range(len(segment_durations))] - - snippets, _ = generate_snippets(sampling_frequency=sampling_frequency, durations=segment_durations) - - NpySnippetsExtractor.write_snippets(snippets, file_path) - npy_snippets = NpySnippetsExtractor(file_path, sampling_frequency=sampling_frequency) - - assert npy_snippets.get_num_segments() == snippets.get_num_segments() - assert npy_snippets.get_num_snippets(1) == snippets.get_num_snippets(1) - assert npy_snippets.snippet_len == snippets.snippet_len - - -if __name__ == "__main__": - test_NpySnippetsExtractor() diff --git a/src/spikeinterface/extractors/extractor_classes.py b/src/spikeinterface/extractors/extractor_classes.py index a2b11726a8..21c383309f 100644 --- a/src/spikeinterface/extractors/extractor_classes.py +++ b/src/spikeinterface/extractors/extractor_classes.py @@ -5,11 +5,9 @@ NumpyRecording, NpzSortingExtractor, NumpySorting, - NpySnippetsExtractor, read_binary, read_zarr, read_npz_sorting, - read_npy_snippets, ) from spikeinterface.core.zarrextractors import ( @@ -72,9 +70,6 @@ read_shybrid_sorting, ) -# snippers -from .waveclussnippetstextractors import WaveClusSnippetsExtractor, read_waveclus_snippets - # misc from .alfsortingextractor import ALFSortingExtractor, read_alf_sorting @@ -154,11 +149,6 @@ # events only from neo _event_extractor_full_dict = neo_event_extractors_dict -_snippets_extractor_full_dict = { - NpySnippetsExtractor: dict(wrapper_string="read_npy_snippets", wrapper_class=read_npy_snippets), - WaveClusSnippetsExtractor: dict(wrapper_string="read_waveclus_snippets", wrapper_class=read_waveclus_snippets), -} - ############################################################################################################ # Organize the possible extractors into a user facing format with keys being extractor names # (e.g. 'intan' , 'kilosort') and values being the appropriate Extractor class returned as its wrapper @@ -181,24 +171,18 @@ event_class.__name__.replace("Event", "").replace("Extractor", "").lower(): event_func["wrapper_class"] for event_class, event_func in _event_extractor_full_dict.items() } -snippets_extractor_full_dict = { - snippets_class.__name__.replace("Snippets", "").replace("Extractor", "").lower(): snippets_func["wrapper_class"] - for snippets_class, snippets_func in _snippets_extractor_full_dict.items() -} # we only do the functions in the init rather than pull in the classes __all__ = [func["wrapper_string"] for func in _recording_extractor_full_dict.values()] __all__ += [func["wrapper_string"] for func in _sorting_extractor_full_dict.values()] __all__ += [func["wrapper_string"] for func in _event_extractor_full_dict.values()] -__all__ += [func["wrapper_string"] for func in _snippets_extractor_full_dict.values()] __all__.extend( [ "read_nwb", # convenience function for multiple nwb formats "recording_extractor_full_dict", "sorting_extractor_full_dict", "event_extractor_full_dict", - "snippets_extractor_full_dict", "read_binary", # convenience function for binary formats "read_zarr", "read_neuroscope", # convenience function for neuroscope diff --git a/src/spikeinterface/extractors/waveclussnippetstextractors.py b/src/spikeinterface/extractors/waveclussnippetstextractors.py deleted file mode 100644 index 72bdbd089a..0000000000 --- a/src/spikeinterface/extractors/waveclussnippetstextractors.py +++ /dev/null @@ -1,154 +0,0 @@ -from pathlib import Path -import numpy as np - -from spikeinterface.core import BaseSnippets, BaseSnippetsSegment -from spikeinterface.core.core_tools import define_function_from_class -from .matlabhelpers import MatlabHelper - - -class WaveClusSnippetsExtractor(MatlabHelper, BaseSnippets): - - def __init__(self, file_path): - file_path = Path(file_path) if isinstance(file_path, str) else file_path - MatlabHelper.__init__(self, file_path) - wc_snippets = self._getfield("spikes") - # handle both types of waveclus results - - # the spikes can be in the times_file - if file_path.name.startswith("times_"): - times = self._getfield("cluster_class")[:, 1] - elif file_path.name.endswith("_spikes.mat"): - times = np.ravel(self._getfield("index")) - else: - raise ("Filename not compatible with waveclus file.") - - sampling_frequency = float(self._getfield("par/sr")) - pre = int(self._getfield("par/w_pre")) - 1 - post = int(self._getfield("par/w_post")) + 1 - - sp_len = pre + post - - nchannels = int(wc_snippets.shape[1] / sp_len) - # waveclus use: #snippets,#concatenated_samples(sample * nchannels) - # spikeinterface use: #snippets,#samples, #nchannels - snippets = np.dstack(np.array_split(wc_snippets, nchannels, 1)) - - BaseSnippets.__init__( - self, - sampling_frequency=sampling_frequency, - nbefore=pre, - snippet_len=sp_len, - dtype=wc_snippets.dtype, - channel_ids=np.arange(nchannels), - ) - - snp_segment = WaveClustSnippetsSegment( - snippets=snippets, spikesframes=np.round(times * (sampling_frequency / 1000)) - ) - self.add_snippets_segment(snp_segment) - - self._kwargs = {"file_path": str(Path(file_path).absolute())} - - @staticmethod - def write_snippets(snippets_extractor, save_file_path): - assert snippets_extractor.is_aligned(), "Waveclus requires aligned spikes" - save_file_path = Path(save_file_path) - assert save_file_path.name.endswith("_spikes.mat"), "Waveclus snippets files should end with _spikes.mat" - frame_to_ms = snippets_extractor.get_sampling_frequency() / 1000 - index = np.concatenate( - [ - snippets_extractor.get_frames(segment_index=sinx) * frame_to_ms - for sinx in range(snippets_extractor.get_num_segments()) - ] - ) - spikes = np.concatenate( - [ - snippets_extractor.get_snippets(segment_index=sinx) * frame_to_ms - for sinx in range(snippets_extractor.get_num_segments()) - ] - ) - spikes = np.swapaxes(spikes, 1, 2).reshape([spikes.shape[0], spikes.shape[1] * spikes.shape[2]], order="C") - par = dict( - sr=snippets_extractor.get_sampling_frequency(), - w_pre=snippets_extractor.nbefore + 1, # waveclus includes the peak in the pre samples - w_post=snippets_extractor.nafter - 1, - ) - MatlabHelper.write_dict_to_mat( - mat_file_path=save_file_path, dict_to_write={"index": index, "spikes": spikes, "par": par} - ) - - -class WaveClustSnippetsSegment(BaseSnippetsSegment): - def __init__(self, snippets, spikesframes): - BaseSnippetsSegment.__init__(self) - self._snippets = snippets - self._spikestimes = spikesframes - - def get_snippets( - self, - indices, - channel_indices: list | None = None, - ) -> np.ndarray: - """ - Return the snippets, optionally for a subset of samples and/or channels - - Parameters - ---------- - indices: list[int] - Indices of the snippets to return - channel_indices: list | None, default: None - Indices of channels to return, or all channels if None - - Returns - ------- - snippets: np.ndarray - Array of snippets, num_snippets x num_samples x num_channels - """ - if indices is None: - return self._snippets[:, :, channel_indices] - return self._snippets[indices, :, channel_indices] - - def get_num_snippets(self): - return self._spikestimes.shape[0] - - def frames_to_indices(self, start_frame: int | None = None, end_frame: int | None = None): - """ - Return the slice of snippets - - Parameters - ---------- - start_frame: int | None, default: None - start sample index, or zero if Non - end_frame: int | None, default: None - end_sample, or number of samples if None - - Returns - ------- - snippets: slice - slice of selected snippets - """ - # must be implemented in subclass - if start_frame is None: - init = 0 - else: - init = np.searchsorted(self._spikestimes, start_frame, side="left") - if end_frame is None: - endi = self._spikestimes.shape[0] - else: - endi = np.searchsorted(self._spikestimes, end_frame, side="left") - return slice(init, endi, 1) - - def get_frames(self, indices=None): - """Returns the frames of the snippets in this segment - - Returns: - SampleIndex: Number of samples in the segment - """ - if indices is None: - return self._spikestimes - return self._spikestimes[indices] - - -read_waveclus_snippets = define_function_from_class( - source_class=WaveClusSnippetsExtractor, name="read_waveclus_snippets" -) diff --git a/src/spikeinterface/generation/__init__.py b/src/spikeinterface/generation/__init__.py index e0fb98dbbd..ea4cccb911 100644 --- a/src/spikeinterface/generation/__init__.py +++ b/src/spikeinterface/generation/__init__.py @@ -34,7 +34,6 @@ from spikeinterface.core.generate import ( generate_recording, generate_sorting, - generate_snippets, generate_templates, generate_recording_by_size, generate_ground_truth_recording, diff --git a/src/spikeinterface/preprocessing/__init__.py b/src/spikeinterface/preprocessing/__init__.py index 957b5a97cc..502c3748da 100644 --- a/src/spikeinterface/preprocessing/__init__.py +++ b/src/spikeinterface/preprocessing/__init__.py @@ -21,6 +21,3 @@ ) from .detect_artifacts import detect_artifact_periods, detect_artifact_periods_by_envelope, detect_saturation_periods - -# for snippets -from .align_snippets import AlignSnippets diff --git a/src/spikeinterface/preprocessing/align_snippets.py b/src/spikeinterface/preprocessing/align_snippets.py deleted file mode 100644 index e6227cc628..0000000000 --- a/src/spikeinterface/preprocessing/align_snippets.py +++ /dev/null @@ -1,113 +0,0 @@ -import numpy as np - -from spikeinterface.core import BaseSnippets, BaseSnippetsSegment - - -class AlignSnippets(BaseSnippets): - installation_mesg = "" # err - - def __init__(self, snippets, new_nbefore, new_nafter, mode="main_peak", interpolate=1, det_sign=0): - assert isinstance(snippets, BaseSnippets), "'snippets' must be a SnippetsExtractor" - assert mode in ("ch_peak", "main_peak"), "mode must be " "ch_peak" " or " "main_peak" "" - - BaseSnippets.__init__( - self, - sampling_frequency=snippets.get_sampling_frequency(), - nbefore=new_nbefore, - snippet_len=new_nbefore + new_nafter, - channel_ids=snippets.channel_ids, - dtype=snippets.get_dtype(), - ) - assert self.snippet_len >= new_nbefore + new_nafter, "snippet_len smaller than new_nbefore+new_nafter" - - snippets.copy_metadata(self, only_main=False, ids=None) - self._parent = snippets - - for i in range(snippets.get_num_segments()): - self.add_snippets_segment( - AlignSnippetsSegment( - snippets._snippets_segments[i], - snippets.snippet_len, - new_nbefore, - new_nafter, - mode, - interpolate, - det_sign, - ) - ) - - self._kwargs = dict( - snippets=snippets, - new_nbefore=new_nbefore, - new_nafter=new_nafter, - mode=mode, - interpolate=interpolate, - det_sign=det_sign, - ) - - -class AlignSnippetsSegment(BaseSnippetsSegment): - def __init__(self, parent_snippets_segment, org_splen, new_nbefore, new_nafter, mode, interpolate, det_sign): - BaseSnippetsSegment.__init__(self) - self.parent_snippets_segment = parent_snippets_segment - self._interpolate = interpolate - self._new_nbefore = new_nbefore - self._new_nafter = new_nafter - self._org_splen = org_splen - self._det_sign = det_sign - start_search = int(new_nbefore * interpolate) - end_search = int((org_splen - new_nafter + 1) * interpolate) - if det_sign == 0: - self._find_peak = lambda x: start_search + np.argmax(np.abs(x[:, start_search:end_search, :]), axis=1) - elif det_sign > 0: - self._find_peak = lambda x: start_search + np.argmax(x[:, start_search:end_search, :], axis=1) - else: - self._find_peak = lambda x: start_search + np.argmin(x[:, start_search:end_search, :], axis=1) - self._mode = mode - if mode == "main_peak": - if det_sign == 0: - self._find_main_ch = lambda x: np.argmax(np.abs(x)) - elif det_sign > 0: - self._find_main_ch = np.argmax - else: - self._find_main_ch = np.argmin - - def get_snippets( - self, - indices=None, - end_frame: int | None = None, - channel_indices: list | None = None, - ) -> np.ndarray: - snippets = self.parent_snippets_segment.get_snippets(indices=indices, channel_indices=channel_indices) - if self._interpolate > 1: - xs = np.arange(0, self._org_splen, 1 / self._interpolate) - from scipy.interpolate import CubicSpline - - snippets = CubicSpline(xs[:: self._interpolate], snippets, axis=1, bc_type="natural")(xs) - - peaks_pos = self._find_peak(snippets) - - aligned_snippets = np.empty( - [snippets.shape[0], self._new_nbefore + self._new_nafter, snippets.shape[2]], - ) - inp = self._interpolate - pres = self._new_nbefore * self._interpolate - posts = self._new_nafter * self._interpolate - if self._mode == "main_peak": - for i, pos in enumerate(peaks_pos): - jpeak = pos[self._find_main_ch([snippets[i, p, ch] for ch, p in enumerate(pos)])] - aligned_snippets[i, :, :] = snippets[i, jpeak - pres : jpeak + posts : inp, :] - else: - for i, pos in enumerate(peaks_pos): - for chi, chpos in enumerate(pos): - aligned_snippets[i, :, chi] = snippets[i, chpos - pres : chpos + posts : inp, chi] - return aligned_snippets - - def get_num_snippets(self): - return self.parent_snippets_segment.get_num_snippets() - - def get_frames(self, indices): - return self.parent_snippets_segment.get_num_snippets(indices) - - def frames_to_indices(self, start_frame: int | None = None, end_frame: int | None = None): - return self.parent_snippets_segment.frames_to_indices(start_frame, end_frame) diff --git a/src/spikeinterface/preprocessing/tests/test_align_snippets.py b/src/spikeinterface/preprocessing/tests/test_align_snippets.py deleted file mode 100644 index 8347fa49ca..0000000000 --- a/src/spikeinterface/preprocessing/tests/test_align_snippets.py +++ /dev/null @@ -1,39 +0,0 @@ -""" -test for BaseSnippets are done with NumpySnippets. -but check only for BaseRecording general methods. -""" - -from spikeinterface.core import generate_snippets -from spikeinterface.preprocessing.align_snippets import AlignSnippets - - -def test_AlignSnippets(): - duration = [4, 3] - num_channels = 3 - nbefore = 20 - nafter = 44 - - snippets, _ = generate_snippets( - durations=duration, num_channels=num_channels, nbefore=nbefore, nafter=nafter, wf_folder=None - ) - # simple realigment test - alined_snippets = AlignSnippets( - snippets, new_nbefore=10, new_nafter=44, mode="main_peak", interpolate=5, det_sign=0 - ) - spikes = alined_snippets.get_snippets(segment_index=0) - - # test interpolation - snippets_ch_peak = AlignSnippets( - snippets, new_nbefore=15, new_nafter=30, mode="ch_peak", interpolate=5, det_sign=-1 - ) - spikes = snippets_ch_peak.get_snippets(segment_index=0) - - # test interpolation - snippets_interpolated = AlignSnippets( - snippets, new_nbefore=15, new_nafter=20, mode="ch_peak", interpolate=5, det_sign=-1 - ) - spikes = snippets_interpolated.get_snippets(segment_index=0) - - -if __name__ == "__main__": - test_AlignSnippets() diff --git a/src/spikeinterface/sorters/external/tests/test_waveclus_snippets.py b/src/spikeinterface/sorters/external/tests/test_waveclus_snippets.py deleted file mode 100644 index 2a39c0290f..0000000000 --- a/src/spikeinterface/sorters/external/tests/test_waveclus_snippets.py +++ /dev/null @@ -1,24 +0,0 @@ -import os -import getpass - -if getpass.getuser() == "localadmin1": - waveclus_path = "F:/GitHub/wave_clus_original/" - os.environ["WAVECLUS_PATH"] = waveclus_path - -import unittest -import pytest - -from spikeinterface.sorters import WaveClusSnippetsSorter -from spikeinterface.sorters.tests.common_tests import SnippetsSorterCommonTestSuite - - -# This run several tests -@pytest.mark.skipif(not WaveClusSnippetsSorter.is_installed(), reason="waveclus not installed") -class WaveClusSnippetsCommonTestSuite(SnippetsSorterCommonTestSuite, unittest.TestCase): - SorterClass = WaveClusSnippetsSorter - - -if __name__ == "__main__": - test = WaveClusSnippetsCommonTestSuite() - test.setUp() - test.test_with_run() diff --git a/src/spikeinterface/sorters/external/waveclus_snippets.py b/src/spikeinterface/sorters/external/waveclus_snippets.py deleted file mode 100644 index 85a5104846..0000000000 --- a/src/spikeinterface/sorters/external/waveclus_snippets.py +++ /dev/null @@ -1,217 +0,0 @@ -from pathlib import Path -import os -import shutil -import sys -import json - -from spikeinterface.sorters.basesorter import BaseSorter -from spikeinterface.sorters.utils import ShellScript - -from spikeinterface.extractors.extractor_classes import WaveClusSortingExtractor -from spikeinterface.extractors.extractor_classes import WaveClusSnippetsExtractor - - -def check_if_installed(waveclus_path: str | None): - if waveclus_path is None: - return False - assert isinstance(waveclus_path, str) - - if waveclus_path.startswith('"'): - waveclus_path = waveclus_path[1:-1] - waveclus_path = str(Path(waveclus_path).absolute()) - - if (Path(waveclus_path) / "wave_clus.m").is_file(): - return True - else: - return False - - -class WaveClusSnippetsSorter(BaseSorter): - """WaveClus Sorter object.""" - - sorter_name: str = "waveclus_snippets" - compiled_name: str = "waveclus_snippets_compiled" - waveclus_path: str | None = os.getenv("WAVECLUS_PATH", None) - requires_locations = False - - _default_params = { - "feature_type": "wav", - "scales": 4, - "min_clus": 20, - "maxtemp": 0.251, - "template_sdnum": 3, - "mintemp": 0, - "stdmax": 50, - "max_spk": 40000, - "keep_good_only": True, - "chunk_memory": "500M", - } - - _params_description = { - "feature_type": "wav (for wavelets) or pca, type of feature extraction applied to the spikes", - "scales": "Levels of the wavelet decomposition used as features", - "min_clus": "Minimum increase of cluster sizes used by the peak selection on the temperature map", - "maxtemp": "Maximum temperature calculated by the SPC method", - "template_sdnum": "Maximum distance (in total variance of the cluster) from the mean waveform to force a " - "spike into a cluster", - "mintemp": "Minimum temperature calculated by the SPC algorithm", - "stdmax": "The events with a value over this number of noise standard deviations will be discarded", - "max_spk": "Maximum number of spikes used by the SPC algorithm", - "keep_good_only": "If True only 'good' units are returned", - "chunk_memory": "Chunk size in Mb to write h5 file (default 500Mb)", - } - - sorter_description = """Wave Clus combines a wavelet-based feature extraction and paramagnetic clustering with a - template-matching approach. It is mainly designed for monotrodes and low-channel count probes. - For more information see https://doi.org/10.1152/jn.00339.2018""" - - installation_mesg = """\nTo use WaveClus run:\n - >>> git clone https://github.com/csn-le/wave_clus - and provide the installation path by setting the WAVECLUS_PATH - environment variables or using WaveClusSorter.set_waveclus_path().\n\n - - More information on WaveClus at: - https://github.com/csn-le/wave_clus/wiki - """ - - @classmethod - def is_installed(cls): - if cls.check_compiled(): - return True - return check_if_installed(cls.waveclus_path) - - @classmethod - def get_sorter_version(cls): - if cls.check_compiled(): - return "compiled" - p = os.getenv("WAVECLUS_PATH", None) - if p is None: - return "unknown" - else: - with open(str(Path(p) / "version.txt"), encoding="utf8") as f: - version = f.readline() - return version - - @classmethod - def set_waveclus_path(cls, waveclus_path: str | Path): - waveclus_path = str(Path(waveclus_path).absolute()) - WaveClusSnippetsSorter.waveclus_path = waveclus_path - try: - print("Setting WAVECLUS_PATH environment variable for subprocess calls to:", waveclus_path) - os.environ["WAVECLUS_PATH"] = waveclus_path - except Exception as e: - print("Could not set WAVECLUS_PATH environment variable:", e) - - @classmethod - def _check_apply_filter_in_params(cls, params): - return False - - @classmethod - def _setup_recording(cls, snippets, sorter_output_folder, params, verbose): - # Generate mat files in the dataset directory - - WaveClusSnippetsExtractor.write_snippets(snippets, sorter_output_folder / "results_spikes.mat") - - if verbose: - num_snippets = snippets.get_total_snippets() - num_channels = snippets.get_num_channels() - print(f"Num. channels = {num_channels}, Num. snippets = {num_snippets}") - - @classmethod - def _run_from_folder(cls, sorter_output_folder, params, verbose): - sorter_output_folder = sorter_output_folder.absolute() - - cls._generate_par_file(params, sorter_output_folder) - if verbose: - print(f"Running waveclus in {sorter_output_folder}..") - - if cls.check_compiled(): - shell_cmd = f""" - #!/bin/bash - {cls.compiled_name} {sorter_output_folder} - """ - else: - source_dir = Path(__file__).parent - shutil.copy(str(source_dir / f"waveclus_snippets_master.m"), str(sorter_output_folder)) - - sorter_path = Path(cls.waveclus_path).absolute() - if "win" in sys.platform and sys.platform != "darwin": - disk_move = str(sorter_output_folder.absolute())[:2] - shell_cmd = f""" - {disk_move} - cd {sorter_output_folder} - matlab -nosplash -wait -log -r "waveclus_snippets_master('{sorter_output_folder}', '{sorter_path}')" - """ - else: - shell_cmd = f""" - #!/bin/bash - cd "{sorter_output_folder}" - matlab -nosplash -nodisplay -log -r "waveclus_snippets_master('{sorter_output_folder}', '{sorter_path}')" - """ - shell_cmd = ShellScript( - shell_cmd, - script_path=sorter_output_folder / f"run_{cls.sorter_name}", - log_path=sorter_output_folder / f"{cls.sorter_name}.log", - verbose=verbose, - ) - shell_cmd.start() - retcode = shell_cmd.wait() - - if retcode != 0: - raise Exception("waveclus returned a non-zero exit code") - - result_fname = sorter_output_folder / "times_results.mat" - if not result_fname.is_file(): - raise Exception(f"Result file does not exist: {result_fname}") - - @classmethod - def _get_result_from_folder(cls, sorter_output_folder): - sorter_output_folder = Path(sorter_output_folder) - result_fname = str(sorter_output_folder / "times_results.mat") - if (sorter_output_folder.parent / "spikeinterface_params.json").is_file(): - params_file = sorter_output_folder.parent / "spikeinterface_params.json" - else: - # back-compatibility - params_file = sorter_output_folder / "spikeinterface_params.json" - with params_file.open("r") as f: - sorter_params = json.load(f)["sorter_params"] - keep_good_only = sorter_params.get("keep_good_only", True) - sorting = WaveClusSortingExtractor(file_path=result_fname, keep_good_only=keep_good_only) - return sorting - - @staticmethod - def _generate_par_file(params, sorter_output_folder): - """ - This function generates parameters data for waveclus and saves as `par_input.mat` - - Loading example in Matlab (shouldn't be assigned to a variable): - >> load('/sorter_output_folder/par_input.mat'); - - Parameters - ---------- - params: dict - Custom parameters dictionary for waveclus - sorter_output_folder: pathlib.Path - Path object to save `par_input.mat` - """ - p = params.copy() - - par_renames = {"feature_type": "features"} - par_input = {} - for key, value in p.items(): - if type(value) == bool: - value = f"{value}".lower() - if key in par_renames: - key = par_renames[key] - par_input[key] = value - - # Converting integer values into float - # matlab interprets numerical fields as double by default - for k, v in par_input.items(): - if isinstance(v, int): - par_input[k] = float(v) - - par_input = {"par_input": par_input} - import scipy.io - - scipy.io.savemat(str(sorter_output_folder / "par_input.mat"), par_input) diff --git a/src/spikeinterface/sorters/external/waveclus_snippets_master.m b/src/spikeinterface/sorters/external/waveclus_snippets_master.m deleted file mode 100644 index ce57f620de..0000000000 --- a/src/spikeinterface/sorters/external/waveclus_snippets_master.m +++ /dev/null @@ -1,78 +0,0 @@ -function waveclus_snippets_master(outputFolder, waveclusPath) - if ~isdeployed - addpath(genpath(waveclusPath)); - end - - % par.mat file should contain a variable called par_input - load(fullfile(outputFolder, 'par_input.mat')); - - % Default parameters - S_par = set_parameters_ss(); - - % Update with custom parameters - S_par = update_parameters(S_par, par_input, 'relevant'); - - try - cd(outputFolder); - vcFile_spikes = fullfile(outputFolder, 'results_spikes.mat'); - vcFile_cluster = fullfile(outputFolder, 'times_results.mat'); - Do_clustering(vcFile_spikes, 'make_plots', false,'save_spikes',false,'par',S_par); - - if ~exist(vcFile_cluster,'file') - load(vcFile_spikes,'par'); - par = update_parameters(S_par, par, 'relevant'); - cluster_class = zeros(0,2); - save(vcFile_cluster,'cluster_class','par'); - end - catch - fprintf('----------------------------------------'); - fprintf(lasterr()); - quit(1); - end - quit(0); -end - - -function par = set_parameters_ss() - - % SPC PARAMETERS - par.mintemp = 0.00; % minimum temperature for SPC - par.maxtemp = 0.251; % maximum temperature for SPC - par.tempstep = 0.01; % temperature steps - par.SWCycles = 100; % SPC iterations for each temperature (default 100) - par.KNearNeighb = 11; % number of nearest neighbors for SPC - par.min_clus = 20; % minimum size of a cluster (default 60) - par.max_clus = 200; % maximum number of clusters allowed (default 200) - par.randomseed = 0; % if 0, random seed is taken as the clock value (default 0) - %par.randomseed = 147; % If not 0, random seed - %par.temp_plot = 'lin'; % temperature plot in linear scale - par.temp_plot = 'log'; % temperature plot in log scale - - par.c_ov = 0.7; % Overlapping coefficient to use for the inclusion criterion. - par.elbow_min = 0.4; %Thr_border parameter for regime border detection. - - % FEATURES PARAMETERS - par.min_inputs = 10; % number of inputs to the clustering - par.max_inputs = 0.75; % number of inputs to the clustering. if < 1 it will the that proportion of the maximum. - par.scales = 4; % number of scales for the wavelet decomposition - par.features = 'wav'; % type of feature ('wav' or 'pca') - %par.features = 'pca' - - - % FORCE MEMBERSHIP PARAMETERS - par.template_sdnum = 3; % max radius of cluster in std devs. - par.template_k = 10; % # of nearest neighbors - par.template_k_min = 10; % min # of nn for vote - %par.template_type = 'mahal'; % nn, center, ml, mahal - par.template_type = 'center'; % nn, center, ml, mahal - par.force_feature = 'spk'; % feature use for forcing (whole spike shape) - %par.force_feature = 'wav'; % feature use for forcing (wavelet coefficients). - par.force_auto = true; %automatically force membership (only for batch scripts). - - % TEMPLATE MATCHING - par.match = 'y'; % for template matching - %par.match = 'n'; % for no template matching - par.max_spk = 40000; % max. # of spikes before starting templ. match. - par.permut = 'y'; % for selection of random 'par.max_spk' spikes before starting templ. match. - % par.permut = 'n'; % for selection of the first 'par.max_spk' spikes before starting templ. match. -end diff --git a/src/spikeinterface/sorters/runsorter.py b/src/spikeinterface/sorters/runsorter.py index dfa96234fe..9b0e6855f5 100644 --- a/src/spikeinterface/sorters/runsorter.py +++ b/src/spikeinterface/sorters/runsorter.py @@ -55,7 +55,6 @@ kilosort2_5="kilosort2_5-compiled", kilosort3="kilosort3-compiled", waveclus="waveclus-compiled", - waveclus_snippets="waveclus-compiled", # archived # klusta="klusta", # yass="yass", diff --git a/src/spikeinterface/sorters/sorterlist.py b/src/spikeinterface/sorters/sorterlist.py index 81279e8ebf..55e853ee5d 100644 --- a/src/spikeinterface/sorters/sorterlist.py +++ b/src/spikeinterface/sorters/sorterlist.py @@ -16,7 +16,6 @@ from .external.spyking_circus import SpykingcircusSorter from .external.tridesclous import TridesclousSorter from .external.waveclus import WaveClusSorter -from .external.waveclus_snippets import WaveClusSnippetsSorter from .external.yass import YassSorter # based on spikeinterface.sortingcomponents @@ -44,7 +43,6 @@ SpykingcircusSorter, TridesclousSorter, WaveClusSorter, - WaveClusSnippetsSorter, YassSorter, # internal Spykingcircus2Sorter, diff --git a/src/spikeinterface/sorters/tests/common_tests.py b/src/spikeinterface/sorters/tests/common_tests.py index 316ea21d40..79e2d5621f 100644 --- a/src/spikeinterface/sorters/tests/common_tests.py +++ b/src/spikeinterface/sorters/tests/common_tests.py @@ -3,7 +3,6 @@ from spikeinterface import generate_ground_truth_recording from spikeinterface.sorters import run_sorter -from spikeinterface.core.snippets_tools import snippets_from_sorting class SorterCommonTestSuite: @@ -64,57 +63,3 @@ def test_with_run(self): def test_get_version(self): version = self.SorterClass.get_sorter_version() print("test_get_versions", self.SorterClass.sorter_name, version) - - -class SnippetsSorterCommonTestSuite: - """ - This class run some basic for a sorter class. - This is the minimal test suite for each sorter class: - * run once - """ - - @pytest.fixture(autouse=True) - def create_cache_folder(self, tmp_path_factory): - self.cache_folder = tmp_path_factory.mktemp("cache_folder") - - SorterClass = None - - def setUp(self): - recording, sorting_gt = generate_ground_truth_recording(num_channels=4, durations=[60], seed=0) - snippets_folder = self.cache_folder / "snippets" - if snippets_folder.is_dir(): - shutil.rmtree(snippets_folder) - - nse = snippets_from_sorting(recording=recording, sorting=sorting_gt, nbefore=20, nafter=44) - - self.snippets = nse.save(folder=snippets_folder, verbose=False, format="npy") - print(self.snippets) - - def test_with_run(self): - # some sorter (TDC, KS, KS2, ...) work by default with the raw binary - # format as input to avoid copy when the recording is already this format - - snippets = self.snippets - - sorter_name = self.SorterClass.sorter_name - - output_folder = self.cache_folder / sorter_name - - sorter_params = self.SorterClass.default_params() - - sorting = run_sorter( - sorter_name, - snippets, - folder=output_folder, - remove_existing_folder=True, - delete_output_folder=False, - verbose=False, - raise_error=True, - **sorter_params, - ) - - del sorting - - def test_get_version(self): - version = self.SorterClass.get_sorter_version() - print("test_get_versions", self.SorterClass.sorter_name, version) From 92b72ca7f9b50120fcde2c5ec2f137da65a50ddd Mon Sep 17 00:00:00 2001 From: Alessio Buccino Date: Wed, 7 Oct 2026 15:00:29 +0200 Subject: [PATCH 4/5] last one standing --- src/spikeinterface/sorters/basesorter.py | 5 +---- 1 file changed, 1 insertion(+), 4 deletions(-) diff --git a/src/spikeinterface/sorters/basesorter.py b/src/spikeinterface/sorters/basesorter.py index d0ec68533d..f25fe33e92 100644 --- a/src/spikeinterface/sorters/basesorter.py +++ b/src/spikeinterface/sorters/basesorter.py @@ -12,7 +12,7 @@ import warnings -from spikeinterface.core import load, BaseRecordingSnippets, BaseRecording +from spikeinterface.core import load, BaseRecording from spikeinterface.core.core_tools import check_json from spikeinterface.core.recording_tools import get_rec_attributes from spikeinterface.core.globals import get_global_job_kwargs @@ -102,9 +102,6 @@ def initialize_folder(cls, recording, output_folder, verbose, remove_existing_fo f"The sorter {cls.sorter_name} is not installed. Please install it with:\n{cls.installation_mesg}" ) - if not isinstance(recording, BaseRecordingSnippets): - raise ValueError("recording must be a Recording or a Snippets!!") - if cls.requires_locations: if not recording.has_probe(): raise RuntimeError( From 236a75521fdaf07c902e1d966e29787cbef67444 Mon Sep 17 00:00:00 2001 From: Alessio Buccino Date: Wed, 7 Oct 2026 15:09:17 +0200 Subject: [PATCH 5/5] fix: warnings in base recording --- src/spikeinterface/core/baserecording.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/spikeinterface/core/baserecording.py b/src/spikeinterface/core/baserecording.py index 9ff382b3b3..e127dc9201 100644 --- a/src/spikeinterface/core/baserecording.py +++ b/src/spikeinterface/core/baserecording.py @@ -339,7 +339,7 @@ def is_filtered(self): def set_channel_groups(self, groups, channel_ids=None): if "probes" in self._annotations: - warn("set_channel_groups() destroys the probe description. Using set_probe() is preferable") + warnings.warn("set_channel_groups() destroys the probe description. Using set_probe() is preferable") self._annotations.pop("probes") self.set_property("group", groups, ids=channel_ids) @@ -621,7 +621,7 @@ def _get_probegroup_based_on_device_channel_indices(self, probegroup: ProbeGroup order = np.argsort(device_channel_indices) probegroup = probegroup.get_slice(order) else: - warn( + warnings.warn( "No connected channels in the probegroup! " "The probegroup will be attached but no channel will be selected." )