Skip to content

Image processing

Apply a Gaussian blur three ways.

Apply a Gaussian blur to an OME-Zarr image three ways with ngio: eagerly on a numpy array, lazily with dask, and through an ngio iterator. Along the way you derive a new container that keeps the metadata of the original and write the processed image into it.

Step 1: set up

Start with a function that applies a Gaussian blur to an image. It takes an image and a sigma value as input, and returns the blurred image.

import numpy as np
import skimage


def gaussian_blur(image: np.ndarray, sigma: float) -> np.ndarray:
    """Apply gaussian blur to an image."""
    original_type = image.dtype
    image = skimage.filters.gaussian(
        image, sigma=sigma, channel_axis=0, preserve_range=True
    )
    # Convert the image back to the original type
    image = image.astype(original_type)
    return image

Step 2: open the OME-Zarr container

from pathlib import Path

from ngio import open_ome_zarr_container
from ngio.utils import download_ome_zarr_dataset

# Download the dataset
download_dir = Path("./data").absolute()

hcs_path = download_ome_zarr_dataset("CardiomyocyteTiny", download_dir=download_dir)
image_path = hcs_path / "B" / "03" / "0"

# Open the OME-Zarr container
ome_zarr = open_ome_zarr_container(image_path)

Step 3: create a new empty OME-Zarr container

ngio can "derive" a new container from an existing one. Use this when you want to apply processing to an image and save the results in a new container that preserves the original metadata and dimensions (unless you change them explicitly when deriving).

# Take the image to read from
image = ome_zarr.get_image()

# Derive a new OME-Zarr container to store the processed image in

blurred_omezarr_path = image_path.parent / "0_blurred"
blurred_omezarr = ome_zarr.derive_image(
    store=blurred_omezarr_path, name="Blurred Image", overwrite=True
)
blurred_image = blurred_omezarr.get_image()

Step 4: apply the Gaussian blur and consolidate the processed image

# `axes_order` sets the order the array comes back in. Here it is
# ["c", "z", "y", "x"], to match the blur function, which expects the
# channel axis first.
image_data = image.get_as_numpy(axes_order=["c", "z", "y", "x"])
# Apply gaussian blur to the image
sigma = 5.0
blurred_image_data = gaussian_blur(image_data, sigma=sigma)

# Write the processed data back to the OME-Zarr image
blurred_image.set_array(patch=blurred_image_data, axes_order=["c", "z", "y", "x"])

# `set_array` wrote to one resolution level only, so the rest of the pyramid is
# still empty. `consolidate` rebuilds the other levels from it.
blurred_image.consolidate()

Plot the results

Finally, visualise the original and blurred images with matplotlib.

original = image.get_as_numpy(c=0, z=1, axes_order=["y", "x"])
blurred = blurred_image.get_as_numpy(c=0, z=1, axes_order=["y", "x"])

# The data does not fill its uint16 range, so window it on percentiles. One window for
# both panels: stretching them separately would misrepresent the difference.
vmin, vmax = np.percentile(original, (1, 99.8))

fig, axs = plt.subplots(2, 1, figsize=(8, 6))
axs[0].set_title("Original image")
axs[0].imshow(original, cmap="gray", vmin=vmin, vmax=vmax)
axs[1].set_title("Blurred image")
axs[1].imshow(blurred, cmap="gray", vmin=vmin, vmax=vmax)
for ax in axs:
    ax.axis("off")
fig.tight_layout()
print(figure_html(fig))
2026-08-11T12:15:59.878433 image/svg+xml Matplotlib v3.11.0, https://matplotlib.org/ Original image Blurred image

Step 5: out-of-memory processing

Some images are larger than memory. In that case, use the dask library to process the image in chunks: with ngio you query the data as a dask array and apply the processing function to it.

from dask import array as da


def dask_gaussian_blur(image: da.Array, sigma: float) -> da.Array:
    """Apply gaussian blur to a dask array."""
    # This introduces edge artefacts at chunk boundaries. In a real application,
    # use map_overlap with a depth chosen from sigma to avoid them.
    return da.map_blocks(gaussian_blur, image, dtype=image.dtype, sigma=sigma)


image_dask = image.get_as_dask(axes_order=["c", "z", "y", "x"])
blurred_image_dask = dask_gaussian_blur(image_dask, sigma=sigma)
print(blurred_image_dask)
dask.array

Step 6: image processing iterators

ngio also processes large images with iterators. This API is not meant to replace dask: it lets you iterate over arbitrary regions, and it supplies default broadcasting behaviour.

from ngio.iterators import ImageProcessingIterator

iterator = ImageProcessingIterator(
    input_image=image,
    output_image=blurred_image,
    axes_order=["c", "z", "y", "x"],
)

# A freshly built iterator covers the entire image in one region.
print(f"Iterator over the whole image: {iterator}")

# Narrow it to an arbitrary ROI table. `product` takes the cartesian product
# of the iterator's regions and the table's.
table = ome_zarr.get_roi_table("FOV_ROI_table")
iterator = iterator.product(table)
print(f"Iterator after product with table: {iterator}")

# Set the broadcasting explicitly. `by_zyx` splits the time axis, so each step
# yields one whole ZYX volume rather than the full time series at once.
iterator = iterator.by_zyx()

# Optionally assert the regions do not overlap each other...
iterator.require_no_regions_overlap()
# ...nor share chunks, which is what makes parallel writes safe.
iterator.require_no_chunks_overlap()

# Map the blur across every region
iterator.map_as_numpy(lambda x: gaussian_blur(x, sigma=sigma))

# No need to consolidate: the iterator does it once every region is processed
Iterator over the whole image: ImageProcessingIterator(regions=1) Iterator after product with table: ImageProcessingIterator(regions=2)

Next steps