Skip to content

Commit 0a865e0

Browse files
jhlegarretajhlegarreta
committed
WIP: ENH: Add watershed with distance map example.
Add watershed with distance map example. Resolves #48.
1 parent 2985782 commit 0a865e0

File tree

11 files changed

+856
-0
lines changed

11 files changed

+856
-0
lines changed

‎src/Segmentation/Watersheds/CMakeLists.txt

Lines changed: 8 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1,4 +1,12 @@
11

2+
add_example(SegmentWithWatershedAndDistanceMap)
3+
compare_to_baseline(EXAMPLE_NAME SegmentWithWatershedAndDistanceMap
4+
TEST_NAME SegmentWithWatershedAndDistanceMapTest01BaselineComparison
5+
BASELINE_PREFIX SegmentWithWatershedAndDistanceMapTest01Baseline
6+
DEPENDS SegmentWithWatershedAndDistanceMapTest01
7+
TEST_IMAGE SegmentWithWatershedAndDistanceMapTest01.png
8+
)
9+
210
add_example(SegmentWithWatershedImageFilter)
311

412
foreach( i RANGE 1 5 )

‎src/Segmentation/Watersheds/SegmentWithWatershedAndDistanceMap/BubbleVolume.csv.sha512

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1 @@
1+
0d5403fda7b3c94d3229bf5e42b523a4a5bb2cd84a45973011520daf5b3001192e74ff3d290b43c0feb0ebd84c17a328f4646a7887408c1f63e52927676f2360

‎src/Segmentation/Watersheds/SegmentWithWatershedAndDistanceMap/CMakeLists.txt

Lines changed: 59 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,59 @@
1+
cmake_minimum_required(VERSION 3.10.2)
2+
3+
project(SegmentWithWatershedAndDistanceMap)
4+
5+
find_package(ITK REQUIRED)
6+
include(${ITK_USE_FILE})
7+
8+
9+
add_executable(SegmentWithWatershedAndDistanceMap Code.cxx)
10+
target_link_libraries(SegmentWithWatershedAndDistanceMap ${ITK_LIBRARIES})
11+
12+
install(TARGETS SegmentWithWatershedAndDistanceMap
13+
DESTINATION bin/ITKExamples/Segmentation/Watersheds
14+
COMPONENT Runtime
15+
)
16+
17+
install(FILES Code.cxx CMakeLists.txt Code.py
18+
DESTINATION share/ITKExamples/Code/Segmentation/Watersheds/SegmentWithWatershedAndDistanceMap/
19+
COMPONENT Code
20+
)
21+
22+
23+
enable_testing()
24+
set(input_image ${CMAKE_CURRENT_BINARY_DIR}/PlateauBorder.tif)
25+
set(binarizingRadius01 2)
26+
set(majorityThreshold01 2)
27+
set(watershedThreshold01 0.01)
28+
set(watershedLevel01 0.5)
29+
set(cleaningStructuringElementRadius01 3)
30+
31+
add_test(NAME SegmentWithWatershedAndDistanceMapTest01
32+
COMMAND ${CMAKE_RUNTIME_OUTPUT_DIRECTORY}/SegmentWithWatershedAndDistanceMap
33+
${input_image}
34+
ReversedInputImageTest01.tif
35+
DistanceMapImageTest01.tif
36+
WatershedImageTest01.tif
37+
SegmentWithWatershedAndDistanceMapTest01.tif
38+
${binarizingRadius01}
39+
${majorityThreshold01}
40+
${watershedThreshold01}
41+
${watershedLevel01}
42+
${cleaningStructuringElementRadius01}
43+
)
44+
45+
if(ITK_WRAP_PYTHON)
46+
add_test(NAME SegmentWithWatershedAndDistanceMapTest01Python
47+
COMMAND ${PYTHON_EXECUTABLE} ${CMAKE_CURRENT_SOURCE_DIR}/Code.py
48+
${input_image}
49+
ReversedInputImageTest01.tif
50+
DistanceMapImageTest01.tif
51+
WatershedImageTest01.tif
52+
SegmentWithWatershedAndDistanceMapTest01Python.tif
53+
${binarizingRadius01}
54+
${majorityThreshold01}
55+
${watershedThreshold01}
56+
${watershedLevel01}
57+
${cleaningStructuringElementRadius01}
58+
)
59+
endif()

‎src/Segmentation/Watersheds/SegmentWithWatershedAndDistanceMap/Code.cxx

Lines changed: 179 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,179 @@
1+
/*=========================================================================
2+
*
3+
* Copyright NumFOCUS
4+
*
5+
* Licensed under the Apache License, Version 2.0 (the "License");
6+
* you may not use this file except in compliance with the License.
7+
* You may obtain a copy of the License at
8+
*
9+
* http://www.apache.org/licenses/LICENSE-2.0.txt
10+
*
11+
* Unless required by applicable law or agreed to in writing, software
12+
* distributed under the License is distributed on an "AS IS" BASIS,
13+
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14+
* See the License for the specific language governing permissions and
15+
* limitations under the License.
16+
*
17+
*=========================================================================*/
18+
19+
#include "itkBinaryBallStructuringElement.h"
20+
#include "itkBinaryMorphologicalOpeningImageFilter.h"
21+
#include "itkImageFileReader.h"
22+
#include "itkImageFileWriter.h"
23+
#include "itkScalarToRGBColormapImageFilter.h"
24+
#include "itkSignedMaurerDistanceMapImageFilter.h"
25+
#include "itkVotingBinaryIterativeHoleFillingImageFilter.h"
26+
#include "itkWatershedImageFilter.h"
27+
28+
// Run with:
29+
// ./SegmentWithWatershedAndDistanceMap <inputImageFile>
30+
// <reversedInputImageFile> <distanceMapOutputImageFile>
31+
// <watershedOutputFileName> <segmentationResultOutputImageFile>
32+
// binarizingRadius majorityThreshold watershedThreshold watershedLevel
33+
// cleaningStructuringElementRadius
34+
// e.g.
35+
// ./SegmentWithWatershedAndDistanceMap PlateauBorder.tif
36+
// reversedInputImage.tif distanceMap.tif watershed.tif segmentationResult.tif
37+
// 2 2 0.01 0.5 3
38+
// (A rule of thumb is to set the threshold to be about 1 / 100 of the level.)
39+
40+
int main( int argc, char *argv[] )
41+
{
42+
if( argc < 10 )
43+
{
44+
std::cerr << "Missing parameters." << std::endl;
45+
std::cerr << "Usage: " << argv[0]
46+
<< " inputImageFile"
47+
<< " reversedInputImageFile"
48+
<< " distanceMapOutputImageFile"
49+
<< " watershedOutputFileName"
50+
<< " segmentationResultOutputImageFile"
51+
<< " binarizingRadius"
52+
<< " majorityThreshold"
53+
<< " watershedThreshold"
54+
<< " watershedLevel"
55+
<< " cleaningStructuringElementRadius"
56+
<< std::endl;
57+
return EXIT_FAILURE;
58+
}
59+
60+
constexpr unsigned int Dimension = 3;
61+
62+
using UnsignedCharPixelType = unsigned char;
63+
using FloatPixelType = float;
64+
65+
using InputImageType = itk::Image< UnsignedCharPixelType, Dimension >;
66+
using FloatImageType = itk::Image< FloatPixelType, Dimension >;
67+
using RGBPixelType = itk::RGBPixel< UnsignedCharPixelType >;
68+
using RGBImageType = itk::Image< RGBPixelType, Dimension >;
69+
using LabeledImageType = itk::Image< itk::IdentifierType, Dimension >;
70+
71+
72+
using FileReaderType = itk::ImageFileReader< InputImageType >;
73+
FileReaderType::Pointer reader = FileReaderType::New();
74+
reader->SetFileName( argv[1] );
75+
reader->Update();
76+
77+
78+
// Create bubble image: get a binarized version of the input image
79+
using VotingBinaryIterativeHoleFillingImageFilterType =
80+
itk::VotingBinaryIterativeHoleFillingImageFilter< InputImageType >;
81+
VotingBinaryIterativeHoleFillingImageFilterType::Pointer votingBinaryHoleFillingImageFilter =
82+
VotingBinaryIterativeHoleFillingImageFilterType::New();
83+
votingBinaryHoleFillingImageFilter->SetInput( reader->GetOutput() );
84+
85+
const unsigned int binarizingRadius = std::stoi( argv[6] );
86+
87+
InputImageType::SizeType indexRadius;
88+
indexRadius.Fill( binarizingRadius );
89+
90+
votingBinaryHoleFillingImageFilter->SetRadius( indexRadius );
91+
92+
votingBinaryHoleFillingImageFilter->SetBackgroundValue( 0 );
93+
votingBinaryHoleFillingImageFilter->SetForegroundValue( 255 );
94+
95+
const unsigned int majorityThreshold = std::stoi( argv[7] );
96+
votingBinaryHoleFillingImageFilter->SetMajorityThreshold( majorityThreshold );
97+
98+
votingBinaryHoleFillingImageFilter->Update();
99+
100+
using FileWriterType = itk::ImageFileWriter< InputImageType >;
101+
FileWriterType::Pointer reversedImageWriter = FileWriterType::New();
102+
reversedImageWriter->SetFileName( argv[2] );
103+
reversedImageWriter->SetInput( votingBinaryHoleFillingImageFilter->GetOutput() );
104+
reversedImageWriter->Update();
105+
106+
107+
// Get the distance map of the input image
108+
using SignedMaurerDistanceMapImageFilterType =
109+
itk::SignedMaurerDistanceMapImageFilter< InputImageType, FloatImageType >;
110+
SignedMaurerDistanceMapImageFilterType::Pointer distanceMapImageFilter =
111+
SignedMaurerDistanceMapImageFilterType::New();
112+
distanceMapImageFilter->SetInput( votingBinaryHoleFillingImageFilter->GetOutput() );
113+
114+
distanceMapImageFilter->SetInsideIsPositive( false );
115+
distanceMapImageFilter->Update();
116+
117+
118+
using DistanceMapFileWriterType = itk::ImageFileWriter< FloatImageType >;
119+
DistanceMapFileWriterType::Pointer distanceMapWriter = DistanceMapFileWriterType::New();
120+
distanceMapWriter->SetFileName( argv[3] );
121+
distanceMapWriter->SetInput( distanceMapImageFilter->GetOutput() );
122+
distanceMapWriter->Update();
123+
124+
125+
// Apply the watershed segmentation
126+
using WatershedFilterType = itk::WatershedImageFilter< FloatImageType >;
127+
WatershedFilterType::Pointer watershed = WatershedFilterType::New();
128+
129+
const float watershedThreshold = std::stod( argv[8] );
130+
const float watershedLevel = std::stod( argv[9] );
131+
132+
watershed->SetThreshold( watershedThreshold );
133+
watershed->SetLevel( watershedLevel );
134+
135+
watershed->SetInput( distanceMapImageFilter->GetOutput() );
136+
watershed->Update();
137+
138+
139+
using RGBFilterType = itk::ScalarToRGBColormapImageFilter< LabeledImageType, RGBImageType>;
140+
RGBFilterType::Pointer colormapImageFilter = RGBFilterType::New();
141+
colormapImageFilter->SetColormap( RGBFilterType::Jet );
142+
colormapImageFilter->SetInput( watershed->GetOutput() );
143+
colormapImageFilter->Update();
144+
145+
using WatershedFileWriterType = itk::ImageFileWriter< RGBImageType >;
146+
WatershedFileWriterType::Pointer watershedWriter = WatershedFileWriterType::New();
147+
watershedWriter->SetFileName( argv[4] );
148+
watershedWriter->SetInput( colormapImageFilter->GetOutput() );
149+
watershedWriter->Update();
150+
151+
152+
// Clean the segmentation image: remove small objects by performing an
153+
// opening morphological operation
154+
using StructuringElementType =
155+
itk::BinaryBallStructuringElement< LabeledImageType::PixelType, LabeledImageType::ImageDimension >;
156+
StructuringElementType structuringElement;
157+
158+
const unsigned int cleaningStructuringElementRadius = std::stoi( argv[10] );
159+
structuringElement.SetRadius( cleaningStructuringElementRadius );
160+
structuringElement.CreateStructuringElement();
161+
162+
using BinaryMorphologicalOpeningImageFilterType =
163+
itk::BinaryMorphologicalOpeningImageFilter< LabeledImageType, LabeledImageType, StructuringElementType >;
164+
BinaryMorphologicalOpeningImageFilterType::Pointer openingFilter =
165+
BinaryMorphologicalOpeningImageFilterType::New();
166+
openingFilter->SetInput( watershed->GetOutput() );
167+
openingFilter->SetKernel( structuringElement );
168+
openingFilter->Update();
169+
170+
171+
using SegmentationFileWriterType = itk::ImageFileWriter< RGBImageType >;
172+
SegmentationFileWriterType::Pointer segmentationWriter = SegmentationFileWriterType::New();
173+
segmentationWriter->SetFileName( argv[5] );
174+
segmentationWriter->SetInput( colormapImageFilter->GetOutput() );
175+
segmentationWriter->Update();
176+
177+
178+
return EXIT_SUCCESS;
179+
}

‎src/Segmentation/Watersheds/SegmentWithWatershedAndDistanceMap/Code.ipynb

Lines changed: 109 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,109 @@
1+
{
2+
"cells": [
3+
{
4+
"cell_type": "code",
5+
"execution_count": 15,
6+
"id": "9e8cfdc6",
7+
"metadata": {},
8+
"outputs": [],
9+
"source": [
10+
"import sys\n",
11+
"\n",
12+
"import itk\n",
13+
"import numpy as np\n",
14+
"\n",
15+
"from itkwidgets import view, cm\n",
16+
"import matplotlib.pyplot as plt"
17+
]
18+
},
19+
{
20+
"cell_type": "code",
21+
"execution_count": 17,
22+
"id": "057c9f2f",
23+
"metadata": {},
24+
"outputs": [
25+
{
26+
"data": {
27+
"application/vnd.jupyter.widget-view+json": {
28+
"model_id": "962ea74edef74bc4b0b026a83f3dde3a",
29+
"version_major": 2,
30+
"version_minor": 0
31+
},
32+
"text/plain": [
33+
"Viewer(cmap=['bone_Matlab'], geometries=[], gradient_opacity=0.22, point_sets=[], rendered_image=<itk.itkImage…"
34+
]
35+
},
36+
"metadata": {},
37+
"output_type": "display_data"
38+
}
39+
],
40+
"source": [
41+
"stack_image = itk.imread('./PlateauBorder.tif')\n",
42+
"view(stack_image, cmap=[cm.bone])"
43+
]
44+
},
45+
{
46+
"cell_type": "code",
47+
"execution_count": 20,
48+
"id": "ab956229",
49+
"metadata": {},
50+
"outputs": [
51+
{
52+
"data": {
53+
"application/vnd.jupyter.widget-view+json": {
54+
"model_id": "c0fc5e027fe149d4b4bc866b054c6e38",
55+
"version_major": 2,
56+
"version_minor": 0
57+
},
58+
"text/plain": [
59+
"Viewer(cmap=['bone_Matlab'], geometries=[], gradient_opacity=0.22, point_sets=[], rendered_image=<itk.itkImage…"
60+
]
61+
},
62+
"metadata": {},
63+
"output_type": "display_data"
64+
}
65+
],
66+
"source": [
67+
"binarizing_radius = 2\n",
68+
"majority_threshold = 2\n",
69+
"\n",
70+
"bubble_image = itk.voting_binary_hole_filling_image_filter(stack_image,\n",
71+
" radius=binarizing_radius,\n",
72+
" majority_threshold=majority_threshold,\n",
73+
" background_value=0,\n",
74+
" foreground_value=255)\n",
75+
"bubble_image = itk.invert_intensity_image_filter(bubble_image)\n",
76+
"view(bubble_image, cmap=[cm.bone,])"
77+
]
78+
},
79+
{
80+
"cell_type": "code",
81+
"execution_count": null,
82+
"id": "bfc02551",
83+
"metadata": {},
84+
"outputs": [],
85+
"source": []
86+
}
87+
],
88+
"metadata": {
89+
"kernelspec": {
90+
"display_name": "Python 3",
91+
"language": "python",
92+
"name": "python3"
93+
},
94+
"language_info": {
95+
"codemirror_mode": {
96+
"name": "ipython",
97+
"version": 3
98+
},
99+
"file_extension": ".py",
100+
"mimetype": "text/x-python",
101+
"name": "python",
102+
"nbconvert_exporter": "python",
103+
"pygments_lexer": "ipython3",
104+
"version": "3.8.6"
105+
}
106+
},
107+
"nbformat": 4,
108+
"nbformat_minor": 5
109+
}

0 commit comments

Comments
 (0)