Merge pull request #4 from catalystneuro/add_raw_fiber_photometry

Add raw fiber photometry

make_env.yml

@@ -4,6 +4,7 @@ channels:
 - defaults
 dependencies:
 - python>=3.9
+- bioformats_jar  # required for BioFormats imaging extractor
 - pip
 - pip:
   - -e .  # This calls the setup and therefore requirements minimal

src/howe_lab_to_nwb/vu2024/extractors/bioformats_utils.py

@@ -0,0 +1,86 @@
+from pathlib import Path
+
+
+import numpy as np
+import aicsimageio
+from aicsimageio.formats import FORMAT_IMPLEMENTATIONS
+from neuroconv.utils import FilePathType
+from ome_types import OME
+
+
+def check_file_format_is_supported(file_path: FilePathType):
+    """
+    Check if the file format is supported by BioformatsReader from aicsimageio.
+
+    Returns ValueError if the file format is not supported.
+
+    Parameters
+    ----------
+    file_path : FilePathType
+        Path to the file.
+    """
+    bioformats_reader = "aicsimageio.readers.bioformats_reader.BioformatsReader"
+    supported_file_suffixes = [
+        suffix_name for suffix_name, reader in FORMAT_IMPLEMENTATIONS.items() if bioformats_reader in reader
+    ]
+
+    file_suffix = Path(file_path).suffix.replace(".", "")
+    if file_suffix not in supported_file_suffixes:
+        raise ValueError(f"File '{file_path}' is not supported by BioformatsReader.")
+
+
+def extract_ome_metadata(
+    file_path: FilePathType,
+) -> OME:
+    """
+    Extract OME metadata from a file using aicsimageio.
+
+    Parameters
+    ----------
+    file_path : FilePathType
+        Path to the file.
+    """
+    check_file_format_is_supported(file_path)
+
+    with aicsimageio.readers.bioformats_reader.BioFile(file_path) as reader:
+        ome_metadata = reader.ome_metadata
+
+    return ome_metadata
+
+
+def parse_ome_metadata(metadata: OME) -> dict:
+    """
+    Parse metadata in OME format to extract relevant information and store it standard keys for ImagingExtractors.
+
+    Currently supports:
+    - num_frames
+    - sampling_frequency
+    - num_channels
+    - num_planes
+    - num_rows (height of the image)
+    - num_columns (width of the image)
+    - dtype
+    - channel_names
+
+    """
+    images_metadata = metadata.images[0]
+    pixels_metadata = images_metadata.pixels
+
+    sampling_frequency = None
+    if pixels_metadata.time_increment is not None:
+        sampling_frequency = 1 / pixels_metadata.time_increment
+
+    channel_names = [channel.id for channel in pixels_metadata.channels]
+
+    metadata_parsed = dict(
+        num_frames=images_metadata.pixels.size_t,
+        sampling_frequency=sampling_frequency,
+        num_channels=images_metadata.pixels.size_c,
+        num_planes=images_metadata.pixels.size_z,
+        num_rows=images_metadata.pixels.size_y,
+        num_columns=images_metadata.pixels.size_x,
+        dtype=np.dtype(pixels_metadata.type.numpy_dtype),
+        channel_names=channel_names,
+    )
+
+    return metadata_parsed

src/howe_lab_to_nwb/vu2024/extractors/cxdimagingextractor.py

@@ -0,0 +1,187 @@
+import os
+from pathlib import Path
+from typing import List, Tuple
+
+import aicsimageio
+import numpy as np
+from neuroconv.utils import FilePathType
+from roiextractors import ImagingExtractor
+from roiextractors.extraction_tools import DtypeType
+
+
+class CxdImagingExtractor(ImagingExtractor):
+    """Imaging extractor for reading Hamamatsu Photonics imaging data from .cxd files."""
+
+    extractor_name = "CxdImaging"
+
+    @classmethod
+    def get_available_channels(cls, file_path) -> List[str]:
+        """Get the available channel names from a CXD file produced by Hamamatsu Photonics.
+
+        Parameters
+        ----------
+        file_path : PathType
+            Path to the Bio-Formats file.
+
+        Returns
+        -------
+        channel_names: list
+            List of channel names.
+        """
+        from .bioformats_utils import extract_ome_metadata, parse_ome_metadata
+
+        ome_metadata = extract_ome_metadata(file_path=file_path)
+        parsed_metadata = parse_ome_metadata(metadata=ome_metadata)
+        channel_names = parsed_metadata["channel_names"]
+        return channel_names
+
+    @classmethod
+    def get_available_planes(cls, file_path):
+        """Get the available plane names from a CXD file produced by Hamamatsu Photonics.
+
+        Parameters
+        ----------
+        file_path : PathType
+            Path to the Bio-Formats file.
+
+        Returns
+        -------
+        plane_names: list
+            List of plane names.
+        """
+        from .bioformats_utils import extract_ome_metadata, parse_ome_metadata
+
+        ome_metadata = extract_ome_metadata(file_path=file_path)
+        parsed_metadata = parse_ome_metadata(metadata=ome_metadata)
+        num_planes = parsed_metadata["num_planes"]
+        plane_names = [f"{i}" for i in range(num_planes)]
+        return plane_names
+
+    def __init__(
+        self,
+        file_path: FilePathType,
+        channel_name: str = None,
+        plane_name: str = None,
+        sampling_frequency: float = None,
+    ):
+        r"""
+        Create a CxdImagingExtractor instance from a CXD file produced by Hamamatsu Photonics.
+
+        This extractor requires `bioformats_jar` to be installed in the environment,
+        and requires the java executable to be available on the path (or via the JAVA_HOME environment variable),
+        along with the mvn executable.
+
+        If you are using conda, you can install with `conda install -c conda-forge bioformats_jar`.
+        Note: you may need to reactivate your conda environment after installing.
+        If you are still getting a JVMNotFoundException, try:
+        # mac and linux:
+        `export JAVA_HOME=$CONDA_PREFIX`
+
+        # windows:
+        `set JAVA_HOME=%CONDA_PREFIX%\\Library`
+
+        Parameters
+        ----------
+        file_path : PathType
+            Path to the CXD file.
+        channel_name : str
+            The name of the channel for this extractor. (default=None)
+        plane_name : str
+            The name of the plane for this extractor. (default=None)
+        sampling_frequency : float
+            The sampling frequency of the imaging data. (default=None)
+            Has to be provided manually if not found in the metadata.
+        """
+        from .bioformats_utils import extract_ome_metadata, parse_ome_metadata
+
+        if ".cxd" not in Path(file_path).suffixes:
+            raise ValueError("The file suffix must be .cxd!")
+
+        if "JAVA_HOME" not in os.environ:
+            conda_home = os.environ.get("CONDA_PREFIX")
+            os.environ["JAVA_HOME"] = conda_home
+
+        self.ome_metadata = extract_ome_metadata(file_path=file_path)
+        parsed_metadata = parse_ome_metadata(metadata=self.ome_metadata)
+
+        self._num_frames = parsed_metadata["num_frames"]
+        self._num_channels = parsed_metadata["num_channels"]
+        self._num_planes = parsed_metadata["num_planes"]
+        self._num_rows = parsed_metadata["num_rows"]
+        self._num_columns = parsed_metadata["num_columns"]
+        self._dtype = parsed_metadata["dtype"]
+        self._sampling_frequency = parsed_metadata["sampling_frequency"]
+        self._channel_names = parsed_metadata["channel_names"]
+        self._plane_names = [f"{i}" for i in range(self._num_planes)]
+
+        if channel_name is None:
+            if self._num_channels > 1:
+                raise ValueError(
+                    "More than one channel is detected! Please specify which channel you wish to load "
+                    "with the `channel_name` argument. To see which channels are available, use "
+                    "`CxdImagingExtractor.get_available_channels(file_path=...)`"
+                )
+            channel_name = self._channel_names[0]
+
+        if channel_name not in self._channel_names:
+            raise ValueError(
+                f"The selected channel '{channel_name}' is not a valid channel name."
+                f" The available channel names are: {self._channel_names}."
+            )
+        self.channel_index = self._channel_names.index(channel_name)
+
+        if plane_name is None:
+            if self._num_planes > 1:
+                raise ValueError(
+                    "More than one plane is detected! Please specify which plane you wish to load "
+                    "with the `plane_name` argument. To see which planes are available, use "
+                    "`CxdImagingExtractor.get_available_planes(file_path=...)`"
+                )
+            plane_name = self._plane_names[0]
+
+        if plane_name not in self._plane_names:
+            raise ValueError(
+                f"The selected plane '{plane_name}' is not a valid plane name."
+                f" The available plane names are: {self._plane_names}."
+            )
+        self.plane_index = self._plane_names.index(plane_name)
+
+        if self._sampling_frequency is None:
+            if sampling_frequency is None:
+                raise ValueError(
+                    "Sampling frequency is not found in the metadata. Please provide it manually with the 'sampling_frequency' argument."
+                )
+            self._sampling_frequency = sampling_frequency
+
+        pixels_metadata = self.ome_metadata.images[0].pixels
+        timestamps = [plane.delta_t for plane in pixels_metadata.planes]
+        if np.any(timestamps):
+            self._times = np.array(timestamps)
+
+        with aicsimageio.readers.bioformats_reader.BioFile(file_path) as reader:
+            self._video = reader.to_dask()
+
+        super().__init__(file_path=file_path, channel_name=channel_name, plane_name=plane_name)
+
+    def get_channel_names(self) -> list:
+        return self._channel_names
+
+    def get_dtype(self) -> DtypeType:
+        return self._dtype
+
+    def get_image_size(self) -> Tuple[int, int]:
+        return self._num_rows, self._num_columns
+
+    def get_num_channels(self) -> int:
+        return self._num_channels
+
+    def get_num_frames(self) -> int:
+        return self._num_frames
+
+    def get_sampling_frequency(self):
+        return self._sampling_frequency
+
+    def get_video(self, start_frame=None, end_frame=None, channel: int = 0) -> np.ndarray:
+        video = self._video[start_frame:end_frame, self.channel_index, self.plane_index, ...]
+
+        return video.compute()

src/howe_lab_to_nwb/vu2024/interfaces/__init__.py

@@ -0,0 +1,2 @@
+from .vu2024_fiberphotometryinterface import Vu2024FiberPhotometryInterface
+from .cxdimaginginterface import CxdImagingInterface

src/howe_lab_to_nwb/vu2024/interfaces/cxdimaginginterface.py

@@ -0,0 +1,72 @@
+from typing import Literal
+
+from neuroconv.datainterfaces.ophys.baseimagingextractorinterface import BaseImagingExtractorInterface
+from neuroconv.utils import FilePathType, DeepDict
+
+from howe_lab_to_nwb.vu2024.extractors.cxdimagingextractor import CxdImagingExtractor
+
+
+class CxdImagingInterface(BaseImagingExtractorInterface):
+    """
+    Interface for reading Hamamatsu Photonics imaging data from .cxd files.
+    """
+
+    display_name = "CXD Imaging"
+    associated_suffixes = (".cxd",)
+    info = "Interface for Hamamatsu Photonics CXD files."
+
+    Extractor = CxdImagingExtractor
+
+    def __init__(
+        self,
+        file_path: FilePathType,
+        channel_name: str = None,
+        plane_name: str = None,
+        sampling_frequency: float = None,
+        verbose: bool = True,
+    ):
+        """
+        DataInterface for reading Hamamatsu Photonics imaging data from .cxd files.
+
+        Parameters
+        ----------
+        file_path : FilePathType
+            Path to the CXD file.
+        channel_name : str, optional
+            The name of the channel for this extractor.
+        plane_name : str, optional
+            The name of the plane for this extractor.
+        sampling_frequency : float, optional
+            The sampling frequency of the data. If None, the sampling frequency will be read from the file.
+            If missing from the file, the sampling frequency must be provided.
+        verbose : bool, default: True
+            controls verbosity.
+        """
+        super().__init__(
+            file_path=file_path,
+            channel_name=channel_name,
+            plane_name=plane_name,
+            sampling_frequency=sampling_frequency,
+            verbose=verbose,
+        )
+
+    def get_metadata(
+        self, photon_series_type: Literal["OnePhotonSeries", "TwoPhotonSeries"] = "TwoPhotonSeries"
+    ) -> DeepDict:
+        metadata = super().get_metadata(photon_series_type=photon_series_type)
+
+        device_name = "HamamatsuMicroscope"
+        metadata["Ophys"]["Device"][0].update(
+            name=device_name,
+            manufacturer="Hamamatsu Photonics",
+        )
+        optical_channel_name = "OpticalChannel"  # TODO: add better channel name
+        imaging_plane_metadata = metadata["Ophys"]["ImagingPlane"][0]
+        optical_channel_metadata = imaging_plane_metadata["optical_channel"][0]
+        optical_channel_metadata.update(name=optical_channel_name)
+        imaging_plane_metadata.update(
+            device=device_name,
+            optical_channel=[optical_channel_metadata],
+        )
+
+        return metadata