Merge pull request #40 from ggirelli/ggirelli-patch-parser-for-stitched-images

Removed empty bytes at channel subchunk ends for issue #31

nd2reader/parser.py

@@ -8,9 +8,9 @@ from pims.base_frames import Frame
 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
+from nd2reader import stitched
 
 
 class Parser(object):
@@ -232,8 +232,7 @@ class Parser(object):
         Returns:
 
         """
-        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"]))
+        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):
@@ -268,6 +267,7 @@ class Parser(object):
         timestamp = struct.unpack("d", data[:8])[0]
         image_group_data = array.array("H", data)
         image_data_start = 4 + channel_offset
+        image_group_data = stitched.remove_parsed_unwanted_bytes(image_group_data, image_data_start, height, width)
 
         # 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
@@ -276,7 +276,8 @@ class Parser(object):
         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))))
+            new_width = len(image_group_data[image_data_start::number_of_true_channels]) // height
+            image_data = np.reshape(image_group_data[image_data_start::number_of_true_channels], (height, new_width))
 
         # 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
@@ -285,11 +286,12 @@ class Parser(object):
         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. 
+        # 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)')
+            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):

nd2reader/stitched.py

+# -*- coding: utf-8 -*-
+import numpy as np  # type: ignore
+import warnings
+
+
+def get_unwanted_bytes_ids(image_group_data, image_data_start, height, width):
+    # Check if the byte array size conforms to the image axes size. If not, check
+    # that the number of unexpected (unwanted) bytes is a multiple of the number of
+    # rows (height), as the same unmber of unwanted bytes is expected to be
+    # appended at the end of each row. Then, returns the indexes of the unwanted
+    # bytes.
+    number_of_true_channels = int(len(image_group_data[4:]) / (height * width))
+    n_unwanted_bytes = (len(image_group_data[image_data_start:])) % (height * width)
+    if not n_unwanted_bytes:
+        return np.arange(0)
+    assert 0 == n_unwanted_bytes % height, (
+        "An unexpected number of extra bytes was encountered based on the expected"
+        + " frame size, therefore the file could not be parsed."
+    )
+    return np.arange(
+        image_data_start + height * number_of_true_channels,
+        len(image_group_data) - n_unwanted_bytes + 1,
+        height * number_of_true_channels,
+    )
+
+
+def remove_bytes_by_id(byte_ids, image_group_data, height):
+    # Remove bytes by ID.
+    bytes_per_row = len(byte_ids) // height
+    warnings.warn(
+        f"{len(byte_ids)} ({bytes_per_row}*{height}) unexpected zero "
+        + "bytes were found in the ND2 file and removed to allow further parsing."
+    )
+    for i in range(len(byte_ids)):
+        del image_group_data[byte_ids[i] : (byte_ids[i] + bytes_per_row)]
+
+
+def remove_parsed_unwanted_bytes(image_group_data, image_data_start, height, width):
+    # Stitched ND2 files have been reported to contain unexpected (according to
+    # image shape) zero bytes at the end of each image data row. This hinders
+    # proper reshaping of the data. Hence, here the unwanted zero bytes are
+    # identified and removed.
+    unwanted_byte_ids = get_unwanted_bytes_ids(
+        image_group_data, image_data_start, height, width
+    )
+    if 0 != len(unwanted_byte_ids):
+        assert np.all(
+            image_group_data[unwanted_byte_ids + np.arange(len(unwanted_byte_ids))] == 0
+        ), (
+            f"{len(unwanted_byte_ids)} unexpected non-zero bytes were found"
+            + " in the ND2 file, the file could not be parsed."
+        )
+        remove_bytes_by_id(unwanted_byte_ids, image_group_data, height)
+    return image_group_data

tests/test_parser.py

@@ -2,8 +2,8 @@ import unittest
 from os import path
 from nd2reader.artificial import ArtificialND2
 from nd2reader.common import check_or_make_dir
-from nd2reader.exceptions import InvalidVersionError
 from nd2reader.parser import Parser
+import urllib.request
 
 
 class TestParser(unittest.TestCase):
@@ -13,15 +13,24 @@ class TestParser(unittest.TestCase):
 
     def setUp(self):
         dir_path = path.dirname(path.realpath(__file__))
-        check_or_make_dir(path.join(dir_path, 'test_data/'))
-        self.test_file = path.join(dir_path, 'test_data/test.nd2')
+        check_or_make_dir(path.join(dir_path, "test_data/"))
+        self.test_file = path.join(dir_path, "test_data/test.nd2")
         self.create_test_nd2()
 
     def test_can_open_test_file(self):
         self.create_test_nd2()
 
-        with open(self.test_file, 'rb') as fh:
+        with open(self.test_file, "rb") as fh:
             parser = Parser(fh)
             self.assertTrue(parser.supported)
 
-
+    def test_get_image(self):
+        stitched_path = "test_data/test_stitched.nd2"
+        if not path.isfile(stitched_path):
+            file_name, header = urllib.request.urlretrieve(
+                "https://downloads.openmicroscopy.org/images/ND2/karl/sample_image.nd2",
+                stitched_path,
+            )
+        with open(stitched_path, "rb") as fh:
+            parser = Parser(fh)
+            parser.get_image(0)