Segmentation module

struct dip::PerObjectEllipsoidFitParameters

#include "diplib/segmentation.h"

Defines the parameters for the PerObjectEllipsoidFit function.

void dip::Watershed(dip::Image const& in, dip::Image const& mask, dip::Image& out, dip::uint connectivity = 1, dip::dfloat maxDepth = 1, dip::uint maxSize = 0, dip::StringSet flags = {})

#include "diplib/morphology.h"

Computes the watershed of in within mask, with on-line merging of regions.

The watershed is a segmentation algorithm that divides the image according to its grey-value ridges.

connectivity determines which pixels are considered neighbors; the default value of 1 leads to vertex-connected watershed lines (i.e. thinnest possible result). See Connectivity for information on the connectivity parameter.

flags determines how the output is computed. There are three options:

• "labels" or "binary": returns either the labels used during processing, with the watershed lines as background (value 0), or a binary image where the watershed lines are set and the regions are not set. "binary" is the default.

• "low first" or "high first": determines the sort order of pixels. The default of "low first" yields the normal watershed, where local minima are origin of the basins, and the watershed lines run along the high ridges in the image. "high first" simply inverts the definition, such that local maxima are at the centers of the basins, and the watershed lines run along the low valleys.

• "fast" or "correct": determines which algorithm is used:

  • "fast" (the default) is an algorithm that takes a few shortcuts, but usually manages to produce good results any way. One shortcut leads to all border pixels being marked as watershed lines. It is possible to extend the image by one pixel before processing to circumvent this. The other shortcut means that plateaus are not handled correctly. A plateau is a region in the image where pixels have exactly the same value. This is usually seen as watershed lines not running through the middle of the plateaus, instead being shifted to one side. Adding a little bit of noise to the image, and setting maxDepth to the range of the noise, usually improves the results in these cases a little bit.

  • "correct" is an algorithm that first finds the local minima through dip::Minima (or maxima if "high first" is set), and then applies dip::SeededWatershed. This always produces correct results, but is significantly slower.

The on-line region merging works as follows: When two regions first meet, a decision is made on whether to keep the regions separate (and thus put a watershed pixel at that point), or to merge the regions. If one of the regions is no deeper than maxDepth (i.e. the intensity difference between the region’s minimum and the pixel where the region meets another), and is no larger than maxSize (i.e. the number of pixels belonging to the region and that have been seen so far), then it can be merged. The merged region is subsequently treated as a single region, and their labels are considered equal. If maxSize is zero, no test for size is done. In this case, the merging is exactly equivalent to applying an H-minima transform to the image before computing the watershed.

Note that for the "fast" algorithm, maxDepth is always at least 0 (negative values will be ignored). That is, two regions without a grey-value difference between them (they are on the same plateau) will always be merged. This is necessary to prevent unexpected results (i.e. a plateau being split into multiple regions). For the "correct" algorithm, any negative value of maxDepth will disable the merging. But note that, due to the way that the region seeds are computed (dip::Minima), setting maxDepth to 0 would lead to the exact same result.

Any pixel that is infinity will be part of the watershed lines, as is any pixel not within mask.

void dip::SeededWatershed(dip::Image const& in, dip::Image const& seeds, dip::Image const& mask, dip::Image& out, dip::uint connectivity = 1, dip::dfloat maxDepth = 1, dip::uint maxSize = 0, dip::StringSet const& flags = {})

#include "diplib/morphology.h"

Computes the watershed of in within mask, starting at seeds, with on-line merging of regions.

seeds is a binary or labeled image (if binary, it is labeled using connectivity). These labels are iteratively expanded in the watershed order (i.e. pixels that have a low value in in go first) until they meet. Pixels where two regions meet are marked as the watershed lines. seeds is commonly used to direct the segmentation, and merging is consequently not necessary. However, this algorithm does include on-line merging. Note that two labeled regions in seeds that do not have a grey-value ridge between them (i.e. they are on the same plateau) will be merged unless merging is disabled (see below). Merged labels will be painted with the label of one of the originating seeds, and the other labels will not be present in the output (only if flags contains "labels").

connectivity determines which pixels are considered neighbors; the default value of 1 leads to vertex-connected watershed lines (i.e. thinnest possible result). See Connectivity for information on the connectivity parameter.

The region merging and the flags parameter work as described in dip::Watershed, with the following additions:

• If maxDepth is negative, regions will never be merged, even if they have no grey-value difference between them.
• The flags values "fast" or "correct" are not allowed.
• flags can contain the string "no gaps", which prevents the formation of watershed lines in between the regions. That is, seeds are grown until they touch. This flag implies the flag "labels", since in a binary image there would be no distinction between initially separate regions. Pixels that have an infinite value in in, or a zero value in mask, will still be excluded from the region growing process.
• flags can contain the string "uphill only", which will limit the region growing to be exclusively uphill (or downhill if "high first" is also given). This means that regions will grow to fill the local catchment basin, but will not grow into neighboring catchment basins that have no seeds. This flag will also disable any merging.

void dip::CompactWatershed(dip::Image const& in, dip::Image const& seeds, dip::Image const& mask, dip::Image& out, dip::uint connectivity = 1, dip::dfloat compactness = 1.0, dip::StringSet const& flags = {})

#include "diplib/morphology.h"

Computes the compact watershed of in within mask, starting at seeds.

seeds is a binary or labeled image (if binary, it is labeled using connectivity). These labels are iteratively expanded in the watershed order (i.e. pixels that have a low value in in go first), modified with a compactness term, until they meet. Pixels where two regions meet are marked as the watershed lines.

The compactness term modifies the normal watershed order by taking into account the distance to the originating seed. This distance, multiplied by compactness, is added to the grey value when determining the processing order. A compactness of 0 leads to the normal seeded watershed, and a very large value for compactness leads to disregarding the pixel values in in, thereby creating a Voronoi diagram.

connectivity determines which pixels are considered neighbors; the default value of 1 leads to vertex-connected watershed lines (i.e. thinnest possible result). See Connectivity for information on the connectivity parameter.

The flags parameter work as described in dip::SeededWatershed, except that "uphill only" is not supported.

void dip::StochasticWatershed(dip::Image const& in, dip::Image& out, dip::Random& random, dip::uint nSeeds = 100, dip::uint nIterations = 50, dip::dfloat noise = 0, dip::String const& seeds = S::HEXAGONAL)

#include "diplib/morphology.h"

Computes the stochastic watershed of in.

The stochastic watershed is computed by applying a watershed with randomly placed seeds nIterations times, and adding the results. The output is an image where each pixel’s value is the likelihood that it belongs to an edge in the image, the values are in the range [0,nIterations]. The input image in should contain high grey values at the edges of the regions to be segmented. Thresholding out at an appropriate value will yield the relevant edges in the image. Alternatively, apply dip::Watershed to the result, with maxDepth set to the appropriate threshold value.

The number of seeds used is given by nSeeds. Actually seeds are chosen with a density of nSeeds / in.NumberOfPixels(), the random process causes the actual number of seeds to differ between runs. Seeds are placed either through a Poisson point process (seeds is "poisson") or a randomly translated and rotated grid (seeds is "rectangular" (any number of dimensions), "hexagonal" (2D only), or "bcc" or "fcc" (3D only)). The output contains counts, in the range [0,nIterations].

If seeds is "exact", or if nIterations is 0, then the exact probabilities are computed (Malmberg and Luengo, 2014). The output contains probabilities, in the range [0,1]. Note that this algorithm requires O(n²) space, and is not suitable for very large images.

The stochastic watershed expects the image to contain roughly equally-sized regions. nSeeds should be approximately equal to the number of expected regions. If there is a strong difference in region sizes, larger regions will be split into smaller ones.

If the image contains regions with different sizes, it is recommended to set noise to a value that is larger than the variation within regions, but smaller than the height of the barrier between regions. Uniform noise will be added to the input image for every iteration of the process, causing non-significant edges to be strongly suppressed (Bernander et al., 2013). In the case of the exact stochastic watershed, the operation is applied three times with random noise added to the input, and the geometric mean of the results is returned (Selig et al., 2015).

in must be real-valued and scalar. out will be of a suitable unsigned integer type (depending on the number of iterations, but typically dip::DT_UINT8), or of type dip::DT_SFLOAT if the exact stochastic watershed is computed.

void dip::WatershedMinima(dip::Image const& in, dip::Image const& mask, dip::Image& out, dip::uint connectivity = 1, dip::dfloat maxDepth = 1, dip::uint maxSize = 0, dip::String const& output = S::BINARY)

#include "diplib/morphology.h"

Marks significant local minima.

This algorithm works exactly like dip::Watershed with the "fast" flag set. All pixels with a value equal to the lowest value within each watershed basin form a local minimum. Note that they can form disconnected regions, use the "labels" flag to recognize such disconnected regions as a single local minimum. See dip::Watershed for a description of all the parameters.

output can be "binary" or "labels", and determines whether the algorithm outputs a binary image or a labeled image.

See Connectivity for information on the connectivity parameter.

void dip::WatershedMaxima(dip::Image const& in, dip::Image const& mask, dip::Image& out, dip::uint connectivity = 1, dip::dfloat maxDepth = 1, dip::uint maxSize = 0, dip::String const& output = S::BINARY)

#include "diplib/morphology.h"

Marks significant local maxima.

This algorithm works exactly like dip::Watershed with the "fast" flag set. All pixels with a value equal to the highest value within each watershed basin form a local maximum. Note that they can form disconnected regions, use the "labels" flag to recognize such disconnected regions as a single local maximum. See dip::Watershed for a description of all the parameters.

output can be "binary" or "labels", and determines whether the algorithm outputs a binary image or a labeled image.

See Connectivity for information on the connectivity parameter.

void dip::Minima(dip::Image const& in, dip::Image& out, dip::uint connectivity = 0, dip::String const& output = S::BINARY)

#include "diplib/morphology.h"

Marks local minima.

This algorithm finds single pixels or plateaus (connected groups of pixels with identical value) that are surrounded by pixels with a higher value. If output is "binary", the result is a binary image where these pixels and plateaus are set. If output is "labels", the result is a labeled image.

See Connectivity for information on the connectivity parameter.

void dip::Maxima(dip::Image const& in, dip::Image& out, dip::uint connectivity = 0, dip::String const& output = S::BINARY)

#include "diplib/morphology.h"

Marks local maxima.

This algorithm finds single pixels or plateaus (connected groups of pixels with identical value) that are surrounded by pixels with a lower value. If output is "binary", the result is a binary image where these pixels and plateaus are set. If output is "labels", the result is a labeled image.

See Connectivity for information on the connectivity parameter.

dip::CoordinateArray dip::KMeansClustering(dip::Image const& in, dip::Image& out, dip::Random& random, dip::uint nClusters = 2)

#include "diplib/segmentation.h"

Applies k-means clustering to an image, yielding nClusters labeled regions.

in is a scalar, real-valued image. nClusters cluster centers are found, centered on regions of high intensity. out is a labeled image with nClusters regions tiling the image. Each region is identified by a different label. Boundaries between regions are the Voronoi tessellation given the identified cluster centers.

Note that this creates a spatial partitioning, not a partitioning of image intensities.

K-means clustering is an iterative process with a random initialization. It is likely to get stuck in local minima. Repeating the clustering several times and picking the best result (e.g. determined by times each cluster center is found) can be necessary.

The returned dip::CoordinateArray contains the cluster centers. Element i in this array corresponds to label i+1.

dip::CoordinateArray dip::KMeansClustering(dip::Image const& in, dip::Image& out, dip::uint nClusters = 2)

#include "diplib/segmentation.h"

Like above, using a default-initialized dip::Random object.

dip::CoordinateArray dip::MinimumVariancePartitioning(dip::Image const& in, dip::Image& out, dip::uint nClusters = 2)

#include "diplib/segmentation.h"

Spatially partitions an image into nClusters partitions iteratively, minimizing the variance of the partitions.

Minimum variance partitioning builds a k-d tree, where, for each node, the orthogonal projection with the largest variance is split using the same logic as Otsu thresholding applies to a histogram. Note that this creates a spatial partitioning, not a partitioning of image intensities. out is a labeled image with nClusters regions tiling the image. Each region is identified by a different label.

Minimum variance partitioning is much faster than k-means clustering, though its result might not be as good. It is also deterministic.

in must be scalar and real-valued.

The returned dip::CoordinateArray contains the centers of gravity for each cluster. Element i in this array corresponds to label i+1.

dip::FloatArray dip::IsodataThreshold(dip::Image const& in, dip::Image const& mask, dip::Image& out, dip::uint nThresholds = 1)

#include "diplib/segmentation.h"

Thresholds the image in using nThresholds thresholds, determined using the Isodata algorithm (k-means clustering), and the histogram of in.

Only those pixels in mask are used to determine the histogram on which the Isodata algorithm is applied, but the threshold is applied to the whole image. in must be scalar and real-valued.

If nThresholds is 1, then out is a binary image. With more thresholds, the output image is labeled.

The output array contains the thresholds used.

See dip::IsodataThreshold for more information on the algorithm used.

dip::dfloat dip::OtsuThreshold(dip::Image const& in, dip::Image const& mask, dip::Image& out)

#include "diplib/segmentation.h"

Thresholds the image in using the maximal inter-class variance method by Otsu, and the histogram of in.

Only those pixels in mask are used to determine the histogram on which the threshold estimation algorithm is applied, but the threshold is applied to the whole image. in must be scalar and real-valued.

Returns the threshold value used.

See dip::OtsuThreshold for more information on the algorithm used.

dip::dfloat dip::MinimumErrorThreshold(dip::Image const& in, dip::Image const& mask, dip::Image& out)

#include "diplib/segmentation.h"

Thresholds the image in using the minimal error method by Kittler and Illingworth, and the histogram of in.

Only those pixels in mask are used to determine the histogram on which the threshold estimation algorithm is applied, but the threshold is applied to the whole image. in must be scalar and real-valued.

Returns the threshold value used.

See dip::MinimumErrorThreshold for more information on the algorithm used.

dip::FloatArray dip::GaussianMixtureModelThreshold(dip::Image const& in, dip::Image const& mask, dip::Image& out, dip::uint nThresholds = 1)

#include "diplib/segmentation.h"

Thresholds the image in using nThresholds thresholds, determined by fitting a Gaussian Mixture Model to the histogram of in.

Only those pixels in mask are used to determine the histogram on which the Gaussian Mixture Model algorithm is applied, but the threshold is applied to the whole image. in must be scalar and real-valued.

If nThresholds is 1, then out is a binary image. With more thresholds, the output image is labeled.

The output array contains the thresholds used.

See dip::GaussianMixtureModelThreshold for more information on the algorithm used.

dip::dfloat dip::TriangleThreshold(dip::Image const& in, dip::Image const& mask, dip::Image& out, dip::dfloat sigma = 4.0)

#include "diplib/segmentation.h"

Thresholds the image in using the chord method (a.k.a. skewed bi-modality, maximum distance to triangle), and the histogram of in.

Only those pixels in mask are used to determine the histogram on which the threshold estimation algorithm is applied, but the threshold is applied to the whole image. in must be scalar and real-valued.

Returns the threshold value used.

See dip::TriangleThreshold for more information on the algorithm used and the sigma parameter.

dip::dfloat dip::BackgroundThreshold(dip::Image const& in, dip::Image const& mask, dip::Image& out, dip::dfloat distance = 2.0, dip::dfloat sigma = 4.0)

#include "diplib/segmentation.h"

Thresholds the image in using the unimodal background-symmetry method, and the histogram of in.

Only those pixels in mask are used to determine the histogram on which the threshold estimation algorithm is applied, but the threshold is applied to the whole image. in must be scalar and real-valued.

Returns the threshold value used.

See dip::BackgroundThreshold for more information on the algorithm used and the sigma parameter.

dip::dfloat dip::VolumeThreshold(dip::Image const& in, dip::Image const& mask, dip::Image& out, dip::dfloat volumeFraction = 0.5)

#include "diplib/segmentation.h"

Thresholds an image such that a fraction volumeFraction of pixels is foreground.

Only pixels within mask are used to determine the threshold value, but the threshold is applied to the whole image. in must be scalar and real-valued.

The return value is the threshold applied.

void dip::FixedThreshold(dip::Image const& in, dip::Image& out, dip::dfloat threshold, dip::dfloat foreground = 1.0, dip::dfloat background = 0.0, dip::String const& output = S::BINARY)

#include "diplib/segmentation.h"

Thresholds an image at the threshold value.

If output is "binary" (the default), FixedThreshold will produce a binary image. Otherwise an image of the same type as the input image is produced, with the pixels set to either foreground or background. In other words, on a pixel-per-pixel basis the following is applied: out = ( in >= threshold ) ? foreground : background.

in must be real-valued, each tensor element is thresholded independently.

Note that, for the “binary” output case, it might be easier to write:

out = in >= threshold;

void dip::RangeThreshold(dip::Image const& in, dip::Image& out, dip::dfloat lowerBound, dip::dfloat upperBound, dip::String const& output = S::BINARY, dip::dfloat foreground = 1.0, dip::dfloat background = 0.0)

#include "diplib/segmentation.h"

Thresholds an image at two values, equivalent to lowerBound <= in && in <= upperBound.

If output is "binary" (the default), RangeThreshold will produce a binary image. If foreground == 0.0, foreground will be set to false and background to true, otherwise the foreground will be true (this is the default).

If output is not "binary", an image of the same type as the input image is produced, with the pixels set to either foreground or background. In other words, on a pixel-per-pixel basis the following is applied:

out = ( lowerBound <= in && in <= upperBound ) ? foreground : background

in must be real-valued, each tensor element is thresholded independently.

void dip::HysteresisThreshold(dip::Image const& in, dip::Image& out, dip::dfloat lowThreshold, dip::dfloat highThreshold)

#include "diplib/segmentation.h"

Hysteresis threshold.

From the binary image in >= lowThreshold only those connected regions are selected for which at least one location also has in >= highThreshold.

The output image will be a binary image with foreground pixels == 1 and background pixels == 0.

in must be scalar and real-valued.

void dip::MultipleThresholds(dip::Image const& in, dip::Image& out, dip::FloatArray const& thresholds)

#include "diplib/segmentation.h"

Thresholds an image at multiple values, yielding a labeled image.

out will be a dip::DT_UINT8, dip::DT_UINT16, dip::DT_UINT32 or dip::DT_UINT64 image, depending on the length of thresholds. All pixels below thresholds[ 0 ] with be assigned the label 0, all pixels greater or equal to thresholds[ 0 ] and smaller than thresholds[ 1 ] will be assigned label 1, etc. Results might not be as expected if thresholds are not sorted.

in must be scalar and real-valued.

dip::dfloat dip::Threshold(dip::Image const& in, dip::Image const& mask, dip::Image& out, dip::String const& method = S::OTSU, dip::dfloat parameter = infinity)

#include "diplib/segmentation.h"

Automated threshold using method.

This function computes an optimal threshold value for in using method, and applies it. Returns the found threshold value. in must be scalar and real-valued. mask can optionally select the pixels used to determine the threshold value. The threshold is applied to the image as a whole, you can combine it with the mask afterwards:

dip::Image bin = dip::Threshold( image, mask, "otsu" );
bin &= mask;

method can be one of:

• "isodata": see dip::IsodataThreshold.
• "otsu": see dip::OtsuThreshold. This is the default method.
• "minerror": see dip::MinimumErrorThreshold”.
• "gmm": see dip::GaussianMixtureModelThreshold.
• "triangle": see dip::TriangleThreshold.
• "background": see dip::BackgroundThreshold.
• "volume": see dip::VolumeThreshold.
• "fixed": see dip::FixedThreshold. The default parameter value is 128.

If parameter is dip::infinity, the default parameter value for the method will be used.

void dip::PerObjectEllipsoidFit(dip::Image const& in, dip::Image& out, dip::PerObjectEllipsoidFitParameters const& parameters)

#include "diplib/segmentation.h"

Finds a per-object threshold such that found objects are maximally ellipsoidal.

This function thresholds the image such that all objects found are approximately ellipsoidal, within the bounds expressed by parameters. Each object is thresholded at a different level, chosen to maximize its fit to an ellipsoid. The measure maximized is the ratio of the object’s size (area or volume) to the size of the fitted ellipsoid. Ellipsoids are fitted by determining the ellipsoid with the same second order central moments as the object at the given threshold level.

in must be scalar, real-valued, and be 2D (TODO: port the 3D version of this function also). out will be binary and of the same sizes as in.

void dip::Canny(dip::Image const& in, dip::Image& out, dip::FloatArray const& sigmas = {1}, dip::dfloat lower = 0.5, dip::dfloat upper = 0.9, dip::String const& selection = S::ALL)

#include "diplib/segmentation.h"

Detect edges in the grey-value image by finding salient ridges in the gradient magnitude

The Canny edge detector finds the ridges in the gradient magnitude of in, which correspond to the edges in the image. The gradient magnitude (see dip::GradientMagnitude) is computed using Gaussian derivatives, with a sigma of sigma. The found ridges are pruned to remove the less salient edges (see dip::NonMaximumSuppression). Next, a threshold t1 is computed such that the 1 - upper fraction of pixels with the highest gradient magnitude are kept. A second threshold, t2 = t1 * lower, is selected that determines the minimal gradient magnitude expected for an edge. All edge pixels equal or larger to t2, and are in the same connected region as at least one pixel that is equal or larger to t1, are selected as the output of this function (see dip::HysteresisThreshold). Finally, a homotopic thinning is applied to reduce the detections to single-pixel–thick lines (see dip::EuclideanSkeleton).

The 1 - upper fraction is computed over all pixels in the image by default. If the image has relatively few edges, this can lead to t1 being equal to 0. If this happens, the hysteresis threshold would select all pixels in the image, and the homotopic thinning will lead to a line across the image that is unrelated to any edges. Instead, t1 will be set to a value slightly larger than 0.

For more control over the thresholds, the selection parameter can be set to "nonzero", in which case the fraction 1 - upper refers to non-zero pixels only; or to "absolute", in which case upper and lower represent absolute threshold values, and t1 will be set to upper and t2 will be set to lower.

in must be scalar, real-valued, and have at least one dimension.

The Canny edge detector was originally described, and typically implemented, for 2D images only. Here we provide an obvious extension to arbitrary dimensions. The final homotopic thinning is only applied in 2D and 3D, since dip::EuclideanSkeleton is not defined for other dimensionalities.

void dip::Superpixels(dip::Image const& in, dip::Image& out, dip::Random& random, dip::dfloat density = 0.005, dip::dfloat compactness = 1.0, dip::String const& method = S::CW, dip::StringSet const& flags = {})

#include "diplib/segmentation.h"

Generates superpixels (oversegmentation)

density indicates how many superpixels, on average, should be created. It is given in superpixels per pixel. That is, 1/density is the average size of the superpixels.

compactness controls the shape of the superpixels. Reducing this value leads to superpixels that more precisely follow image contours, but also are more varied in size and shape. Increasing this value leads to more isotropic superpixels and less variation in size.

method controls the method used to generate superpixels. Currently only "CW" is supported. This is the compact watershed superpixel segmentation (Neubert and Protzel, 2014).

flags can contain the following flags:

• "rectangular" (default) or "hexagonal": controls the basic shape of the superpixels (the shape they tend towards as compactness increases). For 3D images, "hexagonal" implies an FCC grid (see dip::FillRandomGrid). For images with more than 3 dimensions, "rectangular" will always be used.

• "no gaps" indicates that the superpixels must cover the whole image. By default a 1-pixel gap is left in between superpixels.

in must be real-valued. If not scalar, the norm of the gradient magnitude for each tensor element is used to determine where edges are located. In the case of a color image, no color space conversion is performed, the image is used as-is. It is recommended to pass an image in an appropriate color space for superpixel segmentation, such as CIE Lab.

void dip::Superpixels(dip::Image const& in, dip::Image& out, dip::dfloat density = 0.005, dip::dfloat compactness = 1.0, dip::String const& method = S::CW, dip::StringSet const& flags = {})

#include "diplib/segmentation.h"

Like above, using a default-initialized dip::Random object.