From 373741c4dca30f541138c202c462a3e5a3eda3f5 Mon Sep 17 00:00:00 2001 From: Greg Knoll Date: Thu, 28 Mar 2024 12:46:09 +0100 Subject: [PATCH] NEW: probe reader for Neuropixels 1.0 in SpikeGadgets .rec file. --- src/probeinterface/__init__.py | 2 + src/probeinterface/io.py | 181 ++++++++++++++++++++++++++++++++- 2 files changed, 181 insertions(+), 2 deletions(-) diff --git a/src/probeinterface/__init__.py b/src/probeinterface/__init__.py index a35ebe0b..57987173 100644 --- a/src/probeinterface/__init__.py +++ b/src/probeinterface/__init__.py @@ -16,6 +16,8 @@ write_imro, read_BIDS_probe, write_BIDS_probe, + read_spikegadgets, + parse_spikegadgets_header, read_spikeglx, parse_spikeglx_meta, get_saved_channel_indices_from_spikeglx_meta, diff --git a/src/probeinterface/io.py b/src/probeinterface/io.py index 7aa2c5bd..8c7a2300 100644 --- a/src/probeinterface/io.py +++ b/src/probeinterface/io.py @@ -19,6 +19,7 @@ from collections import OrderedDict from packaging.version import Version, parse import numpy as np +from xml.etree import ElementTree from . import __version__ from .probe import Probe @@ -1263,12 +1264,188 @@ def write_imro(file: str | Path, probe: Probe): f.write("".join(ret)) +def read_spikegadgets(file: str | Path) -> ProbeGroup: + """ + Find active channels of the given Neuropixels probe from a SpikeGadgets .rec file. + SpikeGadgets headstages support up to three Neuropixels 1.0 probes (as of March 28, 2024), + and information for all probes will be returned in a ProbeGroup object. + + Within the SpikeConfiguration header element (sconf), there is a SpikeNTrode element + for each electrophysiology channel that contains information relevant to scaling and + otherwise displaying the information from that channel, as well as the id of the electrode + from which it is recording ('id'). + + Nested within each SpikeNTrode element is a SpikeChannel element with information about + the electrode dynamically connected to that channel. This contains information relevant + for spike sorting, i.e., its spatial location along the probe shank and the hardware channel + to which it is connected. + + Excerpt of a sample SpikeConfiguration element: + + + + + + ... + + + Parameters + ---------- + file : Path or str + The .rec file path + + Returns + ------- + probe_group : ProbeGroup object + + """ + # ------------------------- # + # Npix 1.0 constants # + # ------------------------- # + TOTAL_NPIX_ELECTRODES = 960 + MAX_ACTIVE_CHANNELS = 384 + CONTACT_WIDTH = 16 # um + CONTACT_HEIGHT = 20 # um + # ------------------------- # + + # Read the header and get Configuration elements + header_txt = parse_spikegadgets_header(file) + root = ElementTree.fromstring(header_txt) + hconf = root.find("HardwareConfiguration") + sconf = root.find("SpikeConfiguration") + + # Get number of probes present (each has its own Device element) + probe_configs = [device for device in hconf if device.attrib['name'] == 'NeuroPixels1'] + n_probes = len(probe_configs) + print(n_probes, "Neuropixels 1.0 probes found:") + + # Container to store Probe objects + probe_group = ProbeGroup() + + for curr_probe in range(1, n_probes + 1): + print(f'Reading probe{curr_probe}...', flush=True, end='') + probe_config = probe_configs[curr_probe - 1] + + # Get number of active channels from probe Device element + active_channel_str = [option for option in probe_config if option.attrib['name'] == 'channelsOn'][0].attrib['data'] + active_channels = [int(ch) for ch in active_channel_str.split(' ') if ch] + n_active_channels = sum(active_channels) + assert len(active_channels) == TOTAL_NPIX_ELECTRODES + assert n_active_channels <= MAX_ACTIVE_CHANNELS + + # Find all channels/electrodes that belong to the current probe + contact_ids = [] + device_channels = [] + positions = np.zeros((n_active_channels, 2)) + + nt_i = 0 # Both probes are in sconf, so need an independent counter of probe electrodes while iterating through + for ntrode in sconf: + electrode_id = ntrode.attrib['id'] + if int(electrode_id[0]) == curr_probe: # first digit of electrode id is probe number + contact_ids.append(electrode_id) + positions[nt_i, :] = (ntrode[0].attrib['coord_ml'], ntrode[0].attrib['coord_dv']) + device_channels.append(ntrode[0].attrib['hwChan']) + nt_i += 1 + assert len(contact_ids) == n_active_channels + + # Construct Probe object + probe = Probe(ndim=2, si_units="um", model_name="Neuropixels 1.0", manufacturer="IMEC") + probe.set_contacts( + contact_ids=contact_ids, + positions=positions, + shapes="square", + shank_ids=None, + shape_params={"width": CONTACT_WIDTH, "height": CONTACT_HEIGHT}, + ) + + # Wire it (i.e., point contact/electrode ids to corresponding hardware/channel ids) + probe.set_device_channel_indices(device_channels) + + # Create a nice polygon background when plotting the probes + x_min = positions[:, 0].min() + x_max = positions[:, 0].max() + x_mid = 0.5 * (x_max + x_min) + y_min = positions[:, 1].min() + y_max = positions[:, 1].max() + polygon_default = [ + (x_min - 20, y_min - CONTACT_HEIGHT / 2), + (x_mid, y_min - 100), + (x_max + 20, y_min - CONTACT_HEIGHT / 2), + (x_max + 20, y_max + 20), + (x_min - 20, y_max + 20), + ] + probe.set_planar_contour(polygon_default) + + # If there are multiple probes, they must be shifted such that they don't occupy the same coordinates. + probe.move([250 * (curr_probe - 1), 0]) + + # Add the probe to the probe container + print(probe.get_contact_count(), 'active channels found.') + probe_group.add_probe(probe) + + return probe_group + + +def parse_spikegadgets_header(file: str | Path) -> str: + """ + Parse file (SpikeGadgets .rec format) into a string until "", + which is the last tag of the header, after which the binary data begins. + """ + header_size = None + with open(file, mode="rb") as f: + while True: + line = f.readline() + if b"" in line: + header_size = f.tell() + break + + if header_size is None: + ValueError("SpikeGadgets: the xml header does not contain ''") + + f.seek(0) + return f.read(header_size).decode("utf8") + + def read_spikeglx(file: str | Path) -> Probe: """ Read probe position for the meta file generated by SpikeGLX See http://billkarsh.github.io/SpikeGLX/#metadata-guides for implementation. - The x_pitch/y_pitch/width are set automatically depending the NP version. + The x_pitch/y_pitch/width are set automatically depending on the NP version. The shape is auto generated as a shank. @@ -1333,7 +1510,7 @@ def read_spikeglx(file: str | Path) -> Probe: def parse_spikeglx_meta(meta_file: str | Path) -> dict: """ Parse the "meta" file from spikeglx into a dict. - All fiields are kept in txt format and must also parsed themself. + All fields are kept in txt format and must also be parsed themselves. """ meta_file = Path(meta_file) with meta_file.open(mode="r") as f: