1 /*=========================================================================
2 *
3 * Copyright Insight Software Consortium
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 "itkFastMarchingImageFilter.h"
20 #include "itkImageRegionIteratorWithIndex.h"
21 #include "itkTextOutput.h"
22 #include "itkCommand.h"
23
24
25 namespace{
26 // The following class is used to support callbacks
27 // on the filter in the pipeline that follows later
28 class ShowProgressObject
29 {
30 public:
ShowProgressObject(itk::ProcessObject * o)31 ShowProgressObject(itk::ProcessObject* o)
32 {m_Process = o;}
ShowProgress()33 void ShowProgress()
34 {std::cout << "Progress " << m_Process->GetProgress() << std::endl;}
35 itk::ProcessObject::Pointer m_Process;
36 };
37 }
38
itkFastMarchingTest2(int,char * [])39 int itkFastMarchingTest2(int, char* [] )
40 {
41
42 itk::OutputWindow::SetInstance(itk::TextOutput::New().GetPointer());
43
44 // create a fastmarching object
45 using PixelType = float;
46 using FloatImage = itk::Image<PixelType,2>;
47 using FloatFMType = itk::FastMarchingImageFilter<FloatImage,FloatImage>;
48
49 FloatFMType::Pointer marcher = FloatFMType::New();
50
51 ShowProgressObject progressWatch(marcher);
52 itk::SimpleMemberCommand<ShowProgressObject>::Pointer command;
53 command = itk::SimpleMemberCommand<ShowProgressObject>::New();
54 command->SetCallbackFunction(&progressWatch,
55 &ShowProgressObject::ShowProgress);
56 marcher->AddObserver( itk::ProgressEvent(), command);
57
58 using NodeType = FloatFMType::NodeType;
59 using NodeContainer = FloatFMType::NodeContainer;
60
61 // setup alive points
62 NodeContainer::Pointer alivePoints = NodeContainer::New();
63
64 NodeType node;
65
66 FloatImage::OffsetType offset0 = {{28,35}};
67
68 itk::Index<2> index;
69 index.Fill(0);
70
71 node.SetValue( 0.0 );
72 node.SetIndex( index + offset0 );
73 alivePoints->InsertElement(0, node);
74
75 node.SetValue( 42.0 );
76 index.Fill( 200 );
77 node.SetIndex( index ); // this node is out of range
78 alivePoints->InsertElement(1, node);
79
80 marcher->SetAlivePoints( alivePoints );
81
82
83 // setup trial points
84 NodeContainer::Pointer trialPoints = NodeContainer::New();
85
86 node.SetValue( 1.0 );
87
88 index.Fill(0);
89 index += offset0;
90
91 index[0] += 1;
92 node.SetIndex( index );
93 trialPoints->InsertElement(0, node);
94
95 index[0] -= 1;
96 index[1] += 1;
97 node.SetIndex( index );
98 trialPoints->InsertElement(1, node);
99
100 index[0] -= 1;
101 index[1] -= 1;
102 node.SetIndex( index );
103 trialPoints->InsertElement(2, node);
104
105 index[0] += 1;
106 index[1] -= 1;
107 node.SetIndex( index );
108 trialPoints->InsertElement(3, node);
109
110 node.SetValue( 42.0 );
111 index.Fill( 300 ); // this node is out of range
112 node.SetIndex( index );
113 trialPoints->InsertElement(4, node);
114
115 marcher->SetTrialPoints( trialPoints );
116
117 // specify the size of the output image
118 FloatImage::SizeType size = {{64,64}};
119 marcher->SetOutputSize( size );
120
121 // setup a speed image of ones
122 FloatImage::Pointer speedImage = FloatImage::New();
123 FloatImage::RegionType region;
124 region.SetSize( size );
125 speedImage->SetLargestPossibleRegion( region );
126 speedImage->SetBufferedRegion( region );
127 speedImage->Allocate();
128
129 // setup a binary mask image in float (to make sure it works with float)
130 FloatImage::Pointer MaskImage = FloatImage::New();
131 MaskImage->SetLargestPossibleRegion( region );
132 MaskImage->SetBufferedRegion( region );
133 MaskImage->Allocate();
134
135 itk::ImageRegionIterator<FloatImage>
136 speedIter( speedImage, speedImage->GetBufferedRegion() );
137 itk::ImageRegionIteratorWithIndex<FloatImage>
138 maskIter( MaskImage, MaskImage->GetBufferedRegion() );
139 while ( !speedIter.IsAtEnd() )
140 {
141 speedIter.Set( 1.0 );
142 FloatImage::IndexType idx = maskIter.GetIndex();
143 if( ( ( idx[0] > 22 ) && ( idx [0] < 42 ) && ( idx[1] > 27 ) && ( idx[1] < 37 ) ) ||
144 ( ( idx[1] > 22 ) && ( idx [1] < 42 ) && ( idx[0] > 27 ) && ( idx[0] < 37 ) ) )
145 {
146 maskIter.Set( 1.0 );
147 }
148 else
149 {
150 maskIter.Set( 0.0 );
151 }
152
153 ++maskIter;
154 ++speedIter;
155 }
156
157 speedImage->Print( std::cout );
158 marcher->SetInput( speedImage );
159 marcher->SetBinaryMask( MaskImage.GetPointer() );
160 marcher->SetStoppingValue( 100.0 );
161
162 // turn on debugging
163 marcher->DebugOn();
164
165
166 // update the marcher
167 marcher->Update();
168
169
170 // check the results
171 FloatImage::Pointer output = marcher->GetOutput();
172 itk::ImageRegionIterator<FloatImage>
173 iterator( output, output->GetBufferedRegion() );
174
175 bool passed = true;
176 for(; !iterator.IsAtEnd(); ++iterator )
177 {
178 FloatImage::IndexType tempIndex = iterator.GetIndex();
179 auto outputValue = (float) iterator.Get();
180
181 if( ( ( tempIndex[0] > 22 ) && ( tempIndex [0] < 42 ) && ( tempIndex[1] > 27 ) && ( tempIndex[1] < 37 ) ) ||
182 ( ( tempIndex[1] > 22 ) && ( tempIndex [1] < 42 ) && ( tempIndex[0] > 27 ) && ( tempIndex[0] < 37 ) ) )
183 {
184 tempIndex -= offset0;
185 double distance = 0.0;
186 for ( int j = 0; j < 2; j++ )
187 {
188 distance += tempIndex[j] * tempIndex[j];
189 }
190 distance = std::sqrt( distance );
191
192 if (distance < itk::NumericTraits<double>::epsilon() )
193 {
194 continue;
195 }
196
197 if ( itk::Math::abs( outputValue ) / distance > 1.42 )
198 {
199 std::cout << iterator.GetIndex() << " ";
200 std::cout << itk::Math::abs( outputValue ) / distance << " ";
201 std::cout << itk::Math::abs( outputValue ) << " " << distance << std::endl;
202 passed = false;
203 }
204 }
205 else
206 {
207 if( outputValue != 0. )
208 {
209 std::cout << iterator.GetIndex() << " ";
210 std::cout << outputValue << " " << 0.;
211 std::cout << std::endl;
212 passed = false;
213 }
214 }
215 }
216
217 // Exercise other member functions
218 std::cout << "SpeedConstant: " << marcher->GetSpeedConstant() << std::endl;
219 std::cout << "StoppingValue: " << marcher->GetStoppingValue() << std::endl;
220 std::cout << "CollectPoints: " << marcher->GetCollectPoints() << std::endl;
221
222 marcher->SetNormalizationFactor( 2.0 );
223 std::cout << "NormalizationFactor: " << marcher->GetNormalizationFactor();
224 std::cout << std::endl;
225
226 std::cout << "SpeedImage: " << marcher->GetInput();
227 std::cout << std::endl;
228
229 marcher->Print( std::cout );
230
231 if ( passed )
232 {
233 std::cout << "Fast Marching test passed" << std::endl;
234 return EXIT_SUCCESS;
235 }
236 else
237 {
238 std::cout << "Fast Marching test failed" << std::endl;
239 return EXIT_FAILURE;
240 }
241
242 }
243