MathOnco / valis

Virtual Alignment of pathoLogy Image Series
https://valis.readthedocs.io/en/latest/
MIT License
121 stars 29 forks source link

'NoneType' object has no attribute 'shape' #76

Open Justin-Krieg opened 1 year ago

Justin-Krieg commented 1 year ago

Hi @cdgatenbee

Not sure if this is an error specific to me or not, but trying to use czi files with JPEG-XR compression with version (1.0.2). I followed #60 and tried different slide readers (bioformats and czijpegxr) and im getting this error:

==== Rigid registration

Detecting features   :   0%|                                                                                                                                                                      | 0/3 [00:00<?, ?image/s] 
C:\Users\Home\AppData\Local\Programs\Python\Python310\lib\site-packages\valis\valtils.py:22: UserWarning: 'NoneType' object has no attribute 'shape'
  warnings.warn(warning_msg, warning_type)
Traceback (most recent call last):
  File "C:\Users\Home\AppData\Local\Programs\Python\Python310\lib\site-packages\valis\registration.py", line 4116, in register
    rigid_registrar = self.rigid_register()
  File "C:\Users\Home\AppData\Local\Programs\Python\Python310\lib\site-packages\valis\registration.py", line 2981, in rigid_register
    rigid_registrar = serial_rigid.register_images(self.processed_dir,
  File "C:\Users\Home\AppData\Local\Programs\Python\Python310\lib\site-packages\valis\serial_rigid.py", line 1536, in register_images
    registrar.generate_img_obj_list(feature_detector, qt_emitter=qt_emitter)
  File "C:\Users\Home\AppData\Local\Programs\Python\Python310\lib\site-packages\valis\serial_rigid.py", line 511, in generate_img_obj_list
    img_obj.kp_pos_xy, img_obj.desc = feature_detector.detect_and_compute(img)
  File "C:\Users\Home\AppData\Local\Programs\Python\Python310\lib\site-packages\valis\feature_detectors.py", line 166, in detect_and_compute
    if desc.shape[0] > MAX_FEATURES:
AttributeError: 'NoneType' object has no attribute 'shape'

I then thought id convert these files to ome.tiff using bfconvert (with the -noflat option) and i'm having the same error. Leaving the slide reader option blank yields the same result. This is the relevant code being used. I've again tried using different slide readers, and tried removing the crop.

registrar = registration.Valis(slide_src_dir,  results_dir, img_list = img_list_path, crop="overlap")

rigid_registrar, non_rigid_registrar, error_df = registrar.register()
registrar.register_micro(max_non_rigid_registration_dim_px=2000)

registrar.warp_and_merge_slides(name_list_Stitched, channel_name_dict=channel_name_dict_defined, drop_duplicates=True, crop=True, compression="jp2k", Q=85)

I am using windows 11 with the natively installed programs rather than docker. I've used the ome.tiff files in previous versions of Valis (rc12)

cdgatenbee commented 1 year ago

Hi @kriegy98, Sorry you've run into this :( Based on the error message, it looks like the images are actually be being read OK, since I would think an error related to reading the slides would be thrown earlier in the pipeline, either when generating thumbnails or pre-processing the images. Instead, it looks like they're being converted and processed, but that there is a problem in the feature detection step. Could you check if there are images in results_dir/processed? If something is there, and they look ok (e.g. DAPI channel only?), then it suggests that the SlideReader is working ok and the problem is elsewhere. If not, then maybe it is the SlideReader. Anyway, please let me know what you see and then we'll work from there.

Best, -Chandler

Justin-Krieg commented 1 year ago

Hi @cdgatenbee So having a look at the masks helped, it seems as though Valis is reading in the slide image preview on these images. I have some files from the same round work fine, i converted them to ome.tiff and compared it to ones that didnt work and i couldnt find any major differences from the bftools showinf tool. They both had series 6 and 7 be the slide label preview and the slide image preview, you can tell as these are RGB and the tile is non consistent. I tried specifying the series=0 (which is the largest image size) in the main registratration.Valis class but that still didnt yield any better results (it seems as though Valis isnt observing it). Is there any method to get valis to ignore these embedded extra images? Here is the output from showinf:

filename = B:\R1R3\JK023-15_R1.ome.tiff
Series count = 8
Series #0 :
        Image count = 3
        RGB = false (1)
        Interleaved = false
        Indexed = false (false color)
        Width = 50058
        Height = 62937
        SizeZ = 1
        SizeT = 1
        SizeC = 3
        Tile size = 1024 x 1024
        Thumbnail size = 101 x 128
        Endianness = intel (little)
        Dimension order = XYCZT (certain)
        Pixel type = uint16
        Valid bits per pixel = 16
        Metadata complete = true
        Thumbnail series = false
        -----
        Plane #0 <=> Z 0, C 0, T 0
        Plane #1 <=> Z 0, C 1, T 0
        Plane #2 <=> Z 0, C 2, T 0

Series #1 :
        Image count = 3
        RGB = false (1)
        Interleaved = false
        Indexed = false (false color)
        Width = 16686
        Height = 20979
        SizeZ = 1
        SizeT = 1
        SizeC = 3
        Tile size = 1024 x 1024
        Thumbnail size = 101 x 128
        Endianness = intel (little)
        Dimension order = XYCZT (certain)
        Pixel type = uint16
        Valid bits per pixel = 16
        Metadata complete = true
        Thumbnail series = false
        -----
        Plane #0 <=> Z 0, C 0, T 0
        Plane #1 <=> Z 0, C 1, T 0
        Plane #2 <=> Z 0, C 2, T 0

Series #2 :
        Image count = 3
        RGB = false (1)
        Interleaved = false
        Indexed = false (false color)
        Width = 5562
        Height = 6993
        SizeZ = 1
        SizeT = 1
        SizeC = 3
        Tile size = 1024 x 1024
        Thumbnail size = 101 x 128
        Endianness = intel (little)
        Dimension order = XYCZT (certain)
        Pixel type = uint16
        Valid bits per pixel = 16
        Metadata complete = true
        Thumbnail series = false
        -----
        Plane #0 <=> Z 0, C 0, T 0
        Plane #1 <=> Z 0, C 1, T 0
        Plane #2 <=> Z 0, C 2, T 0

Series #3 :
        Image count = 3
        RGB = false (1)
        Interleaved = false
        Indexed = false (false color)
        Width = 1854
        Height = 2331
        SizeZ = 1
        SizeT = 1
        SizeC = 3
        Tile size = 1854 x 1
        Thumbnail size = 101 x 128
        Endianness = intel (little)
        Dimension order = XYCZT (certain)
        Pixel type = uint16
        Valid bits per pixel = 16
        Metadata complete = true
        Thumbnail series = false
        -----
        Plane #0 <=> Z 0, C 0, T 0
        Plane #1 <=> Z 0, C 1, T 0
        Plane #2 <=> Z 0, C 2, T 0

Series #4 :
        Image count = 3
        RGB = false (1)
        Interleaved = false
        Indexed = false (false color)
        Width = 618
        Height = 777
        SizeZ = 1
        SizeT = 1
        SizeC = 3
        Tile size = 618 x 1
        Thumbnail size = 101 x 128
        Endianness = intel (little)
        Dimension order = XYCZT (certain)
        Pixel type = uint16
        Valid bits per pixel = 16
        Metadata complete = true
        Thumbnail series = false
        -----
        Plane #0 <=> Z 0, C 0, T 0
        Plane #1 <=> Z 0, C 1, T 0
        Plane #2 <=> Z 0, C 2, T 0

Series #5 :
        Image count = 3
        RGB = false (1)
        Interleaved = false
        Indexed = false (false color)
        Width = 206
        Height = 259
        SizeZ = 1
        SizeT = 1
        SizeC = 3
        Tile size = 206 x 1
        Thumbnail size = 101 x 128
        Endianness = intel (little)
        Dimension order = XYCZT (certain)
        Pixel type = uint16
        Valid bits per pixel = 16
        Metadata complete = true
        Thumbnail series = false
        -----
        Plane #0 <=> Z 0, C 0, T 0
        Plane #1 <=> Z 0, C 1, T 0
        Plane #2 <=> Z 0, C 2, T 0

Series #6 :
        Image count = 1
        RGB = true (3)
        Interleaved = false
        Indexed = false (false color)
        Width = 701
        Height = 530
        SizeZ = 1
        SizeT = 1
        SizeC = 3 (effectively 1)
        Tile size = 701 x 1
        Thumbnail size = 128 x 96
        Endianness = intel (little)
        Dimension order = XYCZT (certain)
        Pixel type = uint16
        Valid bits per pixel = 12
        Metadata complete = true
        Thumbnail series = false
        -----
        Plane #0 <=> Z 0, C 0, T 0

Series #7 :
        Image count = 1
        RGB = true (3)
        Interleaved = false
        Indexed = false (false color)
        Width = 1571
        Height = 691
        SizeZ = 1
        SizeT = 1
        SizeC = 3 (effectively 1)
        Tile size = 1571 x 1
        Thumbnail size = 128 x 56
        Endianness = intel (little)
        Dimension order = XYCZT (certain)
        Pixel type = uint16
        Valid bits per pixel = 12
        Metadata complete = true
        Thumbnail series = false
        -----
        Plane #0 <=> Z 0, C 0, T 0

and here is the mask that Valis outputted. Note it did make a processed image, but its just black. JK023-15

Link to file here

cdgatenbee commented 1 year ago

Hi @kriegy98, Thanks for all the info and file, it helps a lot! I'll dig into this to see if I can figure out why valis seems to be reading the preview image instead of the whole slide.

Best, -Chandler

Justin-Krieg commented 1 year ago

Hi @cdgatenbee Any updates with this? It might be worthwhile trying to get Valis to clean the image before reading it in stripping the image of any small slide labels and preview scans.

cdgatenbee commented 1 year ago

Hi @kriegy98, Sorry, nothing yet :( We're having a big workshop this week, so unfortunately I don't think I will be able to work on this issue until next week. I'll keep you updated though.

Best, -Chandler

Justin-Krieg commented 12 months ago

Hi @cdgatenbee Have you had a chance to work on this yet?

cdgatenbee commented 11 months ago

Hi @kriegy98, Sorry, I was a conference most of last week and so didn't get a chance to work on this. Hoping to get this and a few other issues (including issue 82) sorted out this week. I'll keep you updated.

Best, -Chandler

cdgatenbee commented 11 months ago

Hi @kriegy98, Just wanted to give you a quick update on this issue. The problem seems to be related to trying to extract the channel names, which ended up being a problem because the structure of the metadata is slightly different when the images are RGB and when they are IF. I hadn't caught this before because the other czi images I have are RGB or single channel images. I just need to test with the other jpegxr compressed czi images I have, but this should be fixed in the next update. Thanks again for your patience.

Best, -Chandler

cdgatenbee commented 11 months ago

Hi @kriegy98, If you'd like to try out the updated czi jpegxr reader, here is the code:

from valis.slide_io import *

class CziJpgxrReader(SlideReader):
    """Read slides and get metadata

    Attributes
    ----------
    slide_f : str
        Path to slide

    metadata : MetaData
        MetaData containing some basic metadata about the slide

    series : int
        Image series

    """

    def __init__(self, src_f, series=None, *args, **kwargs):
        """
        Parameters
        -----------
        src_f : str
            Path to slide

        series : int
            The series to be read. If `series` is None, the the `series`
            will be set to the series associated with the largest image.

        """

        try:
            from aicspylibczi import CziFile
        except Exception as e:
            msg = "Please install aicspylibczi"
            print(e)
            valtils.print_warning(msg)

        czi_reader = CziFile(src_f)
        self.original_meta_dict = valtils.etree_to_dict(czi_reader.meta)

        self.is_bgr = False
        self.meta_list = [None]

        super().__init__(src_f=src_f, *args, **kwargs)

        try:
            self.meta_list = self.create_metadata()
        except Exception as e:
            print(e)

        self.n_series = len(self.meta_list)
        if series is None:
            img_areas = [np.multiply(*meta.slide_dimensions[0]) for meta in self.meta_list]
            series = np.argmax(img_areas)
            msg = (f"No series provided. "
                   f"Selecting series with largest image, "
                   f"which is series {series}")

            valtils.print_warning(msg, warning_type=None, rgb=valtils.Fore.GREEN)

        self._series = series
        self.series = series

    def _set_series(self, series):
        self._series = series
        self.metadata = self.meta_list[series]

    def _get_series(self):
        return self._series

    series = property(fget=_get_series,
                      fset=_set_series,
                      doc="Slide scene")

    def slide2vips(self, level=0, xywh=None, *args, **kwargs):
        czi_reader = CziFile(self.src_f)

        out_shape_wh = self.metadata.slide_dimensions[0]
        tile_bboxes = czi_reader.get_all_mosaic_tile_bounding_boxes(C=0)

        vips_img = pyvips.Image.black(*out_shape_wh, bands=self.metadata.n_channels)
        print(f"Building CZI mosaic for {valtils.get_name(self.src_f)}")
        for tile_info, tile_bbox in tqdm(tile_bboxes.items()):
            m = tile_info.m_index
            x = tile_bbox.x
            y = tile_bbox.y

            np_tile, tile_dims = czi_reader.read_image(M=m)

            slice_dims = [v - 1 for k, v in tile_dims if k not in ["Y", "X", "A"]]

            np_tile = np_tile[(*slice_dims, ...)]
            if self.is_bgr:
                np_tile = np_tile[..., ::-1]

            vips_tile = warp_tools.numpy2vips(np_tile)
            vips_img = vips_img.insert(vips_tile, *(x, y))

        if xywh is not None:
            vips_img = vips_img.extract_area(*xywh)

        if level != 0:
            scaling = self.metadata._zoom_levels[level]
            vips_img = warp_tools.rescale_img(vips_img, scaling)
        if self.is_bgr:
            vips_img = vips_img.copy(interpretation="srgb")

        np_type = slide_tools.CZI_FORMAT_TO_BF_FORMAT[czi_reader.pixel_type]
        vips_type = slide_tools.NUMPY_FORMAT_VIPS_DTYPE[np_type]
        vips_img = vips_img.cast(vips_type)

        return vips_img

    def slide2image(self, level, xywh=None, *args, **kwargs):
        """Convert slide to image

        Parameters
        -----------
        level : int
            Pyramid level

        xywh : tuple of int, optional
            The region to be sliced from the slide. If None,
            then the entire slide will be converted. Otherwise
            xywh is the (top left x, top left y, width, height) of
            the region to be sliced.

        Returns
        -------
        img : ndarray
            An image of the slide or the region defined by xywh

        """

        vips_img = self.slide2vips(level=level, xywh=xywh, *args, **kwargs)
        np_img = warp_tools.vips2numpy(vips_img)

        return np_img

    def create_metadata(self):
        """ Create and fill in a MetaData object

        Returns
        -------
        metadata : MetaData
            MetaData object containing metadata about slide

        """

        czi_reader = CziFile(self.src_f)
        dims_dict = czi_reader.get_dims_shape()

        n_scenes = len(dims_dict)
        meta_list = [None] * n_scenes
        phys_size = self._get_pixel_physical_size()

        with valtils.HiddenPrints():
            bf_reader = BioFormatsSlideReader(self.src_f)
            original_xml = bf_reader.metadata.original_xml

        for i in range(n_scenes):

            temp_name = f"{os.path.split(self.src_f)[1]}".strip("_")
            full_name = f"{temp_name}_Scene_{i}"

            series_meta = MetaData(full_name, "aicspylibczi", series=i)

            series_meta.is_rgb = self._check_rgb()

            if series_meta.is_rgb:
                n_channels = dims_dict[i]["A"][1]
                series_meta.pyvips_interpretation = 'srgb'
            else:
                n_channels = dims_dict[i]["C"][1]
                if n_channels == 1:
                    series_meta.pyvips_interpretation = 'b-w'
                else:
                    series_meta.pyvips_interpretation = 'multiband'

            series_meta.n_channels = n_channels
            series_meta.slide_dimensions = self._get_slide_dimensions(i)
            series_meta.bf_datatype = slide_tools.CZI_FORMAT_TO_BF_FORMAT[czi_reader.pixel_type]
            series_meta.channel_names = self._get_channel_names(meta=series_meta)
            series_meta.pixel_physical_size_xyu = phys_size
            series_meta.original_xml = original_xml
            series_meta._zoom_levels = self._get_zoom_levels(i)

            meta_list[i] = series_meta

        return meta_list

    def _get_img_meta_dict(self):
        return self.original_meta_dict["ImageDocument"]["Metadata"]["Information"]["Image"]

    def _check_rgb(self, *args, **kwargs):
        """Determine if image is RGB

        Returns
        -------
        is_rgb : bool
            Whether or not the image is RGB

        """
        czi_reader = CziFile(self.src_f)
        self.is_bgr = czi_reader.pixel_type.startswith("bgr")
        _is_rgb = czi_reader.pixel_type.startswith("rgb")
        is_rgb  =_is_rgb or self.is_bgr

        return is_rgb

    def _get_channel_names(self, meta, *args, **kwargs):
        """Get names of each channel

        Get list of channel names

        Returns
        -------
        channel_names : list
            List of channel names

        """
        if meta.is_rgb:
            return None

        img_dict = self._get_img_meta_dict()
        channels = img_dict["Dimensions"]["Channels"]
        if "Channel" in channels:
            channels = channels["Channel"]

        if isinstance(channels, dict):
            channels = list(channels.values())

        try:
            all_channel_ids = [eval(x["@Id"].split(":")[1]) for x in channels]
            max_c = max([eval(img_dict["SizeC"]), max(all_channel_ids)+1])
            channel_names = [None] * max_c
        except:
            channel_names = [None] * eval(img_dict["SizeC"])

        for chnl_attr in channels:
            chnl_name = chnl_attr["@Name"]
            chnl_idx = eval(chnl_attr["@Id"].split(":")[1])
            channel_names[chnl_idx] = chnl_name

        channel_names = [x for x in channel_names if x is not None]

        return channel_names

    def _get_slide_dimensions(self, scene=0, *args, **kwargs):
        """Get dimensions of slide at all pyramid levels

        Returns
        -------
        slide_dims : ndarray
            Dimensions of all images in the pyramid (width, height).

        """
        zoom_levels = self._get_zoom_levels(scene)

        czi_reader = CziFile(self.src_f)
        scene_bbox = czi_reader.get_all_scene_bounding_boxes()[scene]
        scence_l0_wh = np.array([scene_bbox.w, scene_bbox.h])
        slide_dimensions = np.round(scence_l0_wh*zoom_levels[..., np.newaxis]).astype(int)

        return slide_dimensions

    def _get_zoom_levels(self, scene=0):

        img_dict = self._get_img_meta_dict()
        if "S" not in img_dict["Dimensions"]:
            # No pyramid levels
            return np.array([1.0])

        scene_dict = img_dict["Dimensions"]["S"]["Scenes"]["Scene"]
        if scene in scene_dict:
            pyramid_info = scene_dict[scene]["PyramidInfo"]
        else:
            pyramid_info = scene_dict["PyramidInfo"]
        n_levels = eval(pyramid_info["PyramidLayersCount"])
        downsampling = eval(pyramid_info["MinificationFactor"])
        zoom_levels = (1/downsampling)**(np.arange(0, n_levels))

        return zoom_levels

    def _get_pixel_physical_size(self, *args, **kwargs):
        """Get resolution of slide

        Returns
        -------
        res_xyu : tuple
            Physical size per pixel and the unit, e.g. u'\u00B5m'

        Notes
        -----
            If physical unit is micron, it must be u'\u00B5m',
            not mu (u'\u03bcm') or u.

        """

        physical_sizes = self.original_meta_dict["ImageDocument"]["Metadata"]["Scaling"]["Items"]["Distance"]

        physical_size_xyu = [None] * 3
        physical_unit = physical_sizes[0]["DefaultUnitFormat"]
        physical_size_xyu[2] = physical_unit

        if physical_unit == u'\u00B5m':
            physical_scaling = 10**6
        elif physical_unit == "mm":
            physical_scaling = 10**3
        elif physical_unit == "cm":
            physical_scaling = 10**2
        else:
            physical_scaling = 1

        for ps in physical_sizes:
            if ps["@Id"] == "X":
                physical_size_xyu[0] = eval(ps["Value"])*physical_scaling
            elif ps["@Id"] == "Y":
                physical_size_xyu[1] = eval(ps["Value"])*physical_scaling

        return tuple(physical_size_xyu)

    def _read_rgb(self, level=0, xywh=None, *args, **kwargs):

        czi_reader = CziFile(self.src_f)
        img_dict = self._get_img_meta_dict()
        bg_hex = img_dict["Dimensions"]["Channels"]["Channel"]["Color"]
        bg_rgba = valtils.hex_to_rgb(bg_hex)[::-1]

        out_shape_wh = self.metadata.slide_dimensions[0]
        tile_bboxes = czi_reader.get_all_mosaic_tile_bounding_boxes(C=0)

        vips_img = pyvips.Image.black(*out_shape_wh, bands=self.metadata.n_channels) + bg_rgba[0:3]
        print(f"Building CZI mosaic for {valtils.get_name(self.src_f)}")
        for tile_info, tile_bbox in tqdm(tile_bboxes.items()):
            m = tile_info.m_index
            x = tile_bbox.x
            y = tile_bbox.y

            np_tile, tile_dims = czi_reader.read_image(M=m)

            slice_dims = [v - 1 for k, v in tile_dims if k not in ["Y", "X", "A"]]

            np_tile = np_tile[(*slice_dims, ...)]
            if self.is_bgr:
                np_tile = np_tile[..., ::-1]

            vips_tile = warp_tools.numpy2vips(np_tile)
            vips_img = vips_img.insert(vips_tile, *(x, y))

        if xywh is not None:
            vips_img = vips_img.extract_area(*xywh)

        if level != 0:
            scaling = self.metadata._zoom_levels[level]
            vips_img = warp_tools.rescale_img(vips_img, scaling)

        vips_img = vips_img.copy(interpretation="srgb")
        np_type = slide_tools.CZI_FORMAT_TO_BF_FORMAT[czi_reader.pixel_type]
        vips_type = slide_tools.NUMPY_FORMAT_VIPS_DTYPE[np_type]
        vips_img = vips_img.cast(vips_type)

        return vips_img

You can use it by passing an uninstantiated object to the reader_cls argument of Valis.register:

registrar = registration.Valis(src_dir, dst_dir)
rigid_registrar, non_rigid_registrar, error_df = registrar.register(reader_cls=CziJpgxrReader)

This is the image that is extracted when I use the above reader on your file:

JK023-15

Hopefully that looks correct, but if not, or if you run into any other issues, please let me know.

Best, -Chandler

Justin-Krieg commented 11 months ago

Hi @cdgatenbee thank you for providing an update, unfortunately i still have an error when i try to load more than one image. If you have a look at the output below it seems to fail on loading the second image. This holds true when i swap the order of images being loaded, suggesting that its not an image issue but perhaps an issue with delinating between images. There is no output with this either making it difficult to diagnose further.

==== Converting images

Converting images:   0%|                                                                                                                                                                                           | 0/2 [00:00<?, ?image/s]SLF4J: Failed to load class "org.slf4j.impl.StaticLoggerBinder".
SLF4J: Defaulting to no-operation (NOP) logger implementation
SLF4J: See http://www.slf4j.org/codes.html#StaticLoggerBinder for further details.
Building CZI mosaic for JK023-15_R3
100%|█████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████| 796/796 [00:43<00:00, 18.35it/s]
Converting images:  50%|█████████████████████████████████████████████████████████████████████████████████████████                                                                                         | 1/2 [01:52<01:52, 112.66s/image]
Building CZI mosaic for JK023-15
100%|█████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████| 765/765 [00:43<00:00, 17.39it/s]
Valis process returned an error███████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████| 765/765 [00:43<00:00, 17.98it/s] 
JVM has been killed. If this was due to an error, then a new Python session will need to be started
PS C:\Slides> 
cdgatenbee commented 11 months ago

Hi @kriegy98, that's so strange. I'll see if I can reproduce it and debug it. I'll keep you updated.

Best, -Chandler

cdgatenbee commented 11 months ago

Hi @kriegy98, I've updated the reader class (see below), and tested it by reading 5 different czi images: 3 RGB (2 mosaic, 1 not mosaic); 1 multichannel non-RGB (mosaic); 1 single channel (mosaic). Three of those 5 have JPEGXR compression. All images were read correctly in the same Python session, so it looks like the reader is successfully delineating the images. Please try out the updated reader and let me know if works for you too.

Best, -Chandler


class CziJpgxrReader(SlideReader):
    """Read slides and get metadata

    Attributes
    ----------
    slide_f : str
        Path to slide

    metadata : MetaData
        MetaData containing some basic metadata about the slide

    series : int
        Image series

    """
    def __init__(self, src_f, series=None, *args, **kwargs):
        """
        Parameters
        -----------
        src_f : str
            Path to slide

        series : int
            The series to be read. If `series` is None, the the `series`
            will be set to the series associated with the largest image.

        """
        try:
            from aicspylibczi import CziFile
        except Exception as e:
            msg = "Please install aicspylibczi"
            print(e)
            valtils.print_warning(msg)

        czi_reader = CziFile(src_f)
        self.original_meta_dict = valtils.etree_to_dict(czi_reader.meta)

        self.is_bgr = False
        self.meta_list = [None]

        super().__init__(src_f=src_f, *args, **kwargs)

        try:
            self.meta_list = self.create_metadata()
        except Exception as e:
            print(e)

        self.n_series = len(self.meta_list)
        if series is None:
            img_areas = [np.multiply(*meta.slide_dimensions[0]) for meta in self.meta_list]
            series = np.argmax(img_areas)
            msg = (f"No series provided. "
                   f"Selecting series with largest image, "
                   f"which is series {series}")

            valtils.print_warning(msg, warning_type=None, rgb=valtils.Fore.GREEN)

        self._series = series
        self.series = series

    def _set_series(self, series):
        self._series = series
        self.metadata = self.meta_list[series]

    def _get_series(self):
        return self._series

    series = property(fget=_get_series,
                      fset=_set_series,
                      doc="Slide scene")

    def _read_whole_img(self, level=0, xywh=None, *args, **kwargs):
        """

        Return
        ------
        img : ndarray

        """

        czi_reader = CziFile(self.src_f)
        img, shp = czi_reader.read_image(S=self.series)
        shp = dict(shp)
        if img.ndim == 4 and shp["T"] == 1:
            img = img[0]

        if self.is_bgr:
            img = img[..., ::-1]

        vips_img = warp_tools.numpy2vips(img)

        return vips_img

    def _read_mosaic(self, level=0, xywh=None, *args, **kwargs):
        czi_reader = CziFile(self.src_f)

        out_shape_wh = self.metadata.slide_dimensions[0]
        tile_bboxes = czi_reader.get_all_mosaic_tile_bounding_boxes(C=0)

        vips_img = pyvips.Image.black(*out_shape_wh, bands=self.metadata.n_channels)
        print(f"Building CZI mosaic for {valtils.get_name(self.src_f)}")
        for tile_info, tile_bbox in tqdm(tile_bboxes.items()):
            m = tile_info.m_index
            x = tile_bbox.x
            y = tile_bbox.y

            np_tile, tile_dims = czi_reader.read_image(S=self.series, M=m)

            slice_dims = [v - 1 for k, v in tile_dims if k not in ["Y", "X", "A"]]

            np_tile = np_tile[(*slice_dims, ...)]
            if self.is_bgr:
                np_tile = np_tile[..., ::-1]

            vips_tile = warp_tools.numpy2vips(np_tile)
            vips_img = vips_img.insert(vips_tile, x, y)

        return vips_img

    def slide2vips(self, level=0, xywh=None, *args, **kwargs):
        try:
            # Image is mosaic
            vips_img = self._read_mosaic(level=level, xywh=xywh,*args, **kwargs)

        except Exception as e:
            print(e)
            print("Reading whole image")
            vips_img = self._read_whole_img(level=level, xywh=xywh,*args, **kwargs)

        czi_reader = CziFile(self.src_f)
        if xywh is not None:
            vips_img = vips_img.extract_area(*xywh)

        if level != 0:
            scaling = self.metadata._zoom_levels[level]
            vips_img = warp_tools.rescale_img(vips_img, scaling)
        if self.is_bgr:
            vips_img = vips_img.copy(interpretation="srgb")

        np_type = slide_tools.CZI_FORMAT_TO_BF_FORMAT[czi_reader.pixel_type]
        vips_type = slide_tools.NUMPY_FORMAT_VIPS_DTYPE[np_type]
        vips_img = vips_img.cast(vips_type)

        return vips_img

    def slide2image(self, level, xywh=None, *args, **kwargs):
        """Convert slide to image

        Parameters
        -----------
        level : int
            Pyramid level

        xywh : tuple of int, optional
            The region to be sliced from the slide. If None,
            then the entire slide will be converted. Otherwise
            xywh is the (top left x, top left y, width, height) of
            the region to be sliced.

        Returns
        -------
        img : ndarray
            An image of the slide or the region defined by xywh

        """

        vips_img = self.slide2vips(level=level, xywh=xywh, *args, **kwargs)
        np_img = warp_tools.vips2numpy(vips_img)

        return np_img

    def create_metadata(self):
        """ Create and fill in a MetaData object

        Returns
        -------
        metadata : MetaData
            MetaData object containing metadata about slide

        """

        czi_reader = CziFile(self.src_f)
        dims_dict = czi_reader.get_dims_shape()

        n_scenes = len(dims_dict)
        meta_list = [None] * n_scenes
        phys_size = self._get_pixel_physical_size()

        with valtils.HiddenPrints():
            bf_reader = BioFormatsSlideReader(self.src_f)
            original_xml = bf_reader.metadata.original_xml

        for i in range(n_scenes):

            temp_name = f"{os.path.split(self.src_f)[1]}".strip("_")
            full_name = f"{temp_name}_Scene_{i}"

            series_meta = MetaData(full_name, "aicspylibczi", series=i)

            series_meta.is_rgb = self._check_rgb()

            if series_meta.is_rgb:
                n_channels = dims_dict[i]["A"][1]
                series_meta.pyvips_interpretation = 'srgb'
            else:
                n_channels = dims_dict[i]["C"][1]
                if n_channels == 1:
                    series_meta.pyvips_interpretation = 'b-w'
                else:
                    series_meta.pyvips_interpretation = 'multiband'

            series_meta.n_channels = n_channels
            series_meta.slide_dimensions = self._get_slide_dimensions(i)
            series_meta.bf_datatype = slide_tools.CZI_FORMAT_TO_BF_FORMAT[czi_reader.pixel_type]
            series_meta.channel_names = self._get_channel_names(meta=series_meta)

            series_meta.pixel_physical_size_xyu = phys_size
            series_meta.original_xml = original_xml
            series_meta._zoom_levels = self._get_zoom_levels(i)

            meta_list[i] = series_meta

        return meta_list

    def _get_img_meta_dict(self):
        return self.original_meta_dict["ImageDocument"]["Metadata"]["Information"]["Image"]

    def _check_rgb(self, *args, **kwargs):
        """Determine if image is RGB

        Returns
        -------
        is_rgb : bool
            Whether or not the image is RGB

        """
        czi_reader = CziFile(self.src_f)
        self.is_bgr = czi_reader.pixel_type.startswith("bgr")
        _is_rgb = czi_reader.pixel_type.startswith("rgb")
        is_rgb  =_is_rgb or self.is_bgr

        return is_rgb

    def _get_channel_names_aics(self, meta, *args, **kwargs):
        """Get names of each channel

        Get list of channel names

        Returns
        -------
        channel_names : list
            List of channel names

        """

        if meta.is_rgb:
            return None

        img_dict = self._get_img_meta_dict()
        channels = img_dict["Dimensions"]["Channels"]
        if "Channel" in channels:
            channels = channels["Channel"]

        if meta.n_channels == 1 and "@Name" in channels:
            channel_names = [channels["@Name"]]

            return channel_names

        if isinstance(channels, dict):
            channels = list(channels.values())

        try:
            all_channel_ids = [x["@Id"].split(":") for x in channels]
            all_channel_ids = [eval(x["@Id"].split(":")[1]) for x in channels]
            max_c = max([eval(img_dict["SizeC"]), max(all_channel_ids)+1])
            channel_names = [None] * max_c
        except:
            channel_names = [None] * eval(img_dict["SizeC"])

        for chnl_attr in channels:
            chnl_name = chnl_attr["@Name"]
            chnl_idx = eval(chnl_attr["@Id"].split(":")[1])
            channel_names[chnl_idx] = chnl_name

        channel_names = [x for x in channel_names if x is not None]

        return channel_names

    def _get_channel_names_bf(self, meta, *args, **kwargs):
        """Get names of each channel

        Get list of channel names

        Returns
        -------
        channel_names : list
            List of channel names

        """
        if meta.is_rgb:
            return None

        with valtils.HiddenPrints():
            bf_reader = BioFormatsSlideReader(self.src_f)

        rdr, bf_meta = bf_reader._get_bf_objects()
        channel_names = bf_reader._get_channel_names(rdr, bf_meta)

        return channel_names

    def _get_channel_names(self, meta, *args, **kwargs):
        channel_names = self._get_channel_names_aics(meta)
        return channel_names

    def _get_slide_dimensions(self, scene=0, *args, **kwargs):
        """Get dimensions of slide at all pyramid levels

        Returns
        -------
        slide_dims : ndarray
            Dimensions of all images in the pyramid (width, height).

        """
        zoom_levels = self._get_zoom_levels(scene)

        czi_reader = CziFile(self.src_f)
        scene_bbox = czi_reader.get_all_scene_bounding_boxes()[scene]
        scence_l0_wh = np.array([scene_bbox.w, scene_bbox.h])
        slide_dimensions = np.round(scence_l0_wh*zoom_levels[..., np.newaxis]).astype(int)

        return slide_dimensions

    def _get_zoom_levels(self, scene=0):

        img_dict = self._get_img_meta_dict()
        if "S" not in img_dict["Dimensions"]:
            # No pyramid levels
            return np.array([1.0])

        scene_dict = img_dict["Dimensions"]["S"]["Scenes"]["Scene"]
        if scene in scene_dict:
            pyramid_info = scene_dict[scene]["PyramidInfo"]
        else:
            pyramid_info = scene_dict["PyramidInfo"]
        n_levels = eval(pyramid_info["PyramidLayersCount"])
        downsampling = eval(pyramid_info["MinificationFactor"])
        zoom_levels = (1/downsampling)**(np.arange(0, n_levels))

        return zoom_levels

    def _get_pixel_physical_size(self, *args, **kwargs):
        """Get resolution of slide

        Returns
        -------
        res_xyu : tuple
            Physical size per pixel and the unit, e.g. u'\u00B5m'

        Notes
        -----
            If physical unit is micron, it must be u'\u00B5m',
            not mu (u'\u03bcm') or u.

        """

        physical_sizes = self.original_meta_dict["ImageDocument"]["Metadata"]["Scaling"]["Items"]["Distance"]

        physical_size_xyu = [None] * 3
        physical_unit = physical_sizes[0]["DefaultUnitFormat"]
        physical_size_xyu[2] = physical_unit

        if physical_unit == u'\u00B5m':
            physical_scaling = 10**6
        elif physical_unit == "mm":
            physical_scaling = 10**3
        elif physical_unit == "cm":
            physical_scaling = 10**2
        else:
            physical_scaling = 1

        for ps in physical_sizes:
            if ps["@Id"] == "X":
                physical_size_xyu[0] = eval(ps["Value"])*physical_scaling
            elif ps["@Id"] == "Y":
                physical_size_xyu[1] = eval(ps["Value"])*physical_scaling

        return tuple(physical_size_xyu)
Justin-Krieg commented 11 months ago

Hi @cdgatenbee unfortunately i'm still having the same error persist, to try and aid with this i've uploaded a second czi file, along with the code i've implemented to run VALIS in the hopes there might be something in there. These files are located here. I'm away for the next two weeks so i wont be able to test any code.

Justin-Krieg commented 10 months ago

Hi @cdgatenbee did you have any luck with this?

cdgatenbee commented 10 months ago

Hi @kriegy98, I was able to open both images in the same Python session within a loop without error. Next will be to try aligning them and/or running the script you shared to see if I can replicate the error that way. I'll let you know how it goes.

Best, -Chandler

Justin-Krieg commented 9 months ago

Hi @cdgatenbee how did you go with this?

cdgatenbee commented 9 months ago

Hi @kriegy98, Good news :) I was able to successfully run your scripts on your images. I've made quite a few updates to the slide_io methods, so it's a little hard to say what the exact fix was, but hopefully this next update (1.0.3) will get things working for you again. I hope to have the update out in the next few days, but I'll post again when it's been pushed.

Best, -Chandler

Justin-Krieg commented 9 months ago

Hi @cdgatenbee, i tried the new release and when i specify to use the rigid_registrar, non_rigid_registrar, error_df = registrar.register(reader_cls=CziJpgxrReader) i get the same crash as before. When omitting the reader class the code runs well! But opening up the output in qupath leads to a bunch of empty channels and no channel naming - see screenshot (channels 1-5 are working, aligned but incorrectly labelled. channels 6-25 are empty and non-functional). Screenshot 2024-01-30 203733 I suspect that its perhaps trying to add all layers in the pyramid? This is just merging the two files i uploaded last year. Also just need to update your pyproject as open cv is set to an old version: opencv-contrib-python-headless = "4.6.0.66". this causes an error: AttributeError: module 'cv2' has no attribute 'xfeatures2d' which i was able to circumvent by updating openCV to the latest version: 4.9.0.80. Otherwise the new code runs faster for me, it does hang when saving the file around 99-100% but otherwise its great to use CZI files directly.

cdgatenbee commented 9 months ago

Hi @kriegy98, I was able to re-create the issue you're having, where the saved image had several extra channels. It looks like the image was being interpreted as RGB, or a stack of RGB images (even though it had 6 channels). The fix seems to be ensuring that the merged image has the "multiband" interpretation in pyvips, and that SamplesPerPixel=1 is set for each channel in the ome-xml. After these changes, the image seems to open as expected in QuPath (although the tumbnail is distorted) (see below). I also updated the opencv version in the pyproject (thanks for the suggestion). These fixes are in the most recent release (1.0.4), so hopefully this version will fix these issues you've been dealing with. Thanks again for your patience.

Best, -Chandler

image