ITK  6.0.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
44namespace itk
45{
46
52inline constexpr double DefaultImageCoordinateTolerance = 1e-6;
53inline constexpr double DefaultImageDirectionTolerance = 1e-6;
54
55
113template <unsigned int VImageDimension = 2>
114class ITK_TEMPLATE_EXPORT ImageBase : public DataObject
115{
116public:
117 ITK_DISALLOW_COPY_AND_MOVE(ImageBase);
118
124
126 itkNewMacro(Self);
127
129 itkOverrideGetNameOfClassMacro(ImageBase);
130
132 using ImageDimensionType = unsigned int;
133
138 static constexpr ImageDimensionType ImageDimension = VImageDimension;
139
143
148
152
155
162
167
172
174 void
175 Initialize() override;
176
178 static constexpr unsigned int
180 {
181 return VImageDimension;
182 }
183
188 itkSetMacro(Origin, PointType);
189 virtual void
190 SetOrigin(const double origin[VImageDimension]);
191 virtual void
192 SetOrigin(const float origin[VImageDimension]);
220 virtual void
221 SetDirection(const DirectionType & direction);
222
226 itkGetConstReferenceMacro(Direction, DirectionType);
227
231 itkGetConstReferenceMacro(InverseDirection, DirectionType);
232
237 itkGetConstReferenceMacro(Spacing, SpacingType);
238
243 itkGetConstReferenceMacro(Origin, PointType);
244
252 virtual void
253 Allocate(bool initialize = false);
254
257 void
259 {
260 return this->Allocate(true);
261 }
262
269 virtual void
271
278 virtual const RegionType &
280 {
281 return m_LargestPossibleRegion;
282 }
283
287 virtual void
289
293 virtual const RegionType &
295 {
296 return m_BufferedRegion;
297 }
298
306 virtual void
308
316 void
317 SetRequestedRegion(const DataObject * data) override;
318
323 virtual const RegionType &
325 {
326 return m_RequestedRegion;
327 }
328
332 virtual void
333 SetRegions(const RegionType & region)
334 {
335 this->SetLargestPossibleRegion(region);
336 this->SetBufferedRegion(region);
337 this->SetRequestedRegion(region);
338 }
341 virtual void
342 SetRegions(const SizeType & size)
343 {
344 RegionType region;
345 region.SetSize(size);
346
347 this->Self::SetRegions(region);
348 }
349
360 const OffsetValueType *
362 {
363 return m_OffsetTable;
364 }
373 inline OffsetValueType
374 ComputeOffset(const IndexType & ind) const
375 {
376 OffsetValueType offset = 0;
377
379 this->GetBufferedRegion().GetIndex(), ind, m_OffsetTable, offset);
380 return offset;
381 /* NON TEMPLATE_META_PROGRAMMING_LOOP_UNROLLING data version
382 * Leaving here for documentation purposes
383 * OffsetValueType ComputeOffset(const IndexType & ind) const
384 * {
385 * // need to add bounds checking for the region/buffer?
386 * OffsetValueType offset = 0;
387 * const IndexType & bufferedRegionIndex = this->GetBufferedRegion().GetIndex();
388 * // data is arranged as [][][][slice][row][col]
389 * // with Index[0] = col, Index[1] = row, Index[2] = slice
390 * for ( int i = VImageDimension - 1; i > 0; i-- )
391 * {
392 * offset += ( ind[i] - bufferedRegionIndex[i] ) * m_OffsetTable[i];
393 * }
394 * offset += ( ind[0] - bufferedRegionIndex[0] );
395 * return offset;
396 * }
397 */
398 }
399
407 inline IndexType
409 {
410 IndexType index;
411 const IndexType & bufferedRegionIndex = this->GetBufferedRegion().GetIndex();
414 ImageHelper<VImageDimension, VImageDimension>::ComputeIndex(bufferedRegionIndex, offset, m_OffsetTable, index);
415 return index;
416 /* NON TEMPLATE_META_PROGRAMMING_LOOP_UNROLLING data version
417 * Leaving here for documentation purposes
418 * IndexType ComputeIndex(OffsetValueType offset) const
419 * {
420 * IndexType index;
421 * const IndexType & bufferedRegionIndex = this->GetBufferedRegion().GetIndex();
422 * for ( int i = VImageDimension - 1; i > 0; i-- )
423 * {
424 * index[i] = static_cast< IndexValueType >( offset / m_OffsetTable[i] );
425 * offset -= ( index[i] * m_OffsetTable[i] );
426 * index[i] += bufferedRegionIndex[i];
427 * }
428 * index[0] = bufferedRegionIndex[0] + static_cast< IndexValueType >( offset );
429 * return index;
430 * }
431 */
432 }
433
440 virtual void
441 SetSpacing(const SpacingType & spacing);
442 virtual void
443 SetSpacing(const double spacing[VImageDimension]);
444 virtual void
445 SetSpacing(const float spacing[VImageDimension]);
463 template <typename TCoordRep>
464 [[nodiscard]] IndexType
466 {
467 IndexType index;
470 for (unsigned int i = 0; i < VImageDimension; ++i)
471 {
472 TCoordRep sum{};
473 for (unsigned int j = 0; j < VImageDimension; ++j)
474 {
475 sum += this->m_PhysicalPointToIndex[i][j] * (point[j] - this->m_Origin[j]);
476 }
477 index[i] = Math::RoundHalfIntegerUp<IndexValueType>(sum);
478 }
479 return index;
480 }
481
490 template <typename TCoordRep>
491 ITK_NODISCARD("Call the overload which has the point as the only parameter and returns the index")
492 bool TransformPhysicalPointToIndex(const Point<TCoordRep, VImageDimension> & point, IndexType & index) const
493 {
494 index = TransformPhysicalPointToIndex(point);
495
496 // Now, check to see if the index is within allowed bounds
497 const bool isInside = this->GetLargestPossibleRegion().IsInside(index);
498 return isInside;
499 }
500
501
516 template <typename TIndexRep, typename TCoordRep>
519 {
522
523 for (unsigned int k = 0; k < VImageDimension; ++k)
524 {
525 cvector[k] = point[k] - this->m_Origin[k];
526 }
527 cvector = m_PhysicalPointToIndex * cvector;
528 for (unsigned int i = 0; i < VImageDimension; ++i)
529 {
530 index[i] = static_cast<TIndexRep>(cvector[i]);
531 }
532 return index;
533 }
534
543 template <typename TCoordRep, typename TIndexRep>
544 ITK_NODISCARD("Call the overload which has the point as the only parameter and returns the index")
545 bool TransformPhysicalPointToContinuousIndex(const Point<TCoordRep, VImageDimension> & point,
546 ContinuousIndex<TIndexRep, VImageDimension> & index) const
547 {
548 index = TransformPhysicalPointToContinuousIndex<TIndexRep>(point);
549
550 // Now, check to see if the index is within allowed bounds
551 const bool isInside = this->GetLargestPossibleRegion().IsInside(index);
552 return isInside;
553 }
554
559 template <typename TCoordRep, typename TIndexRep>
560 void
563 {
564 for (unsigned int r = 0; r < VImageDimension; ++r)
565 {
566 TCoordRep sum{};
567 for (unsigned int c = 0; c < VImageDimension; ++c)
568 {
569 sum += this->m_IndexToPhysicalPoint(r, c) * index[c];
570 }
571 point[r] = sum + this->m_Origin[r];
572 }
573 }
580 template <typename TCoordRep, typename TIndexRep>
583 {
585 TransformContinuousIndexToPhysicalPoint(index, point);
586 return point;
587 }
595 template <typename TCoordRep>
596 void
598 {
599 for (unsigned int i = 0; i < VImageDimension; ++i)
600 {
601 point[i] = this->m_Origin[i];
602 for (unsigned int j = 0; j < VImageDimension; ++j)
603 {
604 point[i] += m_IndexToPhysicalPoint[i][j] * index[j];
605 }
606 }
607 }
615 template <typename TCoordRep>
618 {
620 TransformIndexToPhysicalPoint(index, point);
621 return point;
622 }
639 template <typename TCoordRep>
640 void
642 FixedArray<TCoordRep, VImageDimension> & outputGradient) const
643 {
644 const DirectionType & direction = this->GetDirection();
645
646 itkAssertInDebugAndIgnoreInReleaseMacro(inputGradient.GetDataPointer() != outputGradient.GetDataPointer());
647
648 for (unsigned int i = 0; i < VImageDimension; ++i)
649 {
650 using CoordSumType = typename NumericTraits<TCoordRep>::AccumulateType;
651 CoordSumType sum{};
652 for (unsigned int j = 0; j < VImageDimension; ++j)
653 {
654 sum += direction[i][j] * inputGradient[j];
655 }
656 outputGradient[i] = static_cast<TCoordRep>(sum);
657 }
658 }
659
666 template <typename TVector>
667 [[nodiscard]] TVector
668 TransformLocalVectorToPhysicalVector(const TVector & inputGradient) const
669 {
670 TVector outputGradient;
671 TransformLocalVectorToPhysicalVector(inputGradient, outputGradient);
672 return outputGradient;
673 }
689 template <typename TCoordRep>
690 void
692 FixedArray<TCoordRep, VImageDimension> & outputGradient) const
693 {
694 const DirectionType & inverseDirection = this->GetInverseDirection();
695
696 itkAssertInDebugAndIgnoreInReleaseMacro(inputGradient.GetDataPointer() != outputGradient.GetDataPointer());
697
698 for (unsigned int i = 0; i < VImageDimension; ++i)
699 {
700 using CoordSumType = typename NumericTraits<TCoordRep>::AccumulateType;
701 CoordSumType sum{};
702 for (unsigned int j = 0; j < VImageDimension; ++j)
703 {
704 sum += inverseDirection[i][j] * inputGradient[j];
705 }
706 outputGradient[i] = static_cast<TCoordRep>(sum);
707 }
708 }
709
716 template <typename TVector>
717 [[nodiscard]] TVector
718 TransformPhysicalVectorToLocalVector(const TVector & inputGradient) const
719 {
720 TVector outputGradient;
721 TransformPhysicalVectorToLocalVector(inputGradient, outputGradient);
722 return outputGradient;
723 }
735 void
736 CopyInformation(const DataObject * data) override;
737
748 virtual void
749 Graft(const Self * image);
750
758 void
760
768 void
770
774 void
776
786 bool
788
797 bool
799
803 bool
804 IsCongruentImageGeometry(const ImageBase * otherImage, double coordinateTolerance, double directionTolerance) const;
805
813 bool
815 double coordinateTolerance = DefaultImageCoordinateTolerance,
816 double directionTolerance = DefaultImageDirectionTolerance) const;
817
836 virtual unsigned int
838 virtual void
842protected:
843 ImageBase() = default;
844 ~ImageBase() override = default;
845 void
846 PrintSelf(std::ostream & os, Indent indent) const override;
847
852 void
854
860 virtual void
862
863protected:
867 SpacingType m_Spacing{ MakeFilled<SpacingType>(1.0) };
868 PointType m_Origin{};
869 DirectionType m_Direction{ DirectionType::GetIdentity() };
870 DirectionType m_InverseDirection{ DirectionType::GetIdentity() };
875 DirectionType m_IndexToPhysicalPoint{ DirectionType::GetIdentity() };
876 DirectionType m_PhysicalPointToIndex{ DirectionType::GetIdentity() };
883 virtual void
885
897 FastComputeOffset(const IndexType & ind) const
898 {
899 OffsetValueType offset = 0;
901 Self::GetBufferedRegion().GetIndex(), ind, m_OffsetTable, offset);
902 return offset;
903 }
919 {
920 IndexType index;
921 const IndexType & bufferedRegionIndex = Self::GetBufferedRegion().GetIndex();
922 ImageHelper<VImageDimension, VImageDimension>::ComputeIndex(bufferedRegionIndex, offset, m_OffsetTable, index);
923 return index;
924 }
927 void
928 Graft(const DataObject * data) override;
929
930private:
931 OffsetValueType m_OffsetTable[VImageDimension + 1]{};
932
933 RegionType m_LargestPossibleRegion{};
934 RegionType m_RequestedRegion{};
935 RegionType m_BufferedRegion{};
936};
937} // end namespace itk
938
939#ifndef ITK_MANUAL_INSTANTIATION
940# include "itkImageBase.hxx"
941#endif
942
943#endif
A templated class holding a point in n-Dimensional image space.
Base class for all data objects in ITK.
Simulate a standard C array with copy semantics.
Definition: itkFixedArray.h:54
ValueType * GetDataPointer()
Base class for templated image classes.
Definition: itkImageBase.h:115
IndexType FastComputeIndex(OffsetValueType offset) const
Definition: itkImageBase.h:918
virtual void SetRequestedRegion(const RegionType &region)
IndexType ComputeIndex(OffsetValueType offset) const
Definition: itkImageBase.h:408
virtual void InitializeBufferedRegion()
void SetRequestedRegionToLargestPossibleRegion() override
bool VerifyRequestedRegion() override
SpacePrecisionType SpacingValueType
Definition: itkImageBase.h:160
void TransformContinuousIndexToPhysicalPoint(const ContinuousIndex< TIndexRep, VImageDimension > &index, Point< TCoordRep, VImageDimension > &point) const
Definition: itkImageBase.h:561
TVector TransformPhysicalVectorToLocalVector(const TVector &inputGradient) const
Definition: itkImageBase.h:718
virtual void SetDirection(const DirectionType &direction)
bool IsSameImageGeometryAs(const ImageBase *otherImage, double coordinateTolerance=DefaultImageCoordinateTolerance, double directionTolerance=DefaultImageDirectionTolerance) const
void AllocateInitialized()
Definition: itkImageBase.h:258
unsigned int ImageDimensionType
Definition: itkImageBase.h:132
void PrintSelf(std::ostream &os, Indent indent) const override
virtual const RegionType & GetBufferedRegion() const
Definition: itkImageBase.h:294
bool RequestedRegionIsOutsideOfTheBufferedRegion() override
virtual void SetRegions(const SizeType &size)
Definition: itkImageBase.h:342
virtual void SetBufferedRegion(const RegionType &region)
ContinuousIndex< TIndexRep, VImageDimension > TransformPhysicalPointToContinuousIndex(const Point< TCoordRep, VImageDimension > &point) const
Returns the continuous index from a physical point.
Definition: itkImageBase.h:518
virtual void SetSpacing(const SpacingType &spacing)
void UpdateOutputData() override
virtual void SetNumberOfComponentsPerPixel(unsigned int)
OffsetValueType FastComputeOffset(const IndexType &ind) const
Definition: itkImageBase.h:897
virtual void Allocate(bool initialize=false)
TVector TransformLocalVectorToPhysicalVector(const TVector &inputGradient) const
Definition: itkImageBase.h:668
virtual void SetOrigin(const float origin[VImageDimension])
void TransformLocalVectorToPhysicalVector(const FixedArray< TCoordRep, VImageDimension > &inputGradient, FixedArray< TCoordRep, VImageDimension > &outputGradient) const
Definition: itkImageBase.h:641
typename OffsetType::OffsetValueType OffsetValueType
Definition: itkImageBase.h:147
typename IndexType::IndexValueType IndexValueType
Definition: itkImageBase.h:142
void UpdateOutputInformation() override
virtual void Graft(const Self *image)
void TransformIndexToPhysicalPoint(const IndexType &index, Point< TCoordRep, VImageDimension > &point) const
Definition: itkImageBase.h:597
virtual unsigned int GetNumberOfComponentsPerPixel() const
virtual void SetSpacing(const float spacing[VImageDimension])
void TransformPhysicalVectorToLocalVector(const FixedArray< TCoordRep, VImageDimension > &inputGradient, FixedArray< TCoordRep, VImageDimension > &outputGradient) const
Definition: itkImageBase.h:691
const OffsetValueType * GetOffsetTable() const
Definition: itkImageBase.h:361
virtual const RegionType & GetRequestedRegion() const
Definition: itkImageBase.h:324
virtual void ComputeIndexToPhysicalPointMatrices()
Point< TCoordRep, VImageDimension > TransformContinuousIndexToPhysicalPoint(const ContinuousIndex< TIndexRep, VImageDimension > &index) const
Definition: itkImageBase.h:582
void SetRequestedRegion(const DataObject *data) override
OffsetValueType ComputeOffset(const IndexType &ind) const
Definition: itkImageBase.h:374
virtual void SetRegions(const RegionType &region)
Definition: itkImageBase.h:333
SpacePrecisionType PointValueType
Definition: itkImageBase.h:165
void CopyInformation(const DataObject *data) override
virtual void SetOrigin(const double origin[VImageDimension])
virtual void SetSpacing(const double spacing[VImageDimension])
ImageBase()=default
virtual void SetLargestPossibleRegion(const RegionType &region)
virtual const RegionType & GetLargestPossibleRegion() const
Definition: itkImageBase.h:279
Point< TCoordRep, VImageDimension > TransformIndexToPhysicalPoint(const IndexType &index) const
Definition: itkImageBase.h:617
void ComputeOffsetTable()
~ImageBase() override=default
IndexType TransformPhysicalPointToIndex(const Point< TCoordRep, VImageDimension > &point) const
Definition: itkImageBase.h:465
void Initialize() override
bool IsCongruentImageGeometry(const ImageBase *otherImage, double coordinateTolerance, double directionTolerance) const
typename SizeType::SizeValueType SizeValueType
Definition: itkImageBase.h:151
void Graft(const DataObject *data) override
static constexpr unsigned int GetImageDimension()
Definition: itkImageBase.h:179
static void ComputeIndex(const IndexType &bufferedRegionIndex, OffsetValueType offset, const OffsetValueType offsetTable[], IndexType &index)
static void ComputeOffset(const IndexType &bufferedRegionIndex, const IndexType &index, const OffsetValueType offsetTable[], OffsetValueType &offset)
An image region represents a structured region of data.
void SetSize(const SizeType &size)
Control indentation during Print() invocation.
Definition: itkIndent.h:50
Base class for most ITK classes.
Definition: itkObject.h:62
static constexpr double e
Definition: itkMath.h:56
The "itk" namespace contains all Insight Segmentation and Registration Toolkit (ITK) classes....
constexpr double DefaultImageDirectionTolerance
Definition: itkImageBase.h:53
class ITK_FORWARD_EXPORT DataObject
Definition: itkDataObject.h:42
constexpr double DefaultImageCoordinateTolerance
Definition: itkImageBase.h:52
*par Constraints *The filter 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
long IndexValueType
Definition: itkIntTypes.h:93
double SpacePrecisionType
Definition: itkFloatTypes.h:30
class ITK_TEMPLATE_EXPORT ImageBase
unsigned long SizeValueType
Definition: itkIntTypes.h:86
long OffsetValueType
Definition: itkIntTypes.h:97
const IndexValueType * GetIndex() const
Definition: itkIndex.h:230
Represent a n-dimensional offset between two n-dimensional indexes of n-dimensional image.
Definition: itkOffset.h:67