diff --git a/neo/rawio/baserawio.py b/neo/rawio/baserawio.py index 4207f3c03..462732a33 100644 --- a/neo/rawio/baserawio.py +++ b/neo/rawio/baserawio.py @@ -1561,6 +1561,11 @@ class BaseRawWithBufferApiIO(BaseRawIO): * np.memmap * hdf5 + * mtscomp + + An mtscomp buffer description contains ``type="mtscomp"``, the compressed + ``file_path`` and JSON ``metadata_path``, plus ``dtype``, ``shape``, and + ``time_axis=0``. In theses cases _get_signal_size and _get_analogsignal_chunk are totaly generic and do not need to be implemented in the class. @@ -1667,6 +1672,30 @@ def _get_analogsignal_chunk( else: raise ValueError(f"time_axis must be 0 or 1, got {time_axis}") + elif buffer_desc["type"] == "mtscomp": + if time_axis != 0: + raise ValueError(f"mtscomp buffers require time_axis=0, got {time_axis}") + + try: + import mtscomp + except ImportError as exc: + raise ImportError( + "Reading mtscomp-compressed signals requires mtscomp. " + 'Install it with `pip install "neo[mtscomp]"` or `pip install mtscomp tqdm`.' + ) from exc + + if not hasattr(self, "_mtscomp_analogsignal_buffers"): + self._mtscomp_analogsignal_buffers = {} + block_readers = self._mtscomp_analogsignal_buffers.setdefault(block_index, {}) + segment_readers = block_readers.setdefault(seg_index, {}) + + if buffer_id not in segment_readers: + reader = mtscomp.Reader() + reader.open(buffer_desc["file_path"], buffer_desc["metadata_path"]) + segment_readers[buffer_id] = reader + + raw_sigs = segment_readers[buffer_id][i_start:i_stop] + elif buffer_desc["type"] == "hdf5": # open files on demand and keep reference to opened file @@ -1726,6 +1755,25 @@ def __del__(self): h5_file.close() del self._hdf5_analogsignal_buffers + self._close_mtscomp_analogsignal_buffers() + + def _close_mtscomp_analogsignal_buffers(self): + """Close cached mtscomp readers, tolerating partially constructed instances.""" + readers_by_block = getattr(self, "_mtscomp_analogsignal_buffers", None) + if readers_by_block is None: + return + + for readers_by_segment in readers_by_block.values(): + for readers_by_buffer in readers_by_segment.values(): + for reader in readers_by_buffer.values(): + try: + reader.close() + except Exception: + # Destructors must remain safe during interpreter shutdown or + # after a Reader failed part-way through opening its files. + pass + del self._mtscomp_analogsignal_buffers + def pprint_vector(vector, lim: int = 8): vector = np.asarray(vector) diff --git a/neo/rawio/openephysbinaryrawio.py b/neo/rawio/openephysbinaryrawio.py index c139a9ae6..aaf417ae6 100644 --- a/neo/rawio/openephysbinaryrawio.py +++ b/neo/rawio/openephysbinaryrawio.py @@ -28,6 +28,105 @@ from .utils import get_memmap_shape +def _resolve_continuous_storage(stream_folder): + """Resolve the signal data backend for one Open Ephys continuous stream.""" + stream_folder = Path(stream_folder) + dat_file = stream_folder / "continuous.dat" + cbin_file = stream_folder / "continuous.cbin" + ch_file = stream_folder / "continuous.ch" + + # Preserve the historical raw-file behavior when both representations exist. + if dat_file.is_file(): + return {"type": "raw", "file_path": dat_file} + + if cbin_file.is_file(): + if not ch_file.is_file(): + raise FileNotFoundError( + f"Found mtscomp-compressed signal data at {cbin_file}, but its required " + f"metadata file {ch_file} is missing." + ) + return {"type": "mtscomp", "file_path": cbin_file, "metadata_path": ch_file} + + if ch_file.is_file(): + raise FileNotFoundError( + f"Found mtscomp metadata at {ch_file}, but no signal data file was found. " + f"Expected {dat_file} or {cbin_file}." + ) + + raise FileNotFoundError(f"No signal data file was found in {stream_folder}. Expected {dat_file} or {cbin_file}.") + + +def _read_mtscomp_metadata(metadata_path): + """Read and validate the metadata needed to describe an mtscomp buffer.""" + metadata_path = Path(metadata_path) + with open(metadata_path, encoding="utf8") as file: + metadata = json.load(file) + + if not isinstance(metadata, dict): + raise ValueError(f"Invalid mtscomp metadata in {metadata_path}: expected a JSON object.") + + required_keys = ("n_channels", "sample_rate", "dtype", "chunk_bounds") + missing_keys = [key for key in required_keys if key not in metadata] + if missing_keys: + raise ValueError( + f"Invalid mtscomp metadata in {metadata_path}: missing required " + f"{'key' if len(missing_keys) == 1 else 'keys'} {', '.join(missing_keys)}." + ) + + n_channels = metadata["n_channels"] + if isinstance(n_channels, bool) or not isinstance(n_channels, int) or n_channels <= 0: + raise ValueError( + f"Invalid mtscomp metadata in {metadata_path}: n_channels must be a positive integer, " + f"observed {n_channels!r}." + ) + + sample_rate = metadata["sample_rate"] + if isinstance(sample_rate, bool): + valid_sample_rate = False + else: + try: + sample_rate = float(sample_rate) + except (TypeError, ValueError): + valid_sample_rate = False + else: + valid_sample_rate = np.isfinite(sample_rate) and sample_rate > 0 + if not valid_sample_rate: + raise ValueError( + f"Invalid mtscomp metadata in {metadata_path}: sample_rate must be finite and positive, " + f"observed {metadata['sample_rate']!r}." + ) + + try: + np.dtype(metadata["dtype"]) + except (TypeError, ValueError) as exc: + raise ValueError( + f"Invalid mtscomp metadata in {metadata_path}: dtype {metadata['dtype']!r} " + f"cannot be interpreted as a NumPy dtype." + ) from exc + + chunk_bounds = metadata["chunk_bounds"] + if not isinstance(chunk_bounds, list) or len(chunk_bounds) == 0: + raise ValueError(f"Invalid mtscomp metadata in {metadata_path}: chunk_bounds must be a non-empty list.") + if any(isinstance(bound, bool) or not isinstance(bound, int) for bound in chunk_bounds): + raise ValueError( + f"Invalid mtscomp metadata in {metadata_path}: chunk_bounds must contain integer sample indices." + ) + if chunk_bounds[0] != 0: + raise ValueError( + f"Invalid mtscomp metadata in {metadata_path}: chunk_bounds must start at 0, " + f"observed {chunk_bounds[0]}." + ) + if chunk_bounds[-1] < 0: + raise ValueError( + f"Invalid mtscomp metadata in {metadata_path}: the final chunk bound must be non-negative, " + f"observed {chunk_bounds[-1]}." + ) + if any(next_bound <= bound for bound, next_bound in zip(chunk_bounds, chunk_bounds[1:])): + raise ValueError(f"Invalid mtscomp metadata in {metadata_path}: chunk_bounds must be monotonically increasing.") + + return metadata + + class OpenEphysBinaryRawIO(BaseRawWithBufferApiIO): """ Handle several Blocks and several Segments. @@ -57,7 +156,9 @@ class OpenEphysBinaryRawIO(BaseRawWithBufferApiIO): │ │ │ ├── structure.oebin # Required: JSON metadata file │ │ │ ├── continuous/ # Signal data streams │ │ │ │ └── AP_band/ # Stream folder (becomes "Record Node 102#AP_band") - │ │ │ │ ├── continuous.dat # Raw binary signal data + │ │ │ │ ├── continuous.dat # Raw binary signal data, or: + │ │ │ │ ├── continuous.cbin # Lossless mtscomp signal data + │ │ │ │ ├── continuous.ch # mtscomp metadata │ │ │ │ ├── timestamps.npy # Sample timestamps (pre-v0.6) │ │ │ │ └── sample_numbers.npy # Sample numbers (v0.6+) │ │ │ └── events/ # Event data streams (optional) @@ -92,6 +193,14 @@ class OpenEphysBinaryRawIO(BaseRawWithBufferApiIO): - Point to specific recording folder: Single recording session only The reader will automatically detect the level and parse accordingly. + Signal data may be stored as ``continuous.dat`` or as an mtscomp + ``continuous.cbin``/``continuous.ch`` pair. If both representations exist, + ``continuous.dat`` takes precedence. Compressed data are decompressed losslessly + on demand without creating a temporary ``.dat`` file. Parsing metadata, events, + and timestamps does not require mtscomp; accessing compressed traces requires the + optional ``mtscomp`` package, available through ``pip install "neo[mtscomp]"`` + or ``pip install mtscomp tqdm``. + # Correspondencies Neo OpenEphys block[n-1] experiment[n] New device start/stop @@ -103,7 +212,7 @@ class OpenEphysBinaryRawIO(BaseRawWithBufferApiIO): """ - extensions = ["xml", "oebin", "txt", "dat", "npy"] + extensions = ["xml", "oebin", "txt", "dat", "cbin", "ch", "npy"] rawmode = "one-dir" def __init__(self, dirname="", experiment_names=None): @@ -119,6 +228,11 @@ def __init__(self, dirname="", experiment_names=None): def _source_name(self): return self.dirname + @property + def uses_mtscomp(self): + """Whether any parsed continuous stream uses mtscomp compression.""" + return getattr(self, "_uses_mtscomp", False) + def _parse_header(self): # Use the static private methods directly folder_structure_dict, possible_experiments = OpenEphysBinaryRawIO._parse_folder_structure( @@ -273,15 +387,66 @@ def _parse_header(self): num_channels = len(info["channels"]) stream_id = str(stream_index) buffer_id = str(stream_index) - shape = get_memmap_shape(info["raw_filename"], info["dtype"], num_channels=num_channels, offset=0) - self._buffer_descriptions[block_index][seg_index][buffer_id] = { - "type": "raw", - "file_path": str(info["raw_filename"]), - "dtype": info["dtype"], - "order": "C", - "file_offset": 0, - "shape": shape, - } + if info["data_format"] == "raw": + shape = get_memmap_shape( + info["raw_filename"], info["dtype"], num_channels=num_channels, offset=0 + ) + buffer_description = { + "type": "raw", + "file_path": str(info["raw_filename"]), + "dtype": info["dtype"], + "order": "C", + "file_offset": 0, + "shape": shape, + } + elif info["data_format"] == "mtscomp": + metadata = _read_mtscomp_metadata(info["compression_metadata_filename"]) + cbin_path = info["data_filename"] + stream_name = info["stream_name"] + + expected_n_channels = len(info["channels"]) + observed_n_channels = metadata["n_channels"] + if observed_n_channels != expected_n_channels: + raise ValueError( + f"Open Ephys stream {stream_name!r} at {cbin_path} declares " + f"{expected_n_channels} channels in structure.oebin, but " + f"continuous.ch declares {observed_n_channels}." + ) + + expected_dtype = np.dtype(info["dtype"]) + observed_dtype = np.dtype(metadata["dtype"]) + if observed_dtype != expected_dtype: + raise ValueError( + f"Open Ephys stream {stream_name!r} at {cbin_path} declares dtype " + f"{expected_dtype} in structure.oebin, but continuous.ch declares " + f"{observed_dtype}." + ) + + expected_sample_rate = float(info["sample_rate"]) + observed_sample_rate = float(metadata["sample_rate"]) + if not np.isclose(observed_sample_rate, expected_sample_rate): + raise ValueError( + f"Open Ephys stream {stream_name!r} at {cbin_path} declares sample rate " + f"{expected_sample_rate} Hz in structure.oebin, but continuous.ch declares " + f"{observed_sample_rate} Hz." + ) + + shape = (metadata["chunk_bounds"][-1], metadata["n_channels"]) + buffer_description = { + "type": "mtscomp", + "file_path": cbin_path, + "metadata_path": info["compression_metadata_filename"], + "dtype": str(observed_dtype), + "shape": shape, + "time_axis": 0, + } + else: + raise ValueError( + f"Unsupported continuous data format {info['data_format']!r} " + f"for Open Ephys stream {info['stream_name']!r}." + ) + + self._buffer_descriptions[block_index][seg_index][buffer_id] = buffer_description has_sync_trace = self._sig_streams[block_index][seg_index][stream_index]["has_sync_trace"] @@ -351,6 +516,13 @@ def _parse_header(self): else: self._stream_buffer_slice[stream_id_non_neural] = slice(num_neural_channels, None) + self._uses_mtscomp = any( + description["type"] == "mtscomp" + for descriptions_by_segment in self._buffer_descriptions.values() + for descriptions in descriptions_by_segment.values() + for description in descriptions.values() + ) + # events zone # channel map: one channel one stream event_channels = [] @@ -725,7 +897,10 @@ def _parse_folder_structure(dirname, experiment_names=None): { "channels": [{"channel_name": "CH1", "bit_volts": 0.195, ...}, ...], "sample_rate": 30000.0, - "raw_filename": "/path/to/continuous.dat", + "data_format": "raw" or "mtscomp", + "data_filename": "/path/to/continuous.dat" or "/path/to/continuous.cbin", + "raw_filename": "/path/to/continuous.dat", # raw streams only; compatibility alias + "compression_metadata_filename": "/path/to/continuous.ch", # mtscomp only "dtype": "int16", "timestamp0": 123456, "t_start": 0.0 @@ -849,7 +1024,8 @@ def _parse_folder_structure(dirname, experiment_names=None): ) continue - raw_filename = recording_folder / "continuous" / info["folder_name"] / "continuous.dat" + stream_folder = recording_folder / "continuous" / info["folder_name"] + storage = _resolve_continuous_storage(stream_folder) # Updates for OpenEphys v0.6: # In new vesion (>=0.6) timestamps.npy is now called sample_numbers.npy @@ -869,7 +1045,15 @@ def _parse_folder_structure(dirname, experiment_names=None): # TODO for later : gap checking signal_stream = info.copy() - signal_stream["raw_filename"] = str(raw_filename) + signal_stream["data_format"] = storage["type"] + signal_stream["data_filename"] = str(storage["file_path"]) + if storage["type"] == "raw": + # Keep the field historically returned by explore_folder(). + # Compressed streams deliberately omit this alias because a + # .cbin file cannot be treated as raw binary data. + signal_stream["raw_filename"] = signal_stream["data_filename"] + else: + signal_stream["compression_metadata_filename"] = str(storage["metadata_path"]) signal_stream["dtype"] = "int16" signal_stream["timestamp0"] = timestamp0 signal_stream["t_start"] = t_start diff --git a/neo/test/rawiotest/test_openephysbinaryrawio.py b/neo/test/rawiotest/test_openephysbinaryrawio.py index 0119a4b53..651c7c83f 100644 --- a/neo/test/rawiotest/test_openephysbinaryrawio.py +++ b/neo/test/rawiotest/test_openephysbinaryrawio.py @@ -1,9 +1,93 @@ +import builtins +import json +import shutil import unittest +from pathlib import Path +from unittest.mock import patch -from neo.rawio.openephysbinaryrawio import OpenEphysBinaryRawIO +from neo.rawio.openephysbinaryrawio import ( + OpenEphysBinaryRawIO, + _read_mtscomp_metadata, + _resolve_continuous_storage, + explore_folder, +) from neo.test.rawiotest.common_rawio_test import BaseTestRawIO import numpy as np +import pytest + + +def _write_synthetic_open_ephys_recording( + dataset_folder, + recordings, + channel_names_by_stream, + sample_rate=100.0, + experiment_index=1, +): + """Create a small raw Open Ephys dataset for compression tests.""" + dataset_folder = Path(dataset_folder) + for recording_index, data_by_stream in enumerate(recordings, start=1): + recording_folder = dataset_folder / f"experiment{experiment_index}" / f"recording{recording_index}" + continuous_metadata = [] + for stream_name, data in data_by_stream.items(): + data = np.asarray(data, dtype=np.int16) + stream_folder = recording_folder / "continuous" / stream_name + stream_folder.mkdir(parents=True) + data.tofile(stream_folder / "continuous.dat") + np.save(stream_folder / "sample_numbers.npy", np.arange(data.shape[0], dtype=np.int64)) + + channels = [ + { + "channel_name": channel_name, + "bit_volts": 0.195, + "units": "uV", + } + for channel_name in channel_names_by_stream[stream_name] + ] + continuous_metadata.append( + { + "folder_name": stream_name, + "sample_rate": sample_rate, + "dtype": "int16", + "channels": channels, + } + ) + + with open(recording_folder / "structure.oebin", "w", encoding="utf8") as file: + json.dump({"continuous": continuous_metadata, "events": []}, file) + + return dataset_folder + + +def _compress_open_ephys_recording(source, destination, remove_dat=True, chunk_duration=1.0): + """Copy an Open Ephys fixture and compress each continuous.dat in the copy.""" + mtscomp = pytest.importorskip("mtscomp") + source = Path(source) + destination = Path(destination) + shutil.copytree(source, destination) + + for structure_path in destination.rglob("structure.oebin"): + with open(structure_path, encoding="utf8") as file: + structure = json.load(file) + for stream_info in structure["continuous"]: + stream_folder = structure_path.parent / "continuous" / stream_info["folder_name"] + dat_path = stream_folder / "continuous.dat" + mtscomp.compress( + dat_path, + stream_folder / "continuous.cbin", + stream_folder / "continuous.ch", + sample_rate=stream_info["sample_rate"], + n_channels=len(stream_info["channels"]), + dtype=np.int16, + chunk_duration=chunk_duration, + n_threads=1, + check_after_compress=False, + quiet=True, + ) + if remove_dat: + dat_path.unlink() + + return destination class TestOpenEphysBinaryRawIO(BaseTestRawIO, unittest.TestCase): @@ -84,5 +168,296 @@ def test_separating_stream_for_non_neural_data(self): assert rawio.header["signal_streams"]["name"].tolist() == ["Rhythm_FPGA-100.0", "Rhythm_FPGA-100.0_ADC"] +def test_resolve_continuous_storage(tmp_path): + dat_path = tmp_path / "continuous.dat" + cbin_path = tmp_path / "continuous.cbin" + ch_path = tmp_path / "continuous.ch" + + dat_path.touch() + assert _resolve_continuous_storage(tmp_path) == {"type": "raw", "file_path": dat_path} + + cbin_path.touch() + ch_path.touch() + # Raw data takes precedence when both representations are present. + assert _resolve_continuous_storage(tmp_path) == {"type": "raw", "file_path": dat_path} + + dat_path.unlink() + assert _resolve_continuous_storage(tmp_path) == { + "type": "mtscomp", + "file_path": cbin_path, + "metadata_path": ch_path, + } + + +@pytest.mark.parametrize( + ("present_filename", "match"), + [ + ("continuous.cbin", "required metadata file"), + ("continuous.ch", "no signal data file"), + (None, "No signal data file"), + ], +) +def test_resolve_continuous_storage_missing_files(tmp_path, present_filename, match): + if present_filename is not None: + (tmp_path / present_filename).touch() + with pytest.raises(FileNotFoundError, match=match): + _resolve_continuous_storage(tmp_path) + + +def test_raw_and_mtscomp_open_ephys_signals_are_equal(tmp_path): + sample_count = 245 + probe_data = (np.arange(sample_count * 3, dtype=np.int16) - 400).reshape(sample_count, 3) + mixed_data = (np.arange(sample_count * 2, dtype=np.int16) + 800).reshape(sample_count, 2) + raw_folder = _write_synthetic_open_ephys_recording( + tmp_path / "raw", + [{"ProbeA-AP": probe_data, "Rhythm": mixed_data}], + { + "ProbeA-AP": ["AP0", "AP1", "Probe_SYNC"], + "Rhythm": ["CH0", "ADC0"], + }, + ) + compressed_folder = _compress_open_ephys_recording(raw_folder, tmp_path / "compressed") + + # Even invalid compressed artifacts must not affect a valid raw recording. + for dat_path in raw_folder.rglob("continuous.dat"): + (dat_path.parent / "continuous.cbin").write_bytes(b"not mtscomp data") + (dat_path.parent / "continuous.ch").write_text("not JSON", encoding="utf8") + + _, raw_streams, _, _, _ = explore_folder(raw_folder) + raw_stream_info = raw_streams[0][0]["continuous"]["ProbeA-AP"] + assert raw_stream_info["raw_filename"] == raw_stream_info["data_filename"] + assert raw_stream_info["raw_filename"].endswith("continuous.dat") + + _, compressed_streams, _, _, _ = explore_folder(compressed_folder) + compressed_stream_info = compressed_streams[0][0]["continuous"]["ProbeA-AP"] + assert "raw_filename" not in compressed_stream_info + assert compressed_stream_info["data_filename"].endswith("continuous.cbin") + + raw_reader = OpenEphysBinaryRawIO(raw_folder) + compressed_reader = OpenEphysBinaryRawIO(compressed_folder) + raw_reader.parse_header() + compressed_reader.parse_header() + + assert raw_reader.uses_mtscomp is False + assert compressed_reader.uses_mtscomp is True + with pytest.raises(AttributeError): + compressed_reader.uses_mtscomp = False + assert not list(compressed_folder.rglob("continuous.dat")) + np.testing.assert_array_equal(raw_reader.header["signal_buffers"], compressed_reader.header["signal_buffers"]) + np.testing.assert_array_equal(raw_reader.header["signal_streams"], compressed_reader.header["signal_streams"]) + np.testing.assert_array_equal(raw_reader.header["signal_channels"], compressed_reader.header["signal_channels"]) + raw_buffer_description = raw_reader.get_analogsignal_buffer_description(0, 0, "0") + assert set(raw_buffer_description) == { + "type", + "file_path", + "dtype", + "order", + "file_offset", + "shape", + } + assert raw_buffer_description["type"] == "raw" + assert raw_buffer_description["file_path"].endswith("continuous.dat") + assert raw_reader.segment_t_start(0, 0) == compressed_reader.segment_t_start(0, 0) + assert raw_reader.segment_t_stop(0, 0) == compressed_reader.segment_t_stop(0, 0) + + for stream_index in range(raw_reader.signal_streams_count()): + signal_size = raw_reader.get_signal_size(0, 0, stream_index) + assert signal_size == compressed_reader.get_signal_size(0, 0, stream_index) + for i_start, i_stop in ((0, 1), (signal_size - 1, signal_size), (95, 105), (35, 190)): + raw_chunk = raw_reader.get_analogsignal_chunk(0, 0, i_start, i_stop, stream_index) + compressed_chunk = compressed_reader.get_analogsignal_chunk(0, 0, i_start, i_stop, stream_index) + np.testing.assert_array_equal(raw_chunk, compressed_chunk) + + stream_names = raw_reader.header["signal_streams"]["name"].tolist() + neural_stream_index = stream_names.index("ProbeA-AP") + mixed_neural_stream_index = stream_names.index("Rhythm") + adc_stream_index = stream_names.index("Rhythm_ADC") + sync_stream_index = stream_names.index("ProbeA-APSYNC") + + neural = compressed_reader.get_analogsignal_chunk(0, 0, 0, sample_count, neural_stream_index) + mixed_neural = compressed_reader.get_analogsignal_chunk(0, 0, 0, sample_count, mixed_neural_stream_index) + adc = compressed_reader.get_analogsignal_chunk(0, 0, 0, sample_count, adc_stream_index) + sync = compressed_reader.get_analogsignal_chunk(0, 0, 0, sample_count, sync_stream_index) + np.testing.assert_array_equal(neural, probe_data[:, :2]) + np.testing.assert_array_equal(mixed_neural, mixed_data[:, :1]) + np.testing.assert_array_equal(adc, mixed_data[:, 1:2]) + np.testing.assert_array_equal(sync, probe_data[:, 2:3]) + + contiguous = compressed_reader.get_analogsignal_chunk( + 0, 0, 20, 40, neural_stream_index, channel_indexes=slice(0, 2) + ) + non_contiguous = compressed_reader.get_analogsignal_chunk(0, 0, 20, 40, neural_stream_index, channel_indexes=[1, 0]) + np.testing.assert_array_equal(contiguous, probe_data[20:40, :2]) + np.testing.assert_array_equal(non_contiguous, probe_data[20:40, [1, 0]]) + assert not list(compressed_folder.rglob("continuous.dat")) + + +def test_mtscomp_readers_are_cached_per_block_segment_and_buffer(tmp_path): + channel_names_by_stream = {"AP": ["AP0", "AP1"], "LFP": ["LFP0"]} + raw_folder = tmp_path / "raw_multi" + expected = {} + for block_index in range(2): + recordings = [] + expected[block_index] = {} + for seg_index, sample_count in enumerate((125 + 10 * block_index, 217 + 10 * block_index)): + offset = 1000 * (2 * block_index + seg_index) + ap = (np.arange(sample_count * 2, dtype=np.int16) + offset).reshape(sample_count, 2) + lfp = (np.arange(sample_count, dtype=np.int16) - offset).reshape(sample_count, 1) + recordings.append({"AP": ap, "LFP": lfp}) + expected[block_index][seg_index] = {"AP": ap, "LFP": lfp} + + _write_synthetic_open_ephys_recording( + raw_folder, + recordings, + channel_names_by_stream, + experiment_index=block_index + 1, + ) + + compressed_folder = _compress_open_ephys_recording(raw_folder, tmp_path / "compressed_multi") + reader = OpenEphysBinaryRawIO(compressed_folder) + reader.parse_header() + + stream_names = reader.header["signal_streams"]["name"].tolist() + for block_index in range(2): + for seg_index in range(2): + for stream_name in ("AP", "LFP"): + stream_index = stream_names.index(stream_name) + first = reader.get_analogsignal_chunk(block_index, seg_index, 5, 25, stream_index) + cached_reader = reader._mtscomp_analogsignal_buffers[block_index][seg_index][str(stream_index)] + second = reader.get_analogsignal_chunk(block_index, seg_index, 10, 30, stream_index) + assert reader._mtscomp_analogsignal_buffers[block_index][seg_index][str(stream_index)] is cached_reader + np.testing.assert_array_equal(first, expected[block_index][seg_index][stream_name][5:25]) + np.testing.assert_array_equal(second, expected[block_index][seg_index][stream_name][10:30]) + + cached_readers = [ + compressed_reader + for block_readers in reader._mtscomp_analogsignal_buffers.values() + for segment_readers in block_readers.values() + for compressed_reader in segment_readers.values() + ] + assert len(cached_readers) == 8 + assert len({id(compressed_reader) for compressed_reader in cached_readers}) == 8 + + +@pytest.mark.parametrize( + ("mutation", "match"), + [ + (lambda metadata: metadata.update(n_channels=3), "declares 2 channels.*declares 3"), + (lambda metadata: metadata.update(dtype="float32"), "declares dtype int16.*declares float32"), + (lambda metadata: metadata.update(sample_rate=123.0), "declares sample rate 100.0.*declares 123.0"), + (lambda metadata: metadata.pop("chunk_bounds"), "missing required key chunk_bounds"), + (lambda metadata: metadata.update(chunk_bounds=[0, 100, 50]), "monotonically increasing"), + (lambda metadata: metadata.update(chunk_bounds=[0, "invalid"]), "integer sample indices"), + ], +) +def test_invalid_mtscomp_metadata_fails_during_header_parsing(tmp_path, mutation, match): + data = np.arange(400, dtype=np.int16).reshape(200, 2) + raw_folder = _write_synthetic_open_ephys_recording( + tmp_path / "raw", + [{"AP": data}], + {"AP": ["AP0", "AP1"]}, + ) + compressed_folder = _compress_open_ephys_recording(raw_folder, tmp_path / "compressed") + metadata_path = next(compressed_folder.rglob("continuous.ch")) + with open(metadata_path, encoding="utf8") as file: + metadata = json.load(file) + mutation(metadata) + with open(metadata_path, "w", encoding="utf8") as file: + json.dump(metadata, file) + + with pytest.raises(ValueError, match=match): + OpenEphysBinaryRawIO(compressed_folder).parse_header() + + +def test_mtscomp_is_only_required_for_trace_access_and_reader_closes(tmp_path): + data = np.arange(400, dtype=np.int16).reshape(200, 2) + raw_folder = _write_synthetic_open_ephys_recording( + tmp_path / "raw", + [{"AP": data}], + {"AP": ["AP0", "AP1"]}, + ) + structure_path = next(raw_folder.rglob("structure.oebin")) + event_folder = structure_path.parent / "events" / "Messages" + event_folder.mkdir(parents=True) + np.save(event_folder / "timestamps.npy", np.array([10, 20], dtype=np.int64)) + np.save(event_folder / "text.npy", np.array([b"start", b"stop"])) + with open(structure_path, encoding="utf8") as file: + structure = json.load(file) + structure["events"] = [ + { + "folder_name": "Messages", + "channel_name": "Messages", + "sample_rate": 100.0, + } + ] + with open(structure_path, "w", encoding="utf8") as file: + json.dump(structure, file) + + compressed_folder = _compress_open_ephys_recording(raw_folder, tmp_path / "compressed") + raw_reader = OpenEphysBinaryRawIO(raw_folder) + reader = OpenEphysBinaryRawIO(compressed_folder) + + original_import = builtins.__import__ + + def import_without_mtscomp(name, *args, **kwargs): + if name == "mtscomp": + raise ImportError("simulated missing optional dependency") + return original_import(name, *args, **kwargs) + + with patch("builtins.__import__", side_effect=import_without_mtscomp): + raw_reader.parse_header() + reader.parse_header() + raw_chunk = raw_reader.get_analogsignal_chunk(0, 0, 0, 10, 0) + np.testing.assert_array_equal(raw_chunk, data[:10]) + + assert raw_reader.event_count(0, 0, 0) == reader.event_count(0, 0, 0) == 2 + raw_events = raw_reader.get_event_timestamps(0, 0, 0) + event_timestamps, event_durations, event_labels = reader.get_event_timestamps(0, 0, 0) + np.testing.assert_array_equal(raw_events[0], event_timestamps) + assert raw_events[1] is event_durations is None + np.testing.assert_array_equal(raw_events[2], event_labels) + np.testing.assert_array_equal(event_timestamps, [10, 20]) + assert event_durations is None + np.testing.assert_array_equal(event_labels, ["start", "stop"]) + with pytest.raises(ImportError, match="pip install mtscomp"): + reader.get_analogsignal_chunk(0, 0, 0, 10, 0) + + reader.get_analogsignal_chunk(0, 0, 0, 10, 0) + mtscomp_reader = reader._mtscomp_analogsignal_buffers[0][0]["0"] + reader.get_analogsignal_chunk(0, 0, 10, 20, 0) + assert reader._mtscomp_analogsignal_buffers[0][0]["0"] is mtscomp_reader + assert not mtscomp_reader.cdata.closed + + reader._close_mtscomp_analogsignal_buffers() + assert mtscomp_reader.cdata.closed + assert not hasattr(reader, "_mtscomp_analogsignal_buffers") + + +def test_read_mtscomp_metadata_rejects_invalid_scalar_fields(tmp_path): + metadata_path = tmp_path / "continuous.ch" + valid_metadata = { + "n_channels": 2, + "sample_rate": 100.0, + "dtype": "int16", + "chunk_bounds": [0, 10], + } + invalid_values = ( + ("n_channels", 0, "positive integer"), + ("n_channels", True, "positive integer"), + ("sample_rate", np.inf, "finite and positive"), + ("sample_rate", 0, "finite and positive"), + ("dtype", "not-a-dtype", "NumPy dtype"), + ("chunk_bounds", [], "non-empty list"), + ("chunk_bounds", [1, 10], "must start at 0"), + ) + for key, value, match in invalid_values: + metadata = valid_metadata.copy() + metadata[key] = value + with open(metadata_path, "w", encoding="utf8") as file: + json.dump(metadata, file) + with pytest.raises(ValueError, match=match): + _read_mtscomp_metadata(metadata_path) + + if __name__ == "__main__": unittest.main() diff --git a/pyproject.toml b/pyproject.toml index d21d4181c..2e55a7b42 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -67,6 +67,7 @@ test = [ "sonpy;python_version<'3.10'", "pynwb", "probeinterface", + "mtscomp>=1.0.2", "zugbruecke>=0.2", "wenv" ] @@ -105,6 +106,7 @@ biocam = ["h5py"] med = ["dhn_med_py<2.0"] # ci failing with 2.0 test future version when stable plexon2 = ["zugbruecke>=0.2; sys_platform!='win32'", "wenv; sys_platform!='win32'"] neuralynx = ["python-dateutil"] +mtscomp = ["mtscomp>=1.0.2", "tqdm"] all = [ "coverage", @@ -116,6 +118,7 @@ all = [ "joblib>=1.0.0", "klusta", "matplotlib", + "mtscomp>=1.0.2", "nixio>=1.5.0", "pillow", "probeinterface",