Channel Arithmetic and Beyond

View on GitHub

This program enables one to specify arithmetic expressions which are used to create new channels. The basic arithmetic operations are supported: +,-,*,/,**. More advanced operations that run short SimpleITK code snippets are also supported. SimpleITK provides hundreds of filters that can be used via this program.

Channels are referenced using square brackets and the channel index, starting at zero. To apply an expression to all channels, use the channel index 'i'.

When creating a single new channel, the arithmetic expression consists of literal channel numbers, one can select a name and color for the new channel. When creating multiple new channels, the arithmetic expression is applied to all channels, the postfix '_modified' is appended to the original channel names and the original color is copied over. Note that for all channels created by the program the channel description will include the arithmetic expression used to create that channel. This transparently supports your efforts to conduct reproducible research.

Because an Imaris image has a specific pixel type (8, 16, 32 bit unsigned integer and 32 bit floating point) all computations are performed using a 32 bit floating point representation and then clamped to the range of the image's pixel type.

The program allows you to use the same expression on multiple files. In this case literal channel values are limited by the number of shared channels. Thus, if working with two files one with three channels and one with four channels, the valid literal channel values are limited to 0, 1, 2. We cannot use 3 as it does not exist in all files. On the other hand, if our autofluorescence channel is one of these channels, e.g. channel 0, we can subtract it from all channels in both files, [i]-[0].

Basic Examples

Multiply channels zero and three:
```
[0]*[3]
```
Multiply channels zero and three and subtract the result from channel four:
```
[4] - ([0]*[3])
```
Duplicate the second channel (indexes are zero based, so second channel index is 1):
```
[1]
```
Duplicate all channels:
```
[i]
```
Subtract channel zero from all channels:
```
[i]-[0]
```
Invert channel intensity (new intensity values are the result of subtracting the current intensity value from the maximal possible value):
```
sitk.InvertIntensity([1])
```
Threshold channel one using a value of 100, resulting image is binary with values in {0,1}:
```
[1]>100
```
Threshold a specific channel to create a binary result using the Otsu filter, resulting image is binary with values in {0,1}:
```
sitk.OtsuThreshold([1], 0, 1)
```

Gaussian blurring (variance specified in metric units, e.g. nm):

# blur channel 3 using a Gaussian with variance of 0.245 for x, y and z
sitk.DiscreteGaussian([3], 0.245)
# blur channel 3 using a Gaussian with different variance per x=0.245, y=0.1, z=0.3
sitk.DiscreteGaussian([3], [0.245, 0.1, 0.3])

Unsharp masking (variance specified in metric units, e.g. nm). This is actually a filter: that emphasizes edges, sharpens the image using a blurred version of the image.
```
# sharpen channel 0 using an unsharp Gaussian mask with variance of 0.245 and
# sharpening increase of 0.5
[0]+([0]-sitk.DiscreteGaussian([0], 0.245))*0.5
```
Median filtering using a radius of x=1, y=2 and z=3 (radius in which to compute the median per dimension, pixel units):
```
sitk.Median([3], [1, 2, 3])
```
Linearly map the intensity values in the interval [windowMinimum, windowMaximum] to the interval [outputMinimum, outputMaximum]. Values lower than windowMinimum are mapped to outputMinimum. Values higher than windowMaximum are mapped to outputMaximum.
```
sitk.IntensityWindowing([1], windowMinimum=20, windowMaximum=200, outputMinimum=0, outputMaximum=255)
```

Advanced Examples

Threshold a specific channel retaining the values above the threshold and setting all values equal or lower to zero:
```
sitk.Cast([1]>100, sitk.sitkFloat32)*[1]
```
Threshold a specific channel retaining the values above the threshold and setting all values equal or lower to the threshold to an arbitrary value (20 in this example):
```
sitk.Cast([1]>100, sitk.sitkFloat32)*[1] + sitk.Cast([1]<=100, sitk.sitkFloat32)*20
```
Threshold a specific channel, get all connected components, then sort the components according to size, discarding those smaller than a minimum size and create a binary mask corresponding to the largest component, which is the first label(second largest component label is 2 etc.):
```
sitk.RelabelComponent(sitk.ConnectedComponent([1]>100), minimumObjectSize = 50)==1
```
Create a binary mask representing the colocalization of two channels, intensity values below 20 are considred noise:
```
([1]>20)*([2]>20)
```
Create a binary mask representing the colocalization of two channels. We are interested in all pixels in channel 2 that have a value above 20 and that are less than 1.0um away from pixels in channel 1 that have a value above 100 (note: this operation yields different results when run using a slice-by-slice approach vs. a volumetric approach):
```
(sitk.Cast([2]>20, sitk.sitkFloat32) *
sitk.Abs(sitk.SignedMaurerDistanceMap([1]>100, insideIsPositive=False, squaredDistance=False, useImageSpacing=True)))<=1.0
```

Create a binary mask using thresholding and then perform morphological closing (dilation followed by erosion) with distance maps, useful for filling holes:

sitk.SignedMaurerDistanceMap(sitk.SignedMaurerDistanceMap([1]>100, insideIsPositive=False, squaredDistance=False, useImageSpacing=True) < 1.0, insideIsPositive=False, squaredDistance=False, useImageSpacing=True)<-1.0

Create a binary mask using thresholding and then perform morphological opening (erosion followed by dilation) with distance maps, useful for removing small islands:

sitk.SignedMaurerDistanceMap(sitk.SignedMaurerDistanceMap([1]>100, insideIsPositive=False, squaredDistance=False, useImageSpacing=True) < -0.2, insideIsPositive=False, squaredDistance=False, useImageSpacing=True)<0.2