Skip to content

Commit 3fbe489

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

File tree

11 files changed

+862
-0
lines changed

11 files changed

+862
-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 )
Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1 @@
1+
0d5403fda7b3c94d3229bf5e42b523a4a5bb2cd84a45973011520daf5b3001192e74ff3d290b43c0feb0ebd84c17a328f4646a7887408c1f63e52927676f2360
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()
Lines changed: 177 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,177 @@
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
41+
main(int argc, char * argv[])
42+
{
43+
if (argc < 10)
44+
{
45+
std::cerr << "Missing parameters." << std::endl;
46+
std::cerr << "Usage: " << argv[0] << " inputImageFile"
47+
<< " reversedInputImageFile"
48+
<< " distanceMapOutputImageFile"
49+
<< " watershedOutputFileName"
50+
<< " segmentationResultOutputImageFile"
51+
<< " binarizingRadius"
52+
<< " majorityThreshold"
53+
<< " watershedThreshold"
54+
<< " watershedLevel"
55+
<< " cleaningStructuringElementRadius" << std::endl;
56+
return EXIT_FAILURE;
57+
}
58+
59+
constexpr unsigned int Dimension = 3;
60+
61+
using UnsignedCharPixelType = unsigned char;
62+
using FloatPixelType = float;
63+
64+
using InputImageType = itk::Image<UnsignedCharPixelType, Dimension>;
65+
using FloatImageType = itk::Image<FloatPixelType, Dimension>;
66+
using RGBPixelType = itk::RGBPixel<UnsignedCharPixelType>;
67+
using RGBImageType = itk::Image<RGBPixelType, Dimension>;
68+
using LabeledImageType = itk::Image<itk::IdentifierType, Dimension>;
69+
70+
71+
using FileReaderType = itk::ImageFileReader<InputImageType>;
72+
FileReaderType::Pointer reader = FileReaderType::New();
73+
reader->SetFileName(argv[1]);
74+
reader->Update();
75+
76+
77+
// Create bubble image: get a binarized version of the input image
78+
using VotingBinaryIterativeHoleFillingImageFilterType =
79+
itk::VotingBinaryIterativeHoleFillingImageFilter<InputImageType>;
80+
VotingBinaryIterativeHoleFillingImageFilterType::Pointer votingBinaryHoleFillingImageFilter =
81+
VotingBinaryIterativeHoleFillingImageFilterType::New();
82+
votingBinaryHoleFillingImageFilter->SetInput(reader->GetOutput());
83+
84+
const unsigned int binarizingRadius = std::stoi(argv[6]);
85+
86+
InputImageType::SizeType indexRadius;
87+
indexRadius.Fill(binarizingRadius);
88+
89+
votingBinaryHoleFillingImageFilter->SetRadius(indexRadius);
90+
91+
votingBinaryHoleFillingImageFilter->SetBackgroundValue(0);
92+
votingBinaryHoleFillingImageFilter->SetForegroundValue(255);
93+
94+
const unsigned int majorityThreshold = std::stoi(argv[7]);
95+
votingBinaryHoleFillingImageFilter->SetMajorityThreshold(majorityThreshold);
96+
97+
votingBinaryHoleFillingImageFilter->Update();
98+
99+
using FileWriterType = itk::ImageFileWriter<InputImageType>;
100+
FileWriterType::Pointer reversedImageWriter = FileWriterType::New();
101+
reversedImageWriter->SetFileName(argv[2]);
102+
reversedImageWriter->SetInput(votingBinaryHoleFillingImageFilter->GetOutput());
103+
reversedImageWriter->Update();
104+
105+
106+
// Get the distance map of the input image
107+
using SignedMaurerDistanceMapImageFilterType =
108+
itk::SignedMaurerDistanceMapImageFilter<InputImageType, FloatImageType>;
109+
SignedMaurerDistanceMapImageFilterType::Pointer distanceMapImageFilter =
110+
SignedMaurerDistanceMapImageFilterType::New();
111+
distanceMapImageFilter->SetInput(votingBinaryHoleFillingImageFilter->GetOutput());
112+
113+
distanceMapImageFilter->SetInsideIsPositive(false);
114+
distanceMapImageFilter->Update();
115+
116+
117+
using DistanceMapFileWriterType = itk::ImageFileWriter<FloatImageType>;
118+
DistanceMapFileWriterType::Pointer distanceMapWriter = DistanceMapFileWriterType::New();
119+
distanceMapWriter->SetFileName(argv[3]);
120+
distanceMapWriter->SetInput(distanceMapImageFilter->GetOutput());
121+
distanceMapWriter->Update();
122+
123+
124+
// Apply the watershed segmentation
125+
using WatershedFilterType = itk::WatershedImageFilter<FloatImageType>;
126+
WatershedFilterType::Pointer watershed = WatershedFilterType::New();
127+
128+
const float watershedThreshold = std::stod(argv[8]);
129+
const float watershedLevel = std::stod(argv[9]);
130+
131+
watershed->SetThreshold(watershedThreshold);
132+
watershed->SetLevel(watershedLevel);
133+
134+
watershed->SetInput(distanceMapImageFilter->GetOutput());
135+
watershed->Update();
136+
137+
138+
using RGBFilterType = itk::ScalarToRGBColormapImageFilter<LabeledImageType, RGBImageType>;
139+
RGBFilterType::Pointer colormapImageFilter = RGBFilterType::New();
140+
colormapImageFilter->SetColormap(RGBFilterType::Jet);
141+
colormapImageFilter->SetInput(watershed->GetOutput());
142+
colormapImageFilter->Update();
143+
144+
using WatershedFileWriterType = itk::ImageFileWriter<RGBImageType>;
145+
WatershedFileWriterType::Pointer watershedWriter = WatershedFileWriterType::New();
146+
watershedWriter->SetFileName(argv[4]);
147+
watershedWriter->SetInput(colormapImageFilter->GetOutput());
148+
watershedWriter->Update();
149+
150+
151+
// Clean the segmentation image: remove small objects by performing an
152+
// opening morphological operation
153+
using StructuringElementType =
154+
itk::BinaryBallStructuringElement<LabeledImageType::PixelType, LabeledImageType::ImageDimension>;
155+
StructuringElementType structuringElement;
156+
157+
const unsigned int cleaningStructuringElementRadius = std::stoi(argv[10]);
158+
structuringElement.SetRadius(cleaningStructuringElementRadius);
159+
structuringElement.CreateStructuringElement();
160+
161+
using BinaryMorphologicalOpeningImageFilterType =
162+
itk::BinaryMorphologicalOpeningImageFilter<LabeledImageType, LabeledImageType, StructuringElementType>;
163+
BinaryMorphologicalOpeningImageFilterType::Pointer openingFilter = BinaryMorphologicalOpeningImageFilterType::New();
164+
openingFilter->SetInput(watershed->GetOutput());
165+
openingFilter->SetKernel(structuringElement);
166+
openingFilter->Update();
167+
168+
169+
using SegmentationFileWriterType = itk::ImageFileWriter<RGBImageType>;
170+
SegmentationFileWriterType::Pointer segmentationWriter = SegmentationFileWriterType::New();
171+
segmentationWriter->SetFileName(argv[5]);
172+
segmentationWriter->SetInput(colormapImageFilter->GetOutput());
173+
segmentationWriter->Update();
174+
175+
176+
return EXIT_SUCCESS;
177+
}
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)