ITK 6.0.0
Insight Toolkit
 
Loading...
Searching...
No Matches
itkLabelGeometryImageFilter.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 itkLabelGeometryImageFilter_h
19#define itkLabelGeometryImageFilter_h
20
22#include "itkNumericTraits.h"
23#include "itkArray.h"
25#include <map>
26#include <vector>
27#include "vnl/algo/vnl_symmetric_eigensystem.h"
28#include "vnl/vnl_det.h"
29#include "itkMath.h"
30
31namespace itk
32{
80template <typename TLabelImage, typename TIntensityImage = TLabelImage>
81class ITK_TEMPLATE_EXPORT
82#if !defined(ITK_LEGACY_SILENT)
83 [[deprecated("This class contains known computational bugs. See class documentation for details.")]]
84#endif
85 LabelGeometryImageFilter : public ImageToImageFilter<TLabelImage, TIntensityImage>
86{
87public:
88 ITK_DISALLOW_COPY_AND_MOVE(LabelGeometryImageFilter);
89
95
97 itkNewMacro(Self);
98
100 itkOverrideGetNameOfClassMacro(LabelGeometryImageFilter);
101
103 using IntensityImageType = TIntensityImage;
104 using InputImagePointer = typename TIntensityImage::Pointer;
105 using RegionType = typename TIntensityImage::RegionType;
106 using SizeType = typename TIntensityImage::SizeType;
107 using IndexType = typename TIntensityImage::IndexType;
108 using PixelType = typename TIntensityImage::PixelType;
109
111 using LabelImageType = TLabelImage;
112 using LabelImagePointer = typename TLabelImage::Pointer;
113 using LabelRegionType = typename TLabelImage::RegionType;
114 using LabelSizeType = typename TLabelImage::SizeType;
115 using LabelIndexType = typename TLabelImage::IndexType;
116 using LabelPixelType = typename TLabelImage::PixelType;
117 using LabelPointType = typename TLabelImage::PointType;
118
120 static constexpr unsigned int ImageDimension = TLabelImage::ImageDimension;
121
124
127
130
134
135 // using BoundingBoxVerticesType = itk::FixedArray<
136 // LabelPointType,std::pow(2.0,Self::ImageDimension)>;
137 using BoundingBoxVerticesType = std::vector<LabelPointType>;
138
141
144
146 using LabelsType = std::vector<LabelPixelType>;
147
149 using LabelIndicesType = std::vector<LabelIndexType>;
150
152 using VectorType = std::vector<double>;
153
155 using MatrixType = vnl_matrix<double>;
156
162 {
163 public:
164 // default constructor
166 : m_Label(0)
167 , m_Sum(RealType{})
168 {
169 // initialized to the default values
170 const unsigned int imageDimension = Self::ImageDimension;
171
172 // m_BoundingBox.resize(imageDimension*2);
173 for (unsigned int i = 0; i < imageDimension * 2; i += 2)
174 {
177 }
178
180 m_BoundingBoxSize.Fill(0);
181 m_PixelIndices.clear();
182 m_Centroid.Fill(0);
183 m_WeightedCentroid.Fill(0);
188 m_Eigenvalues.clear();
190 m_Eigenvectors.fill(0);
191 m_AxesLength.Fill(0);
192 m_Eccentricity = 1;
193 m_Elongation = 1;
194 m_Orientation = 0;
195 LabelPointType emptyPoint{};
196 unsigned int numberOfVertices = 1 << ImageDimension;
197 m_OrientedBoundingBoxVertices.resize(numberOfVertices, emptyPoint);
200 m_OrientedLabelImage = LabelImageType::New();
201 m_OrientedIntensityImage = IntensityImageType::New();
204 m_RotationMatrix.fill(0.0);
205
208 for (unsigned int i = 0; i < ImageDimension; ++i)
209 {
210 for (unsigned int j = 0; j < ImageDimension; ++j)
211 {
212 m_SecondOrderRawMoments(i, j) = 0;
214 }
215 }
216 }
217
242 typename LabelImageType::Pointer m_OrientedLabelImage;
243 typename IntensityImageType::Pointer m_OrientedIntensityImage;
246 };
247
249 // Map from the label to the class storing all of the geometry information.
250 using MapType = std::map<LabelPixelType, LabelGeometry>;
251 using MapIterator = typename std::map<LabelPixelType, LabelGeometry>::iterator;
252 using MapConstIterator = typename std::map<LabelPixelType, LabelGeometry>::const_iterator;
253
254 // Macros for enabling the calculation of additional features.
255 itkGetMacro(CalculatePixelIndices, bool);
256 itkBooleanMacro(CalculatePixelIndices);
257 void
258 SetCalculatePixelIndices(const bool value)
259 {
260 // CalculateOrientedBoundingBox, CalculateOrientedLabelImage, and
261 // CalculateOrientedIntensityImage all need CalculatePixelIndices to be
262 // turned
263 // on if they are turned on. So, CalculatePixelIndices cannot be
264 // turned off if any of these flags are turned on.
265 if (value == false)
266 {
269 {
270 // We cannot change the value, so return.
271 return;
272 }
273 }
274
275 if (this->m_CalculatePixelIndices != value)
276 {
277 this->m_CalculatePixelIndices = value;
278 this->Modified();
279 }
280 }
281
282 itkGetMacro(CalculateOrientedBoundingBox, bool);
283 itkBooleanMacro(CalculateOrientedBoundingBox);
284 void
286 {
287 if (this->m_CalculateOrientedBoundingBox != value)
288 {
289 this->m_CalculateOrientedBoundingBox = value;
290 this->Modified();
291 }
292
293 // CalculateOrientedBoundingBox needs
294 // CalculatePixelIndices to be turned on.
295 if (value)
296 {
297 this->SetCalculatePixelIndices(true);
298 }
299 }
300
301 itkGetMacro(CalculateOrientedLabelRegions, bool);
302 itkBooleanMacro(CalculateOrientedLabelRegions);
303 void
305 {
306 if (this->m_CalculateOrientedLabelRegions != value)
307 {
309 this->Modified();
310
311 // CalculateOrientedLabelImage needs
312 // CalculateOrientedBoundingBox to be turned on.
313 if (value)
314 {
316 }
317 }
318 }
319
320 itkGetMacro(CalculateOrientedIntensityRegions, bool);
321 itkBooleanMacro(CalculateOrientedIntensityRegions);
322 void
324 {
325 if (this->m_CalculateOrientedIntensityRegions != value)
326 {
328 this->Modified();
329
330 // CalculateOrientedIntensityImage needs
331 // CalculateOrientedBoundingBox to be turned on.
332 if (value)
333 {
335 }
336 }
337 }
338
340 void
341 SetIntensityInput(const TIntensityImage * input)
342 {
343 // Process object is not const-correct so the const casting is required.
344 this->SetNthInput(1, const_cast<TIntensityImage *>(input));
345 }
346
348 const TIntensityImage *
350 {
351 return static_cast<TIntensityImage *>(const_cast<DataObject *>(this->ProcessObject::GetInput(1)));
352 }
353
356 bool
358 {
359 return m_LabelGeometryMapper.find(label) != m_LabelGeometryMapper.end();
360 }
361
363 [[nodiscard]] SizeValueType
365 {
366 return m_LabelGeometryMapper.size();
367 }
368
369 [[nodiscard]] SizeValueType
371 {
372 return this->GetNumberOfObjects();
373 }
374
376 std::vector<LabelPixelType>
377 GetLabels() const
378 {
379 return m_AllLabels;
380 }
381
383 LabelIndicesType
385
390
394
398
402
406
410
414
419
424
428
433
437
443
447
451
460
464
468
472
477
481
483 TLabelImage *
485
488 TIntensityImage *
490
491 itkConceptMacro(InputHasNumericTraitsCheck, (Concept::HasNumericTraits<PixelType>));
492
493protected:
495 ~LabelGeometryImageFilter() override = default;
496 void
497 PrintSelf(std::ostream & os, Indent indent) const override;
498
499 void
500 GenerateData() override;
501
502private:
503 bool
504 CalculateOrientedBoundingBoxVertices(vnl_symmetric_eigensystem<double> eig, LabelGeometry & m_LabelGeometry);
505
510
514}; // end of class
515
516} // end namespace itk
517
518#ifndef ITK_MANUAL_INSTANTIATION
519# include "itkLabelGeometryImageFilter.hxx"
520#endif
521
522#endif
Base class for all data objects in ITK.
SmartPointer< Self > Pointer
Simulate a standard C array with copy semantics.
Control indentation during Print() invocation.
Definition itkIndent.h:50
FixedArray< float, Self::ImageDimension > m_AxesLength
typename TLabelImage::RegionType LabelRegionType
void SetCalculateOrientedIntensityRegions(const bool value)
RealType GetMajorAxisLength(LabelPixelType label) const
const TIntensityImage * GetIntensityInput() const
RealType GetElongation(LabelPixelType label) const
SmartPointer< const Self > ConstPointer
typename std::map< LabelPixelType, LabelGeometry >::iterator MapIterator
LabelSizeType GetBoundingBoxSize(LabelPixelType label) const
VectorType GetEigenvalues(LabelPixelType label) const
typename TLabelImage::SizeType LabelSizeType
typename TIntensityImage::SizeType SizeType
typename TIntensityImage::Pointer InputImagePointer
typename DataObject::Pointer DataObjectPointer
ImageToImageFilter< TLabelImage, TIntensityImage > Superclass
RealType GetEccentricity(LabelPixelType label) const
SimpleDataObjectDecorator< RealType > RealObjectType
void SetCalculateOrientedBoundingBox(const bool value)
itk::FixedArray< typename LabelIndexType::IndexValueType, Self::ImageDimension > IndexArrayType
typename TIntensityImage::PixelType PixelType
typename TLabelImage::Pointer LabelImagePointer
LabelPointType GetCentroid(LabelPixelType label) const
void PrintSelf(std::ostream &os, Indent indent) const override
BoundingBoxVerticesType GetOrientedBoundingBoxVertices(LabelPixelType label) const
std::vector< LabelPixelType > LabelsType
LabelIndicesType GetPixelIndices(LabelPixelType label) const
typename TLabelImage::PixelType LabelPixelType
bool CalculateOrientedBoundingBoxVertices(vnl_symmetric_eigensystem< double > eig, LabelGeometry &m_LabelGeometry)
std::vector< LabelPixelType > GetLabels() const
std::vector< LabelPointType > BoundingBoxVerticesType
RealType GetIntegratedIntensity(LabelPixelType label) const
SizeValueType GetVolume(LabelPixelType label) const
MatrixType GetRotationMatrix(LabelPixelType label) const
AxesLengthType GetAxesLength(LabelPixelType label) const
typename TIntensityImage::IndexType IndexType
itk::FixedArray< float, Self::ImageDimension *2 > BoundingBoxFloatType
RealType GetBoundingBoxVolume(LabelPixelType label) const
TLabelImage * GetOrientedLabelImage(LabelPixelType label) const
bool HasLabel(LabelPixelType label) const
LabelPointType GetWeightedCentroid(LabelPixelType label) const
itk::FixedArray< typename LabelIndexType::IndexValueType, Self::ImageDimension *2 > BoundingBoxType
TIntensityImage * GetOrientedIntensityImage(LabelPixelType label) const
typename TIntensityImage::RegionType RegionType
RealType GetOrientedBoundingBoxVolume(LabelPixelType label) const
RealType GetMinorAxisLength(LabelPixelType label) const
MatrixType GetEigenvectors(LabelPixelType label) const
itk::FixedArray< RealType, Self::ImageDimension > AxesLengthType
BoundingBoxType GetBoundingBox(LabelPixelType label) const
typename TLabelImage::PointType LabelPointType
static constexpr unsigned int ImageDimension
LabelPointType GetOrientedBoundingBoxSize(LabelPixelType label) const
std::vector< LabelIndexType > LabelIndicesType
typename TLabelImage::IndexType LabelIndexType
typename std::map< LabelPixelType, LabelGeometry >::const_iterator MapConstIterator
std::map< LabelPixelType, LabelGeometry > MapType
void SetIntensityInput(const TIntensityImage *input)
~LabelGeometryImageFilter() override=default
LabelPointType GetOrientedBoundingBoxOrigin(LabelPixelType label) const
void SetCalculateOrientedLabelRegions(const bool value)
RealType GetOrientation(LabelPixelType label) const
typename NumericTraits< PixelType >::RealType RealType
RegionType GetRegion(LabelPixelType label) const
static constexpr T NonpositiveMin()
static constexpr T max(const T &)
virtual void Modified() const
virtual void SetNthInput(DataObjectPointerArraySizeType idx, DataObject *input)
DataObject * GetInput(const DataObjectIdentifierType &key)
Return an input.
Decorates any "simple" data type (data types without smart pointers) with a DataObject API.
Implements transparent reference counting.
#define itkConceptMacro(name, concept)
The "itk" namespace contains all Insight Segmentation and Registration Toolkit (ITK) classes....
unsigned long SizeValueType
Definition itkIntTypes.h:86