ITK  6.0.0
Insight Toolkit
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;
108 using PixelType = typename TIntensityImage::PixelType;
109
111 using LabelImageType = TLabelImage;
116 using LabelPixelType = typename TLabelImage::PixelType;
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 {
167 // initialized to the default values
168 this->m_Label = 0;
169 this->m_Sum = RealType{};
170
171 const unsigned int imageDimension = Self::ImageDimension;
172
173 // m_BoundingBox.resize(imageDimension*2);
174 for (unsigned int i = 0; i < imageDimension * 2; i += 2)
175 {
178 }
179
180 m_BoundingBoxVolume = 0;
181 m_BoundingBoxSize.Fill(0);
182 m_PixelIndices.clear();
183 m_Centroid.Fill(0);
184 m_WeightedCentroid.Fill(0);
185 m_ZeroOrderMoment = 0;
186 m_FirstOrderRawMoments.Fill(0);
187 m_FirstOrderWeightedRawMoments.Fill(0);
188 m_Eigenvalues.resize(ImageDimension);
189 m_Eigenvalues.clear();
190 m_Eigenvectors.set_size(ImageDimension, ImageDimension);
191 m_Eigenvectors.fill(0);
192 m_AxesLength.Fill(0);
193 m_Eccentricity = 1;
194 m_Elongation = 1;
195 m_Orientation = 0;
196 LabelPointType emptyPoint{};
197 unsigned int numberOfVertices = 1 << ImageDimension;
198 m_OrientedBoundingBoxVertices.resize(numberOfVertices, emptyPoint);
199 m_OrientedBoundingBoxVolume = 0;
200 m_OrientedBoundingBoxSize.Fill(0);
201 m_OrientedLabelImage = LabelImageType::New();
202 m_OrientedIntensityImage = IntensityImageType::New();
203 m_OrientedBoundingBoxOrigin.Fill(0);
204 m_RotationMatrix.set_size(ImageDimension, ImageDimension);
205 m_RotationMatrix.fill(0.0);
206
207 m_SecondOrderRawMoments.set_size(ImageDimension, ImageDimension);
208 m_SecondOrderCentralMoments.set_size(ImageDimension, ImageDimension);
209 for (unsigned int i = 0; i < ImageDimension; ++i)
210 {
211 for (unsigned int j = 0; j < ImageDimension; ++j)
212 {
213 m_SecondOrderRawMoments(i, j) = 0;
214 m_SecondOrderCentralMoments(i, j) = 0;
215 }
216 }
217 }
218
247 };
248
250 // Map from the label to the class storing all of the geometry information.
251 using MapType = std::map<LabelPixelType, LabelGeometry>;
252 using MapIterator = typename std::map<LabelPixelType, LabelGeometry>::iterator;
253 using MapConstIterator = typename std::map<LabelPixelType, LabelGeometry>::const_iterator;
254
255 // Macros for enabling the calculation of additional features.
256 itkGetMacro(CalculatePixelIndices, bool);
257 itkBooleanMacro(CalculatePixelIndices);
258 void 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 {
267 if ((this->m_CalculateOrientedBoundingBox) || (this->m_CalculateOrientedLabelRegions) ||
268 (this->m_CalculateOrientedIntensityRegions))
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 SetCalculateOrientedBoundingBox(const bool value)
285 {
286 if (this->m_CalculateOrientedBoundingBox != value)
287 {
288 this->m_CalculateOrientedBoundingBox = value;
289 this->Modified();
290 }
291
292 // CalculateOrientedBoundingBox needs
293 // CalculatePixelIndices to be turned on.
294 if (value)
295 {
296 this->SetCalculatePixelIndices(true);
297 }
298 }
299
300 itkGetMacro(CalculateOrientedLabelRegions, bool);
301 itkBooleanMacro(CalculateOrientedLabelRegions);
302 void SetCalculateOrientedLabelRegions(const bool value)
303 {
304 if (this->m_CalculateOrientedLabelRegions != value)
305 {
306 this->m_CalculateOrientedLabelRegions = value;
307 this->Modified();
308
309 // CalculateOrientedLabelImage needs
310 // CalculateOrientedBoundingBox to be turned on.
311 if (value)
312 {
313 SetCalculateOrientedBoundingBox(true);
314 }
315 }
316 }
317
318 itkGetMacro(CalculateOrientedIntensityRegions, bool);
319 itkBooleanMacro(CalculateOrientedIntensityRegions);
321 {
322 if (this->m_CalculateOrientedIntensityRegions != value)
323 {
324 this->m_CalculateOrientedIntensityRegions = value;
325 this->Modified();
326
327 // CalculateOrientedIntensityImage needs
328 // CalculateOrientedBoundingBox to be turned on.
329 if (value)
330 {
331 this->SetCalculateOrientedBoundingBox(true);
332 }
333 }
334 }
335
337 void SetIntensityInput(const TIntensityImage * input)
338 {
339 // Process object is not const-correct so the const casting is required.
340 this->SetNthInput(1, const_cast<TIntensityImage *>(input));
341 }
342
344 const TIntensityImage * GetIntensityInput() const
345 {
346 return static_cast<TIntensityImage *>(const_cast<DataObject *>(this->ProcessObject::GetInput(1)));
347 }
348
351 bool HasLabel(LabelPixelType label) const { return m_LabelGeometryMapper.find(label) != m_LabelGeometryMapper.end(); }
352
354 SizeValueType GetNumberOfObjects() const { return m_LabelGeometryMapper.size(); }
355
356 SizeValueType GetNumberOfLabels() const { return this->GetNumberOfObjects(); }
357
359 std::vector<LabelPixelType> GetLabels() const { return m_AllLabels; }
360
363
367
370
373
376
379
382
385
389
393
396
400
403
408
411
414
422
425
428
431
435
438
440 TLabelImage * GetOrientedLabelImage(LabelPixelType label) const;
441
444 TIntensityImage * GetOrientedIntensityImage(LabelPixelType label) const;
445
446#ifdef ITK_USE_CONCEPT_CHECKING
447 // Begin concept checking
448 itkConceptMacro(InputHasNumericTraitsCheck, (Concept::HasNumericTraits<PixelType>));
449 // End concept checking
450#endif
451
452protected:
454 ~LabelGeometryImageFilter() override = default;
455 void PrintSelf(std::ostream & os, Indent indent) const override;
456
457 void GenerateData() override;
458
459private:
460 bool CalculateOrientedBoundingBoxVertices(vnl_symmetric_eigensystem<double> eig, LabelGeometry & m_LabelGeometry);
461
462 bool m_CalculatePixelIndices{};
463 bool m_CalculateOrientedBoundingBox{};
464 bool m_CalculateOrientedLabelRegions{};
465 bool m_CalculateOrientedIntensityRegions{};
466
467 MapType m_LabelGeometryMapper{};
468 LabelGeometry m_LabelGeometry{};
469 LabelsType m_AllLabels{};
470}; // end of class
471
472} // end namespace itk
473
474#ifndef ITK_MANUAL_INSTANTIATION
475# include "itkLabelGeometryImageFilter.hxx"
476#endif
477
478#endif
Base class for all data objects in ITK.
SmartPointer< Self > Pointer
Base class for filters that take an image as input and produce an image as output.
Control indentation during Print() invocation.
Definition: itkIndent.h:50
FixedArray< float, Self::ImageDimension > m_AxesLength
Given a label map and an optional intensity image, compute geometric features.
typename TLabelImage::RegionType LabelRegionType
void SetCalculateOrientedIntensityRegions(const bool value)
RealType GetMajorAxisLength(LabelPixelType label) const
const TIntensityImage * GetIntensityInput() const
RealType GetElongation(LabelPixelType label) const
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
RealType GetEccentricity(LabelPixelType label) const
void SetCalculateOrientedBoundingBox(const bool value)
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
void SetCalculatePixelIndices(const bool value)
SizeValueType GetVolume(LabelPixelType label) const
MatrixType GetRotationMatrix(LabelPixelType label) const
AxesLengthType GetAxesLength(LabelPixelType label) const
typename TIntensityImage::IndexType IndexType
RealType GetBoundingBoxVolume(LabelPixelType label) const
TLabelImage * GetOrientedLabelImage(LabelPixelType label) const
bool HasLabel(LabelPixelType label) const
LabelPointType GetWeightedCentroid(LabelPixelType label) const
TIntensityImage * GetOrientedIntensityImage(LabelPixelType label) const
typename TIntensityImage::RegionType RegionType
RealType GetOrientedBoundingBoxVolume(LabelPixelType label) const
RealType GetMinorAxisLength(LabelPixelType label) const
MatrixType GetEigenvectors(LabelPixelType label) const
BoundingBoxType GetBoundingBox(LabelPixelType label) const
typename TLabelImage::PointType LabelPointType
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
Light weight base class for most itk classes.
Define additional traits for native types such as int or float.
static constexpr T NonpositiveMin()
static constexpr T max(const T &)
DataObject * GetInput(const DataObjectIdentifierType &key)
Return an input.
Decorates any "simple" data type (data types without smart pointers) with a DataObject API.
static Pointer New()
SmartPointer< Self > Pointer
#define itkConceptMacro(name, concept)
The "itk" namespace contains all Insight Segmentation and Registration Toolkit (ITK) classes....
unsigned long SizeValueType
Definition: itkIntTypes.h:86