ITK  5.4.0
Insight Toolkit
itkImageMomentsCalculator.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 itkImageMomentsCalculator_h
19 #define itkImageMomentsCalculator_h
20 
21 #include "itkAffineTransform.h"
22 #include "itkImage.h"
23 #include "itkSpatialObject.h"
24 
25 #include "vnl/vnl_vector_fixed.h"
26 #include "vnl/vnl_matrix_fixed.h"
27 #include "vnl/vnl_diag_matrix.h"
28 
29 namespace itk
30 {
59 template <typename TImage>
60 class ITK_TEMPLATE_EXPORT ImageMomentsCalculator : public Object
61 {
62 public:
63  ITK_DISALLOW_COPY_AND_MOVE(ImageMomentsCalculator);
64 
67  using Superclass = Object;
70 
72  itkNewMacro(Self);
73 
75  itkOverrideGetNameOfClassMacro(ImageMomentsCalculator);
76 
78  static constexpr unsigned int ImageDimension = TImage::ImageDimension;
79 
81  using ScalarType = double;
82 
85 
88 
92 
95 
97  using ImageType = TImage;
98 
102 
106 
108  virtual void
109  SetImage(const ImageType * image)
110  {
111  if (m_Image != image)
112  {
113  m_Image = image;
114  this->Modified();
115  m_Valid = false;
116  }
117  }
121  virtual void
123  {
124  if (m_SpatialObjectMask != so)
125  {
126  m_SpatialObjectMask = so;
127  this->Modified();
128  m_Valid = false;
129  }
130  }
138  void
139  Compute();
140 
145  ScalarType
146  GetTotalMass() const;
147 
153  VectorType
154  GetFirstMoments() const;
155 
161  MatrixType
162  GetSecondMoments() const;
163 
168  VectorType
169  GetCenterOfGravity() const;
170 
175  MatrixType
176  GetCentralMoments() const;
177 
184  VectorType
185  GetPrincipalMoments() const;
186 
199  MatrixType
200  GetPrincipalAxes() const;
206  AffineTransformPointer
207  GetPrincipalAxesToPhysicalAxesTransform() const;
208 
213  AffineTransformPointer
214  GetPhysicalAxesToPrincipalAxesTransform() const;
215 
216 protected:
218  ~ImageMomentsCalculator() override = default;
219  void
220  PrintSelf(std::ostream & os, Indent indent) const override;
221 
222 private:
223  bool m_Valid{}; // Have moments been computed yet?
224  ScalarType m_M0{}; // Zeroth moment
225  VectorType m_M1{}; // First moments about origin
226  MatrixType m_M2{}; // Second moments about origin
227  VectorType m_Cg{}; // Center of gravity (physical units)
228  MatrixType m_Cm{}; // Second central moments (physical)
229  VectorType m_Pm{}; // Principal moments (physical)
230  MatrixType m_Pa{}; // Principal axes (physical)
231 
232  ImageConstPointer m_Image{};
233  SpatialObjectConstPointer m_SpatialObjectMask{};
234 }; // class ImageMomentsCalculator
235 } // end namespace itk
236 
237 #ifndef ITK_MANUAL_INSTANTIATION
238 # include "itkImageMomentsCalculator.hxx"
239 #endif
240 
241 #endif /* itkImageMomentsCalculator_h */
Pointer
SmartPointer< Self > Pointer
Definition: itkAddImageFilter.h:93
itk::ImageMomentsCalculator::ImageConstPointer
typename ImageType::ConstPointer ImageConstPointer
Definition: itkImageMomentsCalculator.h:101
ConstPointer
SmartPointer< const Self > ConstPointer
Definition: itkAddImageFilter.h:94
itk::ImageMomentsCalculator::ImagePointer
typename ImageType::Pointer ImagePointer
Definition: itkImageMomentsCalculator.h:100
itkSpatialObject.h
itk::GTest::TypedefsAndConstructors::Dimension2::VectorType
ImageBaseType::SpacingType VectorType
Definition: itkGTestTypedefsAndConstructors.h:53
itk::Vector< ScalarType, Self::ImageDimension >
itkImage.h
itk::SmartPointer< Self >
itk::Indent
Control indentation during Print() invocation.
Definition: itkIndent.h:49
itkAffineTransform.h
itk::AffineTransform
Definition: itkAffineTransform.h:101
itk::ImageMomentsCalculator::SpatialObjectPointer
typename SpatialObjectType::Pointer SpatialObjectPointer
Definition: itkImageMomentsCalculator.h:90
itk::ImageMomentsCalculator::SpatialObjectConstPointer
typename SpatialObjectType::ConstPointer SpatialObjectConstPointer
Definition: itkImageMomentsCalculator.h:91
itk::ImageMomentsCalculator::SetSpatialObjectMask
virtual void SetSpatialObjectMask(const SpatialObject< Self::ImageDimension > *so)
Definition: itkImageMomentsCalculator.h:122
itk::LightObject
Light weight base class for most itk classes.
Definition: itkLightObject.h:55
itk::ImageMomentsCalculator::ScalarType
double ScalarType
Definition: itkImageMomentsCalculator.h:81
itk::SpatialObject
Implementation of the composite pattern.
Definition: itkSpatialObject.h:58
itk::Matrix< ScalarType, Self::ImageDimension, Self::ImageDimension >
itk::ImageMomentsCalculator::SetImage
virtual void SetImage(const ImageType *image)
Definition: itkImageMomentsCalculator.h:109
itk::ImageMomentsCalculator::AffineTransformPointer
typename AffineTransformType::Pointer AffineTransformPointer
Definition: itkImageMomentsCalculator.h:105
itk
The "itk" namespace contains all Insight Segmentation and Registration Toolkit (ITK) classes....
Definition: itkAnnulusOperator.h:24
itk::ImageMomentsCalculator
Compute moments of an n-dimensional image.
Definition: itkImageMomentsCalculator.h:60
itk::Object
Base class for most ITK classes.
Definition: itkObject.h:61
itk::ImageMomentsCalculator::ImageType
TImage ImageType
Definition: itkImageMomentsCalculator.h:97