zolfa
/
zolfa-nd2reader
1
0
Fork 0
Code Issues 0 Pull Requests 0 Projects 0 Releases 0 Wiki Activity
Browse Source

Regenerated Reader

feature/load_slices
Gabriele Girelli 4 years ago
committed by GitHub
parent
77016ccbaf
commit
b84b52d93d
1 changed files with 135 additions and 243 deletions
Split View
Diff Options
Show Stats Download Patch File Download Diff File
  1. +135
    -243
      nd2reader/reader.py

+ 135
- 243
nd2reader/reader.py View File

@ -1,315 +1,207 @@
# -*- coding: utf-8 -*-
import struct
from pims import Frame
from pims.base_frames import FramesSequenceND
import array
import six
import warnings
from pims.base_frames import Frame
from nd2reader.exceptions import EmptyFileError, InvalidFileType
from nd2reader.parser import Parser
import numpy as np
from nd2reader.common import get_version, read_chunk
from nd2reader.exceptions import InvalidVersionError
from nd2reader.label_map import LabelMap
from nd2reader.raw_metadata import RawMetadata
class ND2Reader(FramesSequenceND):
"""PIMS wrapper for the ND2 parser.
This is the main class: use this to process your .nd2 files.
"""
class Parser(object):
"""Parses ND2 files and creates a Metadata and driver object.
class_priority = 12
"""
CHUNK_HEADER = 0xabeceda
CHUNK_MAP_START = six.b("ND2 FILEMAP SIGNATURE NAME 0001!")
CHUNK_MAP_END = six.b("ND2 CHUNK MAP SIGNATURE 0000001!")
def __init__(self, filename):
super(ND2Reader, self).__init__()
supported_file_versions = {(3, None): True}
if not filename.endswith(".nd2"):
raise InvalidFileType("The file %s you want to read with nd2reader does not have extension .nd2." % filename)
def __init__(self, fh):
self._fh = fh
self._label_map = None
self._raw_metadata = None
self.metadata = None
self.filename = filename
# First check the file version
self.supported = self._check_version_supported()
# first use the parser to parse the file
self._fh = open(filename, "rb")
self._parser = Parser(self._fh)
# Parse the metadata
self._parse_metadata()
# Setup metadata
self.metadata = self._parser.metadata
def calculate_image_properties(self, index):
"""Calculate FOV, channels and z_levels
# Set data type
self._dtype = self._parser.get_dtype_from_metadata()
Args:
index(int): the index (which is simply the order in which the image was acquired)
# Setup the axes
self._setup_axes()
Returns:
tuple: tuple of the field of view, the channel and the z level
# Other properties
self._timesteps = None
"""
field_of_view = self._calculate_field_of_view(index)
channel = self._calculate_channel(index)
z_level = self._calculate_z_level(index)
return field_of_view, channel, z_level
@classmethod
def class_exts(cls):
"""Let PIMS open function use this reader for opening .nd2 files
def get_image(self, index):
"""
Creates an Image object and adds its metadata, based on the index (which is simply the order in which the image
was acquired). May return None if the ND2 contains multiple channels and not all were taken in each cycle (for
example, if you take bright field images every minute, and GFP images every five minutes, there will be some
indexes that do not contain an image. The reason for this is complicated, but suffice it to say that we hope to
eliminate this possibility in future releases. For now, you'll need to check if your image is None if you're
doing anything out of the ordinary.
return {'nd2'} | super(ND2Reader, cls).class_exts()
Args:
index(int): the index (which is simply the order in which the image was acquired)
Returns:
Frame: the image
def close(self):
"""Correctly close the file handle
"""
field_of_view, channel, z_level = self.calculate_image_properties(index)
channel_offset = index % len(self.metadata["channels"])
image_group_number = int(index / len(self.metadata["channels"]))
frame_number = self._calculate_frame_number(image_group_number, field_of_view, z_level)
try:
timestamp, image = self._get_raw_image_data(image_group_number, channel_offset, self.metadata["height"],
self.metadata["width"])
except (TypeError):
return Frame([], frame_no=frame_number, metadata=self._get_frame_metadata())
else:
return Frame(image, frame_no=frame_number, metadata=self._get_frame_metadata())
if self._fh is not None:
self._fh.close()
def get_image_by_attributes(self, frame_number, field_of_view, channel, z_level, height, width):
"""Gets an image based on its attributes alone
def _get_default(self, coord):
try:
return self.default_coords[coord]
except KeyError:
return 0
def get_frame_2D(self, c=0, t=0, z=0, x=0, y=0, v=0):
"""Gets a given frame using the parser
Args:
frame_number: the frame number
field_of_view: the field of view
channel_name: the color channel name
z_level: the z level
height: the height of the image
width: the width of the image
x: The x-index (pims expects this)
y: The y-index (pims expects this)
c: The color channel number
t: The frame number
z: The z stack number
v: The field of view index
Returns:
Frame: the requested image
pims.Frame: The requested frame
"""
frame_number = 0 if frame_number is None else frame_number
field_of_view = 0 if field_of_view is None else field_of_view
channel = 0 if channel is None else channel
z_level = 0 if z_level is None else z_level
# This needs to be set to width/height to return an image
x = self.metadata["width"]
y = self.metadata["height"]
image_group_number = self._calculate_image_group_number(frame_number, field_of_view, z_level)
try:
timestamp, raw_image_data = self._get_raw_image_data(image_group_number, channel,
height, width)
except (TypeError):
return Frame([], frame_no=frame_number, metadata=self._get_frame_metadata())
else:
return Frame(raw_image_data, frame_no=frame_number, metadata=self._get_frame_metadata())
@staticmethod
def get_dtype_from_metadata():
"""Determine the data type from the metadata.
For now, always use float64 to prevent unexpected overflow errors when manipulating the data (calculating sums/
means/etc.)
return self._parser.get_image_by_attributes(t, v, c, z, y, x)
@property
def parser(self):
"""
return np.float64
def _check_version_supported(self):
"""Checks if the ND2 file version is supported by this reader.
Returns the parser object.
Returns:
bool: True on supported
Parser: the parser object
"""
major_version, minor_version = get_version(self._fh)
supported = self.supported_file_versions.get(
(major_version, minor_version)) or self.supported_file_versions.get((major_version, None))
if not supported:
print("Warning: No parser is available for your current ND2 version (%d.%d). " % (
major_version, minor_version) + "This might lead to unexpected behaviour.")
return self._parser
return supported
@property
def pixel_type(self):
"""Return the pixel data type
def _parse_metadata(self):
"""Reads all metadata and instantiates the Metadata object.
Returns:
dtype: the pixel data type
"""
# Retrieve raw metadata from the label mapping
self._label_map = self._build_label_map()
self._raw_metadata = RawMetadata(self._fh, self._label_map)
self.metadata = self._raw_metadata.__dict__
self.acquisition_times = self._raw_metadata.acquisition_times
return self._dtype
def _build_label_map(self):
"""
Every label ends with an exclamation point, however, we can't directly search for those to find all the labels
as some of the bytes contain the value 33, which is the ASCII code for "!". So we iteratively find each label,
grab the subsequent data (always 16 bytes long), advance to the next label and repeat.
@property
def timesteps(self):
"""Get the timesteps of the experiment
Returns:
LabelMap: the computed label map
np.ndarray: an array of times in milliseconds.
"""
# go 8 bytes back from file end
self._fh.seek(-8, 2)
chunk_map_start_location = struct.unpack("Q", self._fh.read(8))[0]
self._fh.seek(chunk_map_start_location)
raw_text = self._fh.read(-1)
return LabelMap(raw_text)
def _calculate_field_of_view(self, index):
"""Determines what field of view was being imaged for a given image.
if self._timesteps is None:
return self.get_timesteps()
return self._timesteps
Args:
index(int): the index (which is simply the order in which the image was acquired)
@property
def events(self):
"""Get the events of the experiment
Returns:
int: the field of view
iterator of events as dict
"""
images_per_cycle = len(self.metadata["z_levels"]) * len(self.metadata["channels"])
return int((index - (index % images_per_cycle)) / images_per_cycle) % len(self.metadata["fields_of_view"])
def _calculate_channel(self, index):
"""Determines what channel a particular image is.
Args:
index(int): the index (which is simply the order in which the image was acquired)
return self._get_metadata_property("events")
@property
def frame_rate(self):
"""The (average) frame rate
Returns:
string: the name of the color channel
float: the (average) frame rate in frames per second
"""
return self.metadata["channels"][index % len(self.metadata["channels"])]
total_duration = 0.0
def _calculate_z_level(self, index):
"""Determines the plane in the z-axis a given image was taken in.
for loop in self.metadata['experiment']['loops']:
total_duration += loop['duration']
In the future, this will be replaced with the actual offset in micrometers.
if total_duration == 0:
total_duration = self.timesteps[-1]
Args:
index(int): the index (which is simply the order in which the image was acquired)
if total_duration == 0:
raise ValueError('Total measurement duration could not be determined from loops')
Returns:
The z level
return self.metadata['num_frames'] / (total_duration/1000.0)
"""
return self.metadata["z_levels"][int(
((index - (index % len(self.metadata["channels"]))) / len(self.metadata["channels"])) % len(
self.metadata["z_levels"]))]
def _get_metadata_property(self, key, default=None):
if self.metadata is None:
return default
def _calculate_image_group_number(self, frame_number, fov, z_level):
"""
Images are grouped together if they share the same time index, field of view, and z-level.
if key not in self.metadata:
return default
Args:
frame_number: the time index
fov: the field of view number
z_level: the z level number
if self.metadata[key] is None:
return default
Returns:
int: the image group number
return self.metadata[key]
def _setup_axes(self):
"""Setup the xyctz axes, iterate over t axis by default
"""
z_length = len(self.metadata['z_levels'])
z_length = z_length if z_length > 0 else 1
fields_of_view = len(self.metadata["fields_of_view"])
fields_of_view = fields_of_view if fields_of_view > 0 else 1
self._init_axis_if_exists('x', self._get_metadata_property("width", default=0))
self._init_axis_if_exists('y', self._get_metadata_property("height", default=0))
self._init_axis_if_exists('c', len(self._get_metadata_property("channels", default=[])), min_size=2)
self._init_axis_if_exists('t', len(self._get_metadata_property("frames", default=[])))
self._init_axis_if_exists('z', len(self._get_metadata_property("z_levels", default=[])), min_size=2)
self._init_axis_if_exists('v', len(self._get_metadata_property("fields_of_view", default=[])), min_size=2)
return frame_number * fields_of_view * z_length + (fov * z_length + z_level)
if len(self.sizes) == 0:
raise EmptyFileError("No axes were found for this .nd2 file.")
def _calculate_frame_number(self, image_group_number, field_of_view, z_level):
"""
Images are in the same frame if they share the same group number and field of view and are taken sequentially.
# provide the default
self.iter_axes = self._guess_default_iter_axis()
Args:
image_group_number: the image group number (see _calculate_image_group_number)
field_of_view: the field of view number
z_level: the z level number
self._register_get_frame(self.get_frame_2D, 'yx')
Returns:
def _init_axis_if_exists(self, axis, size, min_size=1):
if size >= min_size:
self._init_axis(axis, size)
def _guess_default_iter_axis(self):
"""
return (image_group_number - (field_of_view * len(self.metadata["z_levels"]) + z_level)) / (
len(self.metadata["fields_of_view"]) * len(self.metadata["z_levels"]))
@property
def _channel_offset(self):
Guesses the default axis to iterate over based on axis sizes.
Returns:
the axis to iterate over
"""
Image data is interleaved for each image set. That is, if there are four images in a set, the first image
will consist of pixels 1, 5, 9, etc, the second will be pixels 2, 6, 10, and so forth.
priority = ['t', 'z', 'c', 'v']
found_axes = []
for axis in priority:
try:
current_size = self.sizes[axis]
except KeyError:
continue
Returns:
dict: the channel offset for each channel
if current_size > 1:
return axis
"""
return {channel: n for n, channel in enumerate(self.metadata["channels"])}
def _remove_unwanted_bytes(self, image_group_data, image_data_start, height, width):
# Remove unwanted 0-bytes that can appear in stitched images
number_of_true_channels = int(len(image_group_data[4:]) / (height * width))
unwanted_bytes_len = (len(image_group_data[image_data_start:]))%(height*width)
if unwanted_bytes_len:
warnings.warn('Identified unwanted bytes in the ND2 file, possibly stitched.')
byte_ids = range(image_data_start+height*number_of_true_channels, len(image_group_data)-unwanted_bytes_len+1, height*number_of_true_channels)
if all([0 == image_group_data[byte_ids[i]+i] for i in range(len(byte_ids))]):
warnings.warn('All unwanted bytes are zero-bytes, correctly removed.')
for i in range(len(byte_ids)):
del image_group_data[byte_ids[i]]
def _get_raw_image_data(self, image_group_number, channel_offset, height, width):
"""Reads the raw bytes and the timestamp of an image.
found_axes.append(axis)
Args:
image_group_number: the image group number (see _calculate_image_group_number)
channel_offset: the number of the color channel
height: the height of the image
width: the width of the image
return found_axes[0]
def get_timesteps(self):
"""Get the timesteps of the experiment
Returns:
np.ndarray: an array of times in milliseconds.
"""
chunk = self._label_map.get_image_data_location(image_group_number)
data = read_chunk(self._fh, chunk)
# All images in the same image group share the same timestamp! So if you have complicated image data,
# your timestamps may not be entirely accurate. Practically speaking though, they'll only be off by a few
# seconds unless you're doing something super weird.
timestamp = struct.unpack("d", data[:8])[0]
image_group_data = array.array("H", data)
image_data_start = 4 + channel_offset
# The images for the various channels are interleaved within the same array. For example, the second image
# of a four image group will be composed of bytes 2, 6, 10, etc. If you understand why someone would design
# a data structure that way, please send the author of this library a message.
number_of_true_channels = int(len(image_group_data[4:]) / (height * width))
self._remove_unwanted_bytes(image_group_data, image_data_start, height, width)
try:
image_data = np.reshape(image_group_data[image_data_start::number_of_true_channels], (height, width))
except ValueError:
image_data = np.reshape(image_group_data[image_data_start::number_of_true_channels], (height, int(round(len(image_group_data[image_data_start::number_of_true_channels])/height))))
# Skip images that are all zeros! This is important, since NIS Elements creates blank "gap" images if you
# don't have the same number of images each cycle. We discovered this because we only took GFP images every
# other cycle to reduce phototoxicity, but NIS Elements still allocated memory as if we were going to take
# them every cycle.
if np.any(image_data):
return timestamp, image_data
# If a blank "gap" image is encountered, generate an array of corresponding height and width to avoid
# errors with ND2-files with missing frames. Array is filled with nan to reflect that data is missing.
else:
empty_frame = np.full((height, width), np.nan)
warnings.warn('ND2 file contains gap frames which are represented by np.nan-filled arrays; to convert to zeros use e.g. np.nan_to_num(array)')
return timestamp, image_data
def _get_frame_metadata(self):
"""Get the metadata for one frame
if self._timesteps is not None and len(self._timesteps) > 0:
return self._timesteps
Returns:
dict: a dictionary containing the parsed metadata
self._timesteps = np.array(list(self._parser._raw_metadata.acquisition_times), dtype=np.float) * 1000.0
"""
return self.metadata
return self._timesteps

Write Preview
Loading…
Cancel
Save