#!/usr/bin/env python
# -*- coding: utf-8 -*-
# ------------------------------------------------------------------------------
# Copyright (c) 2022 Univ. Lyon, ENS de Lyon, UCB Lyon 1, CNRS, INRAe, Inria
# All rights reserved.
# This file is part of the TimageTK library, and is released under the "GPLv3"
# license. Please see the LICENSE.md file that should have been included as
# part of this package.
# ------------------------------------------------------------------------------
"""Image IO (input/output) module."""
import tempfile
from os.path import exists
from os.path import join
from os.path import sep
from pathlib import Path
import numpy as np
import scipy.ndimage as nd
from tqdm.autonotebook import tqdm
from timagetk import TissueImage2D
from timagetk import TissueImage3D
from timagetk.bin.logger import get_logger
from timagetk.components.labelled_image import LabelledImage
from timagetk.components.metadata import _check_class_def
from timagetk.components.metadata import _check_fname_def
from timagetk.components.metadata import _check_physical_ppty
from timagetk.components.multi_channel import MultiChannelImage
from timagetk.components.spatial_image import SpatialImage
from timagetk.io.util import LSM_EXT
from timagetk.io.util import TIF_EXT
from timagetk.io.util import get_pne
from timagetk.third_party.vt_converter import vt_to_spatial_image
from timagetk.third_party.vt_image import vtImage
from timagetk.util import check_type
from timagetk.util import clean_type
from timagetk.util import now
log = get_logger(__name__)
def _loaded_url(url, tmp_fname=None):
"""Check if the given URL has already been loaded.
Parameters
----------
url : str
The URL to check.
tmp_fname : str, optional
Temporary file name created by `_url_image`, may be replaced by previously loaded file name if any.
Returns
-------
str
File path to temporary file to load.
Notes
-----
Uses a temporary file named `timagetk_image_url.json` to store that info.
Examples
--------
>>> from timagetk.io.image import _loaded_url
>>> url = 'https://zenodo.org/record/3737630/files/sphere_membrane_t0.inr.gz'
>>> tmp_fname = _loaded_url(url)
>>> print(tmp_fname)
/tmp/tmprwxl8e4o.inr.gz
>>> url = "https://zenodo.org/record/3737795/files/qDII-CLV3-PIN1-PI-E35-LD-SAM4.czi"
>>> _loaded_url(url)
>>> print(tmp_fname)
/tmp/tmpgbzdj9yd.czi
"""
import json
if tmp_fname is None:
try:
tmp_fname = Path(url).name
except:
# Create a temporary filename:
tmp_fname = tempfile.NamedTemporaryFile().name
tmpdir = tempfile.gettempdir()
url_json = join(tmpdir, 'timagetk_image_url.json')
try:
# Check if the JSON url-register file exists:
assert exists(url_json)
except AssertionError:
# Create the JSON url-register file and add the {url: tmp_file} dict
with open(url_json, 'w') as json_file:
json.dump({url: tmp_fname}, json_file)
else:
# If exists, load the JSON url-register to a dictionary:
with open(url_json, 'r') as json_file:
url_file = json.load(json_file)
# Replace given temporary file name to the one in JSON if it still exists:
tmp = url_file.get(url, tmp_fname)
if exists(tmp):
tmp_fname = tmp
# Update the url-register and dump it to a JSON file:
url_file.update({url: tmp_fname})
with open(url_json, 'w') as json_file:
json.dump(url_file, json_file)
return join(tmpdir, tmp_fname)
def _image_from_url(url, hash_value=None, hash_method='md5'):
"""Get an image using URL address.
Parameters
----------
url : str
A valid URL pointing toward an image.
hash_value : str, optional
The hash value to use for comparison.
hash_method: {'SHA1', 'SHA224', 'SHA256', 'SHA384', 'SHA512', 'MD5'}, optional
The hash method to use, default is 'md5', case insentitive.
Returns
-------
str
Path to the temporary file containing the retreived image.
Raises
------
ValueError
If the given `hash` is wrong.
See Also
--------
_loaded_url
util.POSS_EXT
Notes
-----
This creates a temporary file and requires a valid image extension.
There is a list of previously loaded temporary files to avoid downloading the same file multiple times.
Examples
--------
>>> from timagetk.io.image import _image_from_url
>>> url = 'https://zenodo.org/record/3737630/files/sphere_membrane_t0.inr.gz'
>>> tmp_file = _image_from_url(url)
>>> print(tmp_file)
/tmp/tmpwv0c8lum.inr.gz
>>> # Read the temporary file image with timagetk ``imread``:
>>> from timagetk.io.image import imread
>>> img = imread(tmp_file)
>>> # Open the stack browser:
>>> from timagetk.visu.stack import stack_browser
>>> stack_browser(img)
>>> # Use secure hash to validate downloaded file
>>> md5_h = '24cc1edfaca092c26486a58fd6979827'
>>> img = _image_from_url(url,hash_value=md5_h,hash_method='md5')
>>> # Use WRONG secure hash to INvalidate downloaded file:
>>> wrong_md5_h = 'd3079693f5d25f0743c19b004e88eaf1'
>>> img = _image_from_url(url,hash_value=md5_h,hash_method='md5')
ValueError: MD5 comparison failed for file /tmp/tmpet2v_mio.inr.gz!
Given hash: d3079693f5d25f0743c19b004e88eaf1
Computed hash: 24cc1edfaca092c26486a58fd6979827
"""
import hashlib
# Uses the JSON file to check if we already downloaded this URL:
tmp_fname = _loaded_url(url)
# Only save the URL to a temporary file if not already defined:
if not exists(tmp_fname):
_save_url_temp_file(url, tmp_fname)
# Performs hash comparison if an `hash_value` is given:
if hash_value is not None:
# Instantiate hash method:
h = hashlib.new(hash_method.lower())
# Compute hash based on file content:
h.update(open(tmp_fname, 'rb').read())
# Compare hash or raise error:
try:
assert hash_value == h.hexdigest()
except AssertionError:
msg = f"{hash_method.upper()} comparison failed for file {tmp_fname}!\n"
msg += f"Given hash: {hash_value}\n"
msg += f"Computed hash: {h.hexdigest()}\n"
raise ValueError(msg)
return tmp_fname
def _save_url_temp_file(url, tmp_fname=None):
"""Save URL to a temporary file.
Parameters
----------
url : str
A valid URL pointing toward a file.
tmp_fname : str, optional
Path to the temporary file, else created.
Returns
-------
str
Path to the temporary file containing the downloaded image.
Examples
--------
>>> from os.path import exists
>>> from timagetk.io.image import _save_url_temp_file
>>> # Example #1 -
>>> url = "https://zenodo.org/record/3737795/files/qDII-CLV3-PIN1-PI-E35-LD-SAM4.czi"
>>> tmp_fname = _save_url_temp_file(url)
>>> tmp_fname
'/tmp/tmp81cyqecs'
>>> url = 'https://zenodo.org/record/3737630/files/sphere_membrane_t0.inr.gz'
>>> tmp_fname = _save_url_temp_file(url, '/tmp/sphere_membrane_t0.inr.gz')
>>> print(tmp_fname)
/tmp/sphere_membrane_t0.inr.gz
>>> exists(tmp_fname)
True
"""
import requests
if tmp_fname is None:
# Fixme: upon reflexion this try/except seems useless...
# TODO: use the response of the GET requests to make sure we have an valid URL pointing to a (image) file (not HTML)?
try:
tmp_fname = Path(url).name
except:
# Create a temporary filename:
tmp_fname = tempfile.NamedTemporaryFile().name
r = requests.get(url, stream=True)
total_size = int(r.headers.get("content-length", 0)) # total size in bytes
block_size = 32 * 1024 # block-size reads
progress = 0 # progress tracker
pbar = tqdm(total=total_size, unit="B", unit_scale=True, unit_divisor=1024)
with open(tmp_fname, "wb") as f:
for chunk in r.iter_content(block_size):
f.write(chunk)
progress = progress + len(chunk)
pbar.update(block_size)
pbar.close()
if total_size != 0 and progress != total_size:
raise IOError(f"Error downloading file {tmp_fname}!")
# Close HTTP(S) connection:
r.close()
return tmp_fname
[docs]
def imread(fname, rtype=None, **kwargs):
"""Read an image (2D/3D +C).
Parameters
----------
fname : str or pathlib.Path
File path to the image, can be a URL.
rtype : type, optional
If defined, force retured type of image to given class.
Other Parameters
----------------
channel_names : list of str
List of channel names, to use with ``MultiChannelImage`` instances to overide those defined in the CZI file.
channel : str
Channel to return, *i.e* only this one will be returned, to use with ``MultiChannelImage`` instances.
hash : str
The hash value to use for comparison.
hash_method: {'SHA1', 'SHA224', 'SHA256', 'SHA384', 'SHA512', 'MD5'}
The hash method to use, default is 'md5', case insentitive.
Returns
-------
timagetk.components.spatial_image.SpatialImage or timagetk.LabelledImage or timagetk.MultiChannelImage
Image and metadata (such as voxel-size, extent, type, etc.)
See Also
--------
timagetk.io.util.POSS_EXT
Example
-------
>>> from timagetk.io.image import imread
>>> # EXAMPLE #1 - Load a file from URL:
>>> img = imread('https://zenodo.org/record/3737630/files/sphere_membrane_t0.inr.gz')
>>> type(img)
timagetk.components.spatial_image.SpatialImage
>>> print((img.filepath, img.filename))
('/tmp', 'sphere_membrane_t0.inr.gz')
>>> # EXAMPLE #2 - Load a file from URL and validate with secure hash:
>>> md5_h = '24cc1edfaca092c26486a58fd6979827'
>>> img = imread('https://zenodo.org/record/3737630/files/sphere_membrane_t0.inr.gz', hash=md5_h, hash_method='md5')
>>> type(img)
timagetk.components.spatial_image.SpatialImage
>>> print((img.filepath, img.filename))
('/tmp', 'sphere_membrane_t0.inr.gz)
>>> # EXAMPLE #3 - Load a file from URL and force returned instance type:
>>> from timagetk import TissueImage3D
>>> url = 'https://zenodo.org/record/3737630/files/sphere_membrane_t0_seg.inr.gz'
>>> img = imread(url, TissueImage3D, background=1, not_a_label=0)
>>> type(img)
timagetk.components.tissue_image.TissueImage3D
>>> # EXAMPLE #4 - Load a file from the shared data:
>>> from timagetk.io.util import shared_dataset
>>> image_path = shared_dataset('sphere')[0]
>>> img = imread(image_path)
>>> type(img)
timagetk.components.spatial_image.SpatialImage
"""
if isinstance(fname, str):
# URL checkpoint:
if fname.startswith('http'):
return imread(_image_from_url(fname), rtype, **kwargs)
else:
fname = Path(fname)
fname = fname.expanduser().absolute()
# Check file exists:
try:
assert fname.exists()
except AssertionError:
log.error(f"Could not find file `{fname}`!")
return None
# Get file path, name, name and extension:
filepath, filename, ext = get_pne(fname, test_format=True)
log.debug("Loading image: {}".format(fname))
if ext == ".czi":
img = read_czi_image(fname, kwargs.get("channel_names", None))
# - Select only one channel:
if 'channel' in kwargs:
img = img.get_channel(kwargs.get('channel'))
elif ext == ".lsm":
img = read_lsm_image(fname, kwargs.get("channel_names", None))
# - Select only one channel:
if 'channel' in kwargs:
img = img.get_channel(kwargs.get('channel'))
elif ext in TIF_EXT:
try:
img = read_tiff_image(fname, kwargs.get("channel_names", None))
except Exception as e:
log.warning(f"TiffFile raised an error! '{e}'")
log.warning(f"Could not load file {fname} using TiffFile, trying with vt...")
img = vt_to_spatial_image(vtImage(str(fname)))
# - Select only one channel:
if isinstance(img, MultiChannelImage) and 'channel' in kwargs:
img = img.get_channel(kwargs.get('channel'))
else:
img = vt_to_spatial_image(vtImage(str(fname)))
img = _check_class_def(img)
img = _check_fname_def(img, filename=filename, filepath=filepath)
if isinstance(img, MultiChannelImage):
channels = img.get_channel_names()
ch = channels[0]
_check_physical_ppty(img.metadata_image[ch].get_dict())
else:
_check_physical_ppty(img.metadata_image.get_dict())
# - Type conversion according to detected class:
md_class = img.metadata.get('timagetk', None)
if md_class == 'LabelledImage':
img = LabelledImage(img, **kwargs)
elif md_class == 'TissueImage2D':
img = TissueImage2D(img, **kwargs)
elif md_class == 'TissueImage3D':
img = TissueImage3D(img, **kwargs)
else:
img.metadata.update({'timagetk': {'class': 'SpatialImage'}})
if rtype is not None:
img = rtype(img, **kwargs)
return img
[docs]
def imsave(fname, sp_img, **kwargs):
"""Save an image (2D/3D).
Parameters
----------
fname : str or pathlib.Path
File path where to save the image.
sp_img : timagetk.components.spatial_image.SpatialImage or timagetk.LabelledImage or timagetk.TissueImage3D or timagetk.MultiChannelImage
Image instance to save on drive.
Other Parameters
----------------
pattern : str
Physical dimension order to use to save the image, "ZCYX" by default.
force : bool
Overwrite any existing file of same `fname`.
See Also
--------
timagetk.io.util.POSS_EXT
Example
-------
>>> from timagetk.array_util import random_spatial_image
>>> import tempfile
>>> import numpy as np
>>> from timagetk.io.image import imsave, imread
>>> from timagetk import SpatialImage
>>> sp_image = random_spatial_image((3, 4, 5), voxelsize=[0.3, 0.4 ,0.5])
>>> # Save the image in a temporary file:
>>> tmp_path = tempfile.NamedTemporaryFile()
>>> imsave(tmp_path.name+'.tif', sp_image)
>>> # Load this file:
>>> sp_image_cp = imread(tmp_path.name+'.tif')
>>> # Compare arrays and voxelizes:
>>> np.testing.assert_array_equal(sp_image, sp_image_cp)
>>> np.testing.assert_array_equal(sp_image.voxelsize, sp_image_cp.voxelsize)
>>> sp_image = random_spatial_image((3, 4), voxelsize=[0.3, 0.4])
>>> # Save the image in a temporary file:
>>> tmp_path = tempfile.NamedTemporaryFile()
>>> imsave(tmp_path.name+'.tif', sp_image)
>>> # Load this file:
>>> sp_image_cp = imread(tmp_path.name+'.tif')
>>> # Compare arrays and voxelizes:
>>> np.testing.assert_array_equal(sp_image, sp_image_cp)
>>> np.testing.assert_array_equal(sp_image.voxelsize, sp_image_cp.voxelsize)
>>> from timagetk import MultiChannelImage
>>> from timagetk.io import imsave
>>> from timagetk.synthetic_data.nuclei_image import example_nuclei_signal_images
>>> n_img, s_img, _, _= example_nuclei_signal_images(n_points=10, extent=20, voxelsize=(1,1,1))
>>> img = MultiChannelImage({'nuclei':n_img, 'signal':s_img})
>>> imsave('~/test.tif', img)
>>> tmp_path = tempfile.NamedTemporaryFile()
>>> imsave(tmp_path.name+'.tif', img)
"""
if isinstance(fname, str):
fname = Path(fname)
# - Assert sp_img is a SpatialImage or a MultiChannelImage instance:
try:
assert isinstance(sp_img, SpatialImage) or isinstance(sp_img, MultiChannelImage)
except AssertionError:
raise TypeError("Parameter 'sp_img' is not a SpatialImage or a MultiChannelImage!")
filepath, filename, ext = get_pne(fname, test_format=True)
# Check file do NOT exist:
try:
assert not fname.exists()
except AssertionError:
if not kwargs.get('force', False):
from timagetk.util import yn_query
log.warning(f"Found existing file `{fname}`!")
if not yn_query(f"Overwrite existing image `{fname}`", default='no'):
log.info("Will NOT overwrite existing image!")
return
else:
kwargs['force'] = True
# - METADATA: Update file path & name before saving:
sp_img = _check_fname_def(sp_img, filename=filename, filepath=filepath)
ctype = clean_type(sp_img)
log.debug("Saving {} under: {}".format(ctype, fname))
if ext in TIF_EXT:
save_tiff_image(fname, sp_img, **kwargs)
else:
# - Assert sp_img is a SpatialImage instance:
try:
assert isinstance(sp_img, SpatialImage)
except AssertionError:
raise TypeError("Parameter 'sp_img' is not a SpatialImage!")
vt_im = sp_img.to_vtimage()
vt_im.write(str(fname))
# FIXME: This does not allow to use metadata!
return
[docs]
def apply_mask(img, mask_filename, masking_value=0, crop_mask=False):
"""Load and apply a z-projection mask (2D) to a SpatialImage (2D/3D).
In case of an intensity image, allows removing unwanted signal intensities.
If the SpatialImage is 3D, it is applied in z-direction.
The mask should contain a distinct "masking value".
Parameters
----------
img : timagetk.components.spatial_image.SpatialImage
Image to modify with the mask.
mask_filename : str
String giving the location of the mask file
masking_value : int, optional
Value (default 0) defining the masked region
crop_mask : bool, optional
If ``True`` (default ``False``), the returned SpatialImage is cropped around the non-masked region
Returns
-------
timagetk.components.spatial_image.SpatialImage
The masked SpatialImage.
"""
try:
from pillow import Image
except ImportError:
from PIL import Image
xsh, ysh, zsh = img.shape
# Read the mask file:
mask_im = Image.open(mask_filename)
# Detect non-masked values positions
mask_im = np.array(mask_im) != masking_value
# Transform mask to 3D by repeating the 2D-mask along the z-axis:
mask_im = np.repeat(mask_im[:, :, np.newaxis], zsh, axis=2)
# Remove masked values from ``img``:
masked_im = img * mask_im
# Crop the mask if required:
if crop_mask:
bbox = nd.find_objects(mask_im)
masked_im = masked_im[bbox]
return masked_im
[docs]
def default_image_attributes(array, log_undef=log.debug, **kwargs):
"""Set the default image attribute values given the array dimensionality.
Parameters
----------
array : numpy.ndarray
The array to use as image.
log_undef : fun
The logging level to use for undefined image attributes.
Other Parameters
----------------
dtype : str or numpy.dtype
A valid ``dtype`` to use with this array, checked against ``AVAIL_TYPES``.
Defaults to ``array.dtype``.
axes_order : str
The physical axes ordering to use with this array.
Defaults to ``default_axes_order(array)``.
origin : list of int
The coordinate of the origin for this array.
Defaults to ``default_origin(array)``.
unit : float
The units associated to the voxel-sizes to use with this array.
Defaults to ``DEFAULT_SIZE_UNIT``.
voxelsize : list of float
The voxel-sizes to use with this array.
Defaults to ``default_voxelsize(array)``.
metadata : dict
The dictionary of metadata to use with this array.
Defaults to ``{}``.
Returns
-------
dict
The image attributes dictionary.
See Also
--------
timagetk.components.spatial_image.AVAIL_TYPES
timagetk.components.spatial_image.DEFAULT_SIZE_UNIT
timagetk.components.spatial_image.default_axes_order
timagetk.components.spatial_image.default_origin
timagetk.components.spatial_image.default_voxelsize
Notes
-----
Any extra keyword argument not listed above will be kept and returned.
Examples
--------
>>> from timagetk.io.image import default_image_attributes
>>> from timagetk.array_util import random_array
>>> # Example 1: set the default image attributes for a 2D array:
>>> img_2D = random_array([7, 5])
>>> attr_2D = default_image_attributes(img_2D)
>>> print(attr_2D)
{'dtype': 'uint8', 'axes_order': 'YX', 'origin': [0, 0], 'unit': 1e-06, 'voxelsize': [1.0, 1.0], 'metadata': {}}
>>> # Example 2: set the default image attributes for a 3D array:
>>> img_3D = random_array([7, 5, 4])
>>> attr_3D = default_image_attributes(img_3D)
>>> print(attr_3D)
{'dtype': 'uint8', 'axes_order': 'ZYX', 'origin': [0, 0, 0], 'unit': 1e-06, 'voxelsize': [1.0, 1.0, 1.0], 'metadata': {}}
>>> # Example 3: set the default image attributes for a 3D array, except voxelsize:
>>> attr_3D = default_image_attributes(img_3D, voxelsize=[0.5, 0.25, 0.25])
>>> print(attr_3D)
{'voxelsize': [0.5, 0.25, 0.25], 'dtype': 'uint8', 'axes_order': 'ZYX', 'origin': [0, 0, 0], 'unit': 1e-06, 'metadata': {}}
"""
from timagetk.components.spatial_image import AVAIL_TYPES
from timagetk.components.spatial_image import default_axes_order
from timagetk.components.spatial_image import default_origin
from timagetk.components.spatial_image import DEFAULT_SIZE_UNIT
from timagetk.components.spatial_image import default_voxelsize
# Get the number of dimension in the array:
ndim = array.ndim
# ----------------------------------------------------------------------
# DTYPE:
# ----------------------------------------------------------------------
dtype = kwargs.get('dtype', None)
if dtype is None:
# - Set it to the array type:
dtype = np.dtype(array.dtype).name
else:
# Try to convert it to a str
try:
dtype = np.dtype(dtype).name
except:
pass
# Check it is a known 'dtype':
try:
assert dtype in AVAIL_TYPES
except AssertionError:
log.error(f"Given ``dtype`` ({dtype}) is not valid!")
log.info(f"Available ``dtype`` are: '{AVAIL_TYPES}'.")
dtype = array.dtype
log.info(f"Falling back to ``array`` dtype: '{dtype}'.")
log.debug(f"Defined parameter ``dtype`` as: {dtype}")
# ----------------------------------------------------------------------
# Physical axes ordering:
# ----------------------------------------------------------------------
axes_order = kwargs.get('axes_order', None)
if axes_order is None:
axes_order = default_axes_order(array)
log_undef(f"Undefined parameter ``axes_order``, using default: {axes_order}")
else:
axes_order = ''.join([ax.upper() for ax in axes_order])
extra_dim = ""
if "B" in axes_order:
# ``BlendImage`` case:
ndim -= 1
extra_dim = "B"
if "C" in axes_order:
# ``MultiChannelImage`` case:
ndim -= 1
extra_dim = "C"
test_axes_order = axes_order.replace(extra_dim, '')
try:
assert len(test_axes_order) == ndim
except AssertionError:
log.error(f"Given ``axes_order`` is of size {len(test_axes_order)} for an array of dimension {ndim}!")
axes_order = default_axes_order(array)
axes_order += extra_dim
log.info(f"Falling back to default value: '{axes_order}'.")
log.debug(f"Defined parameter ``axes_order`` as: {axes_order}")
# ----------------------------------------------------------------------
# Origin:
# ----------------------------------------------------------------------
origin = kwargs.get('origin', None)
if origin is None:
origin = default_origin(array)
log_undef(f"Undefined parameter ``origin``, using default: {origin}")
else:
origin = list(map(int, origin))
try:
assert len(origin) == ndim
except AssertionError:
log.error(f"Given ``origin`` is of size {len(origin)} for an array of dimension {ndim}!")
origin = default_origin(array)
log.info(f"Falling back to default value: '{origin}'.")
log.debug(f"Defined parameter ``origin`` as: {origin}")
# ----------------------------------------------------------------------
# Voxel-size:
# ----------------------------------------------------------------------
voxelsize = kwargs.get('voxelsize', None)
if voxelsize is None:
voxelsize = default_voxelsize(array)
log_undef(f"Undefined parameter ``voxelsize``, using default: {voxelsize}")
else:
voxelsize = list(map(float, voxelsize))
try:
assert len(voxelsize) == ndim
except AssertionError:
log.error(f"Given ``voxelsize`` is of size {len(voxelsize)} for an array of dimension {ndim}!")
voxelsize = default_voxelsize(array)
log.info(f"Falling back to default value: '{voxelsize}'.")
log.debug(f"Defined parameter ``voxelsize`` as: {voxelsize}")
# ----------------------------------------------------------------------
# Size unit:
# ----------------------------------------------------------------------
unit = kwargs.get('unit', None)
if unit is None:
unit = DEFAULT_SIZE_UNIT
log_undef(f"Undefined parameter ``unit``, using default: {unit}")
else:
log.debug(f"Defined parameter ``unit`` as: {unit}")
# ----------------------------------------------------------------------
# Metadata:
# ----------------------------------------------------------------------
metadata = kwargs.get('metadata', None)
# - Initialize (empty) or check type of given metadata dictionary:
if metadata == None:
metadata = {}
log.debug("No metadata dictionary provided!")
else:
check_type(metadata, 'metadata', dict)
kwargs['dtype'] = dtype
kwargs['axes_order'] = axes_order
kwargs['origin'] = origin
kwargs['unit'] = unit
kwargs['voxelsize'] = voxelsize
kwargs['metadata'] = metadata
return kwargs
[docs]
def read_czi_image(czi_file, channel_names=None, **kwargs):
"""Read Carl Zeiss Image (CZI) file.
Parameters
----------
czi_file : str
File path to the image to load.
channel_names : list of str, optional
Defines channel names to use.
By default, try to load channel names from the image metadata.
Other Parameters
----------------
pattern : str, optional
Physical dimension order to use to load the image, ".C.ZXY." by default.
By default, try to load the pattern from the image metadata.
return_order : str
Physical dimension order to use for the returned array.
Defaults to "ZYXC".
Returns
-------
timagetk.SpatialImage or timagetk.MultiChannelImage
If a single channel and time-point is found, returns a ``SpatialImage`` instance.
If multiple channels are found but with a single time-point, returns a ``MultiChannelImage`` instance.
Notes
-----
If more than one channel is found, return a ``MultiChannelImage``.
For the `pattern` variable, use:
- ``C`` as channel;
- ``T`` as time;
- ``X`` as x-axis (rows);
- ``Y`` as y-axis (planes);
- ``Z`` as z-axis (columns);
Examples
--------
>>> from timagetk.io.image import read_czi_image
>>> from timagetk.io.image import _image_from_url
>>> im_url = "https://zenodo.org/record/3737795/files/qDII-CLV3-PIN1-PI-E35-LD-SAM4.czi"
>>> im_fname = _image_from_url(im_url) # download the CZI image
>>> ch_names = ["DII-VENUS-N7", "pPIN1:PIN1-GFP", "Propidium Iodide", "pRPS5a:TagBFP-SV40", "pCLV3:mCherry-N7"]
>>> image = read_czi_image(im_fname, ch_names)
>>> print(image.get_channel_names())
['DII-VENUS-N7', 'pPIN1:PIN1-GFP', 'Propidium Iodide', 'pRPS5a:TagBFP-SV40', 'pCLV3:mCherry-N7']
>>> print(image)
>>> im_url = "https://zenodo.org/record/3737795/files/qDII-CLV3-PIN1-PI-E35-LD-SAM4.czi"
>>> im_fname = _image_from_url(im_url) # download the CZI image
>>> from czifile import CziFile
>>> czi_img = CziFile(im_fname)
>>> import xmltodict
>>> md = xmltodict.parse(czi_img.metadata())['ImageDocument']['Metadata']
"""
import xmltodict as xmltodict
from czifile import CziFile
return_order = kwargs.get('return_order', 'ZYXC')
pattern = kwargs.get('pattern', '.C.ZXY.')
vxs_order = list(return_order.replace('C', ''))
# - Load the file with third-party reader:
log.debug(f"Using `CziFile` to load: {czi_file}")
czi_img = CziFile(czi_file)
log.debug("Done!")
# -- METADATA:
md = xmltodict.parse(czi_img.metadata())['ImageDocument']['Metadata']
# Note: no other dict entries except 'ImageDocument'/'Metadata' so we have all the metadata!
# list(md.keys()) => 'Experiment', 'HardwareSetting', 'ImageScaling', 'Information', 'Scaling', 'DisplaySetting', 'Layers'
# - Get the name of the Microscope:
if 'Microscopes' in md['Information']['Instrument']:
micro_md = md['Information']['Instrument']['Microscopes']['Microscope']
if '@Name' in micro_md and 'Type' in micro_md:
m_name = micro_md['@Name']
m_type = micro_md['Type']
m_system = micro_md['System']
microscope = f"{m_name} - {m_type} {m_system}"
else:
microscope = micro_md['System']
else:
microscope = "Unknown microscope"
log.debug(f"Found metadata from {microscope}")
# - Get the name of the image:
if 'Name' in md['Information']['Document']:
name = md['Information']['Document']['Name']
log.debug(f"Found metadata 'Name': {name}")
else:
name = None
# - Get the acquisition date:
if 'CreationDate' in md['Information']['Document']:
creation_date = md['Information']['Document']['CreationDate']
log.debug(f"Found metadata 'CreationDate': {creation_date}")
else:
creation_date = "Unknown"
# - Get the voxel-size dictionary!
md_spacing = md['Scaling']['Items']['Distance']
md_vxs = {md_spacing[ax]['@Id']: float(md_spacing[ax]['Value']) for ax in range(3)}
log.debug(f"Found metadata 'voxelsize': {md_vxs}")
# - Get the channel & dye names:
md_channels = md['DisplaySetting']['Channels']['Channel']
if isinstance(md_channels, list):
# Multi channel case
n_md_ch = len(md_channels)
md_ch_names = [md_channels[ch]['@Name'] for ch in range(n_md_ch)]
if 'DyeName' in md_channels[0].keys():
md_dye_names = [md_channels[ch]['DyeName'] for ch in range(n_md_ch)]
else:
# Single channel case
n_md_ch = 1
md_ch_names = md_channels['@Name']
if 'DyeName' in md_channels.keys():
md_dye_names = md_channels['DyeName']
log.debug(f"Found {n_md_ch} channels in the metadata...")
# TODO: Some channels may have color info under ``channels[ch]['Color']``
# TODO: Other channels should have a "palette" info under ``channels[ch]['PaletteName']``
# true_ch = md['Experiment']['ExperimentBlocks']['AcquisitionBlock']['MultiTrackSetup']['TrackSetup']
true_ch = md['Information']['Image']['Dimensions']['Channels']['Channel'] # md['Version'] == '1.0'
if isinstance(true_ch, list):
n_true_ch = len(true_ch)
true_ch_names = [true_ch[ch]['@Name'] for ch in range(n_true_ch)]
else:
n_true_ch = 1
true_ch_names = true_ch['@Name']
# true_ich_names = [true_ch[ch]['Detectors']['Detector']['ImageChannelName'] for ch in range(n_true_ch)]
# true_dye_names = [true_ch[ch]['Detectors']['Detector']['Dye'] for ch in range(n_true_ch)]
if len(true_ch_names) == n_md_ch:
log.info(f"Found {n_true_ch} channel names in metadata: {' ,'.join(true_ch_names)}")
if channel_names is None:
channel_names = true_ch_names
# - Indicate if is a time-series: 'true'/'false'
if 'AcquisitionModeSetup' in md['Experiment']['ExperimentBlocks']['AcquisitionBlock']:
time_series = md['Experiment']['ExperimentBlocks']['AcquisitionBlock']['AcquisitionModeSetup']['TimeSeries']
md_pattern = None
for s in czi_img.segments():
if s.SID == "ZISRAWSUBBLOCK":
metadata = str(s).split('\n')
for i, row in enumerate(metadata):
if 'DirectoryEntry' in row:
next_md = metadata[i + 1]
md_pattern = next_md.split(' ')[4]
if md_pattern[:2] == 'b\'': # decoding issue?
md_pattern = md_pattern[2:-1]
# -- voxel-size:
voxelsize = tuple([md_vxs[dim] for dim in vxs_order])
unit = 1
if all(vxs < 1e-5 for vxs in voxelsize):
voxelsize = [vxs * 1e6 for vxs in voxelsize]
unit = 1e-6
# -- ARRAYS:
czi_channels = czi_img.asarray()
czi_ndim = czi_channels.ndim
log.debug(f"Loaded {czi_ndim}D CZI file of shape: {czi_channels.shape}")
# - Replace given pattern by metadata pattern!
ok_pattern_nb = len(md_pattern) == czi_ndim
ok_pattern_dims = np.all([md_pattern.find(d) != -1 for d in return_order])
if md_pattern is not None and ok_pattern_nb and ok_pattern_dims:
log.debug(f"Found dimension pattern in metadata: {md_pattern}...")
pattern = md_pattern
log.warning(f"Found dimension pattern in metadata: {md_pattern}, using it in place of given pattern!")
# - Re-order the array according to ``return_order``
log.debug(f"Array shape BEFORE transpose & reshape: {czi_channels.shape}")
# Axes id order for axes to returns ('X', 'Y', 'Z' & 'C'):
axis = [pattern.find(c) for c in return_order]
# Axes id order for axes NOT to returns (others like 'T'):
extra_axis = [i for i in range(czi_ndim) if pattern[i] not in return_order]
# Transpose axes to get array with returned axes as the 4 first
czi_channels = np.transpose(czi_channels, tuple(axis + extra_axis))
# Clip the array to the 4 first axes:
czi_channels = czi_channels.reshape(czi_channels.shape[:4])
log.debug(f"Array shape AFTER transpose & reshape: {czi_channels.shape}")
ch_axis_id = return_order.find('C')
n_channels = czi_channels.shape[ch_axis_id]
log.debug(f"Found {n_channels} channels in the array...")
# Defines file name & path:
fname = Path(czi_file).name
fpath = Path(czi_file).parent
log.info(f"Loaded file '{fname}' with {n_channels} channels!")
if n_channels > 1:
if channel_names is None:
# channel_names = [f"{n}_{f}" for n, f in zip(md_ch_names, md_dye_names)]
channel_names = [f"CH_{n}" for n in md_ch_names]
channel_images = []
for i_ch, _ in enumerate(channel_names):
array = np.take(czi_channels, i_ch, axis=ch_axis_id)
im_kwargs = default_image_attributes(array, voxelsize=voxelsize, unit=unit)
channel_images.append(SpatialImage(array, **im_kwargs))
metadata = {'acquisition_date': creation_date, 'scan_info': md}
img = MultiChannelImage(channel_images, channel_names=channel_names, filename=str(fname), filepath=str(fpath),
name=name if name is not None else fname, metadata=metadata)
log.info(f"Channel names: {', '.join(img.get_channel_names())}")
else:
array = czi_channels[:, :, :, 0]
metadata = {'acquisition_date': creation_date, 'scan_info': md}
im_kwargs = default_image_attributes(array, voxelsize=voxelsize, unit=unit, metadata=metadata)
img = SpatialImage(array, **im_kwargs)
return img
[docs]
def read_lsm_image(lsm_file, channel_names=None, **kwargs):
"""Read LSM file.
Parameters
----------
lsm_file : str
File path to the image to load.
channel_names : list of str, optional
Defines channel names to use.
By default, try to load channel names from the image metadata.
Other Parameters
----------------
pattern : str, optional
Physical dimension order to use to load the image, "TZCYX" by default.
By default, try to load the pattern from the image metadata.
return_order : str
Physical dimension order to use for the returned array.
Defaults to "ZYXC".
Returns
-------
timagetk.SpatialImage or timagetk.MultiChannelImage
If a single channel and time-point is found, returns a ``SpatialImage`` instance.
If multiple channels are found but with a single time-point, returns a ``MultiChannelImage`` instance.
Notes
-----
If more than one channel is found, return a ``MultiChannelImage``.
For the `pattern` variable, use:
- ``C`` as channel;
- ``T`` as time;
- ``X`` as x-axis (rows);
- ``Y`` as y-axis (planes);
- ``Z`` as z-axis (columns);
Examples
--------
>>> from timagetk.io.image import read_lsm_image
>>> from timagetk.io.util import shared_data
>>> im = read_lsm_image(shared_data("090223-p58-flo-top.lsm", "p58/raw"))
>>> im.get_voxelsize()
[1.0, 0.20031953747828882, 0.20031953747828882]
>>> im.get_shape()
[50, 460, 460]
>>> im.unit
1e-06
>>> im.axes_order
'ZYX'
"""
from tifffile.tifffile import imread as tiffread
from tifffile.tifffile import TiffFile
return_order = kwargs.get('return_order', 'ZYXC')
pattern = kwargs.get('pattern', 'TZCYX')
# Get absolute path and check extension is valid:
if isinstance(lsm_file, str):
lsm_file = Path(lsm_file)
lsm_file = lsm_file.expanduser().absolute()
filename = lsm_file.name
filepath = str(lsm_file.parent)
ext = ''.join(lsm_file.suffixes)
try:
assert ext in LSM_EXT
except AssertionError:
raise NotImplementedError(f"Unknown file extension '{ext}', should be in {TIF_EXT}.")
# - Load the file with third-party reader:
log.debug(f"Using `tifffile.tifffile.TiffFile` to load metadata of {lsm_file}")
tif = TiffFile(lsm_file)
log.debug("Done!")
series = tif.series[0] # get shape and dtype of the first image series
lsm_pattern = series.axes
if lsm_pattern is not None or lsm_pattern != "":
log.debug(f"Found LSM dimensions pattern: {lsm_pattern}") # pattern can be TZCYX or TCZYX...
pattern = lsm_pattern
# Remove "channel" symbol from physical dimension order to use
if 'c' not in pattern.lower():
return_order = return_order.replace('C', '')
lsm_info = series.pages[0].tags[34412]
md_vxs = {"Z": lsm_info.value['VoxelSizeZ'], 'Y': lsm_info.value['VoxelSizeY'], 'X': lsm_info.value['VoxelSizeX']}
log.debug(f"Found metadata 'voxelsize': {md_vxs}")
md_ori = {"Z": lsm_info.value['OriginZ'], 'Y': lsm_info.value['OriginY'], 'X': lsm_info.value['OriginX']}
log.debug(f"Found metadata 'origin': {md_ori}")
md_n_ch = lsm_info.value['DimensionChannels']
md_n_time = lsm_info.value['DimensionTime']
try:
assert md_n_time == 1
except AssertionError:
raise NotImplementedError("Found more than one time point in metadata!")
lsm_md = lsm_info.value['ScanInformation']
name = lsm_md['Name']
# md_ch_names = lsm_info.value['ScanInformation']['Tracks']['Name']
# -- ORIGIN:
origin = tuple([md_ori[dim] for dim in return_order])
# -- VOXELSIZE:
voxelsize = tuple([md_vxs[dim] for dim in return_order])
unit = 1
if all(vxs < 1e-5 for vxs in voxelsize):
voxelsize = [vxs * 1e6 for vxs in voxelsize]
unit = 1e-6
# -- ARRAYS:
log.debug(f"Using `tifffile.tifffile.imread` to load array of {lsm_file}")
lsm_arr = tiffread(lsm_file)
log.debug(f"Loaded {lsm_arr.ndim}D LSM file of shape: {lsm_arr.shape}")
# Defines number of channels in array
ch_axis_id = pattern.find('C') # a value of `-1` means "not found it the str", so no channel defined!
if ch_axis_id != -1:
n_channels = lsm_arr.shape[ch_axis_id]
log.debug(f"Found {n_channels} channels in the array...")
try:
assert n_channels == md_n_ch
except AssertionError:
log.error(f"Loaded an array of shape: {lsm_arr.shape}")
log.error(f"Found a dimension pattern in the metadata: {pattern}")
raise ValueError(f"Metadata ({md_n_ch}) and array ({n_channels}) mismatch on channel number!")
else:
n_channels = 0
if n_channels > 1 and 'C' not in return_order:
return_order += 'C' # add 'channel' to output axis ordering pattern before re-ordering
# - Re-order the array according to ``return_order``
log.debug(f"Array shape BEFORE transpose & reshape: {lsm_arr.shape}")
# Axes id order for axes to returns ('X', 'Y', 'Z' & 'C'):
reorder_axis = [pattern.find(c) for c in return_order]
log.debug(f"Detected axes order: {dict(zip(return_order, reorder_axis))}")
# Axes id order for axes NOT to returns (others like 'T'):
extra_axis = [i for i in range(lsm_arr.ndim) if pattern[i] not in return_order]
# Transpose axes to get array with returned axes as the 4 first
lsm_arr = np.transpose(lsm_arr, tuple(reorder_axis + extra_axis))
# Clip the array to the required axes:
lsm_arr = lsm_arr.reshape(lsm_arr.shape[:len(return_order)])
log.debug(f"Array shape AFTER transpose & reshape: {lsm_arr.shape}")
# Defines file name & path:
if isinstance(lsm_file, str):
fname = lsm_file.split('/')[-1]
fpath = sep.join(lsm_file.split('/')[:-1])
elif isinstance(lsm_file, Path):
fname = lsm_file.name
fpath = lsm_file.parent
else:
fname = ""
fpath = ""
log.debug(f"Loaded file '{fname}' with {n_channels} channels!")
if n_channels > 1:
ch_axis_id = return_order.find('C') # update `ch_axis_id` after reordering
if channel_names is None:
channel_names = [f"CH_{i}" for i in range(n_channels)]
channel_images = []
for i_ch, _ in enumerate(channel_names):
array = np.take(lsm_arr, i_ch, axis=ch_axis_id)
axes_order = return_order[:ch_axis_id:] # remove the "C"
im_kwargs = default_image_attributes(array, axes_order=axes_order, voxelsize=voxelsize, unit=unit,
origin=origin, metadata=lsm_md)
channel_images.append(SpatialImage(array, filepath=filepath, filename=filename, **im_kwargs))
img = MultiChannelImage(channel_images, channel_names=channel_names, filename=fname, filepath=fpath, name=name)
else:
im_kwargs = default_image_attributes(lsm_arr, axes_order=return_order, voxelsize=voxelsize, unit=unit,
origin=origin, metadata=lsm_md)
img = SpatialImage(lsm_arr, filepath=filepath, filename=filename, **im_kwargs)
return img
[docs]
def read_tiff_image(tiff_file, channel_names=None, **kwargs):
"""Read (ImageJ) TIFF image file.
Parameters
----------
tiff_file : str
File path to the image to load.
channel_names : list of str, optional
Defines channel names to use.
By default, try to load channel names from the image metadata.
Other Parameters
----------------
pattern : str, optional
Physical dimension order to use to load the image, "ZCYX" by default.
By default, try to load the pattern from the image metadata.
return_order : str
Physical dimension order to use for the returned array.
Defaults to "ZYXC".
Returns
-------
timagetk.SpatialImage or timagetk.MultiChannelImage
If a single channel and time-point is found, returns a ``SpatialImage`` instance.
If multiple channels are found but with a single time-point, returns a ``MultiChannelImage`` instance.
Notes
-----
If more than one channel is found, return a ``MultiChannelImage``.
For the `pattern` variable, use:
- ``C`` as channel;
- ``X`` as x-axis (rows);
- ``Y`` as y-axis (planes);
- ``Z`` as z-axis (columns);
Examples
--------
>>> from timagetk.io.image import read_tiff_image
>>> from timagetk.io.util import shared_data
>>> im = read_tiff_image(shared_data("p58_t0_fused.tif", "p58/fusion"))
>>> im.get_voxelsize()
[0.2003195434808731, 0.20031954352751363, 0.20031954352751363]
>>> im.get_shape()
[59, 460, 460]
>>> im.axes_order
'ZYX'
"""
from tifffile.tifffile import TiffFile
return_order = kwargs.get('return_order', 'ZYXC')
pattern = kwargs.get('pattern', 'ZCYX')
vxs_order = return_order.replace('C', '')
# Get absolute path and check extension is valid:
if isinstance(tiff_file, str):
tiff_file = Path(tiff_file)
tiff_file = tiff_file.expanduser().absolute()
filename = tiff_file.name
filepath = str(tiff_file.parent)
ext = tiff_file.suffix # should be '.tif' or '.tiff'
try:
assert ext in TIF_EXT
except AssertionError:
raise NotImplementedError(f"Unknown file extension '{ext}', should be in {TIF_EXT}.")
# - Load the file with third-party reader:
tiff_img = TiffFile(tiff_file)
# - Get the array dtype
array_dtype = tiff_img.pages[0].dtype # get dtype of the image in the first page
# - Get the image array:
tiff_channels = tiff_img.asarray().astype(array_dtype)
# - If the first physical dimension is of size one, remove it:
if tiff_channels.shape[0] == 1:
tiff_channels = tiff_channels[0]
# - Gather some metadata:
metadata_dict = {}
int_tags = ['XResolution', 'YResolution', 'ImageWidth', 'ImageLength', 'ImageDescription', 'Location']
try:
page = tiff_img.pages.pages[0]
except:
page = tiff_img.pages[0]
for tag in page.tags.values():
if tag.name in int_tags:
if 'resolution' in tag.name:
if (isinstance(tag.value, tuple) and len(tag.value) == 2):
res = float(tag.value[0] / tag.value[1])
metadata_dict[tag.name] = res
elif isinstance(tag.value, float):
metadata_dict[tag.name] = res
else:
metadata_dict[str(tag.name)] = tag.value
# - Process the 'metadata'
if 'ImageDescription' in metadata_dict:
description = metadata_dict['ImageDescription'].split('\n')
md = dict(x.split('=') for x in description if len(x.split('=')) == 2)
parsed_md = {}
for k, v in md.items():
try:
parsed_md[k] = eval(v)
except:
parsed_md[k] = v
metadata_dict.update(parsed_md)
_ = metadata_dict.pop('ImageDescription')
if 'frames' in parsed_md:
log.warning("Temporal stack!!!")
# TODO: write a temporal stack parser!
# From '/data/Laure/220816_B82-5_bud03/raw/C1-220816_B82-5-allbis.czi - 220816_B82-5-allbis.czi #3.tif':
# {'ImageJ': '1.53q',
# 'images': 4200,
# 'slices': 300,
# 'frames': 14,
# 'hyperstack': 'true',
# 'unit': 'micron',
# 'finterval': 7183.146230769231,
# 'spacing': 0.59,
# 'loop': 'false',
# 'min': 0.0,
# 'max': 65535.0}
# - Process voxel size info and defines the list of values
md_vxs = {'X': 1., 'Y': 1., 'Z': 1.}
if 'XResolution' in metadata_dict:
md_vxs['X'] = float(metadata_dict['XResolution'][1] / metadata_dict['XResolution'][0])
if 'YResolution' in metadata_dict:
md_vxs['Y'] = float(metadata_dict['YResolution'][1] / metadata_dict['YResolution'][0])
if 'spacing' in metadata_dict:
md_vxs['Z'] = float(metadata_dict['spacing'])
else:
vxs_order = vxs_order.replace('Z', '')
voxelsize = tuple([md_vxs[dim] for dim in vxs_order])
# - Defines the number of channels
try:
n_channels = metadata_dict["channels"]
except KeyError:
n_channels = 1 if tiff_channels.ndim <= 3 else tiff_channels.shape[pattern.find('C')]
ch_axis_id = return_order.find('C') # where the channel axes will be after reordering
if n_channels > 1:
if tiff_channels.ndim <= 3:
pattern = pattern.replace('Z', '')
return_order = return_order.replace('Z', '')
voxelsize = voxelsize[1:]
tiff_channels = np.transpose(tiff_channels, tuple([pattern.find(c) for c in return_order]))
if "channel names" in metadata_dict:
try:
channel_names = metadata_dict["channel names"]
assert isinstance(channel_names, list) and len(channel_names) == n_channels
except:
pass
if channel_names is None:
channel_names = ["CH" + str(i) for i in range(n_channels)]
img_channels = []
for i_ch in range(n_channels):
array = np.take(tiff_channels, i_ch, axis=ch_axis_id)
im_kwargs = default_image_attributes(array, voxelsize=voxelsize, metadata=metadata_dict)
img_channels.append(SpatialImage(array, **im_kwargs))
img = MultiChannelImage(img_channels, channel_names=channel_names, filepath=filepath, filename=filename)
else:
pattern = pattern.replace('C', '')
return_order = return_order.replace('C', '')
if tiff_channels.ndim > 3:
tiff_channels = tiff_channels[tuple(
[slice(0, tiff_channels.shape[d]) if (pattern[d] in 'XYZ') else 0 for d in range(tiff_channels.ndim)])]
if tiff_channels.ndim <= 2:
pattern = pattern.replace('Z', '')
return_order = return_order.replace('Z', '')
if len(voxelsize) == 3:
voxelsize = voxelsize[1:]
array = np.transpose(tiff_channels, tuple([pattern.find(c) for c in return_order]))
im_kwargs = default_image_attributes(array, voxelsize=voxelsize, metadata=metadata_dict)
img = SpatialImage(array, filepath=filepath, filename=filename, **im_kwargs)
return img
[docs]
def save_tiff_image(tiff_file, img, **kwargs):
"""Save an image as a TIFF file.
Parameters
----------
tiff_file : str or pathlib.Path
File path where to save the image.
img : timagetk.SpatialImage or timagetk.LabelledImage or timagetk.TissueImage3D or timagetk.MultiChannelImage
Image instance to save on drive.
Other Parameters
----------------
pattern : str, optional
Physical dimension order to use to load the image, "ZCYX" by default.
By default, try to load the pattern from the image metadata.
Notes
-----
For the `pattern` variable, use:
- ``C`` as channel;
- ``X`` as x-axis (rows);
- ``Y`` as y-axis (planes);
- ``Z`` as z-axis (columns);
Example
-------
>>> import tempfile
>>> from timagetk import MultiChannelImage
>>> from timagetk.io.image import save_tiff_image
>>> # EXAMPLE#1 - Save a SpatialImage as a TIFF file:
>>> from timagetk.array_util import random_spatial_image
>>> sp_image = random_spatial_image((3, 4, 5), voxelsize=[0.3, 0.4 ,0.5])
>>> # Save the image in a temporary file:
>>> tmp_path = tempfile.NamedTemporaryFile()
>>> save_tiff_image(tmp_path.name+'.tif', sp_image)
>>> # EXAMPLE#2 - Save a MultiChannel as a TIFF file:
>>> from timagetk.synthetic_data.nuclei_image import example_nuclei_signal_images
>>> n_img, s_img, _, _= example_nuclei_signal_images(n_points=10, extent=20, voxelsize=(1,1,1))
>>> img = MultiChannelImage({'nuclei':n_img, 'signal':s_img})
>>> tmp_path = tempfile.NamedTemporaryFile()
>>> save_tiff_image(tmp_path.name+'.tif', img)
"""
from tifffile.tifffile import TiffWriter
pattern = kwargs.get('pattern', 'ZCYX')
try:
assert isinstance(img, SpatialImage) or isinstance(img, MultiChannelImage)
except AssertionError:
raise TypeError("Parameter 'img' is not a SpatialImage or a MultiChannelImage!")
if isinstance(tiff_file, str):
tiff_file = Path(tiff_file)
tiff_file = tiff_file.expanduser().absolute()
ext = ''.join(tiff_file.suffixes)
try:
assert ext in TIF_EXT
except AssertionError:
raise NotImplementedError(f"Unknown file extension '{ext}', should be in {TIF_EXT}.")
# Check file do NOT exist:
try:
assert not tiff_file.exists()
except AssertionError:
if not kwargs.get('force', False):
from timagetk.util import yn_query
log.warning(f"Found existing file `{tiff_file}`!")
if not yn_query(f"Overwrite existing image `{tiff_file}`", default='no'):
log.info("Will NOT overwrite existing image!")
return
metadata_dict = {}
with TiffWriter(tiff_file, bigtiff=False, imagej=True) as tif:
if isinstance(img, MultiChannelImage):
channel_names = img.channel_names
data = np.array([img.get_channel(c).get_array() for c in channel_names])
n_channels = len(channel_names)
metadata_dict = {'channel names': channel_names}
else:
data = img.get_array()[np.newaxis]
n_channels = 1
order = 'CZYX'
voxelsize = img.voxelsize
n_slices = img.shape[0] if img.is3D() else 1
if img.is2D():
voxelsize = [1.] + list(voxelsize)
pattern = pattern.replace('Z', '')
order = order.replace('Z', '')
data = np.transpose(data, tuple([pattern.find(c) for c in order]))
metadata_dict |= {'slices': n_slices, 'channels': n_channels}
if img.is3D():
metadata_dict |= {'spacing': voxelsize[0]}
tif.write(data, compression=0, resolution=(1.0 / voxelsize[2], 1.0 / voxelsize[1]),
software='timagetk', datetime=now(), metadata=metadata_dict)
tif.close()
log.info(f"Saved image under '{tiff_file}'.")
return
def _brute_search_tiff_metadata(tiff_img):
"""Performs a brute force search of metadata in TIFF image.
Parameters
----------
tiff_img : tifffile.tifffile.TiffFile
The TIFF file to search for metadata.
Returns
-------
dict
A dictionary with parsed metadata.
"""
# List all attributes related to metadata:
md_attrs = [attr for attr in dir(tiff_img) if 'metadata' in attr]
md_dict = {}
for md_attr in md_attrs:
md = getattr(tiff_img, md_attr)
if md is not None:
md_dict[md_attr] = md
if len(md_dict) == 0:
log.error(f"No metadata found!")
return
elif len(md_dict) == 1:
md_attr = list(md_dict.keys())[0]
log.info(f"Found metadata in '{md_attr}' TIFF image attribute!")
return md_attr, md_dict[md_attr]
else:
log.warning(f"Found metadata in more than one TIFF image attribute!")
return md_dict
def _ome_metadata_parser(xml_input):
"""Parse metadata from OME-XML metadata.
Parameters
----------
xml_input : str
The OME-XML string to parse.
Returns
-------
dict
A dictionary with parsed metadata.
"""
import xmltodict
xml_dict = xmltodict.parse(xml_input)
pixel_md = xml_dict['OME']['Image']['Pixels']
out_md = {}
out_md['channels'] = pixel_md['@SizeC']
out_md['frames'] = pixel_md['@SizeT']
out_md['pattern'] = pixel_md['@DimensionOrder']
return out_md
