Working with Images

In the previous tutorial, you already saw that images are loaded as imfusion.SharedImageSet. As the name suggests, an imfusion.SharedImageSet is basically a list of images. You might think of it as a movie clip where each image in the list is a different frame.

Each image in an imfusion.SharedImageSet is an imfusion.SharedImage. The imfusion.SharedImage stores the image’s pixel data in both the CPU and the GPU memory, enabling shared access of the same image across these locations. The framework takes care that they always stay synchronized. In addition, an imfusion.SharedImage stores the position of the image in world coordinates. For more details, see Coordinate Systems.

An imfusion.SharedImageSet behaves similarly to a Python list: you can iterate over it and index individual elements:

>>> imageset, *_ = imfusion.load(us_sweep_dcm)
>>> imageset
imfusion.SharedImageSet(size: 20, [
    imfusion.SharedImage(UBYTE width: 164 height: 552 spacing: 0.228659x0.0724638x1 mm),
    imfusion.SharedImage(UBYTE width: 164 height: 552 spacing: 0.228659x0.0724638x1 mm),
    ...
    imfusion.SharedImage(UBYTE width: 164 height: 552 spacing: 0.228659x0.0724638x1 mm)
])
>>> image = imageset[0]
>>> image
imfusion.SharedImage(UBYTE width: 164 height: 552 spacing: 0.228659x0.0724638x1 mm)

Occasionally, you might come across an imfusion.MemImage, which is a low-level container for the pixel data stored in CPU memory. Therefore, each imfusion.SharedImage contains an imfusion.MemImage.

When working with algorithms and other framework components, you will typically deal with an imfusion.SharedImageSet, as it is the most generic image container. It can hold various types of medical image data, ranging from single 3D CT volumes to lists of 2D ultrasound images.

Accessing image data with numpy

Note

Since numpy arrays can lead to direct image memory access, some of the conversions provided below require knowledge about how we internally store pixel values. In particular, shift and scale are efficient mechanisms that allow modification of the range of pixel values without actually updating all the pixels that are contained in the images. Please refer to Storage and Original values for further details about shift and scale.

The best way to work with the pixel data of an image is by converting it to a numpy.array. imfusion.MemImage, imfusion.SharedImage, and imfusion.SharedImageSet offer three slightly different ways to convert them into a numpy array:

  1. image_or_imageset.numpy(): Creates a new copy of the data, whereby shift and scale are automatically applied to the returned numpy array. The numpy array’s dtype is set such that all image values can be represented without information loss.

  2. numpy.array(image_or_imageset): Creates a new copy of the data, whereby shift and scale are not automatically applied to the returned numpy array. Thus, if needed, shift and scale must be applied manually. The numpy array’s dtype can be specified by passing it as an additional argument. Otherwise, it will be assigned by numpy’s internal rules.

  3. numpy.array(image_or_imageset, copy=False): Does not create a new copy of the data, and shift and scale are not applied to the returned numpy array. This means that any changes to the numpy array will be reflected in the memory pointed to by image_or_imageset. If needed, shift and scale must be applied manually. The numpy array’s dtype can be specified by passing it as an additional argument. Otherwise, it will be assigned by numpy’s internal rules.

As a first example, let’s take a closer look at an imfusion.SharedImage:

>>> import numpy as np
>>> (imageset,) = imfusion.load(ct_image_png)
>>> imageset
imfusion.SharedImageSet(size: 1, [imfusion.SharedImage(USHORT width: 512 height: 512 spacing: 0.661813x0.661813x1 mm)])
>>> image = imageset[0]
>>> arr = np.array(image)
>>> arr.shape
(512, 512, 1)

The number of dimensions can vary between 2 and 4. The shape should be read from right to left: the last dimension is always the number of channels, limited to values between 1 and 4. This is preceded by the width, height, and the number of slices. Dimensions of 1 on the left side are omitted. In this example, the number of slices is omitted as the given 2D image only has one slice.

Unlike the Python Imaging Library (PIL), the color dimension is also present on grayscale images (like in the example above). This way you can always differentiate between a 2D color image (ndim of 3) and a grayscale 3D volume (ndim of 4). To interact with PIL or matplotlib, you might have to remove the color dimension if it is 1:

>>> import matplotlib.pyplot as plt
>>> plt.imshow(arr[:,:,0], cmap='gray')
<matplotlib.image.AxesImage object at 0x...>
_images/ct_image_matplotlib.png

Accordingly, you have to add the color dimension if you want to load images with PIL into an imfusion.SharedImage:

>>> from PIL import Image
>>> arr = np.array(Image.open(ct_image_png))
>>> arr.shape
(512, 512)
>>> imfusion.SharedImage(arr)
imfusion.SharedImage(USHORT width: 512 channels: 512)
>>> arr = np.expand_dims(arr, axis=-1)  # append color dimension
>>> arr.shape
(512, 512, 1)
>>> imfusion.SharedImage(arr)
imfusion.SharedImage(USHORT width: 512 height: 512)

Alternatively, you can also use the optional greyscale argument to tell the imfusion.SharedImage to use an implicit color dimension of 1:

>>> imfusion.SharedImage(arr, greyscale=True)
imfusion.SharedImage(USHORT width: 512 height: 512)

Modifying image data

Using numpy, we can now access the pixel data, but how do we modify it? You will notice that changes made to the numpy array are not reflected in the corresponding imfusion.SharedImage. This is because numpy creates a copy by default. Therefore, we have to copy the altered data back into the imfusion.SharedImage using imfusion.SharedImage.assign_array() (see Direct memory access for how to avoid copies in the first place). Let’s assume we want to apply a simple threshold:

>>> arr = np.array(image)
>>> threshold = 2000
>>> arr[arr < threshold] = 0
>>> arr[arr >= threshold] = 1
>>> image.assign_array(arr)

Note that this way you cannot resize or reshape the image. Neither can you change the type of the image.

SharedImageSet

Since an imfusion.SharedImageSet is essentially a list of imfusion.SharedImage, it has an additional dimension. Let’s load an ultrasound sweep, i.e., a time series of 2D ultrasound frames:

>>> imageset, *_ = imfusion.load(us_sweep_dcm)
>>> imageset
imfusion.SharedImageSet(size: 20, [
    imfusion.SharedImage(UBYTE width: 164 height: 552 spacing: 0.228659x0.0724638x1 mm),
    imfusion.SharedImage(UBYTE width: 164 height: 552 spacing: 0.228659x0.0724638x1 mm),
    ...
    imfusion.SharedImage(UBYTE width: 164 height: 552 spacing: 0.228659x0.0724638x1 mm)
])
>>> arr = np.array(imageset)
>>> arr.shape
(20, 552, 164, 1)

The first dimension of an array, which was created from an imfusion.SharedImageSet, will always be the number of frames. It is important to remember that the number of slices and number of frames will be omitted if they are 1. Consequently, a set of 3D volumes will have 5 dimensions, a set of 2D images 4, and a set of 1D images 3. The following image shows how each dimension in a numpy.shape relates to the corresponding dimension in an imfusion.SharedImageSet:

_images/dimensions.png

Aside from the additional dimension, an imfusion.SharedImageSet can be used similarly to an imfusion.SharedImage:

>>> imageset = imfusion.SharedImageSet(np.zeros([10, 600, 800, 1], dtype='uint8'))
>>> imageset
imfusion.SharedImageSet(size: 10, [
    imfusion.SharedImage(UBYTE width: 800 height: 600),
    imfusion.SharedImage(UBYTE width: 800 height: 600),
    ...
    imfusion.SharedImage(UBYTE width: 800 height: 600)
])
>>> imageset.assign_array(np.ones([10, 600, 800, 1], dtype='uint8'))

Direct memory access

numpy.array() always creates a copy of the input array as long as the optional copy argument is not set (default: copy=True). While it is possible and generally supported to avoid the copy, you have to be very careful with variables referring to such arrays. Since numpy and the C++ backend are accessing the same memory block, your program will crash if the memory gets deleted while numpy is still accessing it. Direct memory access is only available on imfusion.SharedImage and not on imfusion.SharedImageSet.

Without the implicit copy, our thresholding example from before now looks like this:

>>> (imageset,) = imfusion.load('ct_image.png')
>>> image = imageset[0]
>>> arr = np.array(image, copy=False)
>>> threshold = 2000
>>> arr[arr < threshold] = 0
>>> arr[arr >= threshold] = 1
>>> image.set_dirty_mem()

Note the last line: this tells the imfusion.SharedImage that the data has changed and other shared memory locations (like the one on the GPU) must be updated. A missing call to imfusion.SharedImage.setDirtyMem() is a common source of bugs, e.g., GPU algorithms working on outdated data or overwriting changes.

A general recommendation is to write your algorithm by copying the data and only avoid the copy if you notice performance issues. This tends to avoid bugs and leads to less headaches when developing.

Storage and Original values

After converting an image to a numpy array, you might have noticed that the values are not in the range you expect them to be. Especially if you loaded some data with negative values, all the values in the array might be positive. The ImFusionLib distinguishes between Storage- and Original pixel values, whereby their relation is defined by the imfusion.SharedImage.shift and imfusion.SharedImage.scale properties.

  • Original: Pixel values are the same as in their original source (e.g., when loaded from a file). They are identical to the storage pixel values if the image’s scale is 1 and the shift is 0.

  • Storage: Pixel values as they are stored in an imfusion.MemImage. The user may decide to apply such a rescaling in order to better use the available limits of the underlying data type.

The following conversion rules apply:

\[ \begin{align}\begin{aligned}original = (storage / scale) - shift\\storage = (original + shift) * scale\end{aligned}\end{align} \]

Since this is a common operation, there are imfusion.SharedImage.apply_shift_and_scale() and imfusion.SharedImageSet.apply_shift_and_scale() that convert a numpy array from storage to original values:

>>> imageset, *_ = imfusion.load(ct_image_png)
>>> imageset[0].shift = 1024.0
>>> arr = np.array(imageset)
>>> arr.dtype
dtype('uint16')
>>> int(np.min(arr))
0
>>> arr = imageset.apply_shift_and_scale(arr)
>>> arr.dtype
dtype('float64')
>>> float(np.min(arr))
-1024.0

Type conversions

imfusion.MemImage, imfusion.SharedImage and imfusion.SharedImageSet allow to convert the image or images into a new format. It works similarly to conversions in numpy:

>>> from imfusion import PixelType
>>> imageset, *_ = imfusion.load(ct_image_png)
>>> image = imageset[0]
>>> image_as_float = image.astype(PixelType.FLOAT)
>>> image.descriptor.pixel_type, image_as_float.descriptor.pixel_type
(<PixelType.USHORT: 5123>, <PixelType.FLOAT: 5126>)

Besides the above example, the astype method supports as input:

  1. a PixelType (e.g., imfusion.PixelType.UINT);

  2. most of the numpy’s dtypes (e.g., np.uint);

  3. Python’s float or int types.

Selections

class imfusion.Selection(*args, **kwargs)

Utility class for describing a selection of elements out of a set. Conceptually, a Selection pairs a list of bools describing selected items with the index of a “focus” item and provides syntactic sugar on top. For instance, the set of selected items could define which ones to show in general while the focus item is additionally highlighted. The class is fully separate from the item set of which it describes the selection. This means for instance that it cannot know the actual number of items in the set and the user/parent class must manually make sure that they match. Also, a Selection only manages indices and offers no way of accessing the underlying elements. In order to iterate over all selected indices, you can do for instance the following:

>>> for index in range(selection.start, selection.stop):
...     if selection[index]:
...             pass

The same effect can also be achieved in a much more terse fashion:

>>> for selected_index in selection.selected_indices:
...     pass

For convenience, the selection can also be converted to a slice object (if the selection has a regular spacing, see below):

>>> selected_subset = container[selection.to_slice()]

Sometimes it can be more convenient to “thin out” a selection by only selecting every N-th element. To this end, the Selection constructor takes the arguments start, stop and step. Setting step to N will only select every N-th element, mimicking the signature of range, slice, etc.

Overloaded function.

  1. __init__(self: imfusion.Selection) -> None

  2. __init__(self: imfusion.Selection, stop: int) -> None

  3. __init__(self: imfusion.Selection, start: int, stop: int, step: int = 1) -> None

  4. __init__(self: imfusion.Selection, indices: list[int]) -> None