ITK  4.9.0
Insight Segmentation and Registration Toolkit
itkConnectedComponentImageFilter.h
Go to the documentation of this file.
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 #ifndef itkConnectedComponentImageFilter_h
19 #define itkConnectedComponentImageFilter_h
20 
21 #include "itkImageToImageFilter.h"
22 #include "itkImage.h"
23 #include <vector>
24 #include <map>
25 #include "itkProgressReporter.h"
26 #include "itkBarrier.h"
27 
28 namespace itk
29 {
58 template< typename TInputImage, typename TOutputImage, typename TMaskImage = TInputImage >
60  public ImageToImageFilter< TInputImage, TOutputImage >
61 {
62 public:
68 
73 
78  typedef typename TOutputImage::PixelType OutputPixelType;
79  typedef typename TOutputImage::InternalPixelType OutputInternalPixelType;
80  typedef typename TInputImage::PixelType InputPixelType;
81  typedef typename TInputImage::InternalPixelType InputInternalPixelType;
82  typedef typename TMaskImage::PixelType MaskPixelType;
83  itkStaticConstMacro(ImageDimension, unsigned int,
84  TOutputImage::ImageDimension);
85  itkStaticConstMacro(OutputImageDimension, unsigned int,
86  TOutputImage::ImageDimension);
87  itkStaticConstMacro(InputImageDimension, unsigned int,
88  TInputImage::ImageDimension);
90 
94  typedef TInputImage InputImageType;
95  typedef TMaskImage MaskImageType;
96  typedef typename TInputImage::IndexType IndexType;
97  typedef typename TInputImage::SizeType SizeType;
98  typedef typename TInputImage::OffsetType OffsetType;
99 
100  typedef TOutputImage OutputImageType;
101  typedef typename TOutputImage::RegionType RegionType;
102  typedef typename TOutputImage::IndexType OutputIndexType;
103  typedef typename TOutputImage::SizeType OutputSizeType;
104  typedef typename TOutputImage::OffsetType OutputOffsetType;
105  typedef typename TOutputImage::PixelType OutputImagePixelType;
106 
107  typedef std::list< IndexType > ListType;
108  typedef typename MaskImageType::Pointer MaskImagePointer;
109 
115 
120 
124  itkNewMacro(Self);
125 
132  itkSetMacro(FullyConnected, bool);
133  itkGetConstReferenceMacro(FullyConnected, bool);
134  itkBooleanMacro(FullyConnected);
136 
138  typedef IdentifierType LabelType;
139 
140  // only set after completion
141  itkGetConstReferenceMacro(ObjectCount, LabelType);
142 
143  // Concept checking -- input and output dimensions must be the same
144  itkConceptMacro( SameDimension,
145  ( Concept::SameDimension< itkGetStaticConstMacro(InputImageDimension),
146  itkGetStaticConstMacro(OutputImageDimension) > ) );
147  itkConceptMacro( OutputImagePixelTypeIsInteger, ( Concept::IsInteger< OutputImagePixelType > ) );
148 
149  void SetMaskImage(TMaskImage *mask)
150  {
151  this->SetNthInput( 1, const_cast< TMaskImage * >( mask ) );
152  }
153 
154  const TMaskImage * GetMaskImage() const
155  {
156  return ( static_cast< const TMaskImage * >( this->ProcessObject::GetInput(1) ) );
157  }
158 
164  itkSetMacro(BackgroundValue, OutputImagePixelType);
165  itkGetConstMacro(BackgroundValue, OutputImagePixelType);
167 
168 protected:
170  {
171  m_FullyConnected = false;
172  m_ObjectCount = 0;
174  }
175 
177  void PrintSelf(std::ostream & os, Indent indent) const ITK_OVERRIDE;
178 
182  void BeforeThreadedGenerateData() ITK_OVERRIDE;
183 
184  void AfterThreadedGenerateData() ITK_OVERRIDE;
185 
186  void ThreadedGenerateData(const RegionType & outputRegionForThread, ThreadIdType threadId) ITK_OVERRIDE;
187 
191  void GenerateInputRequestedRegion() ITK_OVERRIDE;
192 
197  void EnlargeOutputRequestedRegion( DataObject * itkNotUsed(output) ) ITK_OVERRIDE;
198 
200 
201 private:
202  ConnectedComponentImageFilter(const Self &) ITK_DELETE_FUNCTION;
203  void operator=(const Self &) ITK_DELETE_FUNCTION;
204 
207 
208  // some additional types
209  typedef typename TOutputImage::RegionType::SizeType OutSizeType;
210 
211  // types to support the run length encoding of lines
212  class runLength
213  {
214 public:
215  // run length information - may be a more type safe way of doing this
217  typename TInputImage::IndexType where; // Index of the start of the run
218  LabelType label; // the initial label of the run
219  };
220 
221  typedef std::vector< runLength > lineEncoding;
222 
223  // the map storing lines
224  typedef std::vector< lineEncoding > LineMapType;
225 
226  typedef std::vector< typename TInputImage::OffsetValueType > OffsetVec;
227 
228  // the types to support union-find operations
229  typedef std::vector< LabelType > UnionFindType;
232 
233  // functions to support union-find operations
235  {
236  m_UnionFind = UnionFindType(size + 1);
237  }
238 
239  void InsertSet(const LabelType label);
240 
241  SizeValueType LookupSet(const LabelType label);
242 
243  void LinkLabels(const LabelType lab1, const LabelType lab2);
244 
246 
248  bool CheckNeighbors(const OutputIndexType & A,
249  const OutputIndexType & B);
250 
251  void CompareLines(lineEncoding & current, const lineEncoding & Neighbour);
252 
253  void FillOutput(const LineMapType & LineMap,
254  ProgressReporter & progress);
255 
256  void SetupLineOffsets(OffsetVec & LineOffsets);
257 
258  void Wait()
259  {
260  // use m_NumberOfLabels.size() to get the number of thread used
261  if ( m_NumberOfLabels.size() > 1 )
262  {
263  m_Barrier->Wait();
264  }
265  }
266 
267  typename std::vector< IdentifierType > m_NumberOfLabels;
268  typename std::vector< IdentifierType > m_FirstLineIdToJoin;
269 
271 
272  typename TInputImage::ConstPointer m_Input;
273 #if !defined( ITK_WRAPPING_PARSER )
275 #endif
276 };
277 } // end namespace itk
278 
279 #ifndef ITK_MANUAL_INSTANTIATION
280 #if !defined( ITK_WRAPPING_PARSER )
281 #include "itkConnectedComponentImageFilter.hxx"
282 #endif
283 #endif
284 
285 #endif
typedef(Concept::SameDimension< itkGetStaticConstMacro(InputImageDimension), itkGetStaticConstMacro(OutputImageDimension) >) SameDimension
void LinkLabels(const LabelType lab1, const LabelType lab2)
void InsertSet(const LabelType label)
signed long OffsetValueType
Definition: itkIntTypes.h:154
std::vector< typename TInputImage::OffsetValueType > OffsetVec
InputImageType::Pointer InputImagePointer
ImageToImageFilter< TInputImage, TOutputImage > Superclass
TOutputImage::RegionType::SizeType OutSizeType
TInputImage::InternalPixelType InputInternalPixelType
Base class for all process objects that output image data.
unsigned long SizeValueType
Definition: itkIntTypes.h:143
void FillOutput(const LineMapType &LineMap, ProgressReporter &progress)
SizeValueType IdentifierType
Definition: itkIntTypes.h:147
void EnlargeOutputRequestedRegion(DataObject *) override
TOutputImage::InternalPixelType OutputInternalPixelType
bool CheckNeighbors(const OutputIndexType &A, const OutputIndexType &B)
OutputImageType::PixelType OutputImagePixelType
void SetupLineOffsets(OffsetVec &LineOffsets)
SizeValueType LookupSet(const LabelType label)
unsigned int ThreadIdType
Definition: itkIntTypes.h:159
Implements progress tracking for a filter.
DataObject * GetInput(const DataObjectIdentifierType &key)
Return an input.
void AfterThreadedGenerateData() override
void GenerateInputRequestedRegion() override
void ThreadedGenerateData(const RegionType &outputRegionForThread, ThreadIdType threadId) override
Base class for filters that take an image as input and produce an image as output.
Control indentation during Print() invocation.
Definition: itkIndent.h:49
void BeforeThreadedGenerateData() override
virtual void SetNthInput(DataObjectPointerArraySizeType num, DataObject *input)
#define itkConceptMacro(name, concept)
void PrintSelf(std::ostream &os, Indent indent) const override
Base class for all data objects in ITK.
void CompareLines(lineEncoding &current, const lineEncoding &Neighbour)
Label the objects in a binary image.