API Reference

Zarr-based image with NIfTI compatibility using NgffImage internally.

This class provides chainable operations on OME-Zarr data while maintaining compatibility with NIfTI workflows. It uses NgffImage objects internally for better multiscale support and metadata preservation.

Attributes:

• ngff_image (NgffImage) –

  The internal NgffImage object containing data and metadata.

• axes_order (str) –

  The order of the axes for NIfTI compatibility ('ZYX' or 'XYZ').

• xyz_orientation (str) –

  The anatomical orientation string in XYZ axes order (e.g., 'RAS', 'LPI').

Constructor with backward compatibility for old signature.

Source code in zarrnii/core.py

def __init__(
    self,
    darr=None,
    affine=None,
    axes_order="ZYX",
    orientation="RAS",
    xyz_orientation=None,
    ngff_image=None,
    _omero=None,
    **kwargs,
):
    """
    Constructor with backward compatibility for old signature.
    """
    # Handle backwards compatibility: if xyz_orientation is provided, use it
    # Otherwise, use orientation for backwards compatibility
    final_orientation = (
        xyz_orientation if xyz_orientation is not None else orientation
    )

    if ngff_image is not None:
        # New signature
        object.__setattr__(self, "ngff_image", ngff_image)
        object.__setattr__(self, "axes_order", axes_order)
        object.__setattr__(self, "xyz_orientation", final_orientation)
        object.__setattr__(self, "_omero", _omero)
    elif darr is not None:
        # Legacy signature - delegate to from_darr
        instance = self.from_darr(
            darr=darr,
            affine=affine,
            axes_order=axes_order,
            orientation=final_orientation,
            **kwargs,
        )
        object.__setattr__(self, "ngff_image", instance.ngff_image)
        object.__setattr__(self, "axes_order", instance.axes_order)
        object.__setattr__(self, "xyz_orientation", instance.xyz_orientation)
        object.__setattr__(self, "_omero", instance._omero)
    else:
        raise ValueError("Must provide either ngff_image or darr")

Attributes

zarrnii.ZarrNii.data property writable

Access the image data (dask array).

zarrnii.ZarrNii.darr property writable

Legacy property name for image data.

zarrnii.ZarrNii.shape property

Shape of the image data.

zarrnii.ZarrNii.dims property

Dimension names.

zarrnii.ZarrNii.scale property

Scale information from NgffImage.

zarrnii.ZarrNii.translation property

Translation information from NgffImage.

zarrnii.ZarrNii.name property

Image name from NgffImage.

zarrnii.ZarrNii.orientation property writable

Legacy property for backward compatibility.

Returns the xyz_orientation attribute to maintain backward compatibility with code that expects the 'orientation' property.

Returns:

• str ( str ) –

  The anatomical orientation string in XYZ axes order

zarrnii.ZarrNii.affine property

Affine transformation matrix derived from NgffImage scale and translation.

Returns:

• AffineTransform ( AffineTransform ) –

  4x4 affine transformation matrix in axes order of self.

zarrnii.ZarrNii.axes property

Axes metadata - derived from NgffImage for compatibility.

zarrnii.ZarrNii.coordinate_transformations property

Coordinate transformations - derived from NgffImage scale/translation.

zarrnii.ZarrNii.omero property

Omero metadata object.

Functions

zarrnii.ZarrNii.get_affine_transform(axes_order=None)

Get AffineTransform object from NgffImage metadata.

Parameters:

• axes_order (str, default: None ) –

  Spatial axes order, defaults to self.axes_order

Returns:

• AffineTransform –

  AffineTransform object

Source code in zarrnii/core.py

def get_affine_transform(self, axes_order: str = None) -> AffineTransform:
    """
    Get AffineTransform object from NgffImage metadata.

    Args:
        axes_order: Spatial axes order, defaults to self.axes_order

    Returns:
        AffineTransform object
    """
    matrix = self.get_affine_matrix(axes_order)
    return AffineTransform.from_array(matrix)

zarrnii.ZarrNii.from_ngff_image(ngff_image, axes_order='ZYX', xyz_orientation='RAS', omero=None) classmethod

Create ZarrNii from an existing NgffImage.

Parameters:

• ngff_image (NgffImage) –

  NgffImage to wrap

• axes_order (str, default: 'ZYX' ) –

  Spatial axes order for NIfTI compatibility

• xyz_orientation (str, default: 'RAS' ) –

  Anatomical orientation string in XYZ axes order

• omero (Optional[object], default: None ) –

  Optional omero metadata object

Returns:

• 'ZarrNii' –

  ZarrNii instance

Source code in zarrnii/core.py

@classmethod
def from_ngff_image(
    cls,
    ngff_image: nz.NgffImage,
    axes_order: str = "ZYX",
    xyz_orientation: str = "RAS",
    omero: Optional[object] = None,
) -> "ZarrNii":
    """
    Create ZarrNii from an existing NgffImage.

    Args:
        ngff_image: NgffImage to wrap
        axes_order: Spatial axes order for NIfTI compatibility
        xyz_orientation: Anatomical orientation string in XYZ axes order
        omero: Optional omero metadata object

    Returns:
        ZarrNii instance
    """
    return cls(
        ngff_image=ngff_image,
        axes_order=axes_order,
        xyz_orientation=xyz_orientation,
        _omero=omero,
    )

zarrnii.ZarrNii.from_darr(darr, affine=None, axes_order='ZYX', orientation='RAS', spacing=(1.0, 1.0, 1.0), origin=(0.0, 0.0, 0.0), name='image', omero=None, **kwargs) classmethod

Create ZarrNii from dask array (legacy compatibility constructor).

Parameters:

• darr (Array) –

  Dask array containing image data

• affine (Optional[AffineTransform], default: None ) –

  Optional affine transformation

• axes_order (str, default: 'ZYX' ) –

  Spatial axes order

• orientation (str, default: 'RAS' ) –

  Anatomical orientation string

• spacing (Tuple[float, float, float], default: (1.0, 1.0, 1.0) ) –

  Voxel spacing (used if no affine provided)

• origin (Tuple[float, float, float], default: (0.0, 0.0, 0.0) ) –

  Origin offset (used if no affine provided)

• name (str, default: 'image' ) –

  Image name

• omero (Optional[object], default: None ) –

  Optional omero metadata

Returns:

• 'ZarrNii' –

  ZarrNii instance

Source code in zarrnii/core.py

@classmethod
def from_darr(
    cls,
    darr: da.Array,
    affine: Optional[AffineTransform] = None,
    axes_order: str = "ZYX",
    orientation: str = "RAS",
    spacing: Tuple[float, float, float] = (1.0, 1.0, 1.0),
    origin: Tuple[float, float, float] = (0.0, 0.0, 0.0),
    name: str = "image",
    omero: Optional[object] = None,
    **kwargs,
) -> "ZarrNii":
    """
    Create ZarrNii from dask array (legacy compatibility constructor).

    Args:
        darr: Dask array containing image data
        affine: Optional affine transformation
        axes_order: Spatial axes order
        orientation: Anatomical orientation string
        spacing: Voxel spacing (used if no affine provided)
        origin: Origin offset (used if no affine provided)
        name: Image name
        omero: Optional omero metadata

    Returns:
        ZarrNii instance
    """
    # Create scale and translation from affine if provided
    if affine is not None:
        # Extract scale and translation from affine matrix
        affine_matrix = affine.matrix
        if axes_order == "ZYX":
            scale = {
                "z": affine_matrix[0, 0],
                "y": affine_matrix[1, 1],
                "x": affine_matrix[2, 2],
            }
            translation = {
                "z": affine_matrix[0, 3],
                "y": affine_matrix[1, 3],
                "x": affine_matrix[2, 3],
            }
        else:  # XYZ
            scale = {
                "x": affine_matrix[0, 0],
                "y": affine_matrix[1, 1],
                "z": affine_matrix[2, 2],
            }
            translation = {
                "x": affine_matrix[0, 3],
                "y": affine_matrix[1, 3],
                "z": affine_matrix[2, 3],
            }
    else:
        # Use spacing and origin
        if axes_order == "ZYX":
            scale = {"z": spacing[0], "y": spacing[1], "x": spacing[2]}
            translation = {"z": origin[0], "y": origin[1], "x": origin[2]}
        else:  # XYZ
            scale = {"x": spacing[0], "y": spacing[1], "z": spacing[2]}
            translation = {"x": origin[0], "y": origin[1], "z": origin[2]}

    # Create dimensions based on data shape after dimension adjustments
    final_ndim = len(darr.shape)
    if final_ndim == 4:
        # 4D: (c, z, y, x) or (c, x, y, z) - standard case
        dims = ["c"] + list(axes_order.lower())
    elif final_ndim == 5:
        # 5D: (t, c, z, y, x) or (t, c, x, y, z) - time dimension included
        dims = ["t", "c"] + list(axes_order.lower())
    else:
        # Fallback for other cases
        dims = ["c"] + list(axes_order.lower())

    # Create NgffImage
    ngff_image = nz.NgffImage(
        data=darr, dims=dims, scale=scale, translation=translation, name=name
    )

    return cls(
        ngff_image=ngff_image,
        axes_order=axes_order,
        xyz_orientation=orientation,
        _omero=omero,
    )

zarrnii.ZarrNii.from_ome_zarr(store_or_path, level=0, channels=None, channel_labels=None, timepoints=None, storage_options=None, axes_order='ZYX', orientation=None, downsample_near_isotropic=False, chunks='auto', rechunk=False) classmethod

Load ZarrNii from OME-Zarr store with flexible options.

Creates a ZarrNii instance from an OME-Zarr store, supporting multiscale pyramids, channel/timepoint selection, and various storage backends. Automatically handles metadata extraction and format conversion.

Parameters:

• store_or_path (Union[str, Any]) –

  Store or path to OME-Zarr file. Supports: - Local file paths - Remote URLs (s3://, http://, etc.) - ZIP files (.zip extension) - Zarr store objects

• level (int, default: 0 ) –

  Pyramid level to load (0 = highest resolution). If level exceeds available levels, applies lazy downsampling

• channels (Optional[List[int]], default: None ) –

  List of channel indices to load (0-based). Mutually exclusive with channel_labels

• channel_labels (Optional[List[str]], default: None ) –

  List of channel names to load by label. Requires OMERO metadata. Mutually exclusive with channels

• timepoints (Optional[List[int]], default: None ) –

  List of timepoint indices to load (0-based). If None, loads all available timepoints

• storage_options (Optional[Dict[str, Any]], default: None ) –

  Additional options for zarr storage backend (e.g., credentials for cloud storage)

• axes_order (str, default: 'ZYX' ) –

  Spatial axis order for NIfTI compatibility. Either "ZYX" or "XYZ"

• orientation (Optional[str], default: None ) –

  Default anatomical orientation if not in metadata. Standard orientations like "RAS", "LPI", etc. This is always interpreted in XYZ axes order for consistency. This setting will override any orientation defined in the OME zarr metadata

• downsample_near_isotropic (bool, default: False ) –

  If True, automatically downsample dimensions with smaller voxel sizes to achieve near-isotropic resolution

• chunks (tuple[int, Ellipsis] | Literal['auto'], default: 'auto' ) –

  chunking strategy, or explicit chunk sizes to use if not automatic

• rechunk (bool, default: False ) –

  If True, rechunks the dataset after lazy loading, based on the chunks parameter

Returns:

• 'ZarrNii' –

  ZarrNii instance with loaded data and metadata

Raises:

• ValueError –

  If both channels and channel_labels are specified, or if invalid level/indices are provided

• FileNotFoundError –

  If store_or_path does not exist

• KeyError –

  If specified channel labels are not found

• IOError –

  If unable to read from the storage backend

Examples:

>>> # Load full resolution data
>>> znii = ZarrNii.from_ome_zarr("/path/to/data.zarr")

>>> # Load specific channels and pyramid level
>>> znii = ZarrNii.from_ome_zarr(
...     "/path/to/data.zarr",
...     level=1,
...     channels=[0, 2],
...     orientation="LPI"
... )

>>> # Load from cloud storage
>>> znii = ZarrNii.from_ome_zarr(
...     "s3://bucket/data.zarr",
...     storage_options={"key": "access_key", "secret": "secret"}
... )

Notes

Orientation Metadata Backwards Compatibility:

This method implements backwards compatibility for orientation metadata:

1. Override: Setting the orientation here will override any orientation defined in the OME Zarr metadata.

2. Zarr Metadata: Checks for 'xyz_orientation' first (new format), then falls back to 'orientation' (legacy format)

3. Legacy Fallback: When only legacy 'orientation' is found, the orientation string is automatically reversed to convert from ZYX-based encoding (legacy) to XYZ-based encoding (current standard)

4. Default Fallback: If no orientation metadata is found, uses RAS orientation as the default.

Examples of the conversion: - Legacy 'orientation'='SAR' (ZYX) → 'xyz_orientation'='RAS' (XYZ) - Legacy 'orientation'='IPL' (ZYX) → 'xyz_orientation'='LPI' (XYZ)

This ensures consistent orientation handling while maintaining backwards compatibility with existing OME-Zarr files that use the legacy format.

Source code in zarrnii/core.py

@classmethod
def from_ome_zarr(
    cls,
    store_or_path: Union[str, Any],
    level: int = 0,
    channels: Optional[List[int]] = None,
    channel_labels: Optional[List[str]] = None,
    timepoints: Optional[List[int]] = None,
    storage_options: Optional[Dict[str, Any]] = None,
    axes_order: str = "ZYX",
    orientation: Optional[str] = None,
    downsample_near_isotropic: bool = False,
    chunks: tuple[int, Ellipsis] | Literal["auto"] = "auto",
    rechunk: bool = False,
) -> "ZarrNii":
    """Load ZarrNii from OME-Zarr store with flexible options.

    Creates a ZarrNii instance from an OME-Zarr store, supporting multiscale
    pyramids, channel/timepoint selection, and various storage backends.
    Automatically handles metadata extraction and format conversion.

    Args:
        store_or_path: Store or path to OME-Zarr file. Supports:
            - Local file paths
            - Remote URLs (s3://, http://, etc.)
            - ZIP files (.zip extension)
            - Zarr store objects
        level: Pyramid level to load (0 = highest resolution). If level
            exceeds available levels, applies lazy downsampling
        channels: List of channel indices to load (0-based). Mutually
            exclusive with channel_labels
        channel_labels: List of channel names to load by label. Requires
            OMERO metadata. Mutually exclusive with channels
        timepoints: List of timepoint indices to load (0-based). If None,
            loads all available timepoints
        storage_options: Additional options for zarr storage backend
            (e.g., credentials for cloud storage)
        axes_order: Spatial axis order for NIfTI compatibility.
            Either "ZYX" or "XYZ"
        orientation: Default anatomical orientation if not in metadata.
            Standard orientations like "RAS", "LPI", etc. This is always
            interpreted in XYZ axes order for consistency. This setting will override
            any orientation defined in the OME zarr metadata
        downsample_near_isotropic: If True, automatically downsample
            dimensions with smaller voxel sizes to achieve near-isotropic
            resolution
        chunks: chunking strategy, or explicit chunk sizes to use if not automatic
        rechunk: If True, rechunks the dataset after lazy loading, based
            on the chunks parameter

    Returns:
        ZarrNii instance with loaded data and metadata

    Raises:
        ValueError: If both channels and channel_labels are specified,
            or if invalid level/indices are provided
        FileNotFoundError: If store_or_path does not exist
        KeyError: If specified channel labels are not found
        IOError: If unable to read from the storage backend

    Examples:
        >>> # Load full resolution data
        >>> znii = ZarrNii.from_ome_zarr("/path/to/data.zarr")

        >>> # Load specific channels and pyramid level
        >>> znii = ZarrNii.from_ome_zarr(
        ...     "/path/to/data.zarr",
        ...     level=1,
        ...     channels=[0, 2],
        ...     orientation="LPI"
        ... )

        >>> # Load from cloud storage
        >>> znii = ZarrNii.from_ome_zarr(
        ...     "s3://bucket/data.zarr",
        ...     storage_options={"key": "access_key", "secret": "secret"}
        ... )

    Notes:
        **Orientation Metadata Backwards Compatibility:**

        This method implements backwards compatibility for orientation metadata:

        1. **Override**: Setting the orientation here will override
           any orientation defined in the OME Zarr metadata.

        2. **Zarr Metadata**: Checks for 'xyz_orientation' first (new format),
           then falls back to 'orientation' (legacy format)

        3. **Legacy Fallback**: When only legacy 'orientation' is found, the
           orientation string is automatically reversed to convert from ZYX-based
           encoding (legacy) to XYZ-based encoding (current standard)

        4. **Default Fallback**: If no orientation metadata is found, uses RAS
           orientation as the default.

        Examples of the conversion:
        - Legacy 'orientation'='SAR' (ZYX) → 'xyz_orientation'='RAS' (XYZ)
        - Legacy 'orientation'='IPL' (ZYX) → 'xyz_orientation'='LPI' (XYZ)

        This ensures consistent orientation handling while maintaining backwards
        compatibility with existing OME-Zarr files that use the legacy format.
    """
    # Validate channel and timepoint selection arguments
    if channels is not None and channel_labels is not None:
        raise ValueError("Cannot specify both 'channels' and 'channel_labels'")

    # Load the multiscales object
    try:
        if isinstance(store_or_path, str):
            # Handle ZIP files by creating a ZipStore
            if store_or_path.endswith(".zip"):
                import zarr

                store = zarr.storage.ZipStore(store_or_path, mode="r")
                multiscales = nz.from_ngff_zarr(
                    store, storage_options=storage_options or {}
                )
                # Note: We'll close the store after extracting metadata
            else:
                multiscales = nz.from_ngff_zarr(
                    store_or_path, storage_options=storage_options or {}
                )
        else:
            multiscales = nz.from_ngff_zarr(store_or_path)
    except Exception as e:
        # Fallback for older zarr/ngff_zarr versions
        if isinstance(store_or_path, str):
            if store_or_path.endswith(".zip"):
                import zarr

                store = zarr.storage.ZipStore(store_or_path, mode="r")
                multiscales = nz.from_ngff_zarr(store)
            else:
                store = fsspec.get_mapper(store_or_path, **storage_options or {})
                multiscales = nz.from_ngff_zarr(store)
        else:
            store = store_or_path
            multiscales = nz.from_ngff_zarr(store)

    # Extract omero metadata if available
    omero_metadata = None
    try:
        import zarr

        if isinstance(store_or_path, str):
            if store_or_path.endswith(".zip"):
                zip_store = zarr.storage.ZipStore(store_or_path, mode="r")
                group = zarr.open_group(zip_store, mode="r")
                # Close zip store after getting group
                zip_store.close()
            else:
                group = zarr.open_group(store_or_path, mode="r")

        else:
            group = zarr.open_group(store_or_path, mode="r")

        if "omero" in group.attrs:
            omero_dict = group.attrs["omero"]

            # Create a simple object to hold omero metadata
            class OmeroMetadata:
                def __init__(self, omero_dict):
                    self.channels = []
                    if "channels" in omero_dict:
                        for ch_dict in omero_dict["channels"]:
                            # Create channel objects
                            class ChannelMetadata:
                                def __init__(self, ch_dict):
                                    self.label = ch_dict.get("label", "")
                                    self.color = ch_dict.get("color", "")
                                    if "window" in ch_dict:

                                        class WindowMetadata:
                                            def __init__(self, win_dict):
                                                self.min = win_dict.get("min", 0.0)
                                                self.max = win_dict.get(
                                                    "max", 65535.0
                                                )
                                                self.start = win_dict.get(
                                                    "start", 0.0
                                                )
                                                self.end = win_dict.get(
                                                    "end", 65535.0
                                                )

                                        self.window = WindowMetadata(
                                            ch_dict["window"]
                                        )
                                    else:
                                        self.window = None

                            self.channels.append(ChannelMetadata(ch_dict))

            omero_metadata = OmeroMetadata(omero_dict)
    except Exception:
        # If we can't load omero metadata, that's okay
        pass

    # Read orientation metadata with backwards compatibility support
    # Priority: xyz_orientation (new) > orientation (legacy, with reversal)
    try:
        import zarr

        if orientation is None:

            if isinstance(store_or_path, str):
                if store_or_path.endswith(".zip"):
                    zip_store = zarr.storage.ZipStore(store_or_path, mode="r")
                    group = zarr.open_group(zip_store, mode="r")
                    # Check for new xyz_orientation first, then fallback to legacy orientation
                    if "xyz_orientation" in group.attrs:
                        orientation = group.attrs["xyz_orientation"]
                    elif "orientation" in group.attrs:
                        # Legacy orientation is ZYX-based, reverse it to get XYZ-based orientation
                        legacy_orientation = group.attrs["orientation"]
                        orientation = reverse_orientation_string(legacy_orientation)
                    # If neither found, use the provided default orientation
                    zip_store.close()
                else:
                    group = zarr.open_group(store_or_path, mode="r")
                    # Check for new xyz_orientation first, then fallback to legacy orientation
                    if "xyz_orientation" in group.attrs:
                        orientation = group.attrs["xyz_orientation"]
                    elif "orientation" in group.attrs:
                        # Legacy orientation is ZYX-based, reverse it to get XYZ-based orientation
                        legacy_orientation = group.attrs["orientation"]
                        orientation = reverse_orientation_string(legacy_orientation)
                    # If neither found, use the provided default orientation
            else:
                group = zarr.open_group(store_or_path, mode="r")
                # Check for new xyz_orientation first, then fallback to legacy orientation
                if "xyz_orientation" in group.attrs:
                    orientation = group.attrs["xyz_orientation"]
                elif "orientation" in group.attrs:
                    # Legacy orientation is ZYX-based, reverse it to get XYZ-based orientation
                    legacy_orientation = group.attrs["orientation"]
                    orientation = reverse_orientation_string(legacy_orientation)
                # If neither found, use the provided default orientation

    except Exception:
        # If we can't read orientation metadata, use the provided default
        pass

    # If orientation is still None, use the fallback default
    if orientation is None:
        orientation = "RAS"
    # Determine the available pyramid levels and handle lazy downsampling
    max_level = len(multiscales.images) - 1
    actual_level = min(level, max_level)
    do_downsample = level > max_level

    # Get the highest available level
    ngff_image = multiscales.images[actual_level]

    # Handle channel and timepoint selection and filter omero metadata accordingly
    filtered_omero = omero_metadata
    if channels is not None or channel_labels is not None or timepoints is not None:
        ngff_image, filtered_omero = _select_dimensions_from_image_with_omero(
            ngff_image,
            multiscales,
            channels,
            channel_labels,
            timepoints,
            omero_metadata,
        )

    # Create ZarrNii instance with xyz_orientation
    znimg = cls(
        ngff_image=ngff_image,
        axes_order=axes_order,
        xyz_orientation=orientation,
        _omero=filtered_omero,
    )

    # Apply lazy downsampling if needed
    if do_downsample:
        level_ds = level - max_level
        downsample_factor = 2**level_ds

        # Get spatial dims based on axes order
        spatial_dims = ["z", "y", "x"] if axes_order == "ZYX" else ["x", "y", "z"]

        # Apply downsampling using the existing method
        znimg = znimg.downsample(
            factors=downsample_factor, spatial_dims=spatial_dims
        )

    # Apply near-isotropic downsampling if requested
    if downsample_near_isotropic:
        znimg = _apply_near_isotropic_downsampling(znimg, axes_order)

    if rechunk:
        znimg.data = znimg.data.rechunk(chunks)

    return znimg

zarrnii.ZarrNii.from_nifti(path, chunks='auto', axes_order='XYZ', name=None, as_ref=False, zooms=None) classmethod

Load ZarrNii from NIfTI file with flexible loading options.

Creates a ZarrNii instance from a NIfTI file, automatically converting the data to dask arrays and extracting spatial transformation information. Supports both full data loading and reference-only loading for memory efficiency.

Parameters:

• path (Union[str, bytes]) –

  File path to NIfTI file (.nii, .nii.gz, .img/.hdr)

• chunks (Union[str, Tuple[int, ...]], default: 'auto' ) –

  Dask array chunking strategy. Can be: - "auto": Automatic chunking based on file size - Tuple of ints: Manual chunk sizes for each dimension - Dict mapping axis to chunk size

• axes_order (str, default: 'XYZ' ) –

  Spatial axis ordering convention. Either: - "XYZ": X=left-right, Y=anterior-posterior, Z=inferior-superior - "ZYX": Z=inferior-superior, Y=anterior-posterior, X=left-right

• name (Optional[str], default: None ) –

  Optional name for the resulting NgffImage. If None, uses filename without extension

• as_ref (bool, default: False ) –

  If True, creates empty dask array with correct shape/metadata without loading actual image data (memory efficient for templates)

• zooms (Optional[Tuple[float, float, float]], default: None ) –

  Target voxel spacing as (x, y, z) in mm. Only valid when as_ref=True. Adjusts shape and affine accordingly

Returns:

• 'ZarrNii' –

  ZarrNii instance containing NIfTI data and spatial metadata

Raises:

• ValueError –

  If zooms specified with as_ref=False, or invalid axes_order

• FileNotFoundError –

  If NIfTI file does not exist

• OSError –

  If unable to read NIfTI file

• ImageFileError –

  If file is not valid NIfTI

Examples:

>>> # Load full NIfTI data
>>> znii = ZarrNii.from_nifti("/path/to/brain.nii.gz")

>>> # Load with custom chunking and axis order
>>> znii = ZarrNii.from_nifti(
...     "/path/to/data.nii",
...     chunks=(64, 64, 64),
...     axes_order="ZYX"
... )

>>> # Create reference with target resolution
>>> znii_ref = ZarrNii.from_nifti(
...     "/path/to/template.nii.gz",
...     as_ref=True,
...     zooms=(2.0, 2.0, 2.0)
... )

Notes

The method automatically handles NIfTI orientation codes and converts them to the specified axes_order for consistency with OME-Zarr workflows.

Source code in zarrnii/core.py

@classmethod
def from_nifti(
    cls,
    path: Union[str, bytes],
    chunks: Union[str, Tuple[int, ...]] = "auto",
    axes_order: str = "XYZ",
    name: Optional[str] = None,
    as_ref: bool = False,
    zooms: Optional[Tuple[float, float, float]] = None,
) -> "ZarrNii":
    """Load ZarrNii from NIfTI file with flexible loading options.

    Creates a ZarrNii instance from a NIfTI file, automatically converting
    the data to dask arrays and extracting spatial transformation information.
    Supports both full data loading and reference-only loading for memory
    efficiency.

    Args:
        path: File path to NIfTI file (.nii, .nii.gz, .img/.hdr)
        chunks: Dask array chunking strategy. Can be:
            - "auto": Automatic chunking based on file size
            - Tuple of ints: Manual chunk sizes for each dimension
            - Dict mapping axis to chunk size
        axes_order: Spatial axis ordering convention. Either:
            - "XYZ": X=left-right, Y=anterior-posterior, Z=inferior-superior
            - "ZYX": Z=inferior-superior, Y=anterior-posterior, X=left-right
        name: Optional name for the resulting NgffImage. If None,
            uses filename without extension
        as_ref: If True, creates empty dask array with correct shape/metadata
            without loading actual image data (memory efficient for templates)
        zooms: Target voxel spacing as (x, y, z) in mm. Only valid when
            as_ref=True. Adjusts shape and affine accordingly

    Returns:
        ZarrNii instance containing NIfTI data and spatial metadata

    Raises:
        ValueError: If zooms specified with as_ref=False, or invalid axes_order
        FileNotFoundError: If NIfTI file does not exist
        OSError: If unable to read NIfTI file
        nibabel.filebasedimages.ImageFileError: If file is not valid NIfTI

    Examples:
        >>> # Load full NIfTI data
        >>> znii = ZarrNii.from_nifti("/path/to/brain.nii.gz")

        >>> # Load with custom chunking and axis order
        >>> znii = ZarrNii.from_nifti(
        ...     "/path/to/data.nii",
        ...     chunks=(64, 64, 64),
        ...     axes_order="ZYX"
        ... )

        >>> # Create reference with target resolution
        >>> znii_ref = ZarrNii.from_nifti(
        ...     "/path/to/template.nii.gz",
        ...     as_ref=True,
        ...     zooms=(2.0, 2.0, 2.0)
        ... )

    Notes:
        The method automatically handles NIfTI orientation codes and converts
        them to the specified axes_order for consistency with OME-Zarr workflows.
    """
    if not as_ref and zooms is not None:
        raise ValueError("`zooms` can only be used when `as_ref=True`.")

    # Load NIfTI file
    nifti_img = nib.load(path)
    shape = nifti_img.header.get_data_shape()
    affine_matrix = nifti_img.affine.copy()

    # infer orientation from the affine
    orientation = _affine_to_orientation(affine_matrix)

    # Adjust shape and affine if zooms are provided
    if zooms is not None:
        in_zooms = np.sqrt(
            (affine_matrix[:3, :3] ** 2).sum(axis=0)
        )  # Current voxel spacing
        scaling_factor = in_zooms / zooms
        new_shape = [
            int(np.floor(shape[0] * scaling_factor[2])),  # Z
            int(np.floor(shape[1] * scaling_factor[1])),  # Y
            int(np.floor(shape[2] * scaling_factor[0])),  # X
        ]
        np.fill_diagonal(affine_matrix[:3, :3], zooms)
    else:
        new_shape = shape

    if as_ref:
        # Create an empty dask array with the adjusted shape
        darr = da.empty((1, *new_shape), chunks=chunks, dtype="float32")
    else:
        # Load the NIfTI data and convert to a dask array
        array = nifti_img.get_fdata()
        darr = da.from_array(array, chunks=chunks)

    # Add channel and time dimensions if not present
    original_ndim = len(darr.shape)

    if original_ndim == 3:
        # 3D data: add channel dimension -> (c, z, y, x) or (c, x, y, z)
        darr = darr[np.newaxis, ...]
    elif original_ndim == 4:
        # 4D data: could be (c, z, y, x) or (t, z, y, x) - assume channel by default
        # User can specify if it's time by using appropriate axes_order
        pass  # Keep as is - 4D is already handled
    elif original_ndim == 5:
        # 5D data: assume (t, z, y, x, c) and handle appropriately
        pass  # Keep as is - 5D is already the target format
    else:
        # For 1D, 2D, or >5D data, add channel dimension and let user handle
        darr = darr[np.newaxis, ...]

    # Create dimensions based on data shape after dimension adjustments
    final_ndim = len(darr.shape)
    if final_ndim == 4:
        # 4D: (c, z, y, x) or (c, x, y, z) - standard case
        dims = ["c"] + list(axes_order.lower())
    elif final_ndim == 5:
        # 5D: (t, c, z, y, x) or (t, c, x, y, z) - time dimension included
        dims = ["t", "c"] + list(axes_order.lower())
    else:
        # Fallback for other cases
        dims = ["c"] + list(axes_order.lower())

    # Extract scale and translation from affine
    scale = {}
    translation = {}
    spatial_dims = ["z", "y", "x"] if axes_order == "ZYX" else ["x", "y", "z"]

    for i, dim in enumerate(spatial_dims):
        scale[dim] = np.sqrt((affine_matrix[i, :3] ** 2).sum())
        translation[dim] = affine_matrix[i, 3]

    # Create NgffImage
    if name is None:
        name = f"nifti_image_{path}"

    ngff_image = nz.NgffImage(
        data=darr, dims=dims, scale=scale, translation=translation, name=name
    )

    return cls(
        ngff_image=ngff_image, axes_order=axes_order, xyz_orientation=orientation
    )

zarrnii.ZarrNii.crop(bbox_min, bbox_max=None, spatial_dims=None, physical_coords=False)

Extract a spatial region or multiple regions from the image.

Crops the image to the specified bounding box coordinates, preserving all metadata and non-spatial dimensions (channels, time). The cropping is performed in voxel coordinates by default, or physical coordinates if specified. Can crop a single region or multiple regions at once.

Parameters:

• bbox_min (Union[Tuple[float, ...], List[Tuple[Tuple[float, ...], Tuple[float, ...]]]]) –

  Either: - Minimum corner coordinates of bounding box as tuple (when bbox_max is provided). Length should match number of spatial dimensions (x, y, z order) - List of (bbox_min, bbox_max) tuples for batch cropping (when bbox_max is None)

• bbox_max (Optional[Tuple[float, ...]], default: None ) –

  Maximum corner coordinates of bounding box as tuple. Length should match number of spatial dimensions (x, y, z order). Should be None when bbox_min is a list of bounding boxes.

• spatial_dims (Optional[List[str]], default: None ) –

  Names of spatial dimensions to crop. If None, automatically derived from axes_order ("z","y","x" for ZYX or "x","y","z" for XYZ)

• physical_coords (bool, default: False ) –

  If True, bbox_min and bbox_max are in physical/world coordinates (mm). If False, they are in voxel coordinates. Default is False.

Returns:

• Union['ZarrNii', List['ZarrNii']] –

  New ZarrNii instance with cropped data (single crop) or list of

• Union['ZarrNii', List['ZarrNii']] –

  ZarrNii instances (batch crop) with updated spatial metadata

Raises:

• ValueError –

  If bbox coordinates are invalid or out of bounds, or if both list and bbox_max are provided

• IndexError –

  If bbox dimensions don't match spatial dimensions

Examples:

>>> # Crop 3D region (voxel coordinates)
>>> cropped = znii.crop((10, 20, 30), (110, 120, 130))

>>> # Crop with physical coordinates
>>> cropped = znii.crop((10.5, 20.5, 30.5), (110.5, 120.5, 130.5),
...                      physical_coords=True)

>>> # Crop with explicit spatial dimensions
>>> cropped = znii.crop(
...     (50, 60, 70), (150, 160, 170),
...     spatial_dims=["x", "y", "z"]
... )

>>> # Batch crop multiple regions
>>> bboxes = [
...     ((10, 20, 30), (60, 70, 80)),
...     ((100, 110, 120), (150, 160, 170))
... ]
>>> cropped_list = znii.crop(bboxes, physical_coords=True)

Notes

• Coordinates are in voxel space (0-based indexing) by default
• Physical coordinates are in RAS orientation (Right-Anterior-Superior)
• The cropped region includes bbox_min but excludes bbox_max
• All non-spatial dimensions (channels, time) are preserved
• Spatial transformations are automatically updated
• When batch cropping, all patches share the same spatial_dims and physical_coords settings

Source code in zarrnii/core.py

def crop(
    self,
    bbox_min: Union[
        Tuple[float, ...], List[Tuple[Tuple[float, ...], Tuple[float, ...]]]
    ],
    bbox_max: Optional[Tuple[float, ...]] = None,
    spatial_dims: Optional[List[str]] = None,
    physical_coords: bool = False,
) -> Union["ZarrNii", List["ZarrNii"]]:
    """Extract a spatial region or multiple regions from the image.

    Crops the image to the specified bounding box coordinates, preserving
    all metadata and non-spatial dimensions (channels, time). The cropping
    is performed in voxel coordinates by default, or physical coordinates
    if specified. Can crop a single region or multiple regions at once.

    Args:
        bbox_min: Either:
            - Minimum corner coordinates of bounding box as tuple
              (when bbox_max is provided). Length should match number of
              spatial dimensions (x, y, z order)
            - List of (bbox_min, bbox_max) tuples for batch cropping
              (when bbox_max is None)
        bbox_max: Maximum corner coordinates of bounding box as tuple.
            Length should match number of spatial dimensions (x, y, z order).
            Should be None when bbox_min is a list of bounding boxes.
        spatial_dims: Names of spatial dimensions to crop. If None,
            automatically derived from axes_order ("z","y","x" for ZYX
            or "x","y","z" for XYZ)
        physical_coords: If True, bbox_min and bbox_max are in physical/world
            coordinates (mm). If False, they are in voxel coordinates.
            Default is False.

    Returns:
        New ZarrNii instance with cropped data (single crop) or list of
        ZarrNii instances (batch crop) with updated spatial metadata

    Raises:
        ValueError: If bbox coordinates are invalid or out of bounds, or
            if both list and bbox_max are provided
        IndexError: If bbox dimensions don't match spatial dimensions

    Examples:
        >>> # Crop 3D region (voxel coordinates)
        >>> cropped = znii.crop((10, 20, 30), (110, 120, 130))

        >>> # Crop with physical coordinates
        >>> cropped = znii.crop((10.5, 20.5, 30.5), (110.5, 120.5, 130.5),
        ...                      physical_coords=True)

        >>> # Crop with explicit spatial dimensions
        >>> cropped = znii.crop(
        ...     (50, 60, 70), (150, 160, 170),
        ...     spatial_dims=["x", "y", "z"]
        ... )

        >>> # Batch crop multiple regions
        >>> bboxes = [
        ...     ((10, 20, 30), (60, 70, 80)),
        ...     ((100, 110, 120), (150, 160, 170))
        ... ]
        >>> cropped_list = znii.crop(bboxes, physical_coords=True)

    Notes:
        - Coordinates are in voxel space (0-based indexing) by default
        - Physical coordinates are in RAS orientation (Right-Anterior-Superior)
        - The cropped region includes bbox_min but excludes bbox_max
        - All non-spatial dimensions (channels, time) are preserved
        - Spatial transformations are automatically updated
        - When batch cropping, all patches share the same spatial_dims and
          physical_coords settings
    """
    # Check if this is batch cropping (list of bounding boxes)
    # A batch crop is a list of (bbox_min, bbox_max) tuples
    # Each element should be a tuple/list of two elements
    is_batch_crop = (
        isinstance(bbox_min, list)
        and len(bbox_min) > 0
        and isinstance(bbox_min[0], (tuple, list))
        and len(bbox_min[0]) == 2
    )

    if is_batch_crop:
        if bbox_max is not None:
            raise ValueError(
                "bbox_max should be None when bbox_min is a list of bounding boxes"
            )
        # Batch crop: recursively call crop for each bounding box
        return [
            self.crop(bmin, bmax, spatial_dims, physical_coords)
            for bmin, bmax in bbox_min
        ]

    # Single crop: original implementation
    if bbox_max is None:
        raise ValueError("bbox_max is required when bbox_min is not a list")

    if spatial_dims is None:
        spatial_dims = (
            ["z", "y", "x"] if self.axes_order == "ZYX" else ["x", "y", "z"]
        )

    # Convert physical coordinates to voxel coordinates if needed
    if physical_coords:
        # Physical coords are always in (x, y, z) order
        # Convert to homogeneous coordinates
        phys_min = np.array(list(bbox_min) + [1.0])
        phys_max = np.array(list(bbox_max) + [1.0])

        # Get inverse affine to convert from physical to voxel
        affine_inv = np.linalg.inv(
            self.get_affine_matrix(axes_order="XYZ")
        )  # TODO: should this always be xyz affine??

        # Transform to voxel coordinates
        voxel_min = affine_inv @ phys_min
        voxel_max = affine_inv @ phys_max

        # Extract voxel coordinates (x, y, z)
        voxel_min_xyz = voxel_min[:3]
        voxel_max_xyz = voxel_max[:3]

        # Round to nearest integer voxel indices
        voxel_min_xyz = np.round(voxel_min_xyz).astype(int)
        voxel_max_xyz = np.round(voxel_max_xyz).astype(int)

        # Ensure max >= min
        voxel_min_xyz = np.minimum(voxel_min_xyz, voxel_max_xyz)
        voxel_max_xyz = np.maximum(
            np.round(affine_inv @ phys_min).astype(int)[:3],
            np.round(affine_inv @ phys_max).astype(int)[:3],
        )

        # Create mapping from x,y,z to voxel coordinates
        bbox_min = voxel_min_xyz
        bbox_max = voxel_max_xyz

    # Create mapping from x,y,z to voxel coordinates
    bbox_vox_min = {
        "x": bbox_min[0],
        "y": bbox_min[1],
        "z": bbox_min[2],
    }
    bbox_vox_max = {
        "x": bbox_max[0],
        "y": bbox_max[1],
        "z": bbox_max[2],
    }

    dim_flips = _axcodes2flips(self.orientation)
    cropped_image = crop_ngff_image(
        self.ngff_image, bbox_vox_min, bbox_vox_max, dim_flips
    )
    return ZarrNii(
        ngff_image=cropped_image,
        axes_order=self.axes_order,
        xyz_orientation=self.xyz_orientation,
        _omero=self._omero,
    )

zarrnii.ZarrNii.crop_with_bounding_box(bbox_min, bbox_max, ras_coords=False)

Legacy method name for crop.

Parameters:

• bbox_min –

  Minimum corner coordinates

• bbox_max –

  Maximum corner coordinates

• ras_coords –

  If True, coordinates are in RAS physical space (deprecated, use physical_coords parameter of crop() instead)

Source code in zarrnii/core.py

def crop_with_bounding_box(self, bbox_min, bbox_max, ras_coords=False):
    """Legacy method name for crop.

    Args:
        bbox_min: Minimum corner coordinates
        bbox_max: Maximum corner coordinates
        ras_coords: If True, coordinates are in RAS physical space (deprecated,
            use physical_coords parameter of crop() instead)
    """
    return self.crop(bbox_min, bbox_max, physical_coords=ras_coords)

zarrnii.ZarrNii.crop_centered(centers, patch_size, spatial_dims=None, fill_value=0.0)

Extract fixed-size patches centered at specified coordinates.

Crops the image to extract patches of a fixed size (in voxels) centered at the given physical coordinates. This is particularly useful for machine learning workflows where training patches must have consistent dimensions. The method can process a single center or multiple centers at once.

Patches that extend beyond image boundaries are padded with the fill_value to ensure all patches have exactly the requested size.

Parameters:

• centers (Union[Tuple[float, float, float], List[Tuple[float, float, float]]]) –

  Either: - Single center coordinate as (x, y, z) tuple in physical space (mm) - List of center coordinates for batch processing

• patch_size (Tuple[int, int, int]) –

  Size of the patch in voxels as (x, y, z) tuple. This defines the dimensions of each cropped region in voxel space. All returned patches will have exactly this size.

• spatial_dims (Optional[List[str]], default: None ) –

  Names of spatial dimensions to crop. If None, automatically derived from axes_order ("z","y","x" for ZYX or "x","y","z" for XYZ). Default is None.

• fill_value (float, default: 0.0 ) –

  Value to use for padding when patches extend beyond image boundaries. Default is 0.0.

Returns:

• Union['ZarrNii', List['ZarrNii']] –

  Single ZarrNii instance (when centers is a single tuple) or list of

• Union['ZarrNii', List['ZarrNii']] –

  ZarrNii instances (when centers is a list) with cropped data and

• Union['ZarrNii', List['ZarrNii']] –

  updated spatial metadata. All patches will have exactly the shape

• Union['ZarrNii', List['ZarrNii']] –

  specified by patch_size (plus any non-spatial dimensions).

Raises:

• ValueError –

  If coordinates/dimensions are invalid

• IndexError –

  If patch_size dimensions don't match spatial dimensions

Examples:

>>> # Extract single 256x256x256 voxel patch at a coordinate
>>> center = (50.0, 60.0, 70.0)  # physical coordinates in mm
>>> patch = znii.crop_centered(center, patch_size=(256, 256, 256))
>>>
>>> # Extract multiple patches for ML training
>>> centers = [
...     (50.0, 60.0, 70.0),
...     (100.0, 110.0, 120.0),
...     (150.0, 160.0, 170.0)
... ]
>>> patches = znii.crop_centered(centers, patch_size=(128, 128, 128))
>>> # Returns list of 3 ZarrNii instances, all with shape (1, 128, 128, 128)
>>>
>>> # Use with atlas sampling for ML training workflow
>>> centers = atlas.sample_region_patches(
...     n_patches=100,
...     region_ids="cortex",
...     seed=42
... )
>>> patches = image.crop_centered(centers, patch_size=(256, 256, 256))
>>>
>>> # Use custom fill value for padding
>>> patch = znii.crop_centered(center, patch_size=(256, 256, 256), fill_value=-1.0)

Notes

• Centers are in physical/world coordinates (mm), always in (x, y, z) order
• patch_size is in voxels, in (x, y, z) order
• The patch is centered at the given coordinate, extending ±patch_size/2
• If patch_size is odd, the center voxel is included
• Patches near boundaries are padded with fill_value to maintain size
• All patches are guaranteed to have exactly the requested size
• Useful for ML training where fixed patch sizes are required
• Coordinates from atlas.sample_region_patches() can be used directly

Source code in zarrnii/core.py

def crop_centered(
    self,
    centers: Union[Tuple[float, float, float], List[Tuple[float, float, float]]],
    patch_size: Tuple[int, int, int],
    spatial_dims: Optional[List[str]] = None,
    fill_value: float = 0.0,
) -> Union["ZarrNii", List["ZarrNii"]]:
    """Extract fixed-size patches centered at specified coordinates.

    Crops the image to extract patches of a fixed size (in voxels) centered
    at the given physical coordinates. This is particularly useful for machine
    learning workflows where training patches must have consistent dimensions.
    The method can process a single center or multiple centers at once.

    Patches that extend beyond image boundaries are padded with the fill_value
    to ensure all patches have exactly the requested size.

    Args:
        centers: Either:
            - Single center coordinate as (x, y, z) tuple in physical space (mm)
            - List of center coordinates for batch processing
        patch_size: Size of the patch in voxels as (x, y, z) tuple.
            This defines the dimensions of each cropped region in voxel space.
            All returned patches will have exactly this size.
        spatial_dims: Names of spatial dimensions to crop. If None,
            automatically derived from axes_order ("z","y","x" for ZYX
            or "x","y","z" for XYZ). Default is None.
        fill_value: Value to use for padding when patches extend beyond
            image boundaries. Default is 0.0.

    Returns:
        Single ZarrNii instance (when centers is a single tuple) or list of
        ZarrNii instances (when centers is a list) with cropped data and
        updated spatial metadata. All patches will have exactly the shape
        specified by patch_size (plus any non-spatial dimensions).

    Raises:
        ValueError: If coordinates/dimensions are invalid
        IndexError: If patch_size dimensions don't match spatial dimensions

    Examples:
        >>> # Extract single 256x256x256 voxel patch at a coordinate
        >>> center = (50.0, 60.0, 70.0)  # physical coordinates in mm
        >>> patch = znii.crop_centered(center, patch_size=(256, 256, 256))
        >>>
        >>> # Extract multiple patches for ML training
        >>> centers = [
        ...     (50.0, 60.0, 70.0),
        ...     (100.0, 110.0, 120.0),
        ...     (150.0, 160.0, 170.0)
        ... ]
        >>> patches = znii.crop_centered(centers, patch_size=(128, 128, 128))
        >>> # Returns list of 3 ZarrNii instances, all with shape (1, 128, 128, 128)
        >>>
        >>> # Use with atlas sampling for ML training workflow
        >>> centers = atlas.sample_region_patches(
        ...     n_patches=100,
        ...     region_ids="cortex",
        ...     seed=42
        ... )
        >>> patches = image.crop_centered(centers, patch_size=(256, 256, 256))
        >>>
        >>> # Use custom fill value for padding
        >>> patch = znii.crop_centered(center, patch_size=(256, 256, 256), fill_value=-1.0)

    Notes:
        - Centers are in physical/world coordinates (mm), always in (x, y, z) order
        - patch_size is in voxels, in (x, y, z) order
        - The patch is centered at the given coordinate, extending ±patch_size/2
        - If patch_size is odd, the center voxel is included
        - Patches near boundaries are padded with fill_value to maintain size
        - All patches are guaranteed to have exactly the requested size
        - Useful for ML training where fixed patch sizes are required
        - Coordinates from atlas.sample_region_patches() can be used directly
    """
    # Check if this is batch processing (list of centers)
    is_batch = isinstance(centers, list)

    if is_batch:
        # Batch processing: recursively call crop_centered for each center
        return [
            self.crop_centered(center, patch_size, spatial_dims, fill_value)
            for center in centers
        ]

    # Single center processing
    if spatial_dims is None:
        spatial_dims = (
            ["z", "y", "x"] if self.axes_order == "ZYX" else ["x", "y", "z"]
        )

    # Convert center from physical to voxel coordinates
    # Centers are always in (x, y, z) order
    center_phys = np.array(list(centers) + [1.0])

    # Get inverse affine to convert from physical to voxel
    affine_inv = np.linalg.inv(self.get_affine_matrix(axes_order="XYZ"))

    # Transform to voxel coordinates
    center_voxel = affine_inv @ center_phys
    center_voxel_xyz = center_voxel[:3]

    # patch_size is in voxels, in (x, y, z) order
    patch_size_np = np.array(patch_size)
    half_patch = patch_size_np / 2.0

    # Calculate desired bounding box in voxel coordinates (may extend beyond image)
    voxel_min_xyz = center_voxel_xyz - half_patch
    voxel_max_xyz = center_voxel_xyz + half_patch

    # Round to nearest integer voxel indices
    voxel_min_xyz = np.round(voxel_min_xyz).astype(int)
    voxel_max_xyz = np.round(voxel_max_xyz).astype(int)

    # Ensure we get exactly the requested patch size
    # Adjust max to ensure patch_size is respected
    voxel_max_xyz = voxel_min_xyz + patch_size_np

    # Get image dimensions in voxel space
    # Map spatial dims to their indices
    spatial_dim_indices = {}
    for i, dim in enumerate(self.ngff_image.dims):
        if dim.lower() in [d.lower() for d in spatial_dims]:
            spatial_dim_indices[dim.lower()] = i

    image_shape_xyz = np.array(
        [
            self.ngff_image.data.shape[spatial_dim_indices["x"]],
            self.ngff_image.data.shape[spatial_dim_indices["y"]],
            self.ngff_image.data.shape[spatial_dim_indices["z"]],
        ]
    )

    # Calculate the actual crop region (clipped to image bounds)
    crop_min_xyz = np.maximum(voxel_min_xyz, 0)
    crop_max_xyz = np.minimum(voxel_max_xyz, image_shape_xyz)

    # Ensure crop_max >= crop_min to avoid empty arrays
    crop_max_xyz = np.maximum(crop_min_xyz, crop_max_xyz)

    # Calculate padding needed on each side
    pad_before_xyz = crop_min_xyz - voxel_min_xyz  # How much we're clipped at start
    pad_after_xyz = voxel_max_xyz - crop_max_xyz  # How much we're clipped at end

    # Check if the entire patch is outside the image bounds
    # This happens when crop_min >= crop_max in any dimension after clipping
    is_completely_outside = np.any(crop_min_xyz >= crop_max_xyz)

    if is_completely_outside:
        # The entire patch is outside the image bounds
        # Create a completely padded array with the fill value
        import dask.array as da

        # Build the full patch shape
        full_shape = []
        spatial_idx = 0
        for dim in self.ngff_image.dims:
            if dim.lower() in [d.lower() for d in spatial_dims]:
                full_shape.append(patch_size_np[spatial_idx])
                spatial_idx += 1
            else:
                # Non-spatial dimension - keep original size
                dim_idx = self.ngff_image.dims.index(dim)
                full_shape.append(self.ngff_image.data.shape[dim_idx])

        # Create array filled with fill_value
        padded_data = da.full(
            tuple(full_shape),
            fill_value,
            dtype=self.ngff_image.data.dtype,
            chunks=self.ngff_image.data.chunksize,
        )

        # Calculate translation for the patch center
        # The translation should be at voxel_min_xyz (the desired start of patch)
        new_translation = {}
        for dim in self.ngff_image.dims:
            if dim.lower() in [d.lower() for d in spatial_dims]:
                dim_lower = dim.lower()
                if dim_lower == "x":
                    voxel_start = voxel_min_xyz[0]
                elif dim_lower == "y":
                    voxel_start = voxel_min_xyz[1]
                elif dim_lower == "z":
                    voxel_start = voxel_min_xyz[2]
                else:
                    voxel_start = 0

                # Translation is voxel_start * scale + original translation
                new_translation[dim] = voxel_start * self.ngff_image.scale.get(
                    dim, 1.0
                ) + self.ngff_image.translation.get(dim, 0.0)
            elif dim in self.ngff_image.translation:
                new_translation[dim] = self.ngff_image.translation[dim]

        # Create NgffImage with the padded data
        padded_image = nz.NgffImage(
            data=padded_data,
            dims=self.ngff_image.dims,
            scale=self.ngff_image.scale.copy(),
            translation=new_translation,
            name=self.ngff_image.name,
        )

        return ZarrNii(
            ngff_image=padded_image,
            axes_order=self.axes_order,
            xyz_orientation=self.xyz_orientation,
            _omero=self._omero,
        )

    # Create mapping from x,y,z to voxel coordinates for cropping
    bbox_vox_min = {
        "x": crop_min_xyz[0],
        "y": crop_min_xyz[1],
        "z": crop_min_xyz[2],
    }
    bbox_vox_max = {
        "x": crop_max_xyz[0],
        "y": crop_max_xyz[1],
        "z": crop_max_xyz[2],
    }

    dim_flips = _axcodes2flips(self.orientation)
    # Crop the actual image data that exists
    cropped_image = crop_ngff_image(
        self.ngff_image, bbox_vox_min, bbox_vox_max, dim_flips
    )

    # Check if padding is needed
    needs_padding = np.any(pad_before_xyz > 0) or np.any(pad_after_xyz > 0)

    if needs_padding:
        # Build padding specification for all dimensions
        pad_width = []
        spatial_idx = 0

        for dim in cropped_image.dims:
            if dim.lower() in [d.lower() for d in spatial_dims]:
                # Spatial dimension - may need padding
                dim_lower = dim.lower()
                if dim_lower == "x":
                    pad_width.append((pad_before_xyz[0], pad_after_xyz[0]))
                elif dim_lower == "y":
                    pad_width.append((pad_before_xyz[1], pad_after_xyz[1]))
                elif dim_lower == "z":
                    pad_width.append((pad_before_xyz[2], pad_after_xyz[2]))
                spatial_idx += 1
            else:
                # Non-spatial dimension - no padding
                pad_width.append((0, 0))

        # Apply padding
        import dask.array as da

        padded_data = da.pad(
            cropped_image.data,
            pad_width=pad_width,
            mode="constant",
            constant_values=fill_value,
        )

        # Adjust translation for the padding
        new_translation = cropped_image.translation.copy()

        for i, dim in enumerate(bbox_vox_min.keys()):
            new_translation[dim] = new_translation[dim] + dim_flips[
                dim
            ] * pad_before_xyz[i] * cropped_image.scale.get(dim, 1.0)

        # Create padded NgffImage
        cropped_image = nz.NgffImage(
            data=padded_data,
            dims=cropped_image.dims,
            scale=cropped_image.scale,
            translation=new_translation,
            name=cropped_image.name,
        )

    return ZarrNii(
        ngff_image=cropped_image,
        axes_order=self.axes_order,
        xyz_orientation=self.xyz_orientation,
        _omero=self._omero,
    )

zarrnii.ZarrNii.downsample(factors=None, along_x=1, along_y=1, along_z=1, level=None, spatial_dims=None)

Reduce image resolution by downsampling.

Performs spatial downsampling by averaging blocks of voxels, effectively reducing image resolution and size. Multiple parameter options provide flexibility for different downsampling strategies.

Parameters:

• factors (Optional[Union[int, List[int]]], default: None ) –

  Downsampling factors for spatial dimensions. Can be: - int: Same factor applied to all spatial dimensions - List[int]: Per-dimension factors matching spatial_dims order - None: Use other parameters to determine factors

• along_x (int, default: 1 ) –

  Downsampling factor for X dimension (legacy parameter)

• along_y (int, default: 1 ) –

  Downsampling factor for Y dimension (legacy parameter)

• along_z (int, default: 1 ) –

  Downsampling factor for Z dimension (legacy parameter)

• level (Optional[int], default: None ) –

  Power-of-2 downsampling level (factors = 2^level). Takes precedence over along_* parameters

• spatial_dims (Optional[List[str]], default: None ) –

  Names of spatial dimensions. If None, derived from axes_order

Returns:

• 'ZarrNii' –

  New ZarrNii instance with downsampled data and updated metadata

Raises:

• ValueError –

  If conflicting parameters provided or invalid factors

Examples:

>>> # Isotropic downsampling by factor of 2
>>> downsampled = znii.downsample(factors=2)

>>> # Anisotropic downsampling
>>> downsampled = znii.downsample(factors=[1, 2, 2])

>>> # Using legacy parameters
>>> downsampled = znii.downsample(along_x=2, along_y=2, along_z=1)

>>> # Power-of-2 downsampling
>>> downsampled = znii.downsample(level=2)  # factors = 4

Notes

• Downsampling uses block averaging for anti-aliasing
• Spatial transformations are automatically scaled
• Non-spatial dimensions (channels, time) are preserved
• Original data remains unchanged (creates new instance)

Source code in zarrnii/core.py

def downsample(
    self,
    factors: Optional[Union[int, List[int]]] = None,
    along_x: int = 1,
    along_y: int = 1,
    along_z: int = 1,
    level: Optional[int] = None,
    spatial_dims: Optional[List[str]] = None,
) -> "ZarrNii":
    """Reduce image resolution by downsampling.

    Performs spatial downsampling by averaging blocks of voxels, effectively
    reducing image resolution and size. Multiple parameter options provide
    flexibility for different downsampling strategies.

    Args:
        factors: Downsampling factors for spatial dimensions. Can be:
            - int: Same factor applied to all spatial dimensions
            - List[int]: Per-dimension factors matching spatial_dims order
            - None: Use other parameters to determine factors
        along_x: Downsampling factor for X dimension (legacy parameter)
        along_y: Downsampling factor for Y dimension (legacy parameter)
        along_z: Downsampling factor for Z dimension (legacy parameter)
        level: Power-of-2 downsampling level (factors = 2^level).
            Takes precedence over along_* parameters
        spatial_dims: Names of spatial dimensions. If None, derived
            from axes_order

    Returns:
        New ZarrNii instance with downsampled data and updated metadata

    Raises:
        ValueError: If conflicting parameters provided or invalid factors

    Examples:
        >>> # Isotropic downsampling by factor of 2
        >>> downsampled = znii.downsample(factors=2)

        >>> # Anisotropic downsampling
        >>> downsampled = znii.downsample(factors=[1, 2, 2])

        >>> # Using legacy parameters
        >>> downsampled = znii.downsample(along_x=2, along_y=2, along_z=1)

        >>> # Power-of-2 downsampling
        >>> downsampled = znii.downsample(level=2)  # factors = 4

    Notes:
        - Downsampling uses block averaging for anti-aliasing
        - Spatial transformations are automatically scaled
        - Non-spatial dimensions (channels, time) are preserved
        - Original data remains unchanged (creates new instance)
    """
    # Handle legacy parameters
    if factors is None:
        if level is not None:
            factors = 2**level
        else:
            factors = (
                [along_z, along_y, along_x]
                if self.axes_order == "ZYX"
                else [along_x, along_y, along_z]
            )

    if spatial_dims is None:
        spatial_dims = (
            ["z", "y", "x"] if self.axes_order == "ZYX" else ["x", "y", "z"]
        )

    downsampled_image = downsample_ngff_image(
        self.ngff_image, factors, spatial_dims
    )
    return ZarrNii(
        ngff_image=downsampled_image,
        axes_order=self.axes_order,
        xyz_orientation=self.xyz_orientation,
        _omero=self._omero,
    )

zarrnii.ZarrNii.upsample(along_x=1, along_y=1, along_z=1, to_shape=None)

Upsamples the ZarrNii instance using scipy.ndimage.zoom.

Parameters:

• along_x (int, default: 1 ) –

  Upsampling factor along the X-axis (default: 1).

• along_y (int, default: 1 ) –

  Upsampling factor along the Y-axis (default: 1).

• along_z (int, default: 1 ) –

  Upsampling factor along the Z-axis (default: 1).

• to_shape (tuple, default: None ) –

  Target shape for upsampling. Should include all dimensions (e.g., (c, z, y, x) for ZYX or (c, x, y, z) for XYZ). If provided, along_x, along_y, and along_z are ignored.

Returns:

• ZarrNii –

  A new ZarrNii instance with the upsampled data and updated affine.

Notes

• This method supports both direct scaling via along_* factors or target shape via to_shape.
• If to_shape is provided, chunk sizes and scaling factors are dynamically calculated.
• The affine matrix is updated to reflect the new voxel size after upsampling.

Example

Upsample with scaling factors

upsampled_znimg = znimg.upsample(along_x=2, along_y=2, along_z=2)

Upsample to a specific shape

upsampled_znimg = znimg.upsample(to_shape=(1, 256, 256, 256))

Source code in zarrnii/core.py

def upsample(self, along_x=1, along_y=1, along_z=1, to_shape=None):
    """
    Upsamples the ZarrNii instance using `scipy.ndimage.zoom`.

    Parameters:
        along_x (int, optional): Upsampling factor along the X-axis (default: 1).
        along_y (int, optional): Upsampling factor along the Y-axis (default: 1).
        along_z (int, optional): Upsampling factor along the Z-axis (default: 1).
        to_shape (tuple, optional): Target shape for upsampling. Should include all dimensions
                                     (e.g., `(c, z, y, x)` for ZYX or `(c, x, y, z)` for XYZ).
                                     If provided, `along_x`, `along_y`, and `along_z` are ignored.

    Returns:
        ZarrNii: A new ZarrNii instance with the upsampled data and updated affine.

    Notes:
        - This method supports both direct scaling via `along_*` factors or target shape via `to_shape`.
        - If `to_shape` is provided, chunk sizes and scaling factors are dynamically calculated.
        - The affine matrix is updated to reflect the new voxel size after upsampling.

    Example:
        # Upsample with scaling factors
        upsampled_znimg = znimg.upsample(along_x=2, along_y=2, along_z=2)

        # Upsample to a specific shape
        upsampled_znimg = znimg.upsample(to_shape=(1, 256, 256, 256))
    """
    # Determine scaling and chunks based on input parameters
    if to_shape is None:
        if self.axes_order == "XYZ":
            scaling = (1, along_x, along_y, along_z)
        else:
            scaling = (1, along_z, along_y, along_x)

        chunks_out = tuple(
            tuple(c * scale for c in chunks_i)
            for chunks_i, scale in zip(self.data.chunks, scaling)
        )
    else:
        chunks_out, scaling = self.__get_upsampled_chunks(to_shape)

    # Define block-wise upsampling function
    def zoom_blocks(x, block_info=None):
        """
        Scales blocks to the desired size using `scipy.ndimage.zoom`.

        Parameters:
            x (np.ndarray): Input block data.
            block_info (dict, optional): Metadata about the current block.

        Returns:
            np.ndarray: The upscaled block.
        """
        # Calculate scaling factors based on input and output chunk shapes
        scaling = tuple(
            out_n / in_n
            for out_n, in_n in zip(block_info[None]["chunk-shape"], x.shape)
        )
        return zoom(x, scaling, order=1, prefilter=False)

    # Perform block-wise upsampling
    darr_scaled = da.map_blocks(
        zoom_blocks, self.data, dtype=self.data.dtype, chunks=chunks_out
    )

    # Update the affine matrix to reflect the new voxel size
    if self.axes_order == "XYZ":
        scaling_matrix = np.diag(
            (1 / scaling[1], 1 / scaling[2], 1 / scaling[3], 1)
        )
    else:
        scaling_matrix = np.diag(
            (1 / scaling[-1], 1 / scaling[-2], 1 / scaling[-3], 1)
        )

    # Create new NgffImage with upsampled data
    dims = self.dims
    if self.axes_order == "XYZ":
        new_scale = {
            dims[1]: self.scale[dims[1]] / scaling[1],
            dims[2]: self.scale[dims[2]] / scaling[2],
            dims[3]: self.scale[dims[3]] / scaling[3],
        }
    else:
        new_scale = {
            dims[1]: self.scale[dims[1]] / scaling[1],
            dims[2]: self.scale[dims[2]] / scaling[2],
            dims[3]: self.scale[dims[3]] / scaling[3],
        }

    upsampled_ngff = nz.to_ngff_image(
        darr_scaled,
        dims=dims,
        scale=new_scale,
        translation=self.translation.copy(),
        name=self.name,
    )

    # Return a new ZarrNii instance with the upsampled data
    return ZarrNii.from_ngff_image(
        upsampled_ngff,
        axes_order=self.axes_order,
        xyz_orientation=self.xyz_orientation,
        omero=self.omero,
    )

zarrnii.ZarrNii.get_bounded_subregion(points)

Extracts a bounded subregion of the dask array containing the specified points, along with the grid points for interpolation.

If the points extend beyond the domain of the dask array, the extent is capped at the boundaries. If all points are outside the domain, the function returns (None, None).

Parameters:

• points (ndarray) –

  Nx3 or Nx4 array of coordinates in the array's space. If Nx4, the last column is assumed to be the homogeneous coordinate and is ignored.

Returns:

• tuple –

  grid_points (tuple): A tuple of three 1D arrays representing the grid points along each axis (X, Y, Z) in the subregion. subvol (np.ndarray or None): The extracted subregion as a NumPy array. Returns None if all points are outside the array domain.

Notes

• The function uses compute() on the dask array to immediately load the subregion, as Dask doesn't support the type of indexing required for interpolation.
• A padding of 1 voxel is applied around the extent of the points.

Example

grid_points, subvol = znimg.get_bounded_subregion(points) if subvol is not None: print("Subvolume shape:", subvol.shape)

Source code in zarrnii/core.py

def get_bounded_subregion(self, points: np.ndarray):
    """
    Extracts a bounded subregion of the dask array containing the specified points,
    along with the grid points for interpolation.

    If the points extend beyond the domain of the dask array, the extent is capped
    at the boundaries. If all points are outside the domain, the function returns
    `(None, None)`.

    Parameters:
        points (np.ndarray): Nx3 or Nx4 array of coordinates in the array's space.
                             If Nx4, the last column is assumed to be the homogeneous
                             coordinate and is ignored.

    Returns:
        tuple:
            grid_points (tuple): A tuple of three 1D arrays representing the grid
                                 points along each axis (X, Y, Z) in the subregion.
            subvol (np.ndarray or None): The extracted subregion as a NumPy array.
                                         Returns `None` if all points are outside
                                         the array domain.

    Notes:
        - The function uses `compute()` on the dask array to immediately load the
          subregion, as Dask doesn't support the type of indexing required for
          interpolation.
        - A padding of 1 voxel is applied around the extent of the points.

    Example:
        grid_points, subvol = znimg.get_bounded_subregion(points)
        if subvol is not None:
            print("Subvolume shape:", subvol.shape)
    """
    pad = 1  # Padding around the extent of the points

    # Compute the extent of the points in the array's coordinate space
    min_extent = np.floor(points.min(axis=1)[:3] - pad).astype("int")
    max_extent = np.ceil(points.max(axis=1)[:3] + pad).astype("int")

    # Clip the extents to ensure they stay within the bounds of the array
    clip_min = np.zeros_like(min_extent)
    clip_max = np.array(self.darr.shape[-3:])  # Z, Y, X dimensions

    min_extent = np.clip(min_extent, clip_min, clip_max)
    max_extent = np.clip(max_extent, clip_min, clip_max)

    # Check if all points are outside the domain
    if np.any(max_extent <= min_extent):
        return None, None

    # Extract the subvolume using the computed extents
    subvol = self.darr[
        :,
        min_extent[0] : max_extent[0],
        min_extent[1] : max_extent[1],
        min_extent[2] : max_extent[2],
    ].compute()

    # Generate grid points for interpolation
    grid_points = (
        np.arange(min_extent[0], max_extent[0]),  # Z
        np.arange(min_extent[1], max_extent[1]),  # Y
        np.arange(min_extent[2], max_extent[2]),  # X
    )

    return grid_points, subvol

zarrnii.ZarrNii.apply_transform(*transforms, ref_znimg, spatial_dims=None)

Apply spatial transformations to image data.

Transforms the image data to align with a reference image space using the provided transformation(s). This enables registration, resampling, and coordinate system conversions.

Parameters:

• *transforms (Transform, default: () ) –

  Variable number of Transform objects to apply sequentially. Supported transform types: - AffineTransform: Linear transformations (rotation, scaling, translation) - DisplacementTransform: Non-linear deformation fields

• ref_znimg ('ZarrNii') –

  Reference ZarrNii image defining the target coordinate system, grid spacing, and field of view for the output

• spatial_dims (Optional[List[str]], default: None ) –

  Names of spatial dimensions for transformation. If None, automatically derived from axes_order

Returns:

• 'ZarrNii' –

  New ZarrNii instance with transformed data in reference space

Raises:

• ValueError –

  If no transforms provided or reference image incompatible

• TypeError –

  If transforms are not valid Transform objects

Examples:

>>> # Apply affine transformation
>>> affine = AffineTransform.from_txt("transform.txt")
>>> transformed = moving.apply_transform(affine, ref_znimg=reference)

>>> # Apply multiple transforms sequentially
>>> affine = AffineTransform.identity()
>>> warp = DisplacementTransform.from_nifti("warp.nii.gz")
>>> result = moving.apply_transform(affine, warp, ref_znimg=reference)

Notes

• Transformations are applied in the order specified
• Output data inherits spatial properties from ref_znimg
• Uses interpolation for non-integer coordinate mappings
• Non-spatial dimensions (channels, time) are preserved

Source code in zarrnii/core.py

def apply_transform(
    self,
    *transforms: Transform,
    ref_znimg: "ZarrNii",
    spatial_dims: Optional[List[str]] = None,
) -> "ZarrNii":
    """Apply spatial transformations to image data.

    Transforms the image data to align with a reference image space using
    the provided transformation(s). This enables registration, resampling,
    and coordinate system conversions.

    Args:
        *transforms: Variable number of Transform objects to apply sequentially.
            Supported transform types:
            - AffineTransform: Linear transformations (rotation, scaling, translation)
            - DisplacementTransform: Non-linear deformation fields
        ref_znimg: Reference ZarrNii image defining the target coordinate system,
            grid spacing, and field of view for the output
        spatial_dims: Names of spatial dimensions for transformation. If None,
            automatically derived from axes_order

    Returns:
        New ZarrNii instance with transformed data in reference space

    Raises:
        ValueError: If no transforms provided or reference image incompatible
        TypeError: If transforms are not valid Transform objects

    Examples:
        >>> # Apply affine transformation
        >>> affine = AffineTransform.from_txt("transform.txt")
        >>> transformed = moving.apply_transform(affine, ref_znimg=reference)

        >>> # Apply multiple transforms sequentially
        >>> affine = AffineTransform.identity()
        >>> warp = DisplacementTransform.from_nifti("warp.nii.gz")
        >>> result = moving.apply_transform(affine, warp, ref_znimg=reference)

    Notes:
        - Transformations are applied in the order specified
        - Output data inherits spatial properties from ref_znimg
        - Uses interpolation for non-integer coordinate mappings
        - Non-spatial dimensions (channels, time) are preserved
    """
    if spatial_dims is None:
        spatial_dims = (
            ["z", "y", "x"] if self.axes_order == "ZYX" else ["x", "y", "z"]
        )

    # Initialize the list of transformations to apply
    tfms_to_apply = [ref_znimg.affine]  # Start with the reference image affine

    # Append all transformations passed as arguments
    tfms_to_apply.extend(transforms)

    # Append the inverse of the current image's affine
    tfms_to_apply.append(self.affine.invert())

    # Create new NgffImage from ref image
    interp_ngff_image = nz.NgffImage(
        data=ref_znimg.data,
        dims=ref_znimg.ngff_image.dims.copy(),
        scale=ref_znimg.ngff_image.scale.copy(),
        translation=ref_znimg.ngff_image.translation.copy(),
        name=f"{self.name}_transformed_to_{ref_znimg.name}",
    )
    # Lazily apply the transformations using dask
    interp_ngff_image.data = da.map_blocks(
        interp_by_block,  # Function to interpolate each block
        ref_znimg.data,  # Reference image data
        dtype=np.float32,  # Output data type
        transforms=tfms_to_apply,  # Transformations to apply
        flo_znimg=self,  # Floating image to align
    )

    return ZarrNii.from_ngff_image(
        interp_ngff_image,
        axes_order=ref_znimg.axes_order,
        xyz_orientation=ref_znimg.xyz_orientation,
        omero=self.omero,
    )

zarrnii.ZarrNii.to_ome_zarr(store_or_path, max_layer=4, scale_factors=None, backend='ome-zarr-py', **kwargs)

Save to OME-Zarr store with multiscale pyramid.

Creates an OME-Zarr dataset with automatic multiscale pyramid generation for efficient visualization and processing at multiple resolutions. Preserves spatial metadata and supports various storage backends.

Parameters:

• store_or_path (Union[str, Any]) –

  Target location for OME-Zarr store. Supports: - Local directory path - Remote URLs (s3://, gs://, etc.) - ZIP files (.zip extension for compressed storage) - Zarr store objects

• max_layer (int, default: 4 ) –

  Maximum number of pyramid levels to create (including level 0). Higher values create more downsampled levels

• scale_factors (Optional[List[int]], default: None ) –

  Custom downsampling factors for each pyramid level. If None, uses powers of 2: [2, 4, 8, 16, ...]

• backend (str, default: 'ome-zarr-py' ) –

  Backend library to use for writing. Options: - 'ngff-zarr': Use ngff-zarr library (default) - 'ome-zarr-py': Use ome-zarr-py library for better dask integration

• **kwargs (Any, default: {} ) –

  Additional arguments passed to the save function. For 'ngff-zarr': passed to to_ngff_zarr function For 'ome-zarr-py': passed to write_image (e.g., scaling_method, compute)

Returns:

• 'ZarrNii' –

  Self for method chaining

Raises:

• OSError –

  If unable to write to target location

• ValueError –

  If invalid scale_factors or backend provided

Examples:

>>> # Save with default pyramid levels
>>> znii.to_ome_zarr("/path/to/output.zarr")

>>> # Save to compressed ZIP with custom pyramid
>>> znii.to_ome_zarr(
...     "/path/to/output.zarr.zip",
...     max_layer=3,
...     scale_factors=[2, 4]
... )

>>> # Use ome-zarr-py backend for better dask performance
>>> znii.to_ome_zarr(
...     "/path/to/output.zarr",
...     backend="ome-zarr-py",
...     scaling_method="gaussian"
... )

>>> # Chain with other operations
>>> result = (znii.downsample(2)
...               .crop((0,0,0), (100,100,100))
...               .to_ome_zarr("processed.zarr"))

Notes

• OME-Zarr files are always saved in ZYX axis order
• Automatic axis reordering if current order is XYZ
• Spatial transformations and metadata are preserved
• Orientation information is stored using the new 'xyz_orientation' metadata key for consistency and future compatibility
• The 'ome-zarr-py' backend provides better performance with dask and dask distributed workflows

Source code in zarrnii/core.py

def to_ome_zarr(
    self,
    store_or_path: Union[str, Any],
    max_layer: int = 4,
    scale_factors: Optional[List[int]] = None,
    backend: str = "ome-zarr-py",
    **kwargs: Any,
) -> "ZarrNii":
    """Save to OME-Zarr store with multiscale pyramid.

    Creates an OME-Zarr dataset with automatic multiscale pyramid generation
    for efficient visualization and processing at multiple resolutions.
    Preserves spatial metadata and supports various storage backends.

    Args:
        store_or_path: Target location for OME-Zarr store. Supports:
            - Local directory path
            - Remote URLs (s3://, gs://, etc.)
            - ZIP files (.zip extension for compressed storage)
            - Zarr store objects
        max_layer: Maximum number of pyramid levels to create (including level 0).
            Higher values create more downsampled levels
        scale_factors: Custom downsampling factors for each pyramid level.
            If None, uses powers of 2: [2, 4, 8, 16, ...]
        backend: Backend library to use for writing. Options:
            - 'ngff-zarr': Use ngff-zarr library (default)
            - 'ome-zarr-py': Use ome-zarr-py library for better dask integration
        **kwargs: Additional arguments passed to the save function.
            For 'ngff-zarr': passed to to_ngff_zarr function
            For 'ome-zarr-py': passed to write_image (e.g., scaling_method, compute)

    Returns:
        Self for method chaining

    Raises:
        OSError: If unable to write to target location
        ValueError: If invalid scale_factors or backend provided

    Examples:
        >>> # Save with default pyramid levels
        >>> znii.to_ome_zarr("/path/to/output.zarr")

        >>> # Save to compressed ZIP with custom pyramid
        >>> znii.to_ome_zarr(
        ...     "/path/to/output.zarr.zip",
        ...     max_layer=3,
        ...     scale_factors=[2, 4]
        ... )

        >>> # Use ome-zarr-py backend for better dask performance
        >>> znii.to_ome_zarr(
        ...     "/path/to/output.zarr",
        ...     backend="ome-zarr-py",
        ...     scaling_method="gaussian"
        ... )

        >>> # Chain with other operations
        >>> result = (znii.downsample(2)
        ...               .crop((0,0,0), (100,100,100))
        ...               .to_ome_zarr("processed.zarr"))

    Notes:
        - OME-Zarr files are always saved in ZYX axis order
        - Automatic axis reordering if current order is XYZ
        - Spatial transformations and metadata are preserved
        - Orientation information is stored using the new 'xyz_orientation'
          metadata key for consistency and future compatibility
        - The 'ome-zarr-py' backend provides better performance with dask
          and dask distributed workflows
    """
    # Validate backend parameter
    if backend not in ["ngff-zarr", "ome-zarr-py"]:
        raise ValueError(
            f"Invalid backend '{backend}'. Must be 'ngff-zarr' or 'ome-zarr-py'"
        )

    # Determine the image to save
    if self.axes_order == "XYZ":
        # Need to reorder data from XYZ to ZYX for OME-Zarr
        ngff_image_to_save = self._create_zyx_ngff_image()
    else:
        # Already in ZYX order
        ngff_image_to_save = self.ngff_image

    # Use the appropriate save function based on backend
    if backend == "ngff-zarr":
        save_ngff_image(
            ngff_image_to_save,
            store_or_path,
            max_layer,
            scale_factors,
            xyz_orientation=(
                self.xyz_orientation if hasattr(self, "xyz_orientation") else None
            ),
            **kwargs,
        )
    elif backend == "ome-zarr-py":
        save_ngff_image_with_ome_zarr(
            ngff_image_to_save,
            store_or_path,
            max_layer,
            scale_factors,
            xyz_orientation=(
                self.xyz_orientation if hasattr(self, "xyz_orientation") else None
            ),
            **kwargs,
        )

    # Add orientation metadata to the zarr store (only for non-ZIP files)
    # For ZIP files, orientation is handled inside save_ngff_image
    if not (isinstance(store_or_path, str) and store_or_path.endswith(".zip")):
        try:
            import zarr

            if isinstance(store_or_path, str):
                group = zarr.open_group(store_or_path, mode="r+")
            else:
                group = zarr.open_group(store_or_path, mode="r+")

            # Add metadata for xyz_orientation (new format)
            if hasattr(self, "xyz_orientation") and self.xyz_orientation:
                group.attrs["xyz_orientation"] = self.xyz_orientation
        except Exception:
            # If we can't write orientation metadata, that's not critical
            pass

    return self

zarrnii.ZarrNii.to_nifti(filename=None)

Convert to NIfTI format with automatic dimension handling.

Converts the ZarrNii image to NIfTI-1 format, handling dimension reordering, singleton dimension removal, and spatial transformation conversion. NIfTI files are always written in XYZ axis order.

Parameters:

• filename (Optional[Union[str, bytes]], default: None ) –

  Output file path for saving. Supported extensions: - .nii: Uncompressed NIfTI - .nii.gz: Compressed NIfTI (recommended) If None, returns nibabel image object without saving

Returns:

• Union[Nifti1Image, str] –

  If filename is None: nibabel.Nifti1Image object

• Union[Nifti1Image, str] –

  If filename provided: path to saved file

Raises:

• ValueError –

  If data has non-singleton time or channel dimensions (NIfTI doesn't support >4D data)

• OSError –

  If unable to write to specified filename

Notes

• Automatically reorders data from ZYX to XYZ if necessary
• Removes singleton time/channel dimensions automatically
• Spatial transformations are converted to NIfTI affine format
• For 5D data (T,C,Z,Y,X), only singleton T/C dimensions are supported

Examples:

>>> # Save to compressed NIfTI file
>>> znii.to_nifti("output.nii.gz")

>>> # Get nibabel object without saving
>>> nifti_img = znii.to_nifti()
>>> print(nifti_img.shape)

>>> # Handle multi-channel data by selecting single channel first
>>> znii.select_channels([0]).to_nifti("channel0.nii.gz")

Source code in zarrnii/core.py

def to_nifti(
    self, filename: Optional[Union[str, bytes]] = None
) -> Union[nib.Nifti1Image, str]:
    """Convert to NIfTI format with automatic dimension handling.

    Converts the ZarrNii image to NIfTI-1 format, handling dimension
    reordering, singleton dimension removal, and spatial transformation
    conversion. NIfTI files are always written in XYZ axis order.

    Args:
        filename: Output file path for saving. Supported extensions:
            - .nii: Uncompressed NIfTI
            - .nii.gz: Compressed NIfTI (recommended)
            If None, returns nibabel image object without saving

    Returns:
        If filename is None: nibabel.Nifti1Image object
        If filename provided: path to saved file

    Raises:
        ValueError: If data has non-singleton time or channel dimensions
            (NIfTI doesn't support >4D data)
        OSError: If unable to write to specified filename

    Notes:
        - Automatically reorders data from ZYX to XYZ if necessary
        - Removes singleton time/channel dimensions automatically
        - Spatial transformations are converted to NIfTI affine format
        - For 5D data (T,C,Z,Y,X), only singleton T/C dimensions are supported

    Examples:
        >>> # Save to compressed NIfTI file
        >>> znii.to_nifti("output.nii.gz")

        >>> # Get nibabel object without saving
        >>> nifti_img = znii.to_nifti()
        >>> print(nifti_img.shape)

        >>> # Handle multi-channel data by selecting single channel first
        >>> znii.select_channels([0]).to_nifti("channel0.nii.gz")

    Warnings:
        Large images will be computed in memory during conversion.
        Consider downsampling or cropping first for very large datasets.
    """
    # Get data and dimensions
    data = self.data.compute()

    dims = self.dims

    # Handle dimensional reduction for NIfTI compatibility
    # NIfTI supports up to 4D, so we need to remove singleton dimensions
    squeeze_axes = []
    remaining_dims = []

    for i, dim in enumerate(dims):
        if dim in ["t", "c"] and data.shape[i] == 1:
            # Remove singleton time or channel dimensions
            squeeze_axes.append(i)
        elif dim in ["t", "c"] and data.shape[i] > 1:
            # Non-singleton time or channel dimensions - NIfTI can't handle this
            raise ValueError(
                f"NIfTI format doesn't support non-singleton {dim} dimension. "
                f"Dimension '{dim}' has size {data.shape[i]}. "
                f"Consider selecting specific timepoints/channels first."
            )
        else:
            remaining_dims.append(dim)

    # Squeeze out singleton dimensions
    if squeeze_axes:
        data = np.squeeze(data, axis=tuple(squeeze_axes))

    # Check final dimensionality
    if data.ndim > 4:
        raise ValueError(
            f"Resulting data has {data.ndim} dimensions, but NIfTI supports maximum 4D"
        )

    # Now handle spatial reordering based on axes_order
    if self.axes_order == "ZYX":
        # Data spatial dimensions are in ZYX order, need to transpose to XYZ
        if data.ndim == 3:
            # Pure spatial data: ZYX -> XYZ
            data = data.transpose(2, 1, 0)
        elif data.ndim == 4:
            # 4D data with one non-spatial dimension remaining
            # Could be (T,Z,Y,X) or (C,Z,Y,X) - spatial part needs ZYX->XYZ
            # The non-spatial dimension stays first
            data = data.transpose(0, 3, 2, 1)

        # Get affine matrix in XYZ order
        affine_matrix = self.get_affine_matrix(axes_order="XYZ")
    else:
        # Data is already in XYZ order
        affine_matrix = self.get_affine_matrix(axes_order="XYZ")

    # Create NIfTI image
    nifti_img = nib.Nifti1Image(data, affine_matrix)

    if filename is not None:
        nib.save(nifti_img, filename)
        return filename
    else:
        return nifti_img

zarrnii.ZarrNii.to_tiff_stack(filename_pattern, channel=None, timepoint=None, compress=True, dtype='uint16', rescale=True)

Save data as a stack of 2D TIFF images.

Saves the image data as a series of 2D TIFF files, with each Z-slice saved as a separate file. This format is useful for compatibility with tools that don't support OME-Zarr or napari plugins that require individual TIFF files.

Parameters:

• filename_pattern (str) –

  Output filename pattern. Should contain '{z:04d}' or similar format specifier for the Z-slice number. Examples: - "output_z{z:04d}.tif" - "data/slice_{z:03d}.tiff" If pattern doesn't contain format specifier, '_{z:04d}' is appended before the extension.

• channel (Optional[int], default: None ) –

  Channel index to save (0-based). If None and data has multiple channels, all channels will be saved as separate channel dimensions in each TIFF file (multi-channel TIFFs).

• timepoint (Optional[int], default: None ) –

  Timepoint index to save (0-based). If None and data has multiple timepoints, raises ValueError (must select single timepoint).

• compress (bool, default: True ) –

  Whether to use LZW compression (default: True)

• dtype (Optional[str], default: 'uint16' ) –

  Output data type for TIFF files. Options: - 'uint8': 8-bit unsigned integer (0-255) - 'uint16': 16-bit unsigned integer (0-65535) [default] - 'int16': 16-bit signed integer (-32768 to 32767) - 'float32': 32-bit float (preserves original data) Default 'uint16' provides good range and compatibility.

• rescale (bool, default: True ) –

  Whether to rescale data to fit the output dtype range. If True, data is linearly scaled from [min, max] to the full range of the output dtype. If False, data is clipped to the output dtype range. Default: True

Returns:

• str –

  Base directory path where files were saved

Raises:

• ValueError –

  If data has multiple timepoints but none selected, or if selected channel/timepoint is out of range, or if dtype is not supported

• OSError –

  If unable to write to specified directory

Examples:

>>> # Save as 16-bit with auto-rescaling (default, recommended)
>>> znii.to_tiff_stack("output_z{z:04d}.tif")

>>> # Save as 8-bit for smaller file sizes
>>> znii.to_tiff_stack("output_z{z:04d}.tif", dtype='uint8')

>>> # Save specific channel without rescaling
>>> znii.to_tiff_stack("channel0_z{z:04d}.tif", channel=0, rescale=False)

>>> # Save as float32 to preserve original precision
>>> znii.to_tiff_stack("precise_z{z:04d}.tif", dtype='float32')

Notes

• Z-dimension becomes the stack (file) dimension
• Time and channel dimensions are handled as specified
• Spatial transformations are not preserved in TIFF format
• For 5D data (T,C,Z,Y,X), you must select a single timepoint
• Multi-channel data can be saved as multi-channel TIFFs or selected
• Data type conversion helps ensure compatibility with analysis tools
• uint16 is recommended for most scientific applications (good range + compatibility)

Source code in zarrnii/core.py

def to_tiff_stack(
    self,
    filename_pattern: str,
    channel: Optional[int] = None,
    timepoint: Optional[int] = None,
    compress: bool = True,
    dtype: Optional[str] = "uint16",
    rescale: bool = True,
) -> str:
    """Save data as a stack of 2D TIFF images.

    Saves the image data as a series of 2D TIFF files, with each Z-slice
    saved as a separate file. This format is useful for compatibility with
    tools that don't support OME-Zarr or napari plugins that require
    individual TIFF files.

    Args:
        filename_pattern: Output filename pattern. Should contain '{z:04d}' or similar
            format specifier for the Z-slice number. Examples:
            - "output_z{z:04d}.tif"
            - "data/slice_{z:03d}.tiff"
            If pattern doesn't contain format specifier, '_{z:04d}' is appended
            before the extension.
        channel: Channel index to save (0-based). If None and data has multiple
            channels, all channels will be saved as separate channel dimensions
            in each TIFF file (multi-channel TIFFs).
        timepoint: Timepoint index to save (0-based). If None and data has multiple
            timepoints, raises ValueError (must select single timepoint).
        compress: Whether to use LZW compression (default: True)
        dtype: Output data type for TIFF files. Options:
            - 'uint8': 8-bit unsigned integer (0-255)
            - 'uint16': 16-bit unsigned integer (0-65535) [default]
            - 'int16': 16-bit signed integer (-32768 to 32767)
            - 'float32': 32-bit float (preserves original data)
            Default 'uint16' provides good range and compatibility.
        rescale: Whether to rescale data to fit the output dtype range.
            If True, data is linearly scaled from [min, max] to the full
            range of the output dtype. If False, data is clipped to the
            output dtype range. Default: True

    Returns:
        Base directory path where files were saved

    Raises:
        ValueError: If data has multiple timepoints but none selected,
            or if selected channel/timepoint is out of range,
            or if dtype is not supported
        OSError: If unable to write to specified directory

    Examples:
        >>> # Save as 16-bit with auto-rescaling (default, recommended)
        >>> znii.to_tiff_stack("output_z{z:04d}.tif")

        >>> # Save as 8-bit for smaller file sizes
        >>> znii.to_tiff_stack("output_z{z:04d}.tif", dtype='uint8')

        >>> # Save specific channel without rescaling
        >>> znii.to_tiff_stack("channel0_z{z:04d}.tif", channel=0, rescale=False)

        >>> # Save as float32 to preserve original precision
        >>> znii.to_tiff_stack("precise_z{z:04d}.tif", dtype='float32')

    Warnings:
        This method loads all data into memory. For large datasets,
        consider cropping or downsampling first to reduce memory usage.
        The cellseg3d napari plugin and similar tools work best with
        cropped regions rather than full-resolution whole-brain images.

    Notes:
        - Z-dimension becomes the stack (file) dimension
        - Time and channel dimensions are handled as specified
        - Spatial transformations are not preserved in TIFF format
        - For 5D data (T,C,Z,Y,X), you must select a single timepoint
        - Multi-channel data can be saved as multi-channel TIFFs or selected
        - Data type conversion helps ensure compatibility with analysis tools
        - uint16 is recommended for most scientific applications (good range + compatibility)
    """
    try:
        import tifffile
    except ImportError:
        raise ImportError(
            "tifffile is required for TIFF stack support. "
            "Install with: pip install tifffile"
        )

    # Get data and dimensions
    data = self.data.compute()
    dims = self.dims

    # Create output directory if needed
    import os

    output_dir = os.path.dirname(filename_pattern)
    if output_dir and not os.path.exists(output_dir):
        os.makedirs(output_dir)

    # Handle dimensional selection and validation
    # Remove singleton dimensions first, similar to to_nifti
    squeeze_axes = []
    remaining_dims = []
    time_dim_size = 1
    channel_dim_size = 1

    for i, dim in enumerate(dims):
        if dim == "t":
            time_dim_size = data.shape[i]
            if data.shape[i] == 1:
                squeeze_axes.append(i)
            elif timepoint is None:
                raise ValueError(
                    f"Data has {data.shape[i]} timepoints. "
                    f"Must specify 'timepoint' parameter to select a single timepoint."
                )
            elif timepoint >= data.shape[i]:
                raise ValueError(
                    f"Timepoint {timepoint} is out of range (data has {data.shape[i]} timepoints)"
                )
            else:
                remaining_dims.append(dim)
        elif dim == "c":
            channel_dim_size = data.shape[i]
            if data.shape[i] == 1:
                squeeze_axes.append(i)
            elif channel is None:
                raise ValueError(
                    f"Data has {data.shape[i]} channels. "
                    f"Must specify 'channel' parameter to select a single channel."
                )
            elif channel >= data.shape[i]:
                raise ValueError(
                    f"Channel {channel} is out of range (data has {data.shape[i]} channels)"
                )
            else:
                remaining_dims.append(dim)
        else:
            remaining_dims.append(dim)

    # Select specific timepoint if needed
    if time_dim_size > 1 and timepoint is not None:
        time_axis = dims.index("t")
        data = np.take(data, timepoint, axis=time_axis)
        # Update dims list
        dims = [d for i, d in enumerate(dims) if i != time_axis]

    # Select specific channel if needed
    if channel_dim_size > 1 and channel is not None:
        channel_axis = dims.index("c")
        data = np.take(data, channel, axis=channel_axis)
        # Update dims list
        dims = [d for i, d in enumerate(dims) if i != channel_axis]

    # Squeeze singleton dimensions
    if squeeze_axes:
        # Recalculate squeeze axes after potential dimension removal
        current_squeeze_axes = []
        for axis in squeeze_axes:
            # Count how many axes were removed before this one
            removed_before = sum(
                1
                for removed_axis in [
                    (
                        dims.index("t")
                        if time_dim_size > 1 and timepoint is not None
                        else -1
                    ),
                    (
                        dims.index("c")
                        if channel_dim_size > 1 and channel is not None
                        else -1
                    ),
                ]
                if removed_axis != -1 and removed_axis < axis
            )
            current_squeeze_axes.append(axis - removed_before)

        data = np.squeeze(data, axis=tuple(current_squeeze_axes))
        dims = [dim for i, dim in enumerate(dims) if i not in current_squeeze_axes]

    # Find Z dimension for stacking
    if "z" not in dims:
        raise ValueError("Data must have a Z dimension for TIFF stack export")

    z_axis = dims.index("z")
    z_size = data.shape[z_axis]

    # Check filename pattern contains format specifier
    if "{z" not in filename_pattern:
        # Add default z format before extension
        name, ext = os.path.splitext(filename_pattern)
        filename_pattern = f"{name}_{{z:04d}}{ext}"

    # Move Z axis to first position for easy iteration
    axes_order = list(range(data.ndim))
    axes_order[0], axes_order[z_axis] = axes_order[z_axis], axes_order[0]
    data = data.transpose(axes_order)

    # Handle data type conversion and rescaling
    supported_dtypes = {
        "uint8": np.uint8,
        "uint16": np.uint16,
        "int16": np.int16,
        "float32": np.float32,
    }

    if dtype not in supported_dtypes:
        raise ValueError(
            f"Unsupported dtype '{dtype}'. Supported types: {list(supported_dtypes.keys())}"
        )

    target_dtype = supported_dtypes[dtype]

    if rescale and dtype != "float32":
        # Get the data range
        data_min = np.min(data)
        data_max = np.max(data)

        if data_min == data_max:
            # Handle constant data case
            data_scaled = np.zeros_like(data, dtype=target_dtype)
        else:
            # Get target range for the dtype
            if dtype == "uint8":
                target_min, target_max = 0, 255
            elif dtype == "uint16":
                target_min, target_max = 0, 65535
            elif dtype == "int16":
                target_min, target_max = -32768, 32767

            # Linear rescaling: new_value = (value - data_min) * (target_max - target_min) / (data_max - data_min) + target_min
            data_scaled = (
                (data - data_min)
                * (target_max - target_min)
                / (data_max - data_min)
                + target_min
            ).astype(target_dtype)

        print(
            f"Rescaled data from [{data_min:.3f}, {data_max:.3f}] to {dtype} range"
        )
    else:
        # No rescaling - just clip and convert
        if dtype == "uint8":
            data_scaled = np.clip(data, 0, 255).astype(target_dtype)
        elif dtype == "uint16":
            data_scaled = np.clip(data, 0, 65535).astype(target_dtype)
        elif dtype == "int16":
            data_scaled = np.clip(data, -32768, 32767).astype(target_dtype)
        else:  # float32
            data_scaled = data.astype(target_dtype)

        if dtype != "float32":
            print(f"Converted data to {dtype} with clipping (no rescaling)")

    data = data_scaled

    # Save each Z-slice as a separate TIFF file
    compression = "lzw" if compress else None
    saved_files = []

    for z_idx in range(z_size):
        slice_data = data[z_idx]

        # Generate filename for this slice
        filename = filename_pattern.format(z=z_idx)

        # Save the 2D slice
        tifffile.imwrite(filename, slice_data, compression=compression)
        saved_files.append(filename)

    print(f"Saved {len(saved_files)} TIFF files to {output_dir or '.'}")
    print(
        f"Files: {os.path.basename(saved_files[0])} ... {os.path.basename(saved_files[-1])}"
    )

    return output_dir or "."

zarrnii.ZarrNii.from_imaris(path, level=0, timepoint=0, channel=0, chunks='auto', axes_order='ZYX', orientation='RAS') classmethod

Load from Imaris (.ims) file format.

Imaris files use HDF5 format with specific dataset structure. This method requires the 'imaris' extra dependency (h5py).

Parameters:

• path (str) –

  Path to Imaris (.ims) file

• level (int, default: 0 ) –

  Resolution level to load (0 = full resolution)

• timepoint (int, default: 0 ) –

  Time point to load (default: 0)

• channel (int, default: 0 ) –

  Channel to load (default: 0)

• chunks (str, default: 'auto' ) –

  Chunking strategy for dask array

• axes_order (str, default: 'ZYX' ) –

  Spatial axes order for compatibility (default: "ZYX")

• orientation (str, default: 'RAS' ) –

  Default orientation (default: "RAS")

Returns:

• 'ZarrNii' –

  ZarrNii instance

Raises:

• ImportError –

  If h5py is not available

• ValueError –

  If the file is not a valid Imaris file

Source code in zarrnii/core.py

@classmethod
def from_imaris(
    cls,
    path: str,
    level: int = 0,
    timepoint: int = 0,
    channel: int = 0,
    chunks: str = "auto",
    axes_order: str = "ZYX",
    orientation: str = "RAS",
) -> "ZarrNii":
    """
    Load from Imaris (.ims) file format.

    Imaris files use HDF5 format with specific dataset structure.
    This method requires the 'imaris' extra dependency (h5py).

    Args:
        path: Path to Imaris (.ims) file
        level: Resolution level to load (0 = full resolution)
        timepoint: Time point to load (default: 0)
        channel: Channel to load (default: 0)
        chunks: Chunking strategy for dask array
        axes_order: Spatial axes order for compatibility (default: "ZYX")
        orientation: Default orientation (default: "RAS")

    Returns:
        ZarrNii instance

    Raises:
        ImportError: If h5py is not available
        ValueError: If the file is not a valid Imaris file
    """
    try:
        import h5py
    except ImportError:
        raise ImportError(
            "h5py is required for Imaris support. "
            "Install with: pip install zarrnii[imaris] or pip install h5py"
        )

    # Open Imaris file
    with h5py.File(path, "r") as f:
        # Verify it's an Imaris file by checking for standard structure
        if "DataSet" not in f:
            raise ValueError(
                f"File {path} does not appear to be a valid Imaris file (missing DataSet group)"
            )

        # Navigate to the specific dataset
        dataset_group = f["DataSet"]

        # Find available resolution levels
        resolution_levels = [
            key for key in dataset_group.keys() if key.startswith("ResolutionLevel")
        ]
        if not resolution_levels:
            raise ValueError("No resolution levels found in Imaris file")

        # Validate level parameter
        if level >= len(resolution_levels):
            raise ValueError(
                f"Level {level} not available. Available levels: 0-{len(resolution_levels)-1}"
            )

        # Navigate to specified resolution level
        res_level_key = f"ResolutionLevel {level}"
        if res_level_key not in dataset_group:
            raise ValueError(f"Resolution level {level} not found")

        res_group = dataset_group[res_level_key]

        # Find available timepoints
        timepoints = [
            key for key in res_group.keys() if key.startswith("TimePoint")
        ]
        if not timepoints:
            raise ValueError("No timepoints found in Imaris file")

        # Validate timepoint parameter
        if timepoint >= len(timepoints):
            raise ValueError(
                f"Timepoint {timepoint} not available. Available timepoints: 0-{len(timepoints)-1}"
            )

        # Navigate to specified timepoint
        time_key = f"TimePoint {timepoint}"
        if time_key not in res_group:
            raise ValueError(f"Timepoint {timepoint} not found")

        time_group = res_group[time_key]

        # Find available channels
        channels = [key for key in time_group.keys() if key.startswith("Channel")]
        if not channels:
            raise ValueError("No channels found in Imaris file")

        # Validate channel parameter
        if channel >= len(channels):
            raise ValueError(
                f"Channel {channel} not available. Available channels: 0-{len(channels)-1}"
            )

        # Navigate to specified channel
        channel_key = f"Channel {channel}"
        if channel_key not in time_group:
            raise ValueError(f"Channel {channel} not found")

        channel_group = time_group[channel_key]

        # Load the actual data
        if "Data" not in channel_group:
            raise ValueError("No Data dataset found in channel group")

        data_dataset = channel_group["Data"]

        # Load data into memory first (necessary because HDF5 file will be closed)
        data_numpy = data_dataset[:]

        # Create dask array from numpy array
        data_array = da.from_array(data_numpy, chunks=chunks)

        # Add channel dimension if not present
        if len(data_array.shape) == 3:
            data_array = data_array[np.newaxis, ...]

        # Extract spatial metadata
        # Try to get spacing information from Imaris metadata
        spacing = [1.0, 1.0, 1.0]  # Default spacing
        origin = [0.0, 0.0, 0.0]  # Default origin

        # Look for ImageSizeX, ImageSizeY, ImageSizeZ attributes
        try:
            # Navigate back to get image info
            if "ImageSizeX" in f.attrs:
                x_size = f.attrs["ImageSizeX"]
                y_size = f.attrs["ImageSizeY"]
                z_size = f.attrs["ImageSizeZ"]

                # Calculate spacing based on physical size and voxel count
                if data_array.shape[-1] > 0:  # X dimension
                    spacing[0] = x_size / data_array.shape[-1]
                if data_array.shape[-2] > 0:  # Y dimension
                    spacing[1] = y_size / data_array.shape[-2]
                if data_array.shape[-3] > 0:  # Z dimension
                    spacing[2] = z_size / data_array.shape[-3]
        except (KeyError, IndexError):
            # Use default spacing if metadata is not available
            pass

        # Create dimensions
        dims = ["c"] + list(axes_order.lower())

        # Create scale and translation dictionaries
        scale_dict = {}
        translation_dict = {}
        spatial_dims = ["z", "y", "x"] if axes_order == "ZYX" else ["x", "y", "z"]

        for i, dim in enumerate(spatial_dims):
            scale_dict[dim] = spacing[i]
            translation_dict[dim] = origin[i]

        # Create NgffImage
        ngff_image = nz.NgffImage(
            data=data_array,
            dims=dims,
            scale=scale_dict,
            translation=translation_dict,
            name=f"imaris_image_{path}_{level}_{timepoint}_{channel}",
        )

    # Create and return ZarrNii instance
    return cls(
        ngff_image=ngff_image,
        axes_order=axes_order,
        xyz_orientation=orientation,
        _omero=None,
    )

zarrnii.ZarrNii.to_imaris(path, compression='gzip', compression_opts=6)

Save to Imaris (.ims) file format using HDF5.

This method creates Imaris files compatible with Imaris software by following the exact HDF5 structure from correctly-formed reference files. All attributes use byte-array encoding as required by Imaris.

Parameters:

• path (str) –

  Output path for Imaris (.ims) file

• compression (str, default: 'gzip' ) –

  HDF5 compression method (default: "gzip")

• compression_opts (int, default: 6 ) –

  Compression level (default: 6)

Returns:

• str ( str ) –

  Path to the saved file

Raises:

• ImportError –

  If h5py is not available

Source code in zarrnii/core.py

def to_imaris(
    self, path: str, compression: str = "gzip", compression_opts: int = 6
) -> str:
    """
    Save to Imaris (.ims) file format using HDF5.

    This method creates Imaris files compatible with Imaris software by
    following the exact HDF5 structure from correctly-formed reference files.
    All attributes use byte-array encoding as required by Imaris.

    Args:
        path: Output path for Imaris (.ims) file
        compression: HDF5 compression method (default: "gzip")
        compression_opts: Compression level (default: 6)

    Returns:
        str: Path to the saved file

    Raises:
        ImportError: If h5py is not available
    """
    try:
        import h5py
    except ImportError:
        raise ImportError(
            "h5py is required for Imaris support. "
            "Install with: pip install zarrnii[imaris] or pip install h5py"
        )

    # Ensure path has .ims extension
    if not path.endswith(".ims"):
        path = path + ".ims"

    def _string_to_byte_array(s: str) -> np.ndarray:
        """Convert string to byte array as required by Imaris."""
        return np.array([c.encode() for c in s])

    # Get data and metadata
    if hasattr(self.darr, "compute"):
        data = self.darr.compute()  # Convert Dask array to numpy array
    else:
        data = np.asarray(self.darr)  # Handle numpy arrays directly

    # Handle dimensions: expect ZYX or CZYX
    if len(data.shape) == 4:
        # CZYX format
        n_channels = data.shape[0]
        z, y, x = data.shape[1:]
    elif len(data.shape) == 3:
        # ZYX format - single channel
        n_channels = 1
        z, y, x = data.shape
        data = data[np.newaxis, ...]  # Add channel dimension
    else:
        raise ValueError(
            f"Unsupported data shape: {data.shape}. Expected 3D (ZYX) or 4D (CZYX)"
        )

    # Create Imaris file structure exactly matching reference file
    with h5py.File(path, "w") as f:
        # Root attributes - use exact byte array format from reference
        f.attrs["DataSetDirectoryName"] = _string_to_byte_array("DataSet")
        f.attrs["DataSetInfoDirectoryName"] = _string_to_byte_array("DataSetInfo")
        f.attrs["ImarisDataSet"] = _string_to_byte_array("ImarisDataSet")
        f.attrs["ImarisVersion"] = _string_to_byte_array("5.5.0")
        f.attrs["NumberOfDataSets"] = np.array([1], dtype=np.uint32)
        f.attrs["ThumbnailDirectoryName"] = _string_to_byte_array("Thumbnail")

        # Create main DataSet group structure
        dataset_group = f.create_group("DataSet")
        res_group = dataset_group.create_group("ResolutionLevel 0")
        time_group = res_group.create_group("TimePoint 0")

        # Create channels with proper attributes
        for c in range(n_channels):
            channel_group = time_group.create_group(f"Channel {c}")
            channel_data = data[c]  # (Z, Y, X)

            # Channel attributes - use byte array format exactly like reference
            channel_group.attrs["ImageSizeX"] = _string_to_byte_array(str(x))
            channel_group.attrs["ImageSizeY"] = _string_to_byte_array(str(y))
            channel_group.attrs["ImageSizeZ"] = _string_to_byte_array(str(z))
            channel_group.attrs["ImageBlockSizeX"] = _string_to_byte_array(str(x))
            channel_group.attrs["ImageBlockSizeY"] = _string_to_byte_array(str(y))
            channel_group.attrs["ImageBlockSizeZ"] = _string_to_byte_array(
                str(min(z, 16))
            )

            # Histogram range attributes
            data_min, data_max = float(channel_data.min()), float(
                channel_data.max()
            )
            channel_group.attrs["HistogramMin"] = _string_to_byte_array(
                f"{data_min:.3f}"
            )
            channel_group.attrs["HistogramMax"] = _string_to_byte_array(
                f"{data_max:.3f}"
            )

            # Create data dataset with proper compression
            # Preserve original data type but ensure it's compatible with Imaris
            if channel_data.dtype == np.float32 or channel_data.dtype == np.float64:
                # Keep float data as is for round-trip compatibility
                data_for_storage = channel_data.astype(np.float32)
            elif channel_data.dtype in [np.uint16, np.int16]:
                # Keep 16-bit data as is
                data_for_storage = channel_data
            else:
                # Convert other types to uint8
                data_for_storage = channel_data.astype(np.uint8)

            channel_group.create_dataset(
                "Data",
                data=data_for_storage,
                compression=compression,
                compression_opts=compression_opts,
                chunks=True,
            )

            # Create histogram
            hist_data, _ = np.histogram(
                channel_data.flatten(), bins=256, range=(data_min, data_max)
            )
            channel_group.create_dataset(
                "Histogram", data=hist_data.astype(np.uint64)
            )

        # Get spacing directly from scale dictionary with proper XYZ order
        try:
            # Extract voxel sizes directly from ngff_image scale dictionary
            # This ensures we get X, Y, Z in the correct order regardless of axes_order
            sx = self.ngff_image.scale.get("x", 1.0)
            sy = self.ngff_image.scale.get("y", 1.0)
            sz = self.ngff_image.scale.get("z", 1.0)
        except:
            sx = sy = sz = 1.0

        # Calculate extents (physical coordinates)
        ext_x = sx * x
        ext_y = sy * y
        ext_z = sz * z

        # Create comprehensive DataSetInfo structure matching reference
        info_group = f.create_group("DataSetInfo")

        # Create channel info groups
        for c in range(n_channels):
            channel_info = info_group.create_group(f"Channel {c}")

            # Essential channel attributes in byte array format
            channel_info.attrs["Color"] = _string_to_byte_array(
                "1.000 0.000 0.000"
                if c == 0
                else f"0.000 {1.0 if c == 1 else 0.0:.3f} {1.0 if c == 2 else 0.0:.3f}"
            )
            channel_info.attrs["Name"] = _string_to_byte_array(f"Channel {c}")
            channel_info.attrs["ColorMode"] = _string_to_byte_array("BaseColor")
            channel_info.attrs["ColorOpacity"] = _string_to_byte_array("1.000")
            channel_info.attrs["ColorRange"] = _string_to_byte_array("0 255")
            channel_info.attrs["GammaCorrection"] = _string_to_byte_array("1.000")
            channel_info.attrs["LSMEmissionWavelength"] = _string_to_byte_array(
                "500"
            )
            channel_info.attrs["LSMExcitationWavelength"] = _string_to_byte_array(
                "500"
            )
            channel_info.attrs["LSMPhotons"] = _string_to_byte_array("1")
            channel_info.attrs["LSMPinhole"] = _string_to_byte_array("0")

            # Add description
            description = f"Channel {c} created by ZarrNii"
            channel_info.attrs["Description"] = _string_to_byte_array(description)

        # Create CRITICAL Image group with voxel size information (this was missing!)
        image_info = info_group.create_group("Image")

        # Add essential image metadata with proper voxel size information
        image_info.attrs["X"] = _string_to_byte_array(str(x))
        image_info.attrs["Y"] = _string_to_byte_array(str(y))
        image_info.attrs["Z"] = _string_to_byte_array(str(z))
        image_info.attrs["Unit"] = _string_to_byte_array("um")
        image_info.attrs["Noc"] = _string_to_byte_array(str(n_channels))

        # CRITICAL: Set proper physical extents that define voxel size
        # Imaris reads voxel size from these extent values
        image_info.attrs["ExtMin0"] = _string_to_byte_array(f"{-ext_x/2:.3f}")
        image_info.attrs["ExtMax0"] = _string_to_byte_array(f"{ext_x/2:.3f}")
        image_info.attrs["ExtMin1"] = _string_to_byte_array(f"{-ext_y/2:.3f}")
        image_info.attrs["ExtMax1"] = _string_to_byte_array(f"{ext_y/2:.3f}")
        image_info.attrs["ExtMin2"] = _string_to_byte_array(f"{-ext_z/2:.3f}")
        image_info.attrs["ExtMax2"] = _string_to_byte_array(f"{ext_z/2:.3f}")

        # Add device/acquisition metadata
        image_info.attrs["ManufactorString"] = _string_to_byte_array("ZarrNii")
        image_info.attrs["ManufactorType"] = _string_to_byte_array("Generic")
        image_info.attrs["LensPower"] = _string_to_byte_array("")
        image_info.attrs["NumericalAperture"] = _string_to_byte_array("")
        image_info.attrs["RecordingDate"] = _string_to_byte_array(
            "2024-01-01 00:00:00.000"
        )
        image_info.attrs["Filename"] = _string_to_byte_array(path.split("/")[-1])
        image_info.attrs["Name"] = _string_to_byte_array("ZarrNii Export")
        image_info.attrs["Compression"] = _string_to_byte_array("5794")

        # Add description
        description = (
            f"Imaris file created by ZarrNii from {self.axes_order} format data. "
            f"Original shape: {self.darr.shape}. Converted to Imaris format "
            f"with {n_channels} channel(s) and dimensions {z}x{y}x{x}. "
            f"Voxel size: {sx:.3f} x {sy:.3f} x {sz:.3f} um."
        )
        image_info.attrs["Description"] = _string_to_byte_array(description)

        # Create Imaris metadata group
        imaris_info = info_group.create_group("Imaris")
        imaris_info.attrs["Version"] = _string_to_byte_array("7.0")
        imaris_info.attrs["ThumbnailMode"] = _string_to_byte_array("thumbnailMIP")
        imaris_info.attrs["ThumbnailSize"] = _string_to_byte_array("256")

        # Create ImarisDataSet metadata
        dataset_info = info_group.create_group("ImarisDataSet")
        dataset_info.attrs["Creator"] = _string_to_byte_array("Imaris")
        dataset_info.attrs["Version"] = _string_to_byte_array("7.0")
        dataset_info.attrs["NumberOfImages"] = _string_to_byte_array("1")

        # Add version-specific groups as seen in reference
        dataset_info_ver = info_group.create_group("ImarisDataSet       0.0.0")
        dataset_info_ver.attrs["NumberOfImages"] = _string_to_byte_array("1")
        dataset_info_ver2 = info_group.create_group("ImarisDataSet      0.0.0")
        dataset_info_ver2.attrs["NumberOfImages"] = _string_to_byte_array("1")

        # Create TimeInfo group
        time_info = info_group.create_group("TimeInfo")
        time_info.attrs["DatasetTimePoints"] = _string_to_byte_array("1")
        time_info.attrs["FileTimePoints"] = _string_to_byte_array("1")
        time_info.attrs["TimePoint1"] = _string_to_byte_array(
            "2024-01-01 00:00:00.000"
        )

        # Create Log group (basic processing log)
        log_group = info_group.create_group("Log")
        log_group.attrs["Entries"] = _string_to_byte_array("1")
        log_group.attrs["Entry0"] = _string_to_byte_array(
            f"<ZarrNiiExport channels=\"{' '.join(['on'] * n_channels)}\"/>"
        )

        # Create thumbnail group with proper multi-channel thumbnail
        thumbnail_group = f.create_group("Thumbnail")

        # Create a combined thumbnail (256x1024 for multi-channel as in reference)
        if n_channels > 1:
            # Multi-channel thumbnail: concatenate channels horizontally
            thumb_width = 256 * n_channels
            thumbnail_data = np.zeros((256, thumb_width), dtype=np.uint8)

            for c in range(n_channels):
                # Downsample each channel to 256x256
                channel_data = data[c]
                # Take MIP (Maximum Intensity Projection) along Z
                mip = np.max(channel_data, axis=0)
                # Resize to 256x256 (simple decimation)
                step_y = max(1, mip.shape[0] // 256)
                step_x = max(1, mip.shape[1] // 256)
                thumb_channel = mip[::step_y, ::step_x]

                # Pad or crop to exactly 256x256
                if thumb_channel.shape[0] < 256 or thumb_channel.shape[1] < 256:
                    padded = np.zeros((256, 256), dtype=thumb_channel.dtype)
                    h, w = thumb_channel.shape
                    padded[:h, :w] = thumb_channel
                    thumb_channel = padded
                else:
                    thumb_channel = thumb_channel[:256, :256]

                # Place in thumbnail
                thumbnail_data[:, c * 256 : (c + 1) * 256] = thumb_channel
        else:
            # Single channel: 256x256 thumbnail
            channel_data = data[0]
            mip = np.max(channel_data, axis=0)
            step_y = max(1, mip.shape[0] // 256)
            step_x = max(1, mip.shape[1] // 256)
            thumbnail_data = mip[::step_y, ::step_x]

            if thumbnail_data.shape[0] < 256 or thumbnail_data.shape[1] < 256:
                padded = np.zeros((256, 256), dtype=thumbnail_data.dtype)
                h, w = thumbnail_data.shape
                padded[:h, :w] = thumbnail_data
                thumbnail_data = padded
            else:
                thumbnail_data = thumbnail_data[:256, :256]

        thumbnail_group.create_dataset("Data", data=thumbnail_data.astype(np.uint8))

    return path

zarrnii.ZarrNii.copy()

Create a copy of this ZarrNii.

Returns:

• 'ZarrNii' –

  New ZarrNii with copied data

Source code in zarrnii/core.py

def copy(self) -> "ZarrNii":
    """
    Create a copy of this ZarrNii.

    Returns:
        New ZarrNii with copied data
    """
    # Create a new NgffImage with the same properties
    copied_image = nz.NgffImage(
        data=self.ngff_image.data,  # Dask arrays are lazy so this is efficient
        dims=self.ngff_image.dims.copy(),
        scale=self.ngff_image.scale.copy(),
        translation=self.ngff_image.translation.copy(),
        name=self.ngff_image.name,
    )
    return ZarrNii(
        ngff_image=copied_image,
        axes_order=self.axes_order,
        xyz_orientation=self.xyz_orientation,
        _omero=self._omero,
    )

zarrnii.ZarrNii.compute()

Compute the dask array and return the underlying NgffImage.

This triggers computation of any lazy operations and returns the NgffImage with computed data.

Returns:

• NgffImage –

  NgffImage with computed data

Source code in zarrnii/core.py

def compute(self) -> nz.NgffImage:
    """
    Compute the dask array and return the underlying NgffImage.

    This triggers computation of any lazy operations and returns
    the NgffImage with computed data.

    Returns:
        NgffImage with computed data
    """
    computed_data = self.ngff_image.data.compute()

    # Create new NgffImage with computed data
    computed_image = nz.NgffImage(
        data=computed_data,
        dims=self.ngff_image.dims,
        scale=self.ngff_image.scale,
        translation=self.ngff_image.translation,
        name=self.ngff_image.name,
    )
    return computed_image

zarrnii.ZarrNii.get_orientation()

Get the anatomical orientation of the dataset.

This function returns the orientation string (e.g., 'RAS', 'LPI') of the dataset.

Returns:

• str ( str ) –

  The orientation string corresponding to the dataset's anatomical orientation.

Source code in zarrnii/core.py

def get_orientation(self) -> str:
    """
    Get the anatomical orientation of the dataset.

    This function returns the orientation string (e.g., 'RAS', 'LPI') of the dataset.

    Returns:
        str: The orientation string corresponding to the dataset's anatomical orientation.
    """
    return self.orientation

zarrnii.ZarrNii.get_zooms(axes_order=None)

Get voxel spacing (zooms) from NgffImage scale.

Parameters:

• axes_order (str, default: None ) –

  Spatial axes order, defaults to self.axes_order

Returns:

• ndarray –

  Array of voxel spacings

Source code in zarrnii/core.py

def get_zooms(self, axes_order: str = None) -> np.ndarray:
    """
    Get voxel spacing (zooms) from NgffImage scale.

    Args:
        axes_order: Spatial axes order, defaults to self.axes_order

    Returns:
        Array of voxel spacings
    """
    if axes_order is None:
        axes_order = self.axes_order

    spatial_dims = ["z", "y", "x"] if axes_order == "ZYX" else ["x", "y", "z"]
    zooms = []

    for dim in spatial_dims:
        if dim in self.ngff_image.scale:
            zooms.append(self.ngff_image.scale[dim])
        else:
            zooms.append(1.0)

    return np.array(zooms)

zarrnii.ZarrNii.get_origin(axes_order=None)

Get origin (translation) from NgffImage.

Parameters:

• axes_order (str, default: None ) –

  Spatial axes order, defaults to self.axes_order

Returns:

• ndarray –

  Array of origin coordinates

Source code in zarrnii/core.py

def get_origin(self, axes_order: str = None) -> np.ndarray:
    """
    Get origin (translation) from NgffImage.

    Args:
        axes_order: Spatial axes order, defaults to self.axes_order

    Returns:
        Array of origin coordinates
    """
    if axes_order is None:
        axes_order = self.axes_order

    spatial_dims = ["z", "y", "x"] if axes_order == "ZYX" else ["x", "y", "z"]
    origin = []

    for dim in spatial_dims:
        if dim in self.ngff_image.translation:
            origin.append(self.ngff_image.translation[dim])
        else:
            origin.append(0.0)

    return np.array(origin)

zarrnii.ZarrNii.get_affine_matrix(axes_order=None)

Construct a 4x4 affine matrix from NGFF metadata (scale/translation), and align it to self.orientation (if provided) using nibabel.orientations.

Parameters:

• axes_order (str, default: None ) –

  Spatial axes order, e.g. 'ZYX' or 'XYZ'. Defaults to 'XYZ'.

Returns:

• ndarray –

  np.ndarray: 4x4 affine matrix.

Source code in zarrnii/core.py

def get_affine_matrix(self, axes_order: str = None) -> np.ndarray:
    """
    Construct a 4x4 affine matrix from NGFF metadata (scale/translation),
    and align it to self.orientation (if provided) using nibabel.orientations.

    Args:
        axes_order: Spatial axes order, e.g. 'ZYX' or 'XYZ'. Defaults to 'XYZ'.

    Returns:
        np.ndarray: 4x4 affine matrix.
    """
    if axes_order is None:
        axes_order = self.axes_order

    if axes_order == "ZYX":
        orientation = reverse_orientation_string(self.orientation)
    else:
        orientation = self.orientation

    # Safely pull scale/translation from metadata (dict-like expected)
    scale_meta = getattr(self.ngff_image, "scale", {}) or {}
    trans_meta = getattr(self.ngff_image, "translation", {}) or {}

    scale = np.ones(
        3,
    )
    trans = np.zeros(
        3,
    )

    for i, dim in enumerate(axes_order):
        s = scale_meta.get(dim.lower())
        if s is not None:
            scale[i] = float(s)

    for i, dim in enumerate(axes_order):
        s = trans_meta.get(dim.lower())
        if s is not None:
            trans[i] = float(s)

    affine = _axcodes2aff(orientation, scale=scale, translate=trans)

    return affine

zarrnii.ZarrNii.apply_transform_ref_to_flo_indices(*transforms, ref_znimg, indices)

Transform indices from reference to floating space.

Source code in zarrnii/core.py

def apply_transform_ref_to_flo_indices(self, *transforms, ref_znimg, indices):
    """Transform indices from reference to floating space."""
    # Placeholder implementation - would need full transform logic
    return indices

zarrnii.ZarrNii.apply_transform_flo_to_ref_indices(*transforms, ref_znimg, indices)

Transform indices from floating to reference space.

Source code in zarrnii/core.py

def apply_transform_flo_to_ref_indices(self, *transforms, ref_znimg, indices):
    """Transform indices from floating to reference space."""
    # Placeholder implementation - would need full transform logic
    return indices

zarrnii.ZarrNii.list_channels()

Get list of available channel labels from OMERO metadata.

Extracts channel labels from OMERO metadata if available, providing human-readable names for multi-channel datasets.

Returns:

• List[str] –

  List of channel label strings. Empty list if no OMERO metadata

• List[str] –

  is available or no channels are defined.

Examples:

>>> # Check available channels
>>> labels = znii.list_channels()
>>> print(f"Available channels: {labels}")
>>> # ['DAPI', 'GFP', 'RFP', 'Cy5']

>>> # Select specific channels by label
>>> selected = znii.select_channels(channel_labels=['DAPI', 'GFP'])

Notes

• Requires OMERO metadata to be present in the dataset
• Returns empty list for datasets without channel metadata
• Labels are extracted from the 'label' field of each channel

Source code in zarrnii/core.py

def list_channels(self) -> List[str]:
    """Get list of available channel labels from OMERO metadata.

    Extracts channel labels from OMERO metadata if available, providing
    human-readable names for multi-channel datasets.

    Returns:
        List of channel label strings. Empty list if no OMERO metadata
        is available or no channels are defined.

    Examples:
        >>> # Check available channels
        >>> labels = znii.list_channels()
        >>> print(f"Available channels: {labels}")
        >>> # ['DAPI', 'GFP', 'RFP', 'Cy5']

        >>> # Select specific channels by label
        >>> selected = znii.select_channels(channel_labels=['DAPI', 'GFP'])

    Notes:
        - Requires OMERO metadata to be present in the dataset
        - Returns empty list for datasets without channel metadata
        - Labels are extracted from the 'label' field of each channel
    """
    if self.omero is None or not hasattr(self.omero, "channels"):
        return []

    return [
        ch.label if hasattr(ch, "label") else ch.get("label", "")
        for ch in self.omero.channels
    ]

zarrnii.ZarrNii.select_channels(channels=None, channel_labels=None)

Select specific channels from multi-channel image data.

Creates a new ZarrNii instance containing only the specified channels, reducing memory usage and focusing analysis on channels of interest. Supports selection by both numeric indices and human-readable labels.

Parameters:

• channels (Optional[List[int]], default: None ) –

  List of 0-based channel indices to select. Mutually exclusive with channel_labels

• channel_labels (Optional[List[str]], default: None ) –

  List of channel names to select by label. Requires OMERO metadata. Mutually exclusive with channels

Returns:

• 'ZarrNii' –

  New ZarrNii instance with selected channels and updated metadata

Raises:

• ValueError –

  If both channels and channel_labels specified, or if channel_labels used without OMERO metadata, or if labels not found

• IndexError –

  If channel indices are out of range

Examples:

>>> # Select channels by index
>>> selected = znii.select_channels(channels=[0, 2])

>>> # Select channels by label (requires OMERO metadata)
>>> selected = znii.select_channels(channel_labels=['DAPI', 'GFP'])

>>> # Check available labels first
>>> available = znii.list_channels()
>>> print(f"Available: {available}")
>>> selected = znii.select_channels(channel_labels=available[:2])

Notes

• Preserves all spatial dimensions and timepoints
• Updates OMERO metadata to reflect selected channels
• Maintains spatial transformations and other metadata
• Channel order in output matches selection order

Source code in zarrnii/core.py

def select_channels(
    self,
    channels: Optional[List[int]] = None,
    channel_labels: Optional[List[str]] = None,
) -> "ZarrNii":
    """Select specific channels from multi-channel image data.

    Creates a new ZarrNii instance containing only the specified channels,
    reducing memory usage and focusing analysis on channels of interest.
    Supports selection by both numeric indices and human-readable labels.

    Args:
        channels: List of 0-based channel indices to select.
            Mutually exclusive with channel_labels
        channel_labels: List of channel names to select by label.
            Requires OMERO metadata. Mutually exclusive with channels

    Returns:
        New ZarrNii instance with selected channels and updated metadata

    Raises:
        ValueError: If both channels and channel_labels specified, or if
            channel_labels used without OMERO metadata, or if labels not found
        IndexError: If channel indices are out of range

    Examples:
        >>> # Select channels by index
        >>> selected = znii.select_channels(channels=[0, 2])

        >>> # Select channels by label (requires OMERO metadata)
        >>> selected = znii.select_channels(channel_labels=['DAPI', 'GFP'])

        >>> # Check available labels first
        >>> available = znii.list_channels()
        >>> print(f"Available: {available}")
        >>> selected = znii.select_channels(channel_labels=available[:2])

    Notes:
        - Preserves all spatial dimensions and timepoints
        - Updates OMERO metadata to reflect selected channels
        - Maintains spatial transformations and other metadata
        - Channel order in output matches selection order
    """
    if channels is not None and channel_labels is not None:
        raise ValueError("Cannot specify both 'channels' and 'channel_labels'")

    if channel_labels is not None:
        if self.omero is None:
            raise ValueError(
                "Channel labels were specified but no omero metadata found"
            )

        available_labels = self.list_channels()
        channel_indices = []
        for label in channel_labels:
            if label not in available_labels:
                raise ValueError(f"Channel label '{label}' not found")
            channel_indices.append(available_labels.index(label))
        channels = channel_indices

    if channels is None:
        # Return a copy with all channels
        return self.copy()

    # Check if channel dimension exists
    if "c" not in self.dims:
        raise ValueError("No channel dimension found in the data")

    # Get channel dimension index
    c_idx = self.dims.index("c")

    # Create slice objects for proper dimension indexing
    slices = [slice(None)] * len(self.data.shape)
    slices[c_idx] = channels

    # Select channels from data using proper dimension indexing
    selected_data = self.data[tuple(slices)]

    # Create new NgffImage with selected data
    new_ngff_image = nz.NgffImage(
        data=selected_data,
        dims=self.dims,
        scale=self.scale,
        translation=self.translation,
        name=self.name,
    )

    # Filter omero metadata to match selected channels
    filtered_omero = None
    if self.omero is not None and hasattr(self.omero, "channels"):

        class FilteredOmero:
            def __init__(self, channels):
                self.channels = channels

        filtered_channels = [self.omero.channels[i] for i in channels]
        filtered_omero = FilteredOmero(filtered_channels)

    return ZarrNii(
        ngff_image=new_ngff_image,
        axes_order=self.axes_order,
        xyz_orientation=self.xyz_orientation,
        _omero=filtered_omero,
    )

zarrnii.ZarrNii.select_timepoints(timepoints=None)

Select timepoints from the image data and return a new ZarrNii instance.

Parameters:

• timepoints (Optional[List[int]], default: None ) –

  Timepoint indices to select

Returns:

• 'ZarrNii' –

  New ZarrNii instance with selected timepoints

Source code in zarrnii/core.py

def select_timepoints(self, timepoints: Optional[List[int]] = None) -> "ZarrNii":
    """
    Select timepoints from the image data and return a new ZarrNii instance.

    Args:
        timepoints: Timepoint indices to select

    Returns:
        New ZarrNii instance with selected timepoints
    """
    if timepoints is None:
        # Return a copy with all timepoints
        return self.copy()

    # Check if time dimension exists
    if "t" not in self.dims:
        raise ValueError("No time dimension found in the data")

    # Get time dimension index
    t_idx = self.dims.index("t")

    # Create slice objects
    slices = [slice(None)] * len(self.data.shape)
    slices[t_idx] = timepoints

    # Select timepoints from data
    selected_data = self.data[tuple(slices)]

    # Create new NgffImage with selected data
    new_ngff_image = nz.NgffImage(
        data=selected_data,
        dims=self.dims,
        scale=self.scale,
        translation=self.translation,
        name=self.name,
    )

    return ZarrNii(
        ngff_image=new_ngff_image,
        axes_order=self.axes_order,
        xyz_orientation=self.xyz_orientation,
        _omero=self._omero,  # Timepoint selection doesn't affect omero metadata
    )

zarrnii.ZarrNii.to_ngff_image(name=None)

Convert to NgffImage object.

Parameters:

• name (str, default: None ) –

  Optional name for the image

Returns:

• NgffImage –

  NgffImage representation

Source code in zarrnii/core.py

def to_ngff_image(self, name: str = None) -> nz.NgffImage:
    """
    Convert to NgffImage object.

    Args:
        name: Optional name for the image

    Returns:
        NgffImage representation
    """
    if name is None:
        name = self.name

    return nz.NgffImage(
        data=self.data,
        dims=self.dims,
        scale=self.scale,
        translation=self.translation,
        name=name,
    )

zarrnii.ZarrNii.segment(plugin, chunk_size=None, **kwargs)

Apply segmentation plugin to the image using blockwise processing.

This method applies a segmentation plugin to the image data using dask's blockwise processing for efficient computation on large datasets.

Parameters:

• plugin –

  Segmentation plugin instance or class to apply

• chunk_size (Optional[Tuple[int, ...]], default: None ) –

  Optional chunk size for dask processing. If None, uses current chunks.

• **kwargs –

  Additional arguments passed to the plugin

Returns:

• 'ZarrNii' –

  New ZarrNii instance with segmented data as labels

Source code in zarrnii/core.py

def segment(
    self, plugin, chunk_size: Optional[Tuple[int, ...]] = None, **kwargs
) -> "ZarrNii":
    """
    Apply segmentation plugin to the image using blockwise processing.

    This method applies a segmentation plugin to the image data using dask's
    blockwise processing for efficient computation on large datasets.

    Args:
        plugin: Segmentation plugin instance or class to apply
        chunk_size: Optional chunk size for dask processing. If None, uses current chunks.
        **kwargs: Additional arguments passed to the plugin

    Returns:
        New ZarrNii instance with segmented data as labels
    """
    from .plugins.segmentation import SegmentationPlugin

    # Handle plugin instance or class
    if isinstance(plugin, type) and issubclass(plugin, SegmentationPlugin):
        plugin = plugin(**kwargs)
    elif not isinstance(plugin, SegmentationPlugin):
        raise TypeError(
            "Plugin must be an instance or subclass of SegmentationPlugin"
        )

    # Prepare chunk size
    if chunk_size is not None:
        # Rechunk the data if different chunk size requested
        data = self.data.rechunk(chunk_size)
    else:
        data = self.data

    # Create metadata dict to pass to plugin
    metadata = {
        "axes_order": self.axes_order,
        "orientation": self.xyz_orientation,
        "shape": self.shape,
        "dims": self.dims,
        "scale": self.scale,
        "translation": self.translation,
    }

    # Create a wrapper function for map_blocks
    def segment_block(block):
        """Wrapper function to apply segmentation to a single block."""
        # Handle single blocks
        return plugin.segment(block, metadata)

    # Apply segmentation using dask map_blocks
    segmented_data = da.map_blocks(
        segment_block,
        data,
        dtype=np.uint8,  # Segmentation results are typically uint8
        meta=np.array([], dtype=np.uint8),  # Provide meta information
    )

    # Create new NgffImage with segmented data
    new_ngff_image = nz.NgffImage(
        data=segmented_data,
        dims=self.dims.copy(),
        scale=self.scale.copy(),
        translation=self.translation.copy(),
        name=f"{self.name}_segmented_{plugin.name.lower().replace(' ', '_')}",
    )

    # Return new ZarrNii instance
    return ZarrNii(
        ngff_image=new_ngff_image,
        axes_order=self.axes_order,
        xyz_orientation=self.xyz_orientation,
    )

zarrnii.ZarrNii.segment_otsu(nbins=256, chunk_size=None)

Apply local Otsu thresholding segmentation to the image.

Convenience method for local Otsu thresholding segmentation. This computes the threshold locally for each processing block.

Parameters:

• nbins (int, default: 256 ) –

  Number of bins for histogram computation (default: 256)

• chunk_size (Optional[Tuple[int, ...]], default: None ) –

  Optional chunk size for dask processing

Returns:

• 'ZarrNii' –

  New ZarrNii instance with binary segmentation

Source code in zarrnii/core.py

def segment_otsu(
    self, nbins: int = 256, chunk_size: Optional[Tuple[int, ...]] = None
) -> "ZarrNii":
    """
    Apply local Otsu thresholding segmentation to the image.

    Convenience method for local Otsu thresholding segmentation.
    This computes the threshold locally for each processing block.

    Args:
        nbins: Number of bins for histogram computation (default: 256)
        chunk_size: Optional chunk size for dask processing

    Returns:
        New ZarrNii instance with binary segmentation
    """
    from .plugins.segmentation import LocalOtsuSegmentation

    plugin = LocalOtsuSegmentation(nbins=nbins)
    return self.segment(plugin, chunk_size=chunk_size)

zarrnii.ZarrNii.segment_threshold(thresholds, inclusive=True, chunk_size=None)

Apply threshold-based segmentation to the image.

Convenience method for threshold-based segmentation using either manual threshold values or computed thresholds.

Parameters:

• thresholds (Union[float, List[float]]) –

  Single threshold value or list of threshold values. For single threshold, creates binary segmentation (0/1). For multiple thresholds, creates multi-class segmentation (0/1/2/...).

• inclusive (bool, default: True ) –

  Whether thresholds are inclusive (default: True). If True, pixels >= threshold are labeled as foreground. If False, pixels > threshold are labeled as foreground.

• chunk_size (Optional[Tuple[int, ...]], default: None ) –

  Optional chunk size for dask processing

Returns:

• 'ZarrNii' –

  New ZarrNii instance with labeled segmentation

Examples:

>>> # Binary threshold segmentation
>>> segmented = znimg.segment_threshold(0.5)
>>>
>>> # Multi-level threshold segmentation
>>> thresholds = znimg.compute_otsu_thresholds(classes=3)
>>> segmented = znimg.segment_threshold(thresholds[1:-1])  # Exclude min/max values

Source code in zarrnii/core.py

def segment_threshold(
    self,
    thresholds: Union[float, List[float]],
    inclusive: bool = True,
    chunk_size: Optional[Tuple[int, ...]] = None,
) -> "ZarrNii":
    """
    Apply threshold-based segmentation to the image.

    Convenience method for threshold-based segmentation using either
    manual threshold values or computed thresholds.

    Args:
        thresholds: Single threshold value or list of threshold values.
            For single threshold, creates binary segmentation (0/1).
            For multiple thresholds, creates multi-class segmentation (0/1/2/...).
        inclusive: Whether thresholds are inclusive (default: True).
            If True, pixels >= threshold are labeled as foreground.
            If False, pixels > threshold are labeled as foreground.
        chunk_size: Optional chunk size for dask processing

    Returns:
        New ZarrNii instance with labeled segmentation

    Examples:
        >>> # Binary threshold segmentation
        >>> segmented = znimg.segment_threshold(0.5)
        >>>
        >>> # Multi-level threshold segmentation
        >>> thresholds = znimg.compute_otsu_thresholds(classes=3)
        >>> segmented = znimg.segment_threshold(thresholds[1:-1])  # Exclude min/max values
    """
    from .plugins.segmentation import ThresholdSegmentation

    plugin = ThresholdSegmentation(thresholds=thresholds, inclusive=inclusive)
    return self.segment(plugin, chunk_size=chunk_size)

zarrnii.ZarrNii.compute_histogram(bins=256, range=None, mask=None, **kwargs)

Compute histogram of the image.

This method computes the histogram of image intensities, optionally using a mask to weight the computation. The histogram is computed using dask for efficient processing of large datasets.

Parameters:

• bins (int, default: 256 ) –

  Number of histogram bins (default: 256)

• range (Optional[Tuple[float, float]], default: None ) –

  Optional tuple (min, max) defining histogram range. If None, uses the full range of the data

• mask (Optional['ZarrNii'], default: None ) –

  Optional ZarrNii mask of same shape as image. Only pixels where mask > 0 are included in histogram computation

• **kwargs (Any, default: {} ) –

  Additional arguments passed to dask.array.histogram

Returns:

• Array –

  Tuple of (histogram_counts, bin_edges) where:

• Array –
  • histogram_counts: dask array of histogram bin counts
• Tuple[Array, Array] –
  • bin_edges: dask array of bin edge values (length = bins + 1)

Examples:

>>> # Compute histogram
>>> hist, bin_edges = znimg.compute_histogram(bins=128)
>>>
>>> # Compute histogram with mask
>>> mask = znimg > 0.5
>>> hist_masked, _ = znimg.compute_histogram(mask=mask)

Source code in zarrnii/core.py

def compute_histogram(
    self,
    bins: int = 256,
    range: Optional[Tuple[float, float]] = None,
    mask: Optional["ZarrNii"] = None,
    **kwargs: Any,
) -> Tuple[da.Array, da.Array]:
    """
    Compute histogram of the image.

    This method computes the histogram of image intensities, optionally using
    a mask to weight the computation. The histogram is computed using dask for
    efficient processing of large datasets.

    Args:
        bins: Number of histogram bins (default: 256)
        range: Optional tuple (min, max) defining histogram range. If None,
            uses the full range of the data
        mask: Optional ZarrNii mask of same shape as image. Only pixels
            where mask > 0 are included in histogram computation
        **kwargs: Additional arguments passed to dask.array.histogram

    Returns:
        Tuple of (histogram_counts, bin_edges) where:
        - histogram_counts: dask array of histogram bin counts
        - bin_edges: dask array of bin edge values (length = bins + 1)

    Examples:
        >>> # Compute histogram
        >>> hist, bin_edges = znimg.compute_histogram(bins=128)
        >>>
        >>> # Compute histogram with mask
        >>> mask = znimg > 0.5
        >>> hist_masked, _ = znimg.compute_histogram(mask=mask)
    """
    from .analysis import compute_histogram

    mask_data = mask.darr if mask is not None else None
    return compute_histogram(
        self.darr, bins=bins, range=range, mask=mask_data, **kwargs
    )

zarrnii.ZarrNii.compute_otsu_thresholds(classes=2, bins=256, range=None, mask=None, return_figure=False)

Compute Otsu multi-level thresholds for the image.

This method first computes the histogram of the image, then uses scikit-image's threshold_multiotsu to compute optimal threshold values.

Parameters:

• classes (int, default: 2 ) –

  Number of classes to separate data into (default: 2). Must be >= 2. For classes=2, returns 1 threshold. For classes=k, returns k-1 thresholds.

• bins (int, default: 256 ) –

  Number of histogram bins (default: 256)

• range (Optional[Tuple[float, float]], default: None ) –

  Optional tuple (min, max) defining histogram range. If None, uses the full range of the data

• mask (Optional['ZarrNii'], default: None ) –

  Optional ZarrNii mask of same shape as image. Only pixels where mask > 0 are included in histogram computation

• return_figure (bool, default: False ) –

  If True, returns a tuple containing thresholds and a matplotlib figure with the histogram and annotated threshold lines (default: False).

Returns:

• Union[List[float], Tuple[List[float], Any]] –

  If return_figure is False (default): List of threshold values. For classes=k, returns k+1 values: [0, threshold1, threshold2, ..., threshold_k-1, max_intensity] where 0 represents the minimum and max_intensity represents the maximum.

• Union[List[float], Tuple[List[float], Any]] –

  If return_figure is True: Tuple of (thresholds, figure) where figure is a matplotlib Figure object showing the histogram with annotated threshold lines.

Examples:

>>> # Compute binary threshold (2 classes)
>>> thresholds = znimg.compute_otsu_thresholds(classes=2)
>>> print(f"Binary thresholds: {thresholds}")
>>>
>>> # Compute multi-level thresholds (3 classes)
>>> thresholds = znimg.compute_otsu_thresholds(classes=3)
>>> print(f"Multi-level thresholds: {thresholds}")
>>>
>>> # Get histogram data along with thresholds
>>> thresholds, (hist, bin_edges) = znimg.compute_otsu_thresholds(
...     classes=2, return_histogram=True
... )
>>>
>>> # Generate a figure with annotated thresholds
>>> thresholds, fig = znimg.compute_otsu_thresholds(
...     classes=2, return_figure=True
... )
>>> fig.savefig('otsu_thresholds.png')

Source code in zarrnii/core.py

def compute_otsu_thresholds(
    self,
    classes: int = 2,
    bins: int = 256,
    range: Optional[Tuple[float, float]] = None,
    mask: Optional["ZarrNii"] = None,
    return_figure: bool = False,
) -> Union[
    List[float],
    Tuple[List[float], Any],
]:
    """
    Compute Otsu multi-level thresholds for the image.

    This method first computes the histogram of the image, then uses
    scikit-image's threshold_multiotsu to compute optimal threshold values.

    Args:
        classes: Number of classes to separate data into (default: 2).
            Must be >= 2. For classes=2, returns 1 threshold. For classes=k,
            returns k-1 thresholds.
        bins: Number of histogram bins (default: 256)
        range: Optional tuple (min, max) defining histogram range. If None,
            uses the full range of the data
        mask: Optional ZarrNii mask of same shape as image. Only pixels
            where mask > 0 are included in histogram computation
        return_figure: If True, returns a tuple containing thresholds and a
            matplotlib figure with the histogram and annotated threshold lines
            (default: False).

    Returns:
        If return_figure is False (default):
            List of threshold values. For classes=k, returns k+1 values:
            [0, threshold1, threshold2, ..., threshold_k-1, max_intensity]
            where 0 represents the minimum and max_intensity represents the maximum.

        If return_figure is True:
            Tuple of (thresholds, figure) where figure is a matplotlib Figure
            object showing the histogram with annotated threshold lines.

    Examples:
        >>> # Compute binary threshold (2 classes)
        >>> thresholds = znimg.compute_otsu_thresholds(classes=2)
        >>> print(f"Binary thresholds: {thresholds}")
        >>>
        >>> # Compute multi-level thresholds (3 classes)
        >>> thresholds = znimg.compute_otsu_thresholds(classes=3)
        >>> print(f"Multi-level thresholds: {thresholds}")
        >>>
        >>> # Get histogram data along with thresholds
        >>> thresholds, (hist, bin_edges) = znimg.compute_otsu_thresholds(
        ...     classes=2, return_histogram=True
        ... )
        >>>
        >>> # Generate a figure with annotated thresholds
        >>> thresholds, fig = znimg.compute_otsu_thresholds(
        ...     classes=2, return_figure=True
        ... )
        >>> fig.savefig('otsu_thresholds.png')
    """
    from .analysis import compute_otsu_thresholds

    # First compute histogram
    hist, bin_edges = self.compute_histogram(bins=bins, range=range, mask=mask)

    # Then compute thresholds with optional returns
    return compute_otsu_thresholds(
        hist,
        classes=classes,
        bin_edges=bin_edges,
        return_figure=return_figure,
    )

zarrnii.ZarrNii.apply_scaled_processing(plugin, downsample_factor=4, chunk_size=None, upsampled_ome_zarr_path=None, **kwargs)

Apply scaled processing plugin using multi-resolution approach.

This method implements a multi-resolution processing pipeline where: 1. The image is downsampled for efficient computation 2. The plugin's lowres_func is applied to the downsampled data 3. The result is upsampled using dask-based upsampling 4. The plugin's highres_func applies the result to full-resolution data

Parameters:

• plugin –

  ScaledProcessingPlugin instance or class to apply

• downsample_factor (int, default: 4 ) –

  Factor for downsampling (default: 4)

• chunk_size (Optional[Tuple[int, ...]], default: None ) –

  Optional chunk size for low-res processing. If None, uses (1, 10, 10, 10).

• upsampled_ome_zarr_path (Optional[str], default: None ) –

  Path to save intermediate OME-Zarr, default saved in system temp directory.

• **kwargs –

  Additional arguments passed to the plugin

Returns:

• 'ZarrNii' –

  New ZarrNii instance with processed data

Source code in zarrnii/core.py

def apply_scaled_processing(
    self,
    plugin,
    downsample_factor: int = 4,
    chunk_size: Optional[Tuple[int, ...]] = None,
    upsampled_ome_zarr_path: Optional[str] = None,
    **kwargs,
) -> "ZarrNii":
    """
    Apply scaled processing plugin using multi-resolution approach.

    This method implements a multi-resolution processing pipeline where:
    1. The image is downsampled for efficient computation
    2. The plugin's lowres_func is applied to the downsampled data
    3. The result is upsampled using dask-based upsampling
    4. The plugin's highres_func applies the result to full-resolution data

    Args:
        plugin: ScaledProcessingPlugin instance or class to apply
        downsample_factor: Factor for downsampling (default: 4)
        chunk_size: Optional chunk size for low-res processing. If None, uses (1, 10, 10, 10).
        upsampled_ome_zarr_path: Path to save intermediate OME-Zarr, default saved in system temp directory.
        **kwargs: Additional arguments passed to the plugin

    Returns:
        New ZarrNii instance with processed data
    """
    from .plugins.scaled_processing import ScaledProcessingPlugin

    # Handle plugin instance or class
    if isinstance(plugin, type) and issubclass(plugin, ScaledProcessingPlugin):
        plugin = plugin(**kwargs)
    elif not isinstance(plugin, ScaledProcessingPlugin):
        raise TypeError(
            "Plugin must be an instance or subclass of ScaledProcessingPlugin"
        )

    # Step 1: Downsample the data for low-resolution processing
    lowres_znimg = self.downsample(level=int(np.log2(downsample_factor)))

    # Convert to numpy array for lowres processing
    lowres_array = lowres_znimg.data.compute()

    # Step 2: Apply low-resolution function and prepare for upsampling
    # Use chunk_size parameter for the low-res processing chunks
    lowres_chunks = chunk_size if chunk_size is not None else (1, 10, 10, 10)
    lowres_znimg.data = da.from_array(
        plugin.lowres_func(lowres_array), chunks=lowres_chunks
    )

    # Use temporary OME-Zarr to break up dask graph for performance
    import tempfile

    if upsampled_ome_zarr_path is None:
        upsampled_ome_zarr_path = tempfile.mkdtemp(suffix="_SPIM.ome.zarr")

    # Step 3: Upsample using dask-based upsampling, save to ome zarr
    lowres_znimg.upsample(to_shape=self.shape).to_ome_zarr(
        upsampled_ome_zarr_path, max_layer=1
    )

    upsampled_znimg = ZarrNii.from_ome_zarr(upsampled_ome_zarr_path)

    corrected_znimg = self.copy()

    # Step 4: Apply high-resolution function
    # rechunk original data to use same chunksize as upsampled_data, before multiplying
    corrected_znimg.data = plugin.highres_func(
        self.data.rechunk(upsampled_znimg.data.chunks), upsampled_znimg.data
    )

    return corrected_znimg

Bases: ZarrNii

Brain atlas with segmentation image and region lookup table.

Represents a brain atlas consisting of a segmentation image (dseg) that assigns integer labels to brain regions, and a lookup table (tsv) that maps these labels to region names and other metadata.

Extension of ZarrNii to support atlas label tables.

Inherits all functionality from ZarrNii and adds support for storing region/label metadata in a pandas DataFrame.

Attributes

labels_df : pandas.DataFrame DataFrame containing label information for the atlas. label_column : str Name of the column in labels_df containing label indices. name_column : str Name of the column in labels_df containing region names. abbrev_column : str Name of the column in labels_df containing region abbreviations.

Source code in zarrnii/core.py

def __init__(
    self,
    darr=None,
    affine=None,
    axes_order="ZYX",
    orientation="RAS",
    xyz_orientation=None,
    ngff_image=None,
    _omero=None,
    **kwargs,
):
    """
    Constructor with backward compatibility for old signature.
    """
    # Handle backwards compatibility: if xyz_orientation is provided, use it
    # Otherwise, use orientation for backwards compatibility
    final_orientation = (
        xyz_orientation if xyz_orientation is not None else orientation
    )

    if ngff_image is not None:
        # New signature
        object.__setattr__(self, "ngff_image", ngff_image)
        object.__setattr__(self, "axes_order", axes_order)
        object.__setattr__(self, "xyz_orientation", final_orientation)
        object.__setattr__(self, "_omero", _omero)
    elif darr is not None:
        # Legacy signature - delegate to from_darr
        instance = self.from_darr(
            darr=darr,
            affine=affine,
            axes_order=axes_order,
            orientation=final_orientation,
            **kwargs,
        )
        object.__setattr__(self, "ngff_image", instance.ngff_image)
        object.__setattr__(self, "axes_order", instance.axes_order)
        object.__setattr__(self, "xyz_orientation", instance.xyz_orientation)
        object.__setattr__(self, "_omero", instance._omero)
    else:
        raise ValueError("Must provide either ngff_image or darr")

Attributes

zarrnii.ZarrNiiAtlas.dseg property

Return self as the segmentation image (for compatibility with API).

Functions

zarrnii.ZarrNiiAtlas.create_from_dseg(dseg, labels_df, **kwargs) classmethod

Create ZarrNiiAtlas from a dseg ZarrNii and labels DataFrame.

Parameters:

• dseg (ZarrNii) –

  ZarrNii segmentation image

• labels_df (DataFrame) –

  DataFrame containing label information

• **kwargs –

  Additional keyword arguments for label/name/abbrev columns

Returns:

• –

  ZarrNiiAtlas instance

Source code in zarrnii/atlas.py

@classmethod
def create_from_dseg(cls, dseg: ZarrNii, labels_df: pd.DataFrame, **kwargs):
    """Create ZarrNiiAtlas from a dseg ZarrNii and labels DataFrame.

    Args:
        dseg: ZarrNii segmentation image
        labels_df: DataFrame containing label information
        **kwargs: Additional keyword arguments for label/name/abbrev columns

    Returns:
        ZarrNiiAtlas instance
    """
    if not isinstance(dseg, ZarrNii):
        raise TypeError(f"dseg must be a ZarrNii instance, got {type(dseg)}")

    # Note: attrs strips leading underscore from _omero in __init__ signature
    # so we pass it as 'omero' instead of '_omero'
    return cls(
        ngff_image=dseg.ngff_image,
        axes_order=dseg.axes_order,
        xyz_orientation=dseg.xyz_orientation,
        omero=getattr(dseg, "_omero", None),
        labels_df=labels_df,
        **kwargs,
    )

zarrnii.ZarrNiiAtlas.from_files(dseg_path, labels_path, **kwargs) classmethod

Load ZarrNiiAtlas from dseg image and labels TSV files.

Parameters:

• dseg_path (Union[str, Path]) –

  Path to segmentation image (NIfTI or OME-Zarr)

• labels_path (Union[str, Path]) –

  Path to labels TSV file

• **kwargs –

  Additional arguments passed to ZarrNii.from_file()

Returns:

• –

  ZarrNiiAtlas instance

Source code in zarrnii/atlas.py

@classmethod
def from_files(
    cls, dseg_path: Union[str, Path], labels_path: Union[str, Path], **kwargs
):
    """Load ZarrNiiAtlas from dseg image and labels TSV files.

    Args:
        dseg_path: Path to segmentation image (NIfTI or OME-Zarr)
        labels_path: Path to labels TSV file
        **kwargs: Additional arguments passed to ZarrNii.from_file()

    Returns:
        ZarrNiiAtlas instance
    """
    # Load segmentation image
    dseg = ZarrNii.from_file(str(dseg_path), **kwargs)

    # Load labels dataframe
    labels_df = pd.read_csv(str(labels_path), sep="\t")

    # Create atlas instance using create_from_dseg
    return cls.create_from_dseg(dseg, labels_df)

zarrnii.ZarrNiiAtlas.from_itksnap_lut(path, lut_path, **kwargs) classmethod

Construct from itksnap lut file.

Source code in zarrnii/atlas.py

@classmethod
def from_itksnap_lut(cls, path, lut_path, **kwargs):
    """
    Construct from itksnap lut file.
    """
    znii = super().from_file(path, **kwargs)
    labels_df = cls._import_itksnap_lut(lut_path)
    print(labels_df)
    return cls(
        ngff_image=znii.ngff_image,
        axes_order=znii.axes_order,
        xyz_orientation=znii.xyz_orientation,
        labels_df=labels_df,
        omero=getattr(znii, "_omero", None),
    )

zarrnii.ZarrNiiAtlas.from_csv_lut(path, lut_path, **kwargs) classmethod

Construct from csv lut file.

Source code in zarrnii/atlas.py

@classmethod
def from_csv_lut(cls, path, lut_path, **kwargs):
    """
    Construct from csv lut file.
    """
    znii = super().from_file(path, **kwargs)
    labels_df = cls._import_csv_lut(lut_path)
    return cls(
        ngff_image=znii.ngff_image,
        axes_order=znii.axes_order,
        xyz_orientation=znii.xyz_orientation,
        labels_df=labels_df,
        omero=getattr(znii, "_omero", None),
    )

zarrnii.ZarrNiiAtlas.from_tsv_lut(path, lut_path, **kwargs) classmethod

Construct from tsv lut file.

Source code in zarrnii/atlas.py

@classmethod
def from_tsv_lut(cls, path, lut_path, **kwargs):
    """
    Construct from tsv lut file.
    """
    znii = super().from_file(path, **kwargs)
    labels_df = cls._import_tsv_lut(lut_path)
    return cls(
        ngff_image=znii.ngff_image,
        axes_order=znii.axes_order,
        xyz_orientation=znii.xyz_orientation,
        labels_df=labels_df,
        omero=getattr(znii, "_omero", None),
    )

zarrnii.ZarrNiiAtlas.from_labelmapper_lut(path, lut_path, **kwargs) classmethod

Construct from labelmapper lut file.

Source code in zarrnii/atlas.py

@classmethod
def from_labelmapper_lut(cls, path, lut_path, **kwargs):
    """
    Construct from labelmapper lut file.
    """
    znii = super().from_file(path, **kwargs)
    labels_df = cls._import_labelmapper_lut(lut_path)
    return cls(
        ngff_image=znii.ngff_image,
        axes_order=znii.axes_order,
        xyz_orientation=znii.xyz_orientation,
        labels_df=labels_df,
        omero=getattr(znii, "_omero", None),
    )

zarrnii.ZarrNiiAtlas.get_region_info(region_id)

Get information about a specific region.

Parameters:

• region_id (Union[int, str]) –

  Region identifier (int label, name, or abbreviation)

Returns:

• Dict[str, Any] –

  Dictionary containing region information

Raises:

• ValueError –

  If region not found in atlas

Source code in zarrnii/atlas.py

def get_region_info(self, region_id: Union[int, str]) -> Dict[str, Any]:
    """Get information about a specific region.

    Args:
        region_id: Region identifier (int label, name, or abbreviation)

    Returns:
        Dictionary containing region information

    Raises:
        ValueError: If region not found in atlas
    """
    label = self._resolve_region_identifier(region_id)

    # Find the region in labels DataFrame
    region_row = self.labels_df[self.labels_df[self.label_column] == label]
    if region_row.empty:
        raise ValueError(f"Region with label {label} not found in atlas")

    return region_row.iloc[0].to_dict()

zarrnii.ZarrNiiAtlas.get_region_mask(region_id)

Create binary mask for a specific region.

Parameters:

• region_id (Union[int, str]) –

  Region identifier (int label, name, or abbreviation)

Returns:

• ZarrNii –

  ZarrNii instance containing binary mask (1 for region, 0 elsewhere)

Raises:

• ValueError –

  If region not found in atlas

Source code in zarrnii/atlas.py

def get_region_mask(self, region_id: Union[int, str]) -> ZarrNii:
    """Create binary mask for a specific region.

    Args:
        region_id: Region identifier (int label, name, or abbreviation)

    Returns:
        ZarrNii instance containing binary mask (1 for region, 0 elsewhere)

    Raises:
        ValueError: If region not found in atlas
    """
    label = self._resolve_region_identifier(region_id)

    # Validate that the region exists in our labels_df
    if not (self.labels_df[self.label_column] == label).any():
        raise ValueError(f"Region with label {label} not found in atlas")

    # Create binary mask
    mask_data = (self.dseg.data == label).astype(np.uint8)

    return ZarrNii.from_darr(
        mask_data, affine=self.dseg.affine, axes_order=self.dseg.axes_order
    )

zarrnii.ZarrNiiAtlas.get_region_volume(region_id)

Calculate volume of a specific region in mm³.

Parameters:

• region_id (Union[int, str]) –

  Region identifier (int label, name, or abbreviation)

Returns:

• float –

  Volume in cubic millimeters

Raises:

• ValueError –

  If region not found in atlas

Source code in zarrnii/atlas.py

def get_region_volume(self, region_id: Union[int, str]) -> float:
    """Calculate volume of a specific region in mm³.

    Args:
        region_id: Region identifier (int label, name, or abbreviation)

    Returns:
        Volume in cubic millimeters

    Raises:
        ValueError: If region not found in atlas
    """
    label = self._resolve_region_identifier(region_id)

    # Count voxels with this label
    dseg_data = self.dseg.data
    if hasattr(dseg_data, "compute"):
        voxel_count = int((dseg_data == label).sum().compute())
    else:
        voxel_count = int((dseg_data == label).sum())

    # Calculate volume using voxel size from affine
    # Volume per voxel = abs(det(affine[:3, :3]))
    voxel_volume = abs(np.linalg.det(self.dseg.affine[:3, :3]))

    return float(voxel_count * voxel_volume)

zarrnii.ZarrNiiAtlas.aggregate_image_by_regions(image, aggregation_func='mean', background_label=0, column_suffix='value')

Aggregate image values by atlas regions.

Parameters:

• image (ZarrNii) –

  Image to aggregate (must be compatible with atlas)

• aggregation_func (str, default: 'mean' ) –

  Aggregation function ('mean', 'sum', 'std', 'median', 'min', 'max')

• background_label (int, default: 0 ) –

  Label value to treat as background (excluded from results)

• column_suffix (str, default: 'value' ) –

  String suffix to append to column name. Name will be {agg_func}_{col_suffix}.

Returns: DataFrame with columns: index, name, {aggregation_func}_{column_suffix}, volume_mm3 (e.g., with defaults: index, name, mean_value, volume_mm3)

Raises:

• ValueError –

  If image and atlas are incompatible

Source code in zarrnii/atlas.py

def aggregate_image_by_regions(
    self,
    image: ZarrNii,
    aggregation_func: str = "mean",
    background_label: int = 0,
    column_suffix: str = "value",
) -> pd.DataFrame:
    """Aggregate image values by atlas regions.

    Args:
        image: Image to aggregate (must be compatible with atlas)
        aggregation_func: Aggregation function ('mean', 'sum', 'std', 'median', 'min', 'max')
        background_label: Label value to treat as background (excluded from results)
        column_suffix: String suffix to append to column name. Name will be {agg_func}_{col_suffix}.
    Returns:
        DataFrame with columns: index, name, {aggregation_func}_{column_suffix}, volume_mm3
        (e.g., with defaults: index, name, mean_value, volume_mm3)

    Raises:
        ValueError: If image and atlas are incompatible
    """
    # Validate image compatibility
    if not np.array_equal(image.shape, self.dseg.shape):
        raise ValueError(
            f"Image shape {image.shape} doesn't match atlas shape {self.dseg.shape}"
        )

    if not np.allclose(image.affine, self.dseg.affine, atol=1e-6):
        warnings.warn(
            "Image and atlas affines don't match exactly. "
            "Results may be spatially inconsistent."
        )

    # Get all unique labels (excluding background)
    dseg_data = self.dseg.data
    if hasattr(dseg_data, "compute"):
        dseg_data = dseg_data.compute()
    unique_labels = np.unique(dseg_data)
    unique_labels = unique_labels[unique_labels != background_label]

    results = []
    for label in unique_labels:
        # Create mask for this region
        mask = self.dseg.data == label

        # Extract image values for this region
        region_values = image.data[mask]

        # Skip if no voxels (shouldn't happen with unique labels)
        if region_values.size == 0:
            continue

        # Compute aggregation
        if hasattr(region_values, "compute"):
            # Dask array - need to compute
            if aggregation_func == "mean":
                agg_value = float(region_values.mean().compute())
            elif aggregation_func == "sum":
                agg_value = float(region_values.sum().compute())
            elif aggregation_func == "std":
                agg_value = float(region_values.std().compute())
            elif aggregation_func == "median":
                agg_value = float(np.median(region_values.compute()))
            elif aggregation_func == "min":
                agg_value = float(region_values.min().compute())
            elif aggregation_func == "max":
                agg_value = float(region_values.max().compute())
            else:
                raise ValueError(
                    f"Unknown aggregation function: {aggregation_func}. "
                    "Supported: mean, sum, std, median, min, max"
                )
        else:
            # NumPy array - direct computation
            if aggregation_func == "mean":
                agg_value = float(region_values.mean())
            elif aggregation_func == "sum":
                agg_value = float(region_values.sum())
            elif aggregation_func == "std":
                agg_value = float(region_values.std())
            elif aggregation_func == "median":
                agg_value = float(np.median(region_values))
            elif aggregation_func == "min":
                agg_value = float(region_values.min())
            elif aggregation_func == "max":
                agg_value = float(region_values.max())
            else:
                raise ValueError(
                    f"Unknown aggregation function: {aggregation_func}. "
                    "Supported: mean, sum, std, median, min, max"
                )

        column_name = (
            aggregation_func
            if column_suffix is None
            else f"{aggregation_func}_{column_suffix}"
        )

        # Get region info
        try:
            region_info = self.get_region_info(int(label))
            region_name = region_info[self.name_column]
        except ValueError:
            region_name = f"Unknown_Region_{label}"

        # Calculate volume
        volume = self.get_region_volume(int(label))

        results.append(
            {
                self.label_column: int(label),
                self.name_column: region_name,
                column_name: agg_value,
                "volume_mm3": volume,
            }
        )

    return pd.DataFrame(results)

zarrnii.ZarrNiiAtlas.create_feature_map(feature_data, feature_column, label_column='index')

Create feature map by assigning values to atlas regions.

Parameters:

• feature_data (DataFrame) –

  DataFrame with region labels and feature values

• feature_column (str) –

  Column name containing feature values to map

• label_column (str, default: 'index' ) –

  Column name containing region labels

Returns:

• ZarrNii –

  ZarrNii instance with feature values mapped to regions

Raises:

• ValueError –

  If required columns are missing

Source code in zarrnii/atlas.py

def create_feature_map(
    self,
    feature_data: pd.DataFrame,
    feature_column: str,
    label_column: str = "index",
) -> ZarrNii:
    """Create feature map by assigning values to atlas regions.

    Args:
        feature_data: DataFrame with region labels and feature values
        feature_column: Column name containing feature values to map
        label_column: Column name containing region labels

    Returns:
        ZarrNii instance with feature values mapped to regions

    Raises:
        ValueError: If required columns are missing
    """
    # Validate input
    required_cols = [label_column, feature_column]
    missing_cols = [col for col in required_cols if col not in feature_data.columns]
    if missing_cols:
        raise ValueError(f"Missing columns in feature_data: {missing_cols}")

    dseg_data = self.dseg.data.astype("int")  # dask array of labels

    # make a dense lookup array
    max_label = int(feature_data[label_column].max())
    lut = np.zeros(max_label + 1, dtype=np.float32)
    lut[feature_data[label_column].to_numpy(dtype=int)] = feature_data[
        feature_column
    ].to_numpy(dtype=float)

    # broadcast the mapping in one go
    feature_map = dseg_data.map_blocks(lambda block: lut[block], dtype=np.float32)

    return ZarrNii.from_darr(
        feature_map, affine=self.dseg.affine, axes_order=self.dseg.axes_order
    )

zarrnii.ZarrNiiAtlas.get_region_bounding_box(region_ids=None, regex=None)

Get bounding box in physical coordinates for selected regions.

This method computes the spatial extents (bounding box) of one or more atlas regions in physical/world coordinates. The returned bounding box can be used directly with the crop method to extract a subvolume containing the selected regions.

Parameters:

• region_ids (Union[int, str, List[Union[int, str]]], default: None ) –

  Region identifier(s) to include in bounding box. Can be: - Single int: label index - Single str: region name or abbreviation - List[int/str]: multiple regions by index, name, or abbreviation - None: use regex parameter instead

• regex (Optional[str], default: None ) –

  Regular expression to match region names. If provided, region_ids must be None. Case-insensitive matching.

Returns:

• Tuple[float, float, float] –

  Tuple of (bbox_min, bbox_max) where each is a tuple of (x, y, z)

• Tuple[float, float, float] –

  coordinates in physical/world space (mm). These can be passed

• Tuple[Tuple[float, float, float], Tuple[float, float, float]] –

  directly to ZarrNii.crop() method with physical_coords=True.

Raises:

• ValueError –

  If no regions match the selection criteria, or if both region_ids and regex are provided/omitted

• TypeError –

  If region_ids contains invalid types

Examples:

>>> # Get bounding box for single region
>>> bbox_min, bbox_max = atlas.get_region_bounding_box("Hippocampus")
>>> cropped = image.crop(bbox_min, bbox_max, physical_coords=True)
>>>
>>> # Get bounding box for multiple regions
>>> bbox_min, bbox_max = atlas.get_region_bounding_box(["Hippocampus", "Amygdala"])
>>>
>>> # Use regex to select regions
>>> bbox_min, bbox_max = atlas.get_region_bounding_box(regex="Hip.*")
>>>
>>> # Crop atlas itself to region
>>> cropped_atlas = atlas.crop(bbox_min, bbox_max, physical_coords=True)

Notes

• Bounding box is in physical coordinates (mm), not voxel indices
• Axes ordering is relative to self.axes_order (e.g. ZYX for ome zarr)
• The bounding box is the union of all selected regions
• Use the returned values with crop(physical_coords=True)

Source code in zarrnii/atlas.py

def get_region_bounding_box(
    self,
    region_ids: Union[int, str, List[Union[int, str]]] = None,
    regex: Optional[str] = None,
) -> Tuple[Tuple[float, float, float], Tuple[float, float, float]]:
    """Get bounding box in physical coordinates for selected regions.

    This method computes the spatial extents (bounding box) of one or more
    atlas regions in physical/world coordinates. The returned bounding box
    can be used directly with the crop method to extract a subvolume
    containing the selected regions.

    Args:
        region_ids: Region identifier(s) to include in bounding box. Can be:
            - Single int: label index
            - Single str: region name or abbreviation
            - List[int/str]: multiple regions by index, name, or abbreviation
            - None: use regex parameter instead
        regex: Regular expression to match region names. If provided,
            region_ids must be None. Case-insensitive matching.

    Returns:
        Tuple of (bbox_min, bbox_max) where each is a tuple of (x, y, z)
        coordinates in physical/world space (mm). These can be passed
        directly to ZarrNii.crop() method with physical_coords=True.

    Raises:
        ValueError: If no regions match the selection criteria, or if both
            region_ids and regex are provided/omitted
        TypeError: If region_ids contains invalid types

    Examples:
        >>> # Get bounding box for single region
        >>> bbox_min, bbox_max = atlas.get_region_bounding_box("Hippocampus")
        >>> cropped = image.crop(bbox_min, bbox_max, physical_coords=True)
        >>>
        >>> # Get bounding box for multiple regions
        >>> bbox_min, bbox_max = atlas.get_region_bounding_box(["Hippocampus", "Amygdala"])
        >>>
        >>> # Use regex to select regions
        >>> bbox_min, bbox_max = atlas.get_region_bounding_box(regex="Hip.*")
        >>>
        >>> # Crop atlas itself to region
        >>> cropped_atlas = atlas.crop(bbox_min, bbox_max, physical_coords=True)

    Notes:
        - Bounding box is in physical coordinates (mm), not voxel indices
        - Axes ordering is relative to self.axes_order (e.g. ZYX for ome zarr)
        - The bounding box is the union of all selected regions
        - Use the returned values with crop(physical_coords=True)
    """
    import re

    import dask.array as da

    # Validate input parameters
    if region_ids is None and regex is None:
        raise ValueError("Must provide either region_ids or regex parameter")
    if region_ids is not None and regex is not None:
        raise ValueError("Cannot provide both region_ids and regex parameters")

    # Determine which labels to include
    selected_labels = []

    if regex is not None:
        # Match regions using regex
        pattern = re.compile(regex, re.IGNORECASE)
        for _, row in self.labels_df.iterrows():
            region_name = str(row[self.name_column])
            if pattern.search(region_name):
                selected_labels.append(int(row[self.label_column]))

        if not selected_labels:
            raise ValueError(f"No regions matched regex pattern: {regex}")
    else:
        # Convert region_ids to list if single value
        if not isinstance(region_ids, list):
            region_ids = [region_ids]

        # Resolve each region identifier to label
        for region_id in region_ids:
            label = self._resolve_region_identifier(region_id)
            selected_labels.append(label)

    # Create union mask of all selected regions
    dseg_data = self.dseg.data
    mask = None
    for label in selected_labels:
        region_mask = dseg_data == label
        if mask is None:
            mask = region_mask
        else:
            mask = mask | region_mask

    # Find voxel coordinates where mask is True
    # da.where returns tuple of arrays (one per dimension in data array)
    indices = da.where(mask)

    # Compute the indices to get actual coordinates
    indices_computed = [idx.compute() for idx in indices]

    # Check if any voxels were found
    if any(idx.size == 0 for idx in indices_computed):
        raise ValueError(f"No voxels found for selected regions: {selected_labels}")

    # Get the spatial dimensions from dims (skip non-spatial like 'c', 't')
    spatial_dims_lower = [d.lower() for d in ["x", "y", "z"]]
    spatial_indices = []
    for i, dim in enumerate(self.dseg.dims):
        if dim.lower() in spatial_dims_lower:
            spatial_indices.append(i)

    # Extract spatial coordinates from indices
    # indices_computed has one array per dimension in data
    voxel_coords = []
    for spatial_idx in spatial_indices:
        voxel_coords.append(indices_computed[spatial_idx])

    # Get min and max for each spatial dimension
    voxel_mins = [int(coords.min()) for coords in voxel_coords]
    voxel_maxs = [
        int(coords.max()) + 1 for coords in voxel_coords
    ]  # +1 for inclusive max

    # Now we have voxel coordinates in the order they appear in dims
    # We need to convert to (x, y, z) order for physical coordinates
    dim_to_voxel_range = {}
    for i, spatial_idx in enumerate(spatial_indices):
        dim_name = self.dseg.dims[spatial_idx].lower()
        dim_to_voxel_range[dim_name] = (voxel_mins[i], voxel_maxs[i])

    # Build voxel coordinates in (x, y, z) order
    voxel_min_xyz = np.array(
        [
            dim_to_voxel_range["x"][0],
            dim_to_voxel_range["y"][0],
            dim_to_voxel_range["z"][0],
            1.0,
        ]
    )
    voxel_max_xyz = np.array(
        [
            dim_to_voxel_range["x"][1],
            dim_to_voxel_range["y"][1],
            dim_to_voxel_range["z"][1],
            1.0,
        ]
    )

    # Transform to physical coordinates using affine
    affine_matrix = self.dseg.affine.matrix
    physical_min = affine_matrix @ voxel_min_xyz
    physical_max = affine_matrix @ voxel_max_xyz

    # Return as tuples of (x, y, z) in physical space
    bbox_min = tuple(physical_min[:3].tolist())
    bbox_max = tuple(physical_max[:3].tolist())

    return bbox_min, bbox_max

zarrnii.ZarrNiiAtlas.sample_region_patches(n_patches, region_ids=None, regex=None, seed=None)

Sample random coordinates (centers) within atlas regions.

This method generates a list of center coordinates by randomly sampling voxels within the selected atlas regions. The returned coordinates are in physical/world space (mm) and can be used with crop_centered() to extract fixed-size patches for machine learning training or other workflows.

Parameters:

• n_patches (int) –

  Number of patch centers to sample

• region_ids (Union[int, str, List[Union[int, str]]], default: None ) –

  Region identifier(s) to sample from. Can be: - Single int: label index - Single str: region name or abbreviation - List[int/str]: multiple regions by index, name, or abbreviation - None: use regex parameter instead

• regex (Optional[str], default: None ) –

  Regular expression to match region names. If provided, region_ids must be None. Case-insensitive matching.

• seed (Optional[int], default: None ) –

  Random seed for reproducibility. If None, patches are sampled randomly each time.

Returns:

• List[Tuple[float, float, float]] –

  List of (x, y, z) coordinates in physical/world space (mm).

• List[Tuple[float, float, float]] –

  Each coordinate represents the center of a potential patch and

• List[Tuple[float, float, float]] –

  can be used with crop_centered() to extract fixed-size regions.

Raises:

• ValueError –

  If no regions match the selection criteria, if both region_ids and regex are provided/omitted, or if n_patches is less than 1

• TypeError –

  If region_ids contains invalid types

Examples:

>>> # Sample 10 patch centers from hippocampus
>>> centers = atlas.sample_region_patches(
...     n_patches=10,
...     region_ids="Hippocampus",
...     seed=42
... )
>>> # Extract 256x256x256 voxel patches at each center
>>> patches = image.crop_centered(centers, patch_size=(256, 256, 256))
>>>
>>> # Sample from multiple regions using list
>>> centers = atlas.sample_region_patches(
...     n_patches=20,
...     region_ids=[1, 2, 3],
...     seed=42
... )
>>>
>>> # Sample using regex pattern
>>> centers = atlas.sample_region_patches(
...     n_patches=5,
...     regex=".*cortex.*",
... )

Notes

• Coordinates are in physical space (mm), not voxel indices
• Centers are sampled uniformly from voxels within selected regions
• Use crop_centered() to extract fixed-size patches around these centers
• For ML training with fixed patch sizes (e.g., 256x256x256 voxels), use a lower-resolution atlas to define masks, then crop at higher resolution using physical coordinates

Source code in zarrnii/atlas.py

def sample_region_patches(
    self,
    n_patches: int,
    region_ids: Union[int, str, List[Union[int, str]]] = None,
    regex: Optional[str] = None,
    seed: Optional[int] = None,
) -> List[Tuple[float, float, float]]:
    """Sample random coordinates (centers) within atlas regions.

    This method generates a list of center coordinates by randomly sampling
    voxels within the selected atlas regions. The returned coordinates are
    in physical/world space (mm) and can be used with crop_centered() to
    extract fixed-size patches for machine learning training or other workflows.

    Args:
        n_patches: Number of patch centers to sample
        region_ids: Region identifier(s) to sample from. Can be:
            - Single int: label index
            - Single str: region name or abbreviation
            - List[int/str]: multiple regions by index, name, or abbreviation
            - None: use regex parameter instead
        regex: Regular expression to match region names. If provided,
            region_ids must be None. Case-insensitive matching.
        seed: Random seed for reproducibility. If None, patches are sampled
            randomly each time.

    Returns:
        List of (x, y, z) coordinates in physical/world space (mm).
        Each coordinate represents the center of a potential patch and
        can be used with crop_centered() to extract fixed-size regions.

    Raises:
        ValueError: If no regions match the selection criteria, if both
            region_ids and regex are provided/omitted, or if n_patches is
            less than 1
        TypeError: If region_ids contains invalid types

    Examples:
        >>> # Sample 10 patch centers from hippocampus
        >>> centers = atlas.sample_region_patches(
        ...     n_patches=10,
        ...     region_ids="Hippocampus",
        ...     seed=42
        ... )
        >>> # Extract 256x256x256 voxel patches at each center
        >>> patches = image.crop_centered(centers, patch_size=(256, 256, 256))
        >>>
        >>> # Sample from multiple regions using list
        >>> centers = atlas.sample_region_patches(
        ...     n_patches=20,
        ...     region_ids=[1, 2, 3],
        ...     seed=42
        ... )
        >>>
        >>> # Sample using regex pattern
        >>> centers = atlas.sample_region_patches(
        ...     n_patches=5,
        ...     regex=".*cortex.*",
        ... )

    Notes:
        - Coordinates are in physical space (mm), not voxel indices
        - Centers are sampled uniformly from voxels within selected regions
        - Use crop_centered() to extract fixed-size patches around these centers
        - For ML training with fixed patch sizes (e.g., 256x256x256 voxels),
          use a lower-resolution atlas to define masks, then crop at higher
          resolution using physical coordinates
    """
    import re

    import dask.array as da

    # Validate input
    if n_patches < 1:
        raise ValueError(f"n_patches must be at least 1, got {n_patches}")

    if region_ids is None and regex is None:
        raise ValueError("Must provide either region_ids or regex parameter")
    if region_ids is not None and regex is not None:
        raise ValueError("Cannot provide both region_ids and regex parameters")

    # Set random seed if provided
    if seed is not None:
        np.random.seed(seed)

    # Determine which labels to include (reuse logic from get_region_bounding_box)
    selected_labels = []

    if regex is not None:
        # Match regions using regex
        pattern = re.compile(regex, re.IGNORECASE)
        for _, row in self.labels_df.iterrows():
            region_name = str(row[self.name_column])
            if pattern.search(region_name):
                selected_labels.append(int(row[self.label_column]))

        if not selected_labels:
            raise ValueError(f"No regions matched regex pattern: {regex}")
    else:
        # Convert region_ids to list if single value
        if not isinstance(region_ids, list):
            region_ids = [region_ids]

        # Resolve each region identifier to label
        for region_id in region_ids:
            label = self._resolve_region_identifier(region_id)
            selected_labels.append(label)

    # Create union mask of all selected regions
    dseg_data = self.dseg.data
    mask = None
    for label in selected_labels:
        region_mask = dseg_data == label
        if mask is None:
            mask = region_mask
        else:
            mask = mask | region_mask

    # Find voxel coordinates where mask is True
    indices = da.where(mask)

    # Compute the indices to get actual coordinates
    indices_computed = [idx.compute() for idx in indices]

    # Check if any voxels were found
    if any(idx.size == 0 for idx in indices_computed):
        raise ValueError(f"No voxels found for selected regions: {selected_labels}")

    # Get number of valid voxels
    n_voxels = indices_computed[0].size

    # Sample random voxels
    # If n_patches > n_voxels, sample with replacement
    replace = n_patches > n_voxels
    sampled_indices = np.random.choice(n_voxels, size=n_patches, replace=replace)

    # Get spatial dimensions (skip non-spatial like 'c', 't')
    spatial_dims_lower = [d.lower() for d in ["x", "y", "z"]]
    spatial_indices = []
    for i, dim in enumerate(self.dseg.dims):
        if dim.lower() in spatial_dims_lower:
            spatial_indices.append(i)

    # Build voxel coordinates for sampled centers
    sampled_coords = []
    for spatial_idx in spatial_indices:
        sampled_coords.append(indices_computed[spatial_idx][sampled_indices])

    # Map to x, y, z order
    dim_to_coords = {}
    for i, spatial_idx in enumerate(spatial_indices):
        dim_name = self.dseg.dims[spatial_idx].lower()
        dim_to_coords[dim_name] = sampled_coords[i]

    # Get affine matrix
    affine_matrix = self.dseg.affine.matrix

    # Generate center coordinates in physical space
    centers = []
    for i in range(n_patches):
        # Get center voxel coordinates in (x, y, z) order
        center_voxel_xyz = np.array(
            [
                dim_to_coords["x"][i],
                dim_to_coords["y"][i],
                dim_to_coords["z"][i],
                1.0,
            ]
        )

        # Transform to physical coordinates
        center_physical = affine_matrix @ center_voxel_xyz
        center_xyz = center_physical[:3]

        # Convert to tuple
        center = tuple(center_xyz.tolist())
        centers.append(center)

    return centers

Bases: Transform

Affine transformation for spatial coordinate mapping.

Represents a 4x4 affine transformation matrix that can be used to transform 3D coordinates between different coordinate systems. Supports various operations including matrix multiplication, inversion, and point transformation.

Attributes:

• matrix (ndarray) –

  4x4 affine transformation matrix

Functions

zarrnii.transform.AffineTransform.from_txt(path, invert=False) classmethod

Create AffineTransform from text file containing matrix.

Parameters:

• path (Union[str, bytes]) –

  Path to text file containing 4x4 affine matrix

• invert (bool, default: False ) –

  Whether to invert the matrix after loading

Returns:

• 'AffineTransform' –

  AffineTransform instance with loaded matrix

Raises:

• OSError –

  If file cannot be read

• ValueError –

  If file does not contain valid 4x4 matrix

Source code in zarrnii/transform.py

@classmethod
def from_txt(
    cls, path: Union[str, bytes], invert: bool = False
) -> "AffineTransform":
    """Create AffineTransform from text file containing matrix.

    Args:
        path: Path to text file containing 4x4 affine matrix
        invert: Whether to invert the matrix after loading

    Returns:
        AffineTransform instance with loaded matrix

    Raises:
        OSError: If file cannot be read
        ValueError: If file does not contain valid 4x4 matrix
    """
    matrix = np.loadtxt(path)
    if invert:
        matrix = np.linalg.inv(matrix)
    return cls(matrix=matrix)

zarrnii.transform.AffineTransform.from_array(matrix, invert=False) classmethod

Create AffineTransform from numpy array.

Parameters:

• matrix (ndarray) –

  4x4 numpy array representing affine transformation

• invert (bool, default: False ) –

  Whether to invert the matrix

Returns:

• 'AffineTransform' –

  AffineTransform instance with the matrix

Raises:

• ValueError –

  If matrix is not 4x4

Source code in zarrnii/transform.py

@classmethod
def from_array(cls, matrix: np.ndarray, invert: bool = False) -> "AffineTransform":
    """Create AffineTransform from numpy array.

    Args:
        matrix: 4x4 numpy array representing affine transformation
        invert: Whether to invert the matrix

    Returns:
        AffineTransform instance with the matrix

    Raises:
        ValueError: If matrix is not 4x4
    """
    if matrix.shape != (4, 4):
        raise ValueError(f"Matrix must be 4x4, got shape {matrix.shape}")

    if invert:
        matrix = np.linalg.inv(matrix)
    return cls(matrix=matrix)

zarrnii.transform.AffineTransform.identity() classmethod

Create identity transformation.

Returns:

• 'AffineTransform' –

  AffineTransform representing identity transformation (no change)

Source code in zarrnii/transform.py

@classmethod
def identity(cls) -> "AffineTransform":
    """Create identity transformation.

    Returns:
        AffineTransform representing identity transformation (no change)
    """
    return cls(matrix=np.eye(4, 4))

zarrnii.transform.AffineTransform.apply_transform(vecs)

Apply transformation to coordinate vectors.

Parameters:

• vecs (ndarray) –

  Input coordinates to transform

Returns:

• ndarray –

  Transformed coordinates

Source code in zarrnii/transform.py

def apply_transform(self, vecs: np.ndarray) -> np.ndarray:
    """Apply transformation to coordinate vectors.

    Args:
        vecs: Input coordinates to transform

    Returns:
        Transformed coordinates
    """
    return self @ vecs

zarrnii.transform.AffineTransform.invert()

Return the inverse of the matrix transformation.

Returns:

• 'AffineTransform' –

  New AffineTransform with inverted matrix

Raises:

• LinAlgError –

  If matrix is singular and cannot be inverted

Source code in zarrnii/transform.py

def invert(self) -> "AffineTransform":
    """Return the inverse of the matrix transformation.

    Returns:
        New AffineTransform with inverted matrix

    Raises:
        np.linalg.LinAlgError: If matrix is singular and cannot be inverted
    """
    return AffineTransform.from_array(np.linalg.inv(self.matrix))

zarrnii.transform.AffineTransform.update_for_orientation(input_orientation, output_orientation)

Update the matrix to map from input orientation to output orientation.

Parameters:

• input_orientation (str) –

  Current anatomical orientation (e.g., 'RPI')

• output_orientation (str) –

  Target anatomical orientation (e.g., 'RAS')

Returns:

• 'AffineTransform' –

  New AffineTransform updated for orientation mapping

Raises:

• ValueError –

  If orientations are invalid or cannot be matched

Source code in zarrnii/transform.py

def update_for_orientation(
    self, input_orientation: str, output_orientation: str
) -> "AffineTransform":
    """Update the matrix to map from input orientation to output orientation.

    Args:
        input_orientation: Current anatomical orientation (e.g., 'RPI')
        output_orientation: Target anatomical orientation (e.g., 'RAS')

    Returns:
        New AffineTransform updated for orientation mapping

    Raises:
        ValueError: If orientations are invalid or cannot be matched
    """

    # Define a mapping of anatomical directions to axis indices and flips
    axis_map = {
        "R": (0, 1),
        "L": (0, -1),
        "A": (1, 1),
        "P": (1, -1),
        "S": (2, 1),
        "I": (2, -1),
    }

    # Parse the input and output orientations
    input_axes = [axis_map[ax] for ax in input_orientation]
    output_axes = [axis_map[ax] for ax in output_orientation]

    # Create a mapping from input to output
    reorder_indices = [None] * 3
    flip_signs = [1] * 3

    for out_idx, (out_axis, out_sign) in enumerate(output_axes):
        for in_idx, (in_axis, in_sign) in enumerate(input_axes):
            if out_axis == in_axis:  # Match axis
                reorder_indices[out_idx] = in_idx
                flip_signs[out_idx] = out_sign * in_sign
                break

    # Reorder and flip the affine matrix
    reordered_matrix = np.zeros_like(self.matrix)
    for i, (reorder_idx, flip_sign) in enumerate(zip(reorder_indices, flip_signs)):
        if reorder_idx is None:
            raise ValueError(
                f"Cannot match all axes from {input_orientation} to {output_orientation}."
            )
        reordered_matrix[i, :3] = flip_sign * self.matrix[reorder_idx, :3]
        reordered_matrix[i, 3] = flip_sign * self.matrix[reorder_idx, 3]
    reordered_matrix[3, :] = self.matrix[3, :]  # Preserve the homogeneous row

    return AffineTransform.from_array(reordered_matrix)

Bases: Transform

Non-linear displacement field transformation.

Represents a displacement field transformation where each point in space has an associated displacement vector. Uses interpolation to compute displacements for arbitrary coordinates.

Attributes:

• disp_xyz (ndarray) –

  Displacement vectors at grid points (4D array: x, y, z, vector_component)

• disp_grid (Tuple[ndarray, ...]) –

  Grid coordinates for displacement field

• disp_affine (AffineTransform) –

  Affine transformation from world to displacement field coordinates

Functions

zarrnii.transform.DisplacementTransform.from_nifti(path) classmethod

Create DisplacementTransform from NIfTI file.

Parameters:

• path (Union[str, bytes]) –

  Path to NIfTI displacement field file

Returns:

• 'DisplacementTransform' –

  DisplacementTransform instance loaded from file

Raises:

• OSError –

  If file cannot be read

• ValueError –

  If file format is invalid

Source code in zarrnii/transform.py

@classmethod
def from_nifti(cls, path: Union[str, bytes]) -> "DisplacementTransform":
    """Create DisplacementTransform from NIfTI file.

    Args:
        path: Path to NIfTI displacement field file

    Returns:
        DisplacementTransform instance loaded from file

    Raises:
        OSError: If file cannot be read
        ValueError: If file format is invalid
    """
    disp_nib = nib.load(path)
    disp_xyz = disp_nib.get_fdata().squeeze()
    disp_affine = AffineTransform.from_array(disp_nib.affine)

    # Convert from ITK transform convention
    # ITK uses opposite sign convention for x and y displacements
    disp_xyz[:, :, :, 0] = -disp_xyz[:, :, :, 0]
    disp_xyz[:, :, :, 1] = -disp_xyz[:, :, :, 1]

    # Create grid coordinates
    disp_grid = (
        np.arange(disp_xyz.shape[0]),
        np.arange(disp_xyz.shape[1]),
        np.arange(disp_xyz.shape[2]),
    )

    return cls(
        disp_xyz=disp_xyz,
        disp_grid=disp_grid,
        disp_affine=disp_affine,
    )

zarrnii.transform.DisplacementTransform.apply_transform(vecs)

Apply displacement transformation to coordinate vectors.

Transforms input coordinates by interpolating displacement vectors from the displacement field and adding them to the input coordinates.

Parameters:

• vecs (ndarray) –

  Input coordinates as numpy array. Shape should be (3, N) for N points or (3,) for single point

Returns:

• ndarray –

  Transformed coordinates with same shape as input

Notes

Points outside the displacement field domain are filled with zero displacement (no transformation).

Source code in zarrnii/transform.py

def apply_transform(self, vecs: np.ndarray) -> np.ndarray:
    """Apply displacement transformation to coordinate vectors.

    Transforms input coordinates by interpolating displacement vectors
    from the displacement field and adding them to the input coordinates.

    Args:
        vecs: Input coordinates as numpy array. Shape should be (3, N) for
            N points or (3,) for single point

    Returns:
        Transformed coordinates with same shape as input

    Notes:
        Points outside the displacement field domain are filled with
        zero displacement (no transformation).
    """
    # Transform points to voxel space of the displacement field
    vox_vecs = self.disp_affine.invert() @ vecs

    # Initialize displacement vectors
    disp_vecs = np.zeros(vox_vecs.shape)

    # Interpolate displacement for each spatial dimension (x, y, z)
    for ax in range(3):
        disp_vecs[ax, :] = interpn(
            self.disp_grid,
            self.disp_xyz[:, :, :, ax].squeeze(),
            vox_vecs[:3, :].T,
            method="linear",
            bounds_error=False,
            fill_value=0,
        )

    # Add displacement to original coordinates
    return vecs + disp_vecs