ITK  5.4.0
Insight Toolkit
itkSLICImageFilter.h
Go to the documentation of this file.
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  * https://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 #ifndef itkSLICImageFilter_h
19 #define itkSLICImageFilter_h
20 
21 #include "itkImageToImageFilter.h"
22 #include <mutex>
23 
24 namespace itk
25 {
26 
60 template <typename TInputImage, typename TOutputImage, typename TDistancePixel = float>
61 class ITK_TEMPLATE_EXPORT SLICImageFilter : public ImageToImageFilter<TInputImage, TOutputImage>
62 {
63 public:
64  ITK_DISALLOW_COPY_AND_MOVE(SLICImageFilter);
65 
71 
73  itkNewMacro(Self);
74 
76  itkOverrideGetNameOfClassMacro(SLICImageFilter);
77 
79  static constexpr unsigned int ImageDimension = TInputImage::ImageDimension;
80 
81 
83  using InputImageType = TInputImage;
84  using InputPixelType = typename InputImageType::PixelType;
85  using OutputImageType = TOutputImage;
86  using OutputPixelType = typename OutputImageType::PixelType;
87  using DistanceType = TDistancePixel;
89 
93 
94  using ClusterComponentType = double;
95  using ClusterType = vnl_vector_ref<ClusterComponentType>;
96 
98 
100 
107  itkSetMacro(SpatialProximityWeight, double);
108  itkGetConstMacro(SpatialProximityWeight, double);
109 
114  itkSetMacro(MaximumNumberOfIterations, unsigned int);
115  itkGetConstMacro(MaximumNumberOfIterations, unsigned int);
116 
124  itkSetMacro(SuperGridSize, SuperGridSizeType);
125  itkGetConstMacro(SuperGridSize, SuperGridSizeType);
126  void
127  SetSuperGridSize(unsigned int factor);
128  void
129  SetSuperGridSize(unsigned int i, unsigned int factor);
130 
138  itkSetMacro(InitializationPerturbation, bool);
139  itkGetMacro(InitializationPerturbation, bool);
140  itkBooleanMacro(InitializationPerturbation);
141 
142 
150  itkSetMacro(EnforceConnectivity, bool);
151  itkGetMacro(EnforceConnectivity, bool);
152  itkBooleanMacro(EnforceConnectivity);
153 
154 
161  itkGetConstMacro(AverageResidual, double);
162 
163 protected:
164  SLICImageFilter();
165  ~SLICImageFilter() override = default;
166 
167  void
168  PrintSelf(std::ostream & os, Indent indent) const override;
169 
171  void
172  EnlargeOutputRequestedRegion(DataObject * output) override;
173 
174  void
175  BeforeThreadedGenerateData() override;
176 
177  void
178  ThreadedUpdateDistanceAndLabel(const OutputImageRegionType & outputRegionForThread);
179 
180  void
181  ThreadedUpdateClusters(const OutputImageRegionType & updateRegionForThread);
182 
183  void
184  ThreadedPerturbClusters(SizeValueType clusterIndex);
185 
186  void
187  ThreadedConnectivity(SizeValueType clusterIndex);
188 
189  void
190  SingleThreadedConnectivity();
191 
192  void
193  GenerateData() override;
194 
195 
196  void
197  AfterThreadedGenerateData() override;
198 
200  Distance(const ClusterType & cluster1, const ClusterType & cluster2);
201 
203  Distance(const ClusterType & cluster, const InputPixelType & _v, const PointType & pt);
204 
205 private:
206  SuperGridSizeType m_SuperGridSize{};
207  unsigned int m_MaximumNumberOfIterations{};
208  double m_SpatialProximityWeight{ 10.0 };
209 
211  std::vector<ClusterComponentType> m_Clusters{};
212  std::vector<ClusterComponentType> m_OldClusters{};
213 
214 
215  void
216  RelabelConnectedRegion(const IndexType & seed,
217  OutputPixelType requiredLabel,
218  OutputPixelType outputLabel,
219  std::vector<IndexType> & indexStack);
220 
222  {
223  size_t count;
224  vnl_vector<ClusterComponentType> cluster;
225  };
226 
227  using UpdateClusterMap = std::map<size_t, UpdateCluster>;
228 
230 
231  std::vector<UpdateClusterMap> m_UpdateClusterPerThread{};
232 
233  typename DistanceImageType::Pointer m_DistanceImage{};
234  typename MarkerImageType::Pointer m_MarkerImage{};
235 
236  bool m_EnforceConnectivity{ true };
237 
238  bool m_InitializationPerturbation{ true };
239 
240  double m_AverageResidual{};
241  std::mutex m_Mutex{};
242 };
243 } // end namespace itk
244 
245 #ifndef ITK_MANUAL_INSTANTIATION
246 # include "itkSLICImageFilter.hxx"
247 #endif
248 
249 #endif // itkSLICImageFilter_h
itk::SLICImageFilter::UpdateCluster
Definition: itkSLICImageFilter.h:221
itk::SLICImageFilter::UpdateCluster::count
vcl_size_t count
Definition: itkSLICImageFilter.h:223
itk::GTest::TypedefsAndConstructors::Dimension2::PointType
ImageBaseType::PointType PointType
Definition: itkGTestTypedefsAndConstructors.h:51
itk::SmartPointer< Self >
itk::Indent
Control indentation during Print() invocation.
Definition: itkIndent.h:49
itk::GTest::TypedefsAndConstructors::Dimension2::IndexType
ImageBaseType::IndexType IndexType
Definition: itkGTestTypedefsAndConstructors.h:50
itk::ImageToImageFilter
Base class for filters that take an image as input and produce an image as output.
Definition: itkImageToImageFilter.h:108
itk::SLICImageFilter::ClusterType
vnl_vector_ref< ClusterComponentType > ClusterType
Definition: itkSLICImageFilter.h:95
itk::ImageSource
Base class for all process objects that output image data.
Definition: itkImageSource.h:67
itk::SLICImageFilter::DistanceType
TDistancePixel DistanceType
Definition: itkSLICImageFilter.h:87
itk::SLICImageFilter::InputPixelType
typename InputImageType::PixelType InputPixelType
Definition: itkSLICImageFilter.h:84
itk::SLICImageFilter::IndexType
typename InputImageType::IndexType IndexType
Definition: itkSLICImageFilter.h:90
itk::GTest::TypedefsAndConstructors::Dimension2::RegionType
ImageBaseType::RegionType RegionType
Definition: itkGTestTypedefsAndConstructors.h:54
itk::SLICImageFilter::OutputPixelType
typename OutputImageType::PixelType OutputPixelType
Definition: itkSLICImageFilter.h:86
itk::SLICImageFilter
Simple Linear Iterative Clustering (SLIC) super-pixel segmentation.
Definition: itkSLICImageFilter.h:61
itk::SLICImageFilter::PointType
typename InputImageType::PointType PointType
Definition: itkSLICImageFilter.h:91
itk::ImageToImageFilter::InputImageType
TInputImage InputImageType
Definition: itkImageToImageFilter.h:129
itkImageToImageFilter.h
itk::FixedArray< unsigned int, ImageDimension >
itk::SLICImageFilter::UpdateClusterMap
std::map< vcl_size_t, UpdateCluster > UpdateClusterMap
Definition: itkSLICImageFilter.h:227
itk::ImageSource::OutputImageRegionType
typename OutputImageType::RegionType OutputImageRegionType
Definition: itkImageSource.h:92
itk
The "itk" namespace contains all Insight Segmentation and Registration Toolkit (ITK) classes....
Definition: itkAnnulusOperator.h:24
itk::ContinuousIndex
A templated class holding a point in n-Dimensional image space.
Definition: itkContinuousIndex.h:46
itk::ProcessObject
The base class for all process objects (source, filters, mappers) in the Insight data processing pipe...
Definition: itkProcessObject.h:139
itk::SLICImageFilter::UpdateCluster::cluster
vnl_vector< ClusterComponentType > cluster
Definition: itkSLICImageFilter.h:224
itk::Image
Templated n-dimensional image class.
Definition: itkImage.h:88
itk::SLICImageFilter::ClusterComponentType
double ClusterComponentType
Definition: itkSLICImageFilter.h:94
itk::SizeValueType
unsigned long SizeValueType
Definition: itkIntTypes.h:83
itk::ImageSource::OutputImageType
TOutputImage OutputImageType
Definition: itkImageSource.h:90
itk::DataObject
Base class for all data objects in ITK.
Definition: itkDataObject.h:293