ITK  4.13.0
Insight Segmentation and Registration Toolkit
itkLabelOverlapMeasuresImageFilter.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 itkLabelOverlapMeasuresImageFilter_h
19 #define itkLabelOverlapMeasuresImageFilter_h
20 
21 #include "itkImageToImageFilter.h"
22 #include "itkNumericTraits.h"
23 
24 #include "itksys/hash_map.hxx"
25 
26 namespace itk {
27 
44 template<typename TLabelImage>
45 class ITK_TEMPLATE_EXPORT LabelOverlapMeasuresImageFilter :
46  public ImageToImageFilter<TLabelImage, TLabelImage>
47 {
48 public:
54 
56  itkNewMacro( Self );
57 
60 
62  typedef TLabelImage LabelImageType;
63  typedef typename TLabelImage::Pointer LabelImagePointer;
64  typedef typename TLabelImage::ConstPointer LabelImageConstPointer;
65 
66  typedef typename TLabelImage::RegionType RegionType;
67  typedef typename TLabelImage::SizeType SizeType;
69 
70  typedef typename TLabelImage::PixelType LabelType;
71 
74 
80  {
81  public:
82  // default constructor
84  {
85  m_Source = 0;
86  m_Target = 0;
87  m_Union = 0;
88  m_Intersection = 0;
89  m_SourceComplement = 0;
90  m_TargetComplement = 0;
91  }
92 
93  // added for completeness
95  {
96  if(this != &l)
97  {
98  m_Source = l.m_Source;
99  m_Target = l.m_Target;
100  m_Union = l.m_Union;
101  m_Intersection = l.m_Intersection;
102  m_SourceComplement = l.m_SourceComplement;
103  m_TargetComplement = l.m_TargetComplement;
104  }
105  return *this;
106  }
107 
108  unsigned long m_Source;
109  unsigned long m_Target;
110  unsigned long m_Union;
111  unsigned long m_Intersection;
112  unsigned long m_SourceComplement;
113  unsigned long m_TargetComplement;
114  };
115 
117  typedef itksys::hash_map<LabelType, LabelSetMeasures> MapType;
118  typedef typename MapType::iterator MapIterator;
119  typedef typename MapType::const_iterator MapConstIterator;
120 
122  itkStaticConstMacro( ImageDimension, unsigned int,
123  TLabelImage::ImageDimension );
124 
126  void SetSourceImage( const LabelImageType * image )
127  { this->SetNthInput( 0, const_cast<LabelImageType *>( image ) ); }
128 
130  void SetTargetImage( const LabelImageType * image )
131  { this->SetNthInput( 1, const_cast<LabelImageType *>( image ) ); }
132 
135  { return this->GetInput( 0 ); }
136 
139  { return this->GetInput( 1 ); }
140 
143  { return this->m_LabelSetMeasures; }
144 
145  // Overlap agreement metrics
146 
148  RealType GetTotalOverlap() const;
149 
151  RealType GetTargetOverlap( LabelType ) const;
152 
154  RealType GetUnionOverlap() const;
156  { return this->GetUnionOverlap(); }
158 
161  RealType GetUnionOverlap( LabelType ) const;
163  { return this->GetUnionOverlap( label ); }
165 
167  RealType GetMeanOverlap() const;
169  { return this->GetMeanOverlap(); }
171 
174  RealType GetMeanOverlap( LabelType ) const;
176  { return this->GetMeanOverlap( label ); }
178 
180  RealType GetVolumeSimilarity() const;
181 
183  RealType GetVolumeSimilarity( LabelType ) const;
184 
185  // Overlap error metrics
186 
188  RealType GetFalseNegativeError() const;
189 
191  RealType GetFalseNegativeError( LabelType ) const;
192 
194  RealType GetFalsePositiveError() const;
195 
197  RealType GetFalsePositiveError( LabelType ) const;
198 
199 #ifdef ITK_USE_CONCEPT_CHECKING
200  // Begin concept checking
201  itkConceptMacro( Input1HasNumericTraitsCheck,
203  // End concept checking
204 #endif
205 
206 protected:
209  void PrintSelf( std::ostream& os, Indent indent ) const ITK_OVERRIDE;
210 
216  void AllocateOutputs() ITK_OVERRIDE;
217 
218  void BeforeThreadedGenerateData() ITK_OVERRIDE;
219 
220  void AfterThreadedGenerateData() ITK_OVERRIDE;
221 
223  void ThreadedGenerateData( const RegionType&, ThreadIdType ) ITK_OVERRIDE;
224 
225  // Override since the filter produces all of its output
226  void EnlargeOutputRequestedRegion( DataObject *data ) ITK_OVERRIDE;
227 
228 private:
229  ITK_DISALLOW_COPY_AND_ASSIGN(LabelOverlapMeasuresImageFilter);
230 
231  std::vector<MapType> m_LabelSetMeasuresPerThread;
232  MapType m_LabelSetMeasures;
233 }; // end of class
234 
235 } // end namespace itk
236 
237 #ifndef ITK_MANUAL_INSTANTIATION
238 #include "itkLabelOverlapMeasuresImageFilter.hxx"
239 #endif
240 
241 #endif
Light weight base class for most itk classes.
NumericTraits< LabelType >::RealType RealType
RealType GetJaccardCoefficient(LabelType label) const
itksys::hash_map< LabelType, LabelSetMeasures > MapType
Computes overlap measures between the set same set of labels of pixels of two images. Background is assumed to be 0.
unsigned int ThreadIdType
Definition: itkIntTypes.h:159
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
#define itkConceptMacro(name, concept)
Base class for all data objects in ITK.
ImageToImageFilter< TLabelImage, TLabelImage > Superclass