Engee documentation
Notebook

Indexing and segmentation of images

This demonstration shows how to work with an indexed image using IndirectArrays, and also breaks down how to segment an image using ImageSegmentation.

Indexing

An indexed image consists of two parts: the indexes and the pixels themselves.

The following is a "direct" representation of the image as an array.

In [ ]: 
Pkg.add(["IndirectArrays", "TestImages"])
In [ ]: 
using Images
img = [
    RGB(0.0, 0.0, 0.0) RGB(1.0, 0.0, 0.0) RGB(0.0, 1.0, 0.0) RGB(0.0, 0.0, 1.0) RGB(1.0, 1.0, 1.0)
    RGB(1.0, 0.0, 0.0) RGB(0.0, 1.0, 0.0) RGB(0.0, 0.0, 1.0) RGB(1.0, 1.0, 1.0) RGB(0.0, 0.0, 0.0)
    RGB(0.0, 1.0, 0.0) RGB(0.0, 0.0, 1.0) RGB(1.0, 1.0, 1.0) RGB(0.0, 0.0, 0.0) RGB(1.0, 0.0, 0.0)
    RGB(0.0, 0.0, 1.0) RGB(1.0, 1.0, 1.0) RGB(0.0, 0.0, 0.0) RGB(1.0, 0.0, 0.0) RGB(0.0, 1.0, 0.0)
    RGB(1.0, 1.0, 1.0) RGB(0.0, 0.0, 0.0) RGB(1.0, 0.0, 0.0) RGB(0.0, 1.0, 0.0) RGB(0.0, 0.0, 1.0)]
Out[0]:
No description has been provided for this image

Alternatively, we could save it as an indexed image:

In [ ]: 
indices = [1 2 3 4 5; 2 3 4 5 1; 3 4 5 1 2; 4 5 1 2 3; 5 1 2 3 4] # `i` records the pixel value `palette[i]`
palette = [RGB(0.0, 0.0, 0.0), RGB(1.0, 0.0, 0.0), RGB(0.0, 1.0, 0.0), RGB(0.0, 0.0, 1.0), RGB(1.0, 1.0, 1.0)]

palette[indices] == img
Out[0]:
true

Let's now analyse the pros and cons of using IndirectArrays.

Normal indexed images may require more memory due to storing two separate fields. If we use unique(img), on the other hand, the size is close to the image's own size. Also, standard indexing requires two operations, so it may be slower than direct image representation (and not amenable to SIMD vectorisation).

Also, using the indexed image format with two separate arrays can be awkward, so IndirectArrays provides an array abstraction to combine the two parts:

In [ ]: 
using IndirectArrays

indexed_img = IndirectArray(indices, palette)
img == indexed_img
Out[0]:
true

Because an IndirectArray is simply an array, common image operations apply to this type, such as:

In [ ]: 
new_img = imresize(indexed_img; ratio=2)
Out[0]:
No description has been provided for this image
In [ ]: 
dump(indexed_img)

IndirectArray{RGB{Float64}, 2, Int64, Matrix{Int64}, Vector{RGB{Float64}}}
  index: Array{Int64}((5, 5)) [1 2 … 4 5; 2 3 … 5 1; … ; 4 5 … 2 3; 5 1 … 3 4]
  values: Array{RGB{Float64}}((5,))
    1: RGB{Float64}
      r: Float64 0.0
      g: Float64 0.0
      b: Float64 0.0
    2: RGB{Float64}
      r: Float64 1.0
      g: Float64 0.0
      b: Float64 0.0
    3: RGB{Float64}
      r: Float64 0.0
      g: Float64 1.0
      b: Float64 0.0
    4: RGB{Float64}
      r: Float64 0.0
      g: Float64 0.0
      b: Float64 1.0
    5: RGB{Float64}
      r: Float64 1.0
      g: Float64 1.0
      b: Float64 1.0
In [ ]: 
dump(new_img)

Array{RGB{Float64}}((10, 10))
  1: RGB{Float64}
    r: Float64 0.0
    g: Float64 0.0
    b: Float64 0.0
  2: RGB{Float64}
    r: Float64 0.25
    g: Float64 0.0
    b: Float64 0.0
  3: RGB{Float64}
    r: Float64 0.75
    g: Float64 0.0
    b: Float64 0.0
  4: RGB{Float64}
    r: Float64 0.75
    g: Float64 0.25
    b: Float64 0.0
  5: RGB{Float64}
    r: Float64 0.25
    g: Float64 0.75
    b: Float64 0.0
  ...
  96: RGB{Float64}
    r: Float64 0.75
    g: Float64 0.25
    b: Float64 0.0
  97: RGB{Float64}
    r: Float64 0.25
    g: Float64 0.75
    b: Float64 0.0
  98: RGB{Float64}
    r: Float64 0.0
    g: Float64 0.75
    b: Float64 0.25
  99: RGB{Float64}
    r: Float64 0.0
    g: Float64 0.25
    b: Float64 0.75
  100: RGB{Float64}
    r: Float64 0.0
    g: Float64 0.0
    b: Float64 1.0

As we can see, after the operation on the indexed array, it has been reduced to the form of a standard array of pixels.

Segmentation

Next, we segment the image using watershed algorithm. This is the basic tool of mathematical morphology for image segmentation.

In [ ]: 
using TestImages
using IndirectArrays
In [ ]: 
img = testimage("blobs")
Out[0]:
No description has been provided for this image
In [ ]: 
img_example = zeros(Gray, 5, 5)
img_example[2:4,2:4] .=  Gray(0.6)
bw = Gray.(img) .> 0.5
bw_example = img_example .> 0.5
Out[0]:
5×5 BitMatrix:
 0  0  0  0  0
 0  1  1  1  0
 0  1  1  1  0
 0  1  1  1  0
 0  0  0  0  0

The feature_transform function allows us to find the binary image's feature transform (bw). It finds the nearest "feature" (positions where bw is true) for each location in bw.

In cases where two or more objects are at the same distance, an arbitrary object is chosen.

For example, the nearest true in bw_example[1,1]exists at CartesianIndex(2, 2), so it is assigned the value CartesianIndex(2, 2).

In [ ]: 
bw_transform = feature_transform(bw)
bw_transform_example = feature_transform(bw_example)
Out[0]:
5×5 Matrix{CartesianIndex{2}}:
 CartesianIndex(2, 2)  CartesianIndex(2, 2)  …  CartesianIndex(2, 4)
 CartesianIndex(2, 2)  CartesianIndex(2, 2)     CartesianIndex(2, 4)
 CartesianIndex(3, 2)  CartesianIndex(3, 2)     CartesianIndex(3, 4)
 CartesianIndex(4, 2)  CartesianIndex(4, 2)     CartesianIndex(4, 4)
 CartesianIndex(4, 2)  CartesianIndex(4, 2)     CartesianIndex(4, 4)

the distance transform function performs transformation by specified labels.

For example, bw_transform[1,1] has CartesianIndex(2, 2), and D[i] for this element will be the distance between CartesianIndex(1, 1) and CartesianIndex(2, 2), which is equal to sqrt(2).

In [ ]: 
dist = 1 .- distance_transform(bw_transform)
dist_example = 1 .- distance_transform(bw_transform_example)
Out[0]:
5×5 Matrix{Float64}:
 -0.414214  0.0  0.0  0.0  -0.414214
  0.0       1.0  1.0  1.0   0.0
  0.0       1.0  1.0  1.0   0.0
  0.0       1.0  1.0  1.0   0.0
 -0.414214  0.0  0.0  0.0  -0.414214

The label_components function finds related components in the dist_trans binary array. You can provide a list indicating which dimensions are used to determine connectedness. For example, region = [1,3] will not be checked for connectedness of neighbours on dimension 2. This only corresponds to nearest neighbours, i.e. 4-connectivity in 2d and 6-connectivity in 3d matrices. Default value region = 1:ndims(A): . The output label is an integer array, where 0 is used for background pixels and each individual region of connected pixels gets its own integer index.

In [ ]: 
dist_trans = dist .< 1
markers = label_components(dist_trans)
markers_example = label_components(dist_example .< 0.5)
Gray.(markers/32.0)
Out[0]:
No description has been provided for this image

Function watershed Method segments the image using a watershed transformation.

Each basin formed by the watershed transformation corresponds to a segment. In the marker matrix, zero means no marker. If two markers have the same index, their regions will be merged into one region.

In [ ]: 
segments = watershed(dist, markers)
segments_example = watershed(dist_example , markers_example)
Out[0]:
Segmented Image with:
  labels map: 5×5 Matrix{Int64}
  number of labels: 1

The figure below shows from left to right the original image, the marker regions and the final segmented image, dividing the image segments by colour.

In [ ]: 
labels = labels_map(segments)
colored_labels = IndirectArray(labels, distinguishable_colors(maximum(labels)))
masked_colored_labels = colored_labels .* (1 .- bw)
mosaic(img, colored_labels, masked_colored_labels; nrow=1)
Out[0]:
No description has been provided for this image

Conclusion

In this demonstration, we have shown how to use image indexing to segment individual image regions. These techniques can be useful in isolating regions of interest from the original image.