bugfix in the reader

npy2bdv/examples.py

@@ -23,9 +23,10 @@ def generate_test_image(dim_yx):
     stack.append(plane)
 stack = np.asarray(stack)
 
-if not os.path.exists("./test"):
-    os.mkdir("./test")
-fname = "./test/ex1_t2_ch2_illum2_angle2.h5"
+examples_dir = "./examples/"
+if not os.path.exists(examples_dir):
+    os.mkdir(examples_dir)
+fname = examples_dir + "ex1_t2_ch2_illum2_angle2.h5"
 bdv_writer = npy2bdv.BdvWriter(fname, nchannels=2, nilluminations=2, nangles=2, subsamp=((1, 1, 1),))
 for t in range(2):
     for i_ch in range(2):
@@ -35,7 +36,7 @@ def generate_test_image(dim_yx):
 
 bdv_writer.write_xml_file(ntimes=2)
 bdv_writer.close()
-print("dataset in " + fname)
+print(f"dataset in {fname}")
 
 #########################
 # 2. Writing speed test #
@@ -46,7 +47,7 @@ def generate_test_image(dim_yx):
 start_time_total = time.time()
 i_stacks = 0
 time_list = []
-fname = "./test/ex2_t20_chan2.h5"
+fname = examples_dir + "ex2_t20_chan2.h5"
 bdv_writer = npy2bdv.BdvWriter(fname, nchannels=2, subsamp=((1, 1, 1),))
 for ichannel in range(nchannels):
     for itime in range(ntimes):
@@ -59,9 +60,9 @@ def generate_test_image(dim_yx):
 bdv_writer.write_xml_file(ntimes=ntimes)
 bdv_writer.close()
 time_per_stack = (time.time() - start_time_total) / i_stacks
-print("H5 mean writing time per stack: {:1.3f}".format(time_per_stack) + " sec.")
-print("H5 mean writing speed: " + str(int(sys.getsizeof(stack) / time_per_stack / 1e6)) + " MB/s")
-print("speed test dataset in " + fname)
+print(f"H5 mean writing time per stack: {time_per_stack:1.3f} sec.")
+print(f"H5 mean writing speed: {int(sys.getsizeof(stack) / time_per_stack / 1e6)} MB/s")
+print(f"speed test dataset in {fname}")
 
 #################################################################
 # 3. Writing with affine transformations defined in XML file ####
@@ -73,21 +74,21 @@ def generate_test_image(dim_yx):
                           (0.0, 1.0, 0.0, 0.0),
                           (0.0, 0.0, 1.0, 0.0)))
 
-fname = "./test/ex3_t1_ch1_unshear.h5"
+fname = examples_dir + "ex3_t1_ch1_unshear.h5"
 bdv_writer = npy2bdv.BdvWriter(fname, nchannels=1, subsamp=((1, 1, 1),))
 bdv_writer.append_view(stack, time=0, channel=0,
                        m_affine=affine_matrix,
                        name_affine="unshearing transformation",
                        calibration=(1, 1, 1))
 bdv_writer.write_xml_file(ntimes=1)
 bdv_writer.close()
-print("unsheared stack in " + fname)
+print(f"unsheared stack in {fname}")
 
 ############################################
 # 4. Writing with experiment metadata ######
 ############################################
 print("Example4: writing 1 time point and 1 channel with voxel size, exposure, camera and microscope properties")
-fname = "./test/ex4_t1_ch1_cam_props.h5"
+fname = examples_dir + "ex4_t1_ch1_cam_props.h5"
 bdv_writer = npy2bdv.BdvWriter(fname, nchannels=1, subsamp=((1, 1, 1),))
 bdv_writer.append_view(stack, time=0, channel=0,
                        voxel_size_xyz=(1, 1, 5), voxel_units='um',
@@ -96,29 +97,29 @@ def generate_test_image(dim_yx):
                           microscope_name='Superscope',
                           user_name='nvladimus')
 bdv_writer.close()
-print("dataset is in " + fname)
+print("dataset is in {fname}")
 
 ################################################
 # 5. Writing with subsampling and compression ##
 ################################################
 print("Example5: 1 time point and 1 channel with 3-level subsampling and compression")
-fname = "./test/ex5_t1_ch1_level3_gzip.h5"
+fname = examples_dir + "ex5_t1_ch1_level3_gzip.h5"
 bdv_writer = npy2bdv.BdvWriter(fname, nchannels=1,
                                subsamp=((1, 1, 1), (2, 4, 4), (4, 16, 16)),
                                blockdim=((64, 64, 64),),
                                compression='gzip')
 bdv_writer.append_view(stack, time=0, channel=0)
 bdv_writer.write_xml_file(ntimes=1)
 bdv_writer.close()
-print("dataset is in " + fname)
+print("dataset is in {fname}")
 
 ##########################################################
 # 6. Writing virtual stacks that are too big to fit RAM ##
 ##########################################################
 print("Example6: 1 time point, 2 channels, HUGE virtual stack, 20 GB!")
 stack_dim_zyx = (250, 3648, 5472)
 image_plane = generate_test_image(stack_dim_zyx[1:])
-fname = "./test/ex6_t1_ch1_huge_virtual.h5"
+fname = examples_dir + "ex6_t1_ch1_huge_virtual.h5"
 bdv_writer = npy2bdv.BdvWriter(fname, nchannels=2,
                                blockdim=((1, int(image_plane.shape[0]/4), int(image_plane.shape[1]/4)),))
 
@@ -129,27 +130,27 @@ def generate_test_image(dim_yx):
 
 bdv_writer.write_xml_file(ntimes=1)
 bdv_writer.close()
-print("virtual stack is in " + fname)
+print("virtual stack is in {fname}")
 
 ############################################
 ## 7. Missing views, normal stack writing ##
 ############################################
 print("Example7: Automatic calculation of missing views.")
-fname = "./test/ex7_missing_views.h5"
+fname = examples_dir + "ex7_missing_views.h5"
 bdv_writer = npy2bdv.BdvWriter(fname, nchannels=2, subsamp=((1, 1, 1),))
 bdv_writer.append_view(stack, time=0, channel=0)
 bdv_writer.append_view(stack, time=1, channel=1)
 bdv_writer.write_xml_file(ntimes=2)
 bdv_writer.close()
-print("dataset with missing views in " + fname)
+print("dataset with missing views in {fname}")
 
 #####################################
 ## 8. Missing views, virtual stack ##
 #####################################
 print("Example8: Automatic calculation of missing views, virtual stack.")
 stack_dim_zyx = (50, 1000, 2000)
 image_plane = generate_test_image(stack_dim_zyx[1:])
-fname = "./test/ex8_virtual_stack_missing_views.h5"
+fname = examples_dir + "ex8_virtual_stack_missing_views.h5"
 bdv_writer = npy2bdv.BdvWriter(fname, nchannels=2, subsamp=((1, 1, 1),))
 bdv_writer.append_view(stack=None, virtual_stack_dim=stack_dim_zyx, time=0, channel=0)
 bdv_writer.append_view(stack=None, virtual_stack_dim=stack_dim_zyx, time=1, channel=1)
@@ -160,4 +161,4 @@ def generate_test_image(dim_yx):
 
 bdv_writer.write_xml_file(ntimes=2)
 bdv_writer.close()
-print("dataset with missing views in " + fname)
+print("dataset with missing views in {fname}")

npy2bdv/npy2bdv.py

@@ -9,6 +9,8 @@
 
 
 class BdvWriter:
+    __version__ = "2020.06"
+
     def __init__(self, filename,
                  subsamp=((1, 1, 1),),
                  blockdim=((4, 256, 256),),
@@ -45,7 +47,7 @@ def __init__(self, filename,
             print("Number of blockdim levels (" + str (len(blockdim)) +
                   ") is less than subsamp levels (" + str(len(subsamp)) + ")\n" +
                   "First-level block size " + str(blockdim[0]) + " will be used for all levels")
-
+        self._fmt = 't{:05d}/s{:02d}/{}'
         self.nsetups = nilluminations * nchannels * ntiles * nangles
         self.nilluminations = nilluminations
         self.nchannels = nchannels
@@ -68,7 +70,6 @@ def __init__(self, filename,
         self._write_setups_header()
         self.virtual_stacks = False
         self.setup_id_present = [[False] * self.nsetups]
-        self.__version__ = "2020.03"
 
     def _write_setups_header(self):
         """Write resolutions and subdivisions for all setups into h5 file."""
@@ -102,9 +103,8 @@ def append_plane(self, plane, plane_index, time=0, illumination=0, channel=0, ti
         assert plane.shape == self.stack_shapes[isetup][1:], "Plane dimensions must match (y,x) size of virtual stack."
         assert plane_index < self.stack_shapes[isetup][0], "Plane index must be less than virtual stack z-dimension."
         assert self.nlevels == 1, "No subsampling currently implemented for virtual stack writing."
-        fmt = 't{:05d}/s{:02d}/{}'
         for ilevel in range(self.nlevels):
-            group_name = fmt.format(time, isetup, ilevel)
+            group_name = self._fmt.format(time, isetup, ilevel)
             dataset = self.file_object[group_name]["cells"]
             dataset[plane_index, :, :] = plane.astype('int16')
 
@@ -143,7 +143,6 @@ def append_view(self, stack, virtual_stack_dim=None,
         """
         assert len(calibration) == 3, "Calibration must be a tuple of 3 elements (x, y, z)."
         assert len(voxel_size_xyz) == 3, "Voxel size must be a tuple of 3 elements (x, y, z)."
-        fmt = 't{:05d}/s{:02d}/{}'
         isetup = self._determine_setup_id(illumination, channel, tile, angle)
         self.update_setup_id_present(isetup, time)
         if stack is not None:
@@ -155,7 +154,7 @@ def append_view(self, stack, virtual_stack_dim=None,
             self.virtual_stacks = True
 
         for ilevel in range(self.nlevels):
-            group_name = fmt.format(time, isetup, ilevel)
+            group_name = self._fmt.format(time, isetup, ilevel)
             if group_name in self.file_object:
                 del self.file_object[group_name]
             grp = self.file_object.create_group(group_name)
@@ -391,42 +390,42 @@ def close(self):
 
 
 class BdvReader:
-    def __init__(self, filename_h5=None):
+    __version__ = "2020.06.25"
+
+    def __init__(self, filename_h5):
         """Class for reading a BigDataViewer/BigStitcher HDF5 file into numpy array.
 
         Constructor parameters
             filename_h5: (string), optional, full path to a HDF5 file (default None).
         """
-        self.__version__ = "2020.01.07"
-        self.__fmt = 't{:05d}/s{:02d}/{}'
-        self.filename_h5 = filename_h5
-        if filename_h5 is not None:
-            assert os.path.exists(filename_h5), "Error: HDF5 file not found"
-            self.__file_object = h5py.File(filename_h5, 'r')
-        else:
-            self.__file_object = None
+        self._fmt = 't{:05d}/s{:02d}/{}'
+        assert os.path.exists(filename_h5), "Error: HDF5 file not found"
+        self._file_object = h5py.File(filename_h5, 'r')
 
     def set_path_h5(self, filename_h5):
         """Set the file path to HDF5 file. If another file was already open, it closes it before proceeding"""
         assert os.path.exists(filename_h5), "Error: HDF5 file not found"
-        if self.__file_object is not None:
-            self.__file_object.close()
-        self.__file_object = h5py.File(filename_h5, 'r')
+        if self._file_object:
+            self._file_object.close()
+        self._file_object = h5py.File(filename_h5, 'r')
 
     def read_view(self, time=0, isetup=0, ilevel=0):
-        """Read a view (stack) specified by its time, setup ID, and downsampling level.
+        """Read a view (stack) specified by its time, setup ID, and downsampling level into numpy array (uint16).
 
         Parameters
             time: (int) index of time point (default 0).
             isetup: (int), index of setup ID (default 0)
             ilevel: (int), level of subsampling, if available (default 0, no subsampling)
 
         Returns
-            dataset: a numpy array (ndim=3)"""
-        group_name = self.__fmt.format(time, isetup, ilevel)
-        dataset = self.__file_object[group_name]["cells"]
-        return dataset
+            dataset: a numpy array (dim=3, dtype=uint16)"""
+        group_name = self._fmt.format(time, isetup, ilevel)
+        if self._file_object:
+            dataset = self._file_object[group_name]["cells"].value.astype('uint16')
+            return dataset
+        else:
+            raise ValueError('File object is None')
 
     def close(self):
         """Close the file object."""
-        self.__file_object.close()
+        self._file_object.close()

npy2bdv/test.py

@@ -0,0 +1,20 @@
+import npy2bdv
+import os
+
+examples_dir = "./examples/"
+assert os.path.exists(examples_dir), 'Please run the examples first to generate the datasets.'
+
+
+def check_range():
+    """"Check if the reader imports full uint16 range correctly"""
+    fname = examples_dir + "ex1_t2_ch2_illum2_angle2.h5"
+    assert os.path.exists(fname), f'File {fname} not found.'
+    reader = npy2bdv.BdvReader(fname)
+    view0 = reader.read_view(time=0, isetup=0)
+    print(f"Array min, max, mean: {view0.min()}, {view0.max()}, {int(view0.mean())}")
+    assert view0.min() > 0, "Min() value incorrect: {view0.min()}"
+    assert view0.max() < 65535, "Max() value incorrect: {view0.max()}"
+    reader.close()
+
+check_range()
+

setup.py

@@ -2,15 +2,15 @@
 
 setup(
         name='npy2bdv',
-        version='2020.03.25',
+        version='2020.06.25',
         description='Package for writing/reading 3d numpy arrays to/from HDF5 files (Fiji/BigDataViewer/BigStitcher format).',
         url='https://github.com/nvladimus/npy2bdv',
         author='Nikita Vladimirov',
         author_email="[email protected]",
         install_requires=[
             'h5py',
             'numpy',
-            'scikit-image',
+            'scikit-image'
         ],
         packages=find_packages()
 )