Segmentation of image based fields

Example a_backup/segmentation: Segmentation of image based fields

This example demonstrates how to segment image based computed fields. It also shows how to use the results of the segmentation to generate data points to describe the boundary of an anatomical feature.

Segmentation refers to the process of partitioning an image into various regions. These regions are identified by sets of pixels with the same intensity values. Segmentation of medical images is typically used to identify the boundaries between different types of tissue and pinpoint the position of objects such as bones and organs.

The human brain is very good at segmenting images but it is quite a difficult task to automate. Many different segmentation algorithms have been developed and for a given image one approach may produce better results than another.

Cmgui implements segmentation via computed fields, using the Insight Tookit (ITK) library code. For background information on segmentation see chapter 9 of the ITK software guide. The latest version of the guide can be found at www.itk.org/ItkSoftwareGuide.pdf

Currently cmgui supports the following segmentation algorithms:

Connected threshold, a region growing algorithm (see section 9.1.1 of the ITK software guide)
Fast marching, a level set segmentation algorithm (see section 9.3.1 of the ITK software guide)

This example demonstrates the use of both methods to segment out different areas from a sample brain image

Connected threshold segmentation

The filter pipeline used for the connected threshold segmentation is as follows:

Curvature

anisotropic

diffusion

filter

Connected

threshold

filter

Input image

Smoothed

Segmented image

The images show the stages in segmenting out the left and right ventricle. Two seed points were used for the connected threshold filter, one within the left ventricle and one within the right.

Fast marching segmentation

The filter pipeline used for the fast marching segmentation is as follows:

Curvature

anisotropic

diffusion

filter

Gradient

magnitude

recursive

gaussian

filter

Sigma

filter

Fast

marching

filter

Binary

threshold

filter

Input image

Smoothed

Gradient

Speed term

Time-crossing map

Segmented image

The images show the stages in segmenting out the right ventricle and the gray matter. Two seed points were used for the fast marching filter, one within the right ventricle and one within the left half of the gray matter.

Here the speed term is calculated using a sigmoid function, but it could also be calculated using the negative exponential exp(-x) or the reciprocal 1/(1+x).

Created by Peter Bier, July 2007.

Screenshot of example a_backup/segmentation

The comfile run by this example is as follows:

# This example demonstrates the use of two methods which can be used to
# segment out different regions from a sample brain image.
# The segmentation algorithms are described in more depth in the
# ITK software guides, www.itk.org/ItkSoftwareGuide.pdf
# Once the image has been segmented, data points can be generated
# on the boundaries of the segmented regions, as shown at the
# end of this example.
#

# CONNECTED THRESHOLD SEGMENTATION

# Connected threshold segmentation is a very simple region growing method.
# One or more seed points are specified for the image and the region grows 
# from those points, adding neighbouring pixels to the region if the pixel
# intensity falls within a specific interval.  The regions grow until no 
# more pixels can be added.

# Noise present in an input image can reduce the effectivness of this filter, 
# so it is often useful to remove noise by first passing the image through
# an edge preserving smoothing filter.

# The filter pipeline is as follows (inputs to the filters shown in brackets):

#     Input image
#        |
#        V
# Curvature anisotropic diffusion filter (iterations, time step, conductance)
#        |
#        V
#     Smoothed image
#        |
#        V
# Connected threshold filter (lower threshold, upper threshold, seed points)
#        |
#        V
#     Output image

# We will use the connected threshold filter to segment out the left and
# right ventricle.

# Read in a mesh to hold image which is currently required for cmiss.
gfx read nodes $example/square.exnode
gfx read elements $example/square.exelem

# read in the image to a texture and define a field based on it
gfx create texture input_image image png:$example/BrainProtonDensitySlice.png nearest width 1 height 1;
gfx create material input_image texture input_image

# Set up a dummy image of the correct size for the output field
# This will be used to write out filtered images
gfx create texture output_image image png:$example/BrainProtonDensitySlice.png specify_format i nearest width 1 height 1;
gfx create material output_image texture output_image

# create a spectrum
gfx create spectrum gray_spectrum clear overwrite_colour
gfx modify spectrum gray_spectrum linear range 0 1 extend_above extend_below monochrome colour_range 0 1 ambient diffuse component 1

# Set up the pipeline for the connected threshold filter

# Smooth the image using anisotropic diffusion with 5 iterations
# Note that if trying to duplicate the itk segmentation tests that the 
# conductance and time step do not require dividing by 255
gfx define field smoothed curvature_anisotropic_diffusion_filter field input_image num_iterations 5 time_step 0.125 conductance 9.0;

# Segment the image using the connected threshold filter
# The seed points are defined in terms of their position within the texture
# described by normalised coordinates.
# In terms of pixel value the seed points are (85,107) and (96,107)
# The image pixel dimensions are 181x217 which gives the normallised seed 
# points as (0.47,0.5) and (0.53,0.5)
# Pixels will be added to the region if the fall within the intensity range
# 0.85 to 1.0
gfx define field connected_threshold_segmentation connected_threshold_filter field smoothed num_seed_points 2 dimension 2 seed_points 0.47 0.5 0.53 0.5 lower_threshold 0.85 upper_Threshold 1.0

# write out the segmented image to a file
gfx modify texture output_image nearest width 1 height 1 evaluate_image field connected_threshold_segmentation spectrum gray_spectrum texture_coord xi element_group square;
gfx write texture output_image file connected_threshold_segmentation.png;

# FAST MARCHING SEGMENTATION

# Fast marching segmentation is a level set segmentation method which can be 
# used if the equation governing the level set evolution is very simple. 
# Level set segmentation is a very powerful technique that uses a numerical 
# method for tracking the evolution of contours and surfaces, based on terms
# derived from an input image.

# In the case of the fast marching segmentation filter, the "speed term"
# required to calculate the segmentation is typically some function of the
# image gradient.  As with the connected threshold example, the input image is
# usually smoothed first, before computing the gradient and then some
# function of the gradient to give the speed image.  The fast marching
# filter then produces a time map, indicating how long it takes to get from
# a seed point to any given point on the image.  By using a threshold filter
# on the time map a segemented image can be obtained.

# The fast marching filter segmentation method demonstrated here uses the same 
# filters outlined in the collaboration diagram in the ITK software guide 
# (See Figure 9.15)
#
# The filter pipeline is as follows (inputs to the filters shown in brackets):

#      Input image
#         |
#         V
# Curvature anisotropic diffusion filter (iterations, time step, conductance)
#         |
#         V
#      Smoothed image
#         |
#         V
# Gradient magnitude recursive gaussian (sigma)
#         |
#         V
#      Gradient image
#         |
#         V
# Sigmoid filter (alpha, beta)
#         |
#         V
#      Speed term
#         |
#         V
# Fast marching filter (seed points, stopping value)
#         |
#         V
#      Time-crossing map
#         |
#         V
# Binary Threshold filter (lower threshold, upper threshold)
#         |
#         V
#      Segmented image


# We will use the fast marching image filter to segment out the white
# matter and the right ventricle.

# Assemble the pipeline for the fast marching image filter segmentation

# Smooth the image using anisotropic diffusion with 5 iterations
# Note that if trying to duplicate the itk segmentation tests that the 
# conductance and time step do not require dividing by 255
gfx define field smoothed curvature_anisotropic_diffusion_filter field input_image num_iterations 5 time_step 0.125 conductance 9.0;


# Calculate the gradient of the image using a sigma value of 1.0/255
gfx define field gradient gradient_magnitude_recursive_gaussian_filter field smoothed sigma 0.0039216;

# Caclulate the speed term using a sigmoid with alpha -0.3/255 and 2.0/255 
gfx define field sigma sigmoid_filter field gradient alpha -0.0011765 beta  0.0078431;

# Calculate a time-crossing map indicating how long it takes to reach
# pixels from the starting seed points, using the fast marching image filter
# The seed points are defined in terms of their position within the texture
# described by normalised coordinates.
# In terms of pixel value the seed points are (56,92) and (96,107)
# The image pixel dimensions are 181x217 which gives the normallised seed 
# points as (0.30939,0.423396) and (0.53,0.5)

gfx define field time_cross_map fast_marching_filter field sigma num_seed_points  2 seed_points 0.31 0.42 0.53 0.5 stopping_value 1000000
 
# Caclulate the segmented image by threshold the image.  
# Pixels with intensities between 0 (black) and 0.95 (almost white) 
# will be included in the regions
gfx define field fast_marching_segmentation binary_threshold_filter field time_cross_map lower_threshold 0 upper_threshold 0.95

# write out the output images for each stage of filtering to check that 
# the filter pipeline has worked

gfx modify texture output_image nearest width 1 height 1 evaluate_image field smoothed spectrum gray_spectrum texture_coord xi element_group square;
gfx write texture output_image file smoothed.png;

gfx modify texture output_image nearest width 1 height 1 evaluate_image field gradient spectrum gray_spectrum texture_coord xi element_group square;
gfx write texture output_image file gradient.png;

gfx modify texture output_image nearest width 1 height 1 evaluate_image field sigma spectrum gray_spectrum texture_coord xi element_group square;
gfx write texture output_image file sigma.png;

gfx modify texture output_image nearest width 1 height 1 evaluate_image field time_cross_map spectrum gray_spectrum texture_coord xi element_group square;
gfx write texture output_image file time_cross_map.png;

gfx modify texture output_image nearest width 1 height 1 evaluate_image field fast_marching_segmentation spectrum gray_spectrum texture_coord xi element_group square;
gfx write texture output_image file fast_marching_segmentation.png;


# generate node points around the segmented image

# Use a fine discretization with for generating the isolines.
gfx modify g_element square general element_discretization "512*512*512"

# Define the density field used to control the number of points generated 
# Here we want a constant number of points per segment length.
gfx define field point_density constant 50;

# Create the isolines
gfx modify g_element square iso_surfaces as isolines select_on iso_scalar fast_marching_segmentation iso_value 0.5 use_elements select_on material red selected_material default_selected render_shaded data point_density;


if (!$TESTING)
  {
	# set up a window to view the output of the fast marching segmentation

	gfx create window
	gfx modify window 1 set perturb_lines;
  gfx modify window 1 layout 2d
  gfx mod g_element square surfaces coordinate coordinates select_on material output_image texture_coordinates xi selected_material default_selected render_shaded

  }

# Generate points on the isolines calculated above.
gfx create ngroup data;
gfx define graphics_filter isolines match_graphic_name isolines;
gfx define scene isolines filter isolines;
gfx convert graphics scene isolines render_surface_node_cloud region data coordinate coordinates;

# Write out these points
gfx write nodes group data data_cloud.exnode;

if (!$TESTING)
  {
	 gfx modify g_element data node_points glyph sphere general size ".05*.05*.05" centre 0,0,0 font default select_on material green selected_material default_selected;
  }

Files used by this example are:

Name                         Modified     Size

segmentation.com             20-Apr-2012  9.6k
BrainProtonDensitySlice.png  20-Apr-2012   33k
COPYRIGHT                    19-Apr-2012  504 
square.exelem                20-Apr-2012  1.5k
square.exnode                20-Apr-2012  400

Download the entire example:

Name                                   Modified     Size

examples_a_backup_segmentation.tar.gz  12-Aug-2014  420k

Testing status by version:

	Status	Tested	Real time (s)
i686-linux
cmgui-wx	Success	Sun Mar 6 00:16:36 2016	4
cmgui-wx-debug	Success	Sun Mar 6 00:19:48 2016	5
cmgui-wx-debug-memorycheck	Success	Sun Mar 6 00:23:48 2016	7
cmgui-wx-debug-valgrind	Failure	Sun Mar 6 01:44:25 2016	189
	last break	Sun Mar 6 01:41:00 2016	189
	last success	Tue Feb 10 01:45:00 2015	170
x86_64-linux
cmgui-wx	Failure	Sun Mar 6 00:01:38 2016	1
	last break	Sun Mar 6 00:01:00 2016	1
	last success	Wed Jun 3 00:06:00 2015	2
cmgui-wx-debug	Failure	Sun Mar 6 00:01:38 2016	1
	last break	Sun Mar 6 00:01:00 2016	1
	last success	Wed Jun 3 00:06:00 2015	3
cmgui-wx-debug-memorycheck	Failure	Sun Mar 6 00:01:49 2016	0
	last break	Sun Mar 6 00:01:00 2016	0
	last success	Wed Jun 3 00:07:00 2015	4
cmgui-wx-debug-valgrind	Failure	Sun Mar 6 00:03:23 2016	9
	last break	Sun Mar 6 00:03:00 2016	9
	last success	Wed Jun 3 01:02:00 2015	131

Testing status by file:

status

x86_64-linux

139 cmgui-wx: error exit status 139.

139 cmgui-wx-debug: error exit status 139.

139 cmgui-wx-debug-memorycheck: error exit status 139.

139 cmgui-wx-debug-valgrind: error exit status 139.

cmiss_test.log

i686-linux
Success	cmgui-wx:	diff test: no differences with generic answer.
Success	cmgui-wx-debug:	diff test: no differences with generic answer.
Success	cmgui-wx-debug-memorycheck:	diff test: no differences with generic answer.
Failure	cmgui-wx-debug-valgrind:	diff test: differences with generic answer; Test output.
x86_64-linux
Failure	cmgui-wx:	diff test: differences with generic answer; Test output.
Failure	cmgui-wx-debug:	diff test: differences with generic answer; Test output.
Failure	cmgui-wx-debug-memorycheck:	diff test: differences with generic answer; Test output.
Failure	cmgui-wx-debug-valgrind:	diff test: differences with generic answer; Test output.

Html last generated: Sun Mar 6 05:51:00 2016

Input last modified: Fri Apr 20 16:01:17 2012

CMISS Help / Examples / a_backup / segmentation

x86_64-linux
139	cmgui-wx:	error exit status 139.
139	cmgui-wx-debug:	error exit status 139.
139	cmgui-wx-debug-memorycheck:	error exit status 139.
139	cmgui-wx-debug-valgrind:	error exit status 139.