ITK 6.0.0
Insight Toolkit
 
Loading...
Searching...
No Matches
itkLabelStatisticsImageFilter.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 itkLabelStatisticsImageFilter_h
19#define itkLabelStatisticsImageFilter_h
20
21#include "itkImageSink.h"
22#include "itkNumericTraits.h"
24#include "itkHistogram.h"
25#include "itkPrintHelper.h"
26#include <mutex>
27#include <unordered_map>
28#include <vector>
29
30namespace itk
31{
62template <typename TInputImage, typename TLabelImage>
63class ITK_TEMPLATE_EXPORT LabelStatisticsImageFilter : public ImageSink<TInputImage>
64{
65public:
66 ITK_DISALLOW_COPY_AND_MOVE(LabelStatisticsImageFilter);
67
73
75 itkNewMacro(Self);
76
78 itkOverrideGetNameOfClassMacro(LabelStatisticsImageFilter);
79
81 using InputImagePointer = typename TInputImage::Pointer;
82 using RegionType = typename TInputImage::RegionType;
83 using SizeType = typename TInputImage::SizeType;
84 using IndexType = typename TInputImage::IndexType;
85 using PixelType = typename TInputImage::PixelType;
86
88 using LabelImageType = TLabelImage;
89 using LabelImagePointer = typename TLabelImage::Pointer;
90 using LabelRegionType = typename TLabelImage::RegionType;
91 using LabelSizeType = typename TLabelImage::SizeType;
92 using LabelIndexType = typename TLabelImage::IndexType;
93 using LabelPixelType = typename TLabelImage::PixelType;
94
96 static constexpr unsigned int ImageDimension = TInputImage::ImageDimension;
97
100
103
106
108 using BoundingBoxType = std::vector<IndexValueType>;
109
113
119 {
120 public:
121 // default constructor
123 // initialized to the default values
126 , m_Maximum(NumericTraits<RealType>::NonpositiveMin())
127 , m_Mean(RealType{})
128 , m_Sum(RealType{})
130 , m_Sigma(RealType{})
132 {
133 // Set such that the first pixel encountered can be compared
134
135 const unsigned int imageDimension = Self::ImageDimension;
136 m_BoundingBox.resize(imageDimension * 2);
137 for (unsigned int i = 0; i < imageDimension * 2; i += 2)
138 {
141 }
142 m_Histogram = nullptr;
143 }
144
145 // constructor with histogram enabled
146 LabelStatistics(int size, RealType lowerBound, RealType upperBound)
149 , m_Maximum(NumericTraits<RealType>::NonpositiveMin())
150 , m_Mean(RealType{})
151 , m_Sum(RealType{})
153 , m_Sigma(RealType{})
155 {
156 // Set such that the first pixel encountered can be compared
157 const unsigned int imageDimension = Self::ImageDimension;
158 m_BoundingBox.resize(imageDimension * 2);
159 for (unsigned int i = 0; i < imageDimension * 2; i += 2)
160 {
163 }
164
165 // Histogram
167 typename HistogramType::SizeType hsize;
170 hsize.SetSize(1);
171 lb.SetSize(1);
172 ub.SetSize(1);
173 m_Histogram->SetMeasurementVectorSize(1);
174 hsize[0] = size;
175 lb[0] = lowerBound;
176 ub[0] = upperBound;
177 m_Histogram->Initialize(hsize, lb, ub);
178 }
179
180 // need copy constructor because of smart pointer to histogram
193
195
196 // added for completeness
199 {
200 if (this != &l)
201 {
202 m_Count = l.m_Count;
205 m_Mean = l.m_Mean;
206 m_Sum = l.m_Sum;
208 m_Sigma = l.m_Sigma;
212 }
213 return *this;
214 }
215
216 friend std::ostream &
217 operator<<(std::ostream & os, const LabelStatistics & labelStatistics)
218 {
219 using namespace print_helper;
220
221 os << "Count: " << static_cast<typename NumericTraits<IdentifierType>::PrintType>(labelStatistics.m_Count)
222 << std::endl;
223 os << "Minimum: " << static_cast<typename NumericTraits<RealType>::PrintType>(labelStatistics.m_Minimum)
224 << std::endl;
225 os << "Maximum: " << static_cast<typename NumericTraits<RealType>::PrintType>(labelStatistics.m_Maximum)
226 << std::endl;
227 os << "Mean: " << static_cast<typename NumericTraits<RealType>::PrintType>(labelStatistics.m_Mean) << std::endl;
228 os << "Sum: " << static_cast<typename NumericTraits<RealType>::PrintType>(labelStatistics.m_Sum) << std::endl;
229 os << "SumOfSquares: " << static_cast<typename NumericTraits<RealType>::PrintType>(labelStatistics.m_SumOfSquares)
230 << std::endl;
231 os << "Sigma: " << static_cast<typename NumericTraits<RealType>::PrintType>(labelStatistics.m_Sigma) << std::endl;
232 os << "Variance: " << static_cast<typename NumericTraits<RealType>::PrintType>(labelStatistics.m_Variance)
233 << std::endl;
234 os << "BoundingBox: " << labelStatistics.m_BoundingBox << std::endl;
235
236 os << "Histogram: ";
237 if (labelStatistics.m_Histogram)
238 {
239 labelStatistics.m_Histogram->Print(os);
240 }
241 else
242 {
243 os << "nullptr" << std::endl;
244 }
245
246 return os;
247 }
248
259 };
260
262 using MapType = std::unordered_map<LabelPixelType, LabelStatistics>;
263 using MapIterator = typename MapType::iterator;
264 using MapConstIterator = typename MapType::const_iterator;
266
268 using ValidLabelValuesContainerType = std::vector<LabelPixelType>;
269
270 // macros for Histogram enables
271 itkSetMacro(UseHistograms, bool);
272 itkGetConstMacro(UseHistograms, bool);
273 itkBooleanMacro(UseHistograms);
274
275
276 virtual const ValidLabelValuesContainerType &
278 {
279 return m_ValidLabelValues;
280 }
281
284 itkSetInputMacro(LabelInput, TLabelImage);
285 itkGetInputMacro(LabelInput, TLabelImage);
289 bool
291 {
292 return m_LabelStatistics.find(label) != m_LabelStatistics.end();
293 }
294
296 [[nodiscard]] MapSizeType
298 {
299 return static_cast<MapSizeType>(m_LabelStatistics.size());
300 }
301
302 [[nodiscard]] MapSizeType
304 {
305 return static_cast<MapSizeType>(this->GetNumberOfObjects());
306 }
307
311
315
319
324
328
332
338
342
345 GetSum(LabelPixelType label) const;
346
350
354
356 void
357 SetHistogramParameters(const int numBins, RealType lowerBound, RealType upperBound);
358
359 // Change the access from protected to public to expose streaming option, a using statement can not be used due to
360 // limitations of wrapping.
361 void
362 SetNumberOfStreamDivisions(const unsigned int n) override
363 {
365 }
366 [[nodiscard]] unsigned int
368 {
370 }
371
372 itkConceptMacro(InputHasNumericTraitsCheck, (Concept::HasNumericTraits<PixelType>));
373
374protected:
376 ~LabelStatisticsImageFilter() override = default;
377 void
378 PrintSelf(std::ostream & os, Indent indent) const override;
379
380 void
382 {
383 this->AllocateOutputs();
384 m_LabelStatistics.clear();
385 }
386
389 void
391
392 void
394
395private:
396 void
398
401
403
405
408
409 std::mutex m_Mutex{};
410
411}; // end of class
412} // end namespace itk
413
414#ifndef ITK_MANUAL_INSTANTIATION
415# include "itkLabelStatisticsImageFilter.hxx"
416#endif
417
418#endif
void SetSize(SizeValueType sz)
SmartPointer< Self > Pointer
virtual void AllocateOutputs()
virtual void SetNumberOfStreamDivisions(unsigned int _arg)
virtual unsigned int GetNumberOfStreamDivisions() const
Control indentation during Print() invocation.
Definition itkIndent.h:50
LabelStatistics(int size, RealType lowerBound, RealType upperBound)
LabelStatistics & operator=(const LabelStatistics &l)
typename DataObject::Pointer DataObjectPointer
void SetNumberOfStreamDivisions(const unsigned int n) override
~LabelStatisticsImageFilter() override=default
HistogramPointer GetHistogram(LabelPixelType label) const
RealType GetMean(LabelPixelType label) const
unsigned int GetNumberOfStreamDivisions() const override
typename TLabelImage::RegionType LabelRegionType
itk::Statistics::Histogram< RealType > HistogramType
typename HistogramType::Pointer HistogramPointer
typename TInputImage::IndexType IndexType
SimpleDataObjectDecorator< RealType > RealObjectType
RealType GetSigma(LabelPixelType label) const
BoundingBoxType GetBoundingBox(LabelPixelType label) const
void MergeMap(MapType &, MapType &) const
MapSizeType GetCount(LabelPixelType label) const
typename TLabelImage::Pointer LabelImagePointer
RealType GetMinimum(LabelPixelType label) const
typename TInputImage::RegionType RegionType
RegionType GetRegion(LabelPixelType label) const
bool HasLabel(LabelPixelType label) const
typename TLabelImage::PixelType LabelPixelType
itkSetInputMacro(LabelInput, TLabelImage)
std::vector< IndexValueType > BoundingBoxType
ValidLabelValuesContainerType m_ValidLabelValues
void AfterStreamedGenerateData() override
std::vector< LabelPixelType > ValidLabelValuesContainerType
typename TLabelImage::IndexType LabelIndexType
itkGetInputMacro(LabelInput, TLabelImage)
static constexpr unsigned int ImageDimension
typename TLabelImage::SizeType LabelSizeType
RealType GetSum(LabelPixelType label) const
typename TInputImage::PixelType PixelType
RealType GetMaximum(LabelPixelType label) const
typename TInputImage::Pointer InputImagePointer
std::unordered_map< LabelPixelType, LabelStatistics > MapType
virtual const ValidLabelValuesContainerType & GetValidLabelValues() const
typename TInputImage::SizeType SizeType
void PrintSelf(std::ostream &os, Indent indent) const override
RealType GetMedian(LabelPixelType label) const
typename NumericTraits< PixelType >::RealType RealType
void SetHistogramParameters(const int numBins, RealType lowerBound, RealType upperBound)
typename MapType::const_iterator MapConstIterator
void ThreadedStreamedGenerateData(const RegionType &) override
RealType GetVariance(LabelPixelType label) const
Define additional traits for native types such as int or float.
static constexpr T NonpositiveMin()
static constexpr T max(const T &)
Decorates any "simple" data type (data types without smart pointers) with a DataObject API.
Implements transparent reference counting.
ObjectType * Print(std::ostream &os) const
This class stores measurement vectors in the context of n-dimensional histogram.
Array< itk::SizeValueType > SizeType
#define itkConceptMacro(name, concept)
std::ostream & operator<<(std::ostream &os, const itk::DOMNode &object)
The "itk" namespace contains all Insight Segmentation and Registration Toolkit (ITK) classes....
SizeValueType IdentifierType
Definition itkIntTypes.h:90