From bdd7a5399cf1cf67f9d5b83805a8c4b8f83a1b85 Mon Sep 17 00:00:00 2001 From: Wim Pomp Date: Mon, 14 Aug 2023 17:01:03 +0200 Subject: [PATCH] - set dtype according to pixel type in file - cziread bugfix - add reader for broken files saved by Fiji - make ndread work for arrays with less dimensions than 5 - relative imports - remove some old functions - make bfread check if it can open the file --- ndbioimage/__init__.py | 72 +++++++++++++--------------------- ndbioimage/jvm.py | 18 +++++++-- ndbioimage/readers/bfread.py | 34 ++++++++++------ ndbioimage/readers/cziread.py | 4 +- ndbioimage/readers/fijiread.py | 59 ++++++++++++++++++++++++++++ ndbioimage/readers/ndread.py | 11 +++--- ndbioimage/readers/seqread.py | 2 +- ndbioimage/readers/tifread.py | 7 ++-- pyproject.toml | 2 +- tests/test_open.py | 12 +++--- 10 files changed, 143 insertions(+), 78 deletions(-) create mode 100644 ndbioimage/readers/fijiread.py diff --git a/ndbioimage/__init__.py b/ndbioimage/__init__.py index e0e0034..2d0b4af 100755 --- a/ndbioimage/__init__.py +++ b/ndbioimage/__init__.py @@ -17,13 +17,11 @@ from parfor import parfor from tiffwrite import IJTiffFile from numbers import Number from argparse import ArgumentParser -from typing import List from pathlib import Path from importlib.metadata import version from traceback import print_exc -from ndbioimage.transforms import Transform, Transforms -from ndbioimage.jvm import JVM - +from .transforms import Transform, Transforms +from .jvm import JVM try: __version__ = version(Path(__file__).parent.name) @@ -42,8 +40,13 @@ ureg.default_format = '~P' set_application_registry(ureg) +class ReaderNotFoundError(Exception): + pass + + class ImTransforms(Transforms): """ Transforms class with methods to calculate channel transforms from bead files etc. """ + def __init__(self, path, cyllens, file=None, transforms=None): super().__init__() self.cyllens = cyllens @@ -179,6 +182,7 @@ class ImTransforms(Transforms): class ImShiftTransforms(Transforms): """ Class to handle drift in xy. The image this is applied to must have a channeltransform already, which is then replaced by this class. """ + def __init__(self, im, shifts=None): """ im: Calculate shifts from channel-transformed images im, t x 2 array Sets shifts from array, one row per frame @@ -242,16 +246,19 @@ class ImShiftTransforms(Transforms): @parfor(range(1, im.shape['t']), (im, im0), desc='Calculating image shifts.') def fun(t, im, im0): return Transform(im0, im[:, 0, t].squeeze().transpose(2, 0, 1), 'translation') + transforms = [Transform()] + fun self.shifts = np.array([t.parameters[4:] for t in transforms]) self.set_transforms(transforms, im.transform) def calulate_shifts(self, im): """ Calculate shifts relative to the previous frame """ + @parfor(range(1, im.shape['t']), (im,), desc='Calculating image shifts.') def fun(t, im): - return Transform(im[:, 0, t-1].squeeze().transpose(2, 0, 1), im[:, 0, t].squeeze().transpose(2, 0, 1), + return Transform(im[:, 0, t - 1].squeeze().transpose(2, 0, 1), im[:, 0, t].squeeze().transpose(2, 0, 1), 'translation') + transforms = [Transform()] + fun self.shifts = np.cumsum([t.parameters[4:] for t in transforms]) self.set_transforms(transforms, im.transform) @@ -285,40 +292,12 @@ class DequeDict(OrderedDict): self.__truncate__() -def tolist(item) -> List: - if hasattr(item, 'items'): - return item - elif isinstance(item, str): - return [item] - try: - iter(item) - return list(item) - except TypeError: - return list((item,)) - - def find(obj, **kwargs): for item in obj: if all([getattr(item, key) == value for key, value in kwargs.items()]): return item -def find_rec(obj, **kwargs): - if isinstance(obj, list): - for item in obj: - ans = find_rec(item, **kwargs) - if ans: - return ans - elif isinstance(obj, ome_types._base_type.OMEType): - if all([hasattr(obj, key) for key in kwargs.keys()]) and all( - [getattr(obj, key) == value for key, value in kwargs.items()]): - return obj - for k, v in obj.__dict__.items(): - ans = find_rec(v, **kwargs) - if ans: - return ans - - def try_default(fun, default, *args, **kwargs): try: return fun(*args, **kwargs) @@ -326,9 +305,10 @@ def try_default(fun, default, *args, **kwargs): return default -def ome_subprocess(path): +def get_ome(path): + from .readers.bfread import jars try: - jvm = JVM() + jvm = JVM(jars) ome_meta = jvm.metadata_tools.createOMEXMLMetadata() reader = jvm.image_reader() reader.setMetadataStore(ome_meta) @@ -450,7 +430,7 @@ class Imread(np.lib.mixins.NDArrayOperatorsMixin, metaclass=ABCMeta): def get_ome(path): """ Use java BioFormats to make an ome metadata structure. """ with multiprocessing.get_context('spawn').Pool(1) as pool: - ome = pool.map(ome_subprocess, (path,))[0] + ome = pool.map(get_ome, (path,))[0] return ome def __new__(cls, path=None, *args, **kwargs): @@ -479,6 +459,7 @@ class Imread(np.lib.mixins.NDArrayOperatorsMixin, metaclass=ABCMeta): subclass.do_copy = set(do_copy).union(set(subclass_do_copy)) return super().__new__(subclass) + raise ReaderNotFoundError(f'No reader found for {path}.') @staticmethod def split_path_series(path): @@ -507,7 +488,6 @@ class Imread(np.lib.mixins.NDArrayOperatorsMixin, metaclass=ABCMeta): self.transform = transform self.drift = drift self.beadfile = beadfile - self.dtype = dtype self.reader = None self.pcf = None @@ -523,14 +503,13 @@ class Imread(np.lib.mixins.NDArrayOperatorsMixin, metaclass=ABCMeta): self.frameoffset = 0, 0 # how far apart the centers of frame and sensor are self.flags = dict(C_CONTIGUOUS=False, F_CONTIGUOUS=False, OWNDATA=False, WRITEABLE=False, ALIGNED=False, WRITEBACKIFCOPY=False, UPDATEIFCOPY=False) - self.open() - # extract some metadata from ome instrument = self.ome.instruments[0] if self.ome.instruments else None image = self.ome.images[0] pixels = image.pixels self.shape = pixels.size_x, pixels.size_y, pixels.size_c, pixels.size_z, pixels.size_t + self.dtype = pixels.type.value if dtype is None else dtype self.pxsize = pixels.physical_size_x_quantity try: self.exposuretime = tuple(find(image.pixels.planes, the_c=c).exposure_time_quantity @@ -556,7 +535,7 @@ class Imread(np.lib.mixins.NDArrayOperatorsMixin, metaclass=ABCMeta): try: self.binning = [int(i) for i in image.pixels.channels[0].detector_settings.binning.value.split('x')] self.pxsize *= self.binning[0] - except (AttributeError, IndexError): + except (AttributeError, IndexError, ValueError): self.binning = None self.cnamelist = [channel.name for channel in image.pixels.channels] try: @@ -718,6 +697,7 @@ class Imread(np.lib.mixins.NDArrayOperatorsMixin, metaclass=ABCMeta): s = [f"path/filename: {self.path}", f"series/pos: {self.series}", f"reader: {self.__class__.__module__.split('.')[-1]}", + f"dtype: {self.dtype}", f"shape ({self.axes}):".ljust(15) + f"{' x '.join(str(i) for i in self.shape)}"] if self.pxsize_um: s.append(f'pixel size: {1000 * self.pxsize_um:.2f} nm') @@ -730,10 +710,10 @@ class Imread(np.lib.mixins.NDArrayOperatorsMixin, metaclass=ABCMeta): if self.binning: s.append('binning: {}x{}'.format(*self.binning)) if self.laserwavelengths: - s.append('laser colors: ' + ' | '.join([' & '.join(len(w)*('{:.0f}',)).format(*w) + s.append('laser colors: ' + ' | '.join([' & '.join(len(w) * ('{:.0f}',)).format(*w) for w in self.laserwavelengths]) + ' nm') if self.laserpowers: - s.append('laser powers: ' + ' | '.join([' & '.join(len(p)*('{:.3g}',)).format(*[100 * i for i in p]) + s.append('laser powers: ' + ' | '.join([' & '.join(len(p) * ('{:.3g}',)).format(*[100 * i for i in p]) for p in self.laserpowers]) + ' %') if self.objective: s.append('objective: {}'.format(self.objective.model)) @@ -819,6 +799,7 @@ class Imread(np.lib.mixins.NDArrayOperatorsMixin, metaclass=ABCMeta): for k in b: if k not in a: yield k + for idx in unique_yield(list(self.cache.keys()), product(range(self.shape['c']), range(self.shape['z']), range(self.shape['t']))): xyczt = (slice(None), slice(None)) + idx @@ -1226,6 +1207,7 @@ class Imread(np.lib.mixins.NDArrayOperatorsMixin, metaclass=ABCMeta): @cached_property def piezoval(self): """ gives the height of the piezo and focus motor, only available when CylLensGUI was used """ + def upack(idx): time = list() val = list() @@ -1291,7 +1273,7 @@ class Imread(np.lib.mixins.NDArrayOperatorsMixin, metaclass=ABCMeta): return return - def save_as_tiff(self, fname=None, c=None, z=None, t=None, split=False, bar=True, pixel_type='uint16'): + def save_as_tiff(self, fname=None, c=None, z=None, t=None, split=False, bar=True, pixel_type='uint16', **kwargs): """ saves the image as a tiff-file split: split channels into different files """ if fname is None: @@ -1319,7 +1301,7 @@ class Imread(np.lib.mixins.NDArrayOperatorsMixin, metaclass=ABCMeta): shape = [len(i) for i in n] at_least_one = False with IJTiffFile(fname.with_suffix('.tif'), shape, pixel_type, - pxsize=self.pxsize_um, deltaz=self.deltaz_um) as tif: + pxsize=self.pxsize_um, deltaz=self.deltaz_um, **kwargs) as tif: for i, m in tqdm(zip(product(*[range(s) for s in shape]), product(*n)), total=np.prod(shape), desc='Saving tiff', disable=not bar): if np.any(self(*m)) or not at_least_one: @@ -1353,4 +1335,4 @@ def main(): im.save_as_tiff(out, args.channel, args.zslice, args.time, args.split) -from ndbioimage.readers import * +from .readers import * diff --git a/ndbioimage/jvm.py b/ndbioimage/jvm.py index 0b85bb3..54fbea6 100644 --- a/ndbioimage/jvm.py +++ b/ndbioimage/jvm.py @@ -1,4 +1,5 @@ from pathlib import Path +from urllib import request try: class JVM: @@ -15,10 +16,15 @@ try: cls._instance = object.__new__(cls) return cls._instance - def __init__(self, classpath=None): + def __init__(self, jars=None): if not self.vm_started and not self.vm_killed: - if classpath is None: - classpath = [str(Path(__file__).parent / 'jars' / '*')] + jarpath = Path(__file__).parent / 'jars' + if jars is None: + jars = {} + for jar, src in jars.items(): + if not (jarpath / jar).exists(): + JVM.download(src, jarpath / jar) + classpath = [str(jarpath / jar) for jar in jars.keys()] import jpype jpype.startJVM(classpath=classpath) @@ -39,6 +45,12 @@ try: if self.vm_killed: raise Exception('The JVM was killed before, and cannot be restarted in this Python process.') + @staticmethod + def download(src, dest): + print(f'Downloading {dest.name} to {dest}.') + dest.parent.mkdir(exist_ok=True) + dest.write_bytes(request.urlopen(src).read()) + def kill_vm(self): if self.vm_started and not self.vm_killed: import jpype diff --git a/ndbioimage/readers/bfread.py b/ndbioimage/readers/bfread.py index a6a325d..83c8054 100644 --- a/ndbioimage/readers/bfread.py +++ b/ndbioimage/readers/bfread.py @@ -1,22 +1,17 @@ import multiprocessing import numpy as np from abc import ABC -from ndbioimage import Imread, JVM from multiprocessing import queues from traceback import print_exc -from urllib import request -from pathlib import Path +from .. import Imread, JVM + + +jars = {'bioformats_package.jar': + 'https://downloads.openmicroscopy.org/bio-formats/latest/artifacts/bioformats_package.jar'} class JVMReader: def __init__(self, path, series): - bf_jar = Path(__file__).parent.parent / 'jars' / 'bioformats_package.jar' - if not bf_jar.exists(): - print('Downloading bioformats_package.jar.') - url = 'https://downloads.openmicroscopy.org/bio-formats/latest/artifacts/bioformats_package.jar' - bf_jar.parent.mkdir(exist_ok=True) - bf_jar.write_bytes(request.urlopen(url).read()) - mp = multiprocessing.get_context('spawn') self.path = path self.series = series @@ -50,7 +45,7 @@ class JVMReader: """ Read planes from the image reader file. adapted from python-bioformats/bioformats/formatreader.py """ - jvm = JVM() + jvm = JVM(jars) reader = jvm.image_reader() ome_meta = jvm.metadata_tools.createOMEXMLMetadata() reader.setMetadataStore(ome_meta) @@ -167,6 +162,18 @@ class JVMReader: jvm.kill_vm() +def can_open(path): + try: + jvm = JVM(jars) + reader = jvm.image_reader() + reader.getFormat(str(path)) + return True + except (Exception,): + return False + finally: + jvm.kill_vm() + + class Reader(Imread, ABC): """ This class is used as a last resort, when we don't have another way to open the file. We don't like it because it requires the java vm. @@ -176,7 +183,10 @@ class Reader(Imread, ABC): @staticmethod def _can_open(path): - return not path.is_dir() + """ Use java BioFormats to make an ome metadata structure. """ + with multiprocessing.get_context('spawn').Pool(1) as pool: + ome = pool.map(can_open, (path,))[0] + return ome def open(self): self.reader = JVMReader(self.path, self.series) diff --git a/ndbioimage/readers/cziread.py b/ndbioimage/readers/cziread.py index bee1bd6..25ed732 100644 --- a/ndbioimage/readers/cziread.py +++ b/ndbioimage/readers/cziread.py @@ -4,10 +4,10 @@ import re from lxml import etree from ome_types import model from abc import ABC -from ndbioimage import Imread from functools import cached_property from itertools import product from pathlib import Path +from .. import Imread class Reader(Imread, ABC): @@ -414,7 +414,7 @@ class Reader(Imread, ABC): return ome def __frame__(self, c=0, z=0, t=0): - f = np.zeros(self.shape['xy'], self.dtype) + f = np.zeros(self.file_shape[:2], self.dtype) directory_entries = self.filedict[c, z, t] x_min = min([f.start[f.axes.index('X')] for f in directory_entries]) y_min = min([f.start[f.axes.index('Y')] for f in directory_entries]) diff --git a/ndbioimage/readers/fijiread.py b/ndbioimage/readers/fijiread.py new file mode 100644 index 0000000..6669d4d --- /dev/null +++ b/ndbioimage/readers/fijiread.py @@ -0,0 +1,59 @@ +from abc import ABC +from tifffile import TiffFile +from functools import cached_property +from itertools import product +from ome_types import model +from pathlib import Path +from struct import unpack +from warnings import warn +import numpy as np +from .. import Imread + + +class Reader(Imread, ABC): + """ Can read some tif files written with Fiji which are broken because Fiji didn't finish writing. """ + priority = 90 + do_not_pickle = 'reader' + + @staticmethod + def _can_open(path): + if isinstance(path, Path) and path.suffix in ('.tif', '.tiff'): + with TiffFile(path) as tif: + return tif.is_imagej and not tif.is_bigtiff + else: + return False + + def __frame__(self, c, z, t): # Override this, return the frame at c, z, t + self.reader.filehandle.seek(self.offset + t * self.count) + return np.reshape(unpack(self.fmt, self.reader.filehandle.read(self.count)), self.shape['yx']) + + def open(self): + warn(f'File {self.path.name} is probably damaged, opening with fijiread.') + self.reader = TiffFile(self.path) + assert self.reader.pages[0].compression == 1, "Can only read uncompressed tiff files." + assert self.reader.pages[0].samplesperpixel == 1, "Can only read 1 sample per pixel." + self.offset = self.reader.pages[0].dataoffsets[0] + self.count = self.reader.pages[0].databytecounts[0] + self.bytes_per_sample = self.reader.pages[0].bitspersample // 8 + self.fmt = self.reader.byteorder + self.count // self.bytes_per_sample * 'BHILQ'[self.bytes_per_sample - 1] + + def close(self): + self.reader.close() + + @cached_property + def ome(self): + size_y, size_x = self.reader.pages[0].shape + size_c, size_z = 1, 1 + size_t = int(np.floor((self.reader.filehandle.size - self.reader.pages[0].dataoffsets[0]) / self.count)) + pixel_type = model.simple_types.PixelType(self.reader.pages[0].dtype.name) + ome = model.OME() + ome.instruments.append(model.Instrument()) + ome.images.append( + model.Image( + pixels=model.Pixels( + size_c=size_c, size_z=size_z, size_t=size_t, size_x=size_x, size_y=size_y, + dimension_order="XYCZT", type=pixel_type, + objective_settings=model.ObjectiveSettings(id="Objective:0")))) + for c, z, t in product(range(size_c), range(size_z), range(size_t)): + ome.images[0].pixels.planes.append(model.Plane(the_c=c, the_z=z, the_t=t, delta_t=0)) + return ome diff --git a/ndbioimage/readers/ndread.py b/ndbioimage/readers/ndread.py index f73db50..5477cb2 100644 --- a/ndbioimage/readers/ndread.py +++ b/ndbioimage/readers/ndread.py @@ -2,14 +2,12 @@ import numpy as np from ome_types import model from functools import cached_property from abc import ABC -from ndbioimage import Imread +from .. import Imread from itertools import product class Reader(Imread, ABC): priority = 20 - do_not_pickle = 'ome' - do_not_copy = 'ome' @staticmethod def _can_open(path): @@ -41,8 +39,11 @@ class Reader(Imread, ABC): return ome def open(self): - self.array = np.array(self.path) - self.path = 'numpy array' + if isinstance(self.path, np.ndarray): + self.array = np.array(self.path) + while self.array.ndim < 5: + self.array = np.expand_dims(self.array, -1) + self.path = 'numpy array' def __frame__(self, c, z, t): # xyczt = (slice(None), slice(None), c, z, t) diff --git a/ndbioimage/readers/seqread.py b/ndbioimage/readers/seqread.py index 380fabd..1d39d12 100644 --- a/ndbioimage/readers/seqread.py +++ b/ndbioimage/readers/seqread.py @@ -1,7 +1,6 @@ import tifffile import yaml import re -from ndbioimage import Imread from pathlib import Path from functools import cached_property from ome_types import model @@ -9,6 +8,7 @@ from ome_types._base_type import quantity_property from itertools import product from datetime import datetime from abc import ABC +from .. import Imread def lazy_property(function, field, *arg_fields): diff --git a/ndbioimage/readers/tifread.py b/ndbioimage/readers/tifread.py index f03a3c4..d07de9a 100644 --- a/ndbioimage/readers/tifread.py +++ b/ndbioimage/readers/tifread.py @@ -1,13 +1,12 @@ -from abc import ABC - -from ndbioimage import Imread import numpy as np import tifffile import yaml +from abc import ABC from functools import cached_property from ome_types import model from pathlib import Path from itertools import product +from .. import Imread class Reader(Imread, ABC): @@ -18,7 +17,7 @@ class Reader(Imread, ABC): def _can_open(path): if isinstance(path, Path) and path.suffix in ('.tif', '.tiff'): with tifffile.TiffFile(path) as tif: - return tif.is_imagej + return tif.is_imagej and tif.pages[-1]._nextifd() == 0 else: return False diff --git a/pyproject.toml b/pyproject.toml index adf0a4b..8dc6bb5 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -1,6 +1,6 @@ [tool.poetry] name = "ndbioimage" -version = "2023.7.3" +version = "2023.7.4" description = "Bio image reading, metadata and some affine registration." authors = ["W. Pomp "] license = "GPLv3" diff --git a/tests/test_open.py b/tests/test_open.py index cd3883b..a8691de 100644 --- a/tests/test_open.py +++ b/tests/test_open.py @@ -1,10 +1,12 @@ import pytest from pathlib import Path -from ndbioimage import Imread +from ndbioimage import Imread, ReaderNotFoundError -@pytest.mark.parametrize("file", - [file for file in (Path(__file__).parent / 'files').iterdir() if file.suffix != '.pzl']) +@pytest.mark.parametrize("file", (Path(__file__).parent / 'files').iterdir()) def test_open(file): - with Imread(file) as im: - print(im[dict(c=0, z=0, t=0)].mean()) + try: + with Imread(file) as im: + print(im[dict(c=0, z=0, t=0)].mean()) + except ReaderNotFoundError: + assert len(Imread.__subclasses__()), "No subclasses for Imread found."