ITK  5.4.0
Insight Toolkit
itkImageBase.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 /*=========================================================================
19  *
20  * Portions of this file are subject to the VTK Toolkit Version 3 copyright.
21  *
22  * Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
23  *
24  * For complete copyright, license and disclaimer of warranty information
25  * please refer to the NOTICE file at the top of the ITK source tree.
26  *
27  *=========================================================================*/
28 #ifndef itkImageBase_h
29 #define itkImageBase_h
30 
31 #include "itkDataObject.h"
32 
33 #include "itkImageRegion.h"
34 #include "itkMatrix.h"
35 #include "itkObjectFactory.h"
36 #include "itkOffset.h"
37 #include "itkFixedArray.h"
38 #include "itkImageHelper.h"
39 #include "itkFloatTypes.h"
40 
41 #include <vxl_version.h>
42 #include "vnl/vnl_matrix_fixed.hxx" // Get the templates
43 
44 namespace itk
45 {
46 
104 template <unsigned int VImageDimension = 2>
105 class ITK_TEMPLATE_EXPORT ImageBase : public DataObject
106 {
107 public:
108  ITK_DISALLOW_COPY_AND_MOVE(ImageBase);
109 
111  using Self = ImageBase;
115 
117  itkNewMacro(Self);
118 
120  itkTypeMacro(ImageBase, DataObject);
121 
123  using ImageDimensionType = unsigned int;
124 
129  static constexpr ImageDimensionType ImageDimension = VImageDimension;
130 
134 
139 
143 
146 
153 
158 
163 
165  void
166  Initialize() override;
167 
169  static unsigned int
171  {
172  return VImageDimension;
173  }
174 
179  itkSetMacro(Origin, PointType);
180  virtual void
181  SetOrigin(const double origin[VImageDimension]);
182  virtual void
183  SetOrigin(const float origin[VImageDimension]);
211  virtual void
212  SetDirection(const DirectionType & direction);
213 
217  itkGetConstReferenceMacro(Direction, DirectionType);
218 
222  itkGetConstReferenceMacro(InverseDirection, DirectionType);
223 
228  itkGetConstReferenceMacro(Spacing, SpacingType);
229 
234  itkGetConstReferenceMacro(Origin, PointType);
235 
243  virtual void
244  Allocate(bool initialize = false);
245 
252  virtual void
253  SetLargestPossibleRegion(const RegionType & region);
254 
261  virtual const RegionType &
263  {
264  return m_LargestPossibleRegion;
265  }
266 
270  virtual void
271  SetBufferedRegion(const RegionType & region);
272 
276  virtual const RegionType &
278  {
279  return m_BufferedRegion;
280  }
281 
289  virtual void
290  SetRequestedRegion(const RegionType & region);
291 
299  void
300  SetRequestedRegion(const DataObject * data) override;
301 
306  virtual const RegionType &
308  {
309  return m_RequestedRegion;
310  }
311 
315  virtual void
316  SetRegions(const RegionType & region)
317  {
318  this->SetLargestPossibleRegion(region);
319  this->SetBufferedRegion(region);
320  this->SetRequestedRegion(region);
321  }
324  virtual void
325  SetRegions(const SizeType & size)
326  {
327  RegionType region;
328  region.SetSize(size);
329 
330  this->Self::SetRegions(region);
331  }
332 
343  const OffsetValueType *
345  {
346  return m_OffsetTable;
347  }
356  inline OffsetValueType
357  ComputeOffset(const IndexType & ind) const
358  {
359  OffsetValueType offset = 0;
360 
362  this->GetBufferedRegion().GetIndex(), ind, m_OffsetTable, offset);
363  return offset;
364  /* NON TEMPLATE_META_PROGRAMMING_LOOP_UNROLLING data version
365  * Leaving here for documentation purposes
366  * OffsetValueType ComputeOffset(const IndexType & ind) const
367  * {
368  * // need to add bounds checking for the region/buffer?
369  * OffsetValueType offset = 0;
370  * const IndexType & bufferedRegionIndex = this->GetBufferedRegion().GetIndex();
371  * // data is arranged as [][][][slice][row][col]
372  * // with Index[0] = col, Index[1] = row, Index[2] = slice
373  * for ( int i = VImageDimension - 1; i > 0; i-- )
374  * {
375  * offset += ( ind[i] - bufferedRegionIndex[i] ) * m_OffsetTable[i];
376  * }
377  * offset += ( ind[0] - bufferedRegionIndex[0] );
378  * return offset;
379  * }
380  */
381  }
382 
390  inline IndexType
392  {
393  IndexType index;
394  const IndexType & bufferedRegionIndex = this->GetBufferedRegion().GetIndex();
397  ImageHelper<VImageDimension, VImageDimension>::ComputeIndex(bufferedRegionIndex, offset, m_OffsetTable, index);
398  return index;
399  /* NON TEMPLATE_META_PROGRAMMING_LOOP_UNROLLING data version
400  * Leaving here for documentation purposes
401  * IndexType ComputeIndex(OffsetValueType offset) const
402  * {
403  * IndexType index;
404  * const IndexType & bufferedRegionIndex = this->GetBufferedRegion().GetIndex();
405  * for ( int i = VImageDimension - 1; i > 0; i-- )
406  * {
407  * index[i] = static_cast< IndexValueType >( offset / m_OffsetTable[i] );
408  * offset -= ( index[i] * m_OffsetTable[i] );
409  * index[i] += bufferedRegionIndex[i];
410  * }
411  * index[0] = bufferedRegionIndex[0] + static_cast< IndexValueType >( offset );
412  * return index;
413  * }
414  */
415  }
416 
423  virtual void
424  SetSpacing(const SpacingType & spacing);
425  virtual void
426  SetSpacing(const double spacing[VImageDimension]);
427  virtual void
428  SetSpacing(const float spacing[VImageDimension]);
446  template <typename TCoordRep>
447  [[nodiscard]] IndexType
449  {
450  IndexType index;
453  for (unsigned int i = 0; i < VImageDimension; ++i)
454  {
455  TCoordRep sum{};
456  for (unsigned int j = 0; j < VImageDimension; ++j)
457  {
458  sum += this->m_PhysicalPointToIndex[i][j] * (point[j] - this->m_Origin[j]);
459  }
460  index[i] = Math::RoundHalfIntegerUp<IndexValueType>(sum);
461  }
462  return index;
463  }
464 
473  template <typename TCoordRep>
474  ITK_NODISCARD("Call the overload which has the point as the only parameter and returns the index")
475  bool TransformPhysicalPointToIndex(const Point<TCoordRep, VImageDimension> & point, IndexType & index) const
476  {
477  index = TransformPhysicalPointToIndex(point);
478 
479  // Now, check to see if the index is within allowed bounds
480  const bool isInside = this->GetLargestPossibleRegion().IsInside(index);
481  return isInside;
482  }
483 
484 
499  template <typename TIndexRep, typename TCoordRep>
502  {
505 
506  for (unsigned int k = 0; k < VImageDimension; ++k)
507  {
508  cvector[k] = point[k] - this->m_Origin[k];
509  }
510  cvector = m_PhysicalPointToIndex * cvector;
511  for (unsigned int i = 0; i < VImageDimension; ++i)
512  {
513  index[i] = static_cast<TIndexRep>(cvector[i]);
514  }
515  return index;
516  }
517 
526  template <typename TCoordRep, typename TIndexRep>
527  ITK_NODISCARD("Call the overload which has the point as the only parameter and returns the index")
528  bool TransformPhysicalPointToContinuousIndex(const Point<TCoordRep, VImageDimension> & point,
529  ContinuousIndex<TIndexRep, VImageDimension> & index) const
530  {
531  index = TransformPhysicalPointToContinuousIndex<TIndexRep>(point);
532 
533  // Now, check to see if the index is within allowed bounds
534  const bool isInside = this->GetLargestPossibleRegion().IsInside(index);
535  return isInside;
536  }
537 
542  template <typename TCoordRep, typename TIndexRep>
543  void
546  {
547  for (unsigned int r = 0; r < VImageDimension; ++r)
548  {
549  TCoordRep sum{};
550  for (unsigned int c = 0; c < VImageDimension; ++c)
551  {
552  sum += this->m_IndexToPhysicalPoint(r, c) * index[c];
553  }
554  point[r] = sum + this->m_Origin[r];
555  }
556  }
563  template <typename TCoordRep, typename TIndexRep>
566  {
568  TransformContinuousIndexToPhysicalPoint(index, point);
569  return point;
570  }
578  template <typename TCoordRep>
579  void
581  {
582  for (unsigned int i = 0; i < VImageDimension; ++i)
583  {
584  point[i] = this->m_Origin[i];
585  for (unsigned int j = 0; j < VImageDimension; ++j)
586  {
587  point[i] += m_IndexToPhysicalPoint[i][j] * index[j];
588  }
589  }
590  }
598  template <typename TCoordRep>
601  {
603  TransformIndexToPhysicalPoint(index, point);
604  return point;
605  }
622  template <typename TCoordRep>
623  void
625  FixedArray<TCoordRep, VImageDimension> & outputGradient) const
626  {
627  const DirectionType & direction = this->GetDirection();
628 
629  itkAssertInDebugAndIgnoreInReleaseMacro(inputGradient.GetDataPointer() != outputGradient.GetDataPointer());
630 
631  for (unsigned int i = 0; i < VImageDimension; ++i)
632  {
633  using CoordSumType = typename NumericTraits<TCoordRep>::AccumulateType;
634  CoordSumType sum{};
635  for (unsigned int j = 0; j < VImageDimension; ++j)
636  {
637  sum += direction[i][j] * inputGradient[j];
638  }
639  outputGradient[i] = static_cast<TCoordRep>(sum);
640  }
641  }
642 
649  template <typename TVector>
650  [[nodiscard]] TVector
651  TransformLocalVectorToPhysicalVector(const TVector & inputGradient) const
652  {
653  TVector outputGradient;
654  TransformLocalVectorToPhysicalVector(inputGradient, outputGradient);
655  return outputGradient;
656  }
672  template <typename TCoordRep>
673  void
675  FixedArray<TCoordRep, VImageDimension> & outputGradient) const
676  {
677  const DirectionType & inverseDirection = this->GetInverseDirection();
678 
679  itkAssertInDebugAndIgnoreInReleaseMacro(inputGradient.GetDataPointer() != outputGradient.GetDataPointer());
680 
681  for (unsigned int i = 0; i < VImageDimension; ++i)
682  {
683  using CoordSumType = typename NumericTraits<TCoordRep>::AccumulateType;
684  CoordSumType sum{};
685  for (unsigned int j = 0; j < VImageDimension; ++j)
686  {
687  sum += inverseDirection[i][j] * inputGradient[j];
688  }
689  outputGradient[i] = static_cast<TCoordRep>(sum);
690  }
691  }
692 
699  template <typename TVector>
700  [[nodiscard]] TVector
701  TransformPhysicalVectorToLocalVector(const TVector & inputGradient) const
702  {
703  TVector outputGradient;
704  TransformPhysicalVectorToLocalVector(inputGradient, outputGradient);
705  return outputGradient;
706  }
718  void
719  CopyInformation(const DataObject * data) override;
720 
731  virtual void
732  Graft(const Self * image);
733 
741  void
742  UpdateOutputInformation() override;
743 
751  void
752  UpdateOutputData() override;
753 
757  void
758  SetRequestedRegionToLargestPossibleRegion() override;
759 
769  bool
770  RequestedRegionIsOutsideOfTheBufferedRegion() override;
771 
780  bool
781  VerifyRequestedRegion() override;
782 
801  virtual unsigned int
802  GetNumberOfComponentsPerPixel() const;
803  virtual void
804  SetNumberOfComponentsPerPixel(unsigned int);
807 protected:
808  ImageBase() = default;
809  ~ImageBase() override = default;
810  void
811  PrintSelf(std::ostream & os, Indent indent) const override;
812 
817  void
818  ComputeOffsetTable();
819 
825  virtual void
826  ComputeIndexToPhysicalPointMatrices();
827 
828 protected:
832  SpacingType m_Spacing{ MakeFilled<SpacingType>(1.0) };
833  PointType m_Origin{};
834  DirectionType m_Direction{ DirectionType::GetIdentity() };
835  DirectionType m_InverseDirection{ DirectionType::GetIdentity() };
840  DirectionType m_IndexToPhysicalPoint{ DirectionType::GetIdentity() };
841  DirectionType m_PhysicalPointToIndex{ DirectionType::GetIdentity() };
848  virtual void
849  InitializeBufferedRegion();
850 
862  FastComputeOffset(const IndexType & ind) const
863  {
864  OffsetValueType offset = 0;
866  Self::GetBufferedRegion().GetIndex(), ind, m_OffsetTable, offset);
867  return offset;
868  }
882  IndexType
884  {
885  IndexType index;
886  const IndexType & bufferedRegionIndex = Self::GetBufferedRegion().GetIndex();
887  ImageHelper<VImageDimension, VImageDimension>::ComputeIndex(bufferedRegionIndex, offset, m_OffsetTable, index);
888  return index;
889  }
892  void
893  Graft(const DataObject * data) override;
894 
895 private:
896  OffsetValueType m_OffsetTable[VImageDimension + 1]{};
897 
898  RegionType m_LargestPossibleRegion{};
899  RegionType m_RequestedRegion{};
900  RegionType m_BufferedRegion{};
901 };
902 } // end namespace itk
903 
904 #ifndef ITK_MANUAL_INSTANTIATION
905 # include "itkImageBase.hxx"
906 #endif
907 
908 #endif
itk::ImageBase::TransformPhysicalPointToContinuousIndex
ContinuousIndex< TIndexRep, VImageDimension > TransformPhysicalPointToContinuousIndex(const Point< TCoordRep, VImageDimension > &point) const
Returns the continuous index from a physical point.
Definition: itkImageBase.h:501
itk::ImageBase< TImage::ImageDimension >::OffsetValueType
typename OffsetType::OffsetValueType OffsetValueType
Definition: itkImageBase.h:138
itk::ImageBase::GetBufferedRegion
virtual const RegionType & GetBufferedRegion() const
Definition: itkImageBase.h:277
itk::ImageBase::GetImageDimension
static unsigned int GetImageDimension()
Definition: itkImageBase.h:170
itkObjectFactory.h
itk::ImageBase::TransformLocalVectorToPhysicalVector
TVector TransformLocalVectorToPhysicalVector(const TVector &inputGradient) const
Definition: itkImageBase.h:651
itk::Index< VImageDimension >
itk::ImageBase::SetRegions
virtual void SetRegions(const SizeType &size)
Definition: itkImageBase.h:325
itk::GTest::TypedefsAndConstructors::Dimension2::DirectionType
ImageBaseType::DirectionType DirectionType
Definition: itkGTestTypedefsAndConstructors.h:52
itk::ImageHelper::ComputeIndex
static void ComputeIndex(const IndexType &bufferedRegionIndex, OffsetValueType offset, [[maybe_unused]] const OffsetValueType offsetTable[], IndexType &index)
Definition: itkImageHelper.h:64
itk::Size
Represent a n-dimensional size (bounds) of a n-dimensional image.
Definition: itkSize.h:71
itkMatrix.h
itk::ImageBase::GetOffsetTable
const OffsetValueType * GetOffsetTable() const
Definition: itkImageBase.h:344
itkOffset.h
itkFloatTypes.h
itk::ImageBase
Base class for templated image classes.
Definition: itkImageBase.h:105
itkImageHelper.h
itk::ImageRegion
An image region represents a structured region of data.
Definition: itkImageRegion.h:69
itk::ImageBase< TImage::ImageDimension >::SizeValueType
typename SizeType::SizeValueType SizeValueType
Definition: itkImageBase.h:142
itk::ImageBase< TImage::ImageDimension >::PointValueType
SpacePrecisionType PointValueType
Definition: itkImageBase.h:156
itk::GTest::TypedefsAndConstructors::Dimension2::PointType
ImageBaseType::PointType PointType
Definition: itkGTestTypedefsAndConstructors.h:51
itk::Vector< SpacingValueType, VImageDimension >
itk::NumericTraits::AccumulateType
double AccumulateType
Definition: itkNumericTraits.h:75
itk::ImageBase::TransformIndexToPhysicalPoint
void TransformIndexToPhysicalPoint(const IndexType &index, Point< TCoordRep, VImageDimension > &point) const
Definition: itkImageBase.h:580
itk::SmartPointer< Self >
itk::ImageBase
class ITK_TEMPLATE_EXPORT ImageBase
Definition: itkImageHelper.h:50
itk::Indent
Control indentation during Print() invocation.
Definition: itkIndent.h:49
itkDataObject.h
itk::ImageBase< TImage::ImageDimension >::IndexValueType
typename IndexType::IndexValueType IndexValueType
Definition: itkImageBase.h:133
itk::FixedArray::GetDataPointer
ValueType * GetDataPointer()
Definition: itkFixedArray.h:295
itk::ImageBase::TransformContinuousIndexToPhysicalPoint
void TransformContinuousIndexToPhysicalPoint(const ContinuousIndex< TIndexRep, VImageDimension > &index, Point< TCoordRep, VImageDimension > &point) const
Definition: itkImageBase.h:544
itkImageRegion.h
itk::IndexValueType
long IndexValueType
Definition: itkIntTypes.h:90
itk::GTest::TypedefsAndConstructors::Dimension2::IndexType
ImageBaseType::IndexType IndexType
Definition: itkGTestTypedefsAndConstructors.h:50
itk::ImageBase::FastComputeIndex
IndexType FastComputeIndex(OffsetValueType offset) const
Definition: itkImageBase.h:883
itk::Index::GetIndex
const IndexValueType * GetIndex() const
Definition: itkIndex.h:232
itk::ImageBase::ComputeOffset
OffsetValueType ComputeOffset(const IndexType &ind) const
Definition: itkImageBase.h:357
itk::point
*par Constraints *The filter requires an image with at least two dimensions and a vector *length of at least The theory supports extension to scalar but *the implementation of the itk vector classes do not **The template parameter TRealType must be floating point(float or double) or *a user-defined "real" numerical type with arithmetic operations defined *sufficient to compute derivatives. **\par Performance *This filter will automatically multithread if run with *SetUsePrincipleComponents
itk::ImageBase::TransformContinuousIndexToPhysicalPoint
Point< TCoordRep, VImageDimension > TransformContinuousIndexToPhysicalPoint(const ContinuousIndex< TIndexRep, VImageDimension > &index) const
Definition: itkImageBase.h:565
itk::GTest::TypedefsAndConstructors::Dimension2::RegionType
ImageBaseType::RegionType RegionType
Definition: itkGTestTypedefsAndConstructors.h:54
itk::SpacePrecisionType
double SpacePrecisionType
Definition: itkFloatTypes.h:30
itkFixedArray.h
itk::ImageBase::TransformPhysicalVectorToLocalVector
TVector TransformPhysicalVectorToLocalVector(const TVector &inputGradient) const
Definition: itkImageBase.h:701
itk::ImageHelper::ComputeOffset
static void ComputeOffset(const IndexType &bufferedRegionIndex, const IndexType &index, [[maybe_unused]] const OffsetValueType offsetTable[], OffsetValueType &offset)
Definition: itkImageHelper.h:92
itk::OffsetValueType
long OffsetValueType
Definition: itkIntTypes.h:94
itk::ImageBase::SetRegions
virtual void SetRegions(const RegionType &region)
Definition: itkImageBase.h:316
itk::ImageBase::TransformPhysicalVectorToLocalVector
void TransformPhysicalVectorToLocalVector(const FixedArray< TCoordRep, VImageDimension > &inputGradient, FixedArray< TCoordRep, VImageDimension > &outputGradient) const
Definition: itkImageBase.h:674
itk::DataObject
class ITK_FORWARD_EXPORT DataObject
Definition: itkDataObject.h:42
itk::FixedArray
Simulate a standard C array with copy semantics.
Definition: itkFixedArray.h:53
itk::Matrix< SpacePrecisionType, VImageDimension, VImageDimension >
itk::ImageBase::GetLargestPossibleRegion
virtual const RegionType & GetLargestPossibleRegion() const
Definition: itkImageBase.h:262
itk::Offset
Represent a n-dimensional offset between two n-dimensional indexes of n-dimensional image.
Definition: itkOffset.h:69
itk::ImageBase::GetRequestedRegion
virtual const RegionType & GetRequestedRegion() const
Definition: itkImageBase.h:307
itk::ImageBase< TImage::ImageDimension >::SpacingValueType
SpacePrecisionType SpacingValueType
Definition: itkImageBase.h:151
itk::ImageBase::TransformPhysicalPointToIndex
IndexType TransformPhysicalPointToIndex(const Point< TCoordRep, VImageDimension > &point) const
Definition: itkImageBase.h:448
itk
The "itk" namespace contains all Insight Segmentation and Registration Toolkit (ITK) classes....
Definition: itkAnnulusOperator.h:24
itk::ContinuousIndex
A templated class holding a point in n-Dimensional image space.
Definition: itkContinuousIndex.h:46
itk::ImageBase::TransformLocalVectorToPhysicalVector
void TransformLocalVectorToPhysicalVector(const FixedArray< TCoordRep, VImageDimension > &inputGradient, FixedArray< TCoordRep, VImageDimension > &outputGradient) const
Definition: itkImageBase.h:624
itk::Object
Base class for most ITK classes.
Definition: itkObject.h:61
itk::ImageBase::TransformIndexToPhysicalPoint
Point< TCoordRep, VImageDimension > TransformIndexToPhysicalPoint(const IndexType &index) const
Definition: itkImageBase.h:600
itk::ImageBase::FastComputeOffset
OffsetValueType FastComputeOffset(const IndexType &ind) const
Definition: itkImageBase.h:862
itk::Point< PointValueType, VImageDimension >
itk::ImageBase::ComputeIndex
IndexType ComputeIndex(OffsetValueType offset) const
Definition: itkImageBase.h:391
AddImageFilter
Definition: itkAddImageFilter.h:81
itk::ImageRegion::SetSize
void SetSize(const SizeType &size)
Definition: itkImageRegion.h:174
itk::SizeValueType
unsigned long SizeValueType
Definition: itkIntTypes.h:83
itk::DataObject
Base class for all data objects in ITK.
Definition: itkDataObject.h:293