SimpleITK  
sitkImage.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  * http://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 sitkImage_h
19 #define sitkImage_h
20 
21 #include "sitkCommon.h"
22 #include "sitkTemplateFunctions.h"
23 #include "sitkDetail.h"
24 #include "sitkPixelIDTokens.h"
25 #include "sitkInterpolator.h"
26 
27 #include <vector>
28 #include <memory>
29 #include <type_traits>
30 
31 namespace itk
32 {
33 
34 // Forward declaration for pointer
35 class DataObject;
36 
37 template <class T>
38 class SmartPointer;
39 
40 namespace simple
41 {
42 
43 // This is the forward declaration of a class used internally to the
44 // Image class, but the actual interface is not exposed to simple
45 // ITK. A pointer to the implementation is used as per the pimple
46 // idiom.
47 class PimpleImageBase;
48 
77 {
78 public:
79  using Self = Image;
80 
81  virtual ~Image();
82 
84  Image();
85 
86  // copy constructor
87  Image(const Image & img);
88 
89  Image &
90  operator=(const Image & img);
91 
92 #ifndef SWIG
93 
99  Image(Image && img) noexcept;
100  Image &
101  operator=(Image && img) noexcept;
102 
122  Image
123  ProxyForInPlaceOperation();
124 #endif
125 
126 
142  Image(unsigned int width, unsigned int height, PixelIDValueEnum valueEnum);
143  Image(unsigned int width, unsigned int height, unsigned int depth, PixelIDValueEnum valueEnum);
144  Image(const std::vector<unsigned int> & size, PixelIDValueEnum valueEnum, unsigned int numberOfComponents = 0);
164  template <typename TImageType>
166  : Image(image.GetPointer())
167  {}
168 
169  template <typename TImageType>
170  explicit Image(TImageType * image)
171  : Image()
172  {
174  const unsigned int dimension = TImageType::ImageDimension;
175 
176  static_assert(type != sitkUnknown, "invalid pixel type");
177  static_assert(dimension >= 2 && dimension <= SITK_MAX_DIMENSION, "Unsupported image dimension.");
178 
179  this->InternalInitialization(type, dimension, image);
180  }
198  GetITKBase();
199  const itk::DataObject *
200  GetITKBase() const;
210  GetPixelID() const;
212  GetPixelIDValue() const;
213 
215  std::string
216  GetPixelIDTypeAsString() const;
217 
225  unsigned int
226  GetDimension() const;
227 
234  unsigned int
235  GetNumberOfComponentsPerPixel() const;
236 
245  uint64_t
246  GetNumberOfPixels() const;
247 
252  unsigned int
253  GetSizeOfPixelComponent() const;
254 
258  std::vector<double>
259  GetOrigin() const;
260  void
261  SetOrigin(const std::vector<double> & origin);
270  std::vector<double>
271  GetSpacing() const;
272  void
273  SetSpacing(const std::vector<double> & spacing);
283  std::vector<double>
284  GetDirection() const;
285  void
286  SetDirection(const std::vector<double> & direction);
290  std::vector<double>
291  TransformIndexToPhysicalPoint(const std::vector<int64_t> & index) const;
292 
294  std::vector<int64_t>
295  TransformPhysicalPointToIndex(const std::vector<double> & point) const;
296 
298  std::vector<double>
299  TransformPhysicalPointToContinuousIndex(const std::vector<double> & point) const;
300 
302  std::vector<double>
303  TransformContinuousIndexToPhysicalPoint(const std::vector<double> & index) const;
304 
319  std::vector<double>
320  EvaluateAtContinuousIndex(const std::vector<double> & index, InterpolatorEnum interp = sitkLinear) const;
321 
335  std::vector<double>
336  EvaluateAtPhysicalPoint(const std::vector<double> & point, InterpolatorEnum interp = sitkLinear) const;
337 
338 
346  bool
347  IsCongruentImageGeometry(const Image & otherImage, double coordinateTolerance, double directionTolerance) const;
348 
356  bool
357  IsSameImageGeometryAs(const Image & otherImage,
359  double = DefaultImageDirectionTolerance) const;
360 
364  std::vector<unsigned int>
365  GetSize() const;
366 
368  unsigned int
369  GetWidth() const;
370 
372  unsigned int
373  GetHeight() const;
374 
377  unsigned int
378  GetDepth() const;
379 
380 
391  void
392  CopyInformation(const Image & srcImage);
393 
400  std::vector<std::string>
401  GetMetaDataKeys() const;
402 
405  bool
406  HasMetaDataKey(const std::string & key) const;
407 
416  std::string
417  GetMetaData(const std::string & key) const;
418 
423  void
424  SetMetaData(const std::string & key, const std::string & value);
425 
431  bool
432  EraseMetaData(const std::string & key);
433 
434  std::string
435  ToString() const;
436 
454  Image
455  ToVectorImage(bool inPlace = true);
456 
475  Image
476  ToScalarImage(bool inPlace = true);
477 
493  int8_t
494  GetPixelAsInt8(const std::vector<uint32_t> & idx) const;
495  uint8_t
496  GetPixelAsUInt8(const std::vector<uint32_t> & idx) const;
497  int16_t
498  GetPixelAsInt16(const std::vector<uint32_t> & idx) const;
499  uint16_t
500  GetPixelAsUInt16(const std::vector<uint32_t> & idx) const;
501  int32_t
502  GetPixelAsInt32(const std::vector<uint32_t> & idx) const;
503  uint32_t
504  GetPixelAsUInt32(const std::vector<uint32_t> & idx) const;
505  int64_t
506  GetPixelAsInt64(const std::vector<uint32_t> & idx) const;
507  uint64_t
508  GetPixelAsUInt64(const std::vector<uint32_t> & idx) const;
509  float
510  GetPixelAsFloat(const std::vector<uint32_t> & idx) const;
511  double
512  GetPixelAsDouble(const std::vector<uint32_t> & idx) const;
513 
514  std::vector<int8_t>
515  GetPixelAsVectorInt8(const std::vector<uint32_t> & idx) const;
516  std::vector<uint8_t>
517  GetPixelAsVectorUInt8(const std::vector<uint32_t> & idx) const;
518  std::vector<int16_t>
519  GetPixelAsVectorInt16(const std::vector<uint32_t> & idx) const;
520  std::vector<uint16_t>
521  GetPixelAsVectorUInt16(const std::vector<uint32_t> & idx) const;
522  std::vector<int32_t>
523  GetPixelAsVectorInt32(const std::vector<uint32_t> & idx) const;
524  std::vector<uint32_t>
525  GetPixelAsVectorUInt32(const std::vector<uint32_t> & idx) const;
526  std::vector<int64_t>
527  GetPixelAsVectorInt64(const std::vector<uint32_t> & idx) const;
528  std::vector<uint64_t>
529  GetPixelAsVectorUInt64(const std::vector<uint32_t> & idx) const;
530  std::vector<float>
531  GetPixelAsVectorFloat32(const std::vector<uint32_t> & idx) const;
532  std::vector<double>
533  GetPixelAsVectorFloat64(const std::vector<uint32_t> & idx) const;
534 
535  std::complex<float>
536  GetPixelAsComplexFloat32(const std::vector<uint32_t> & idx) const;
537  std::complex<double>
538  GetPixelAsComplexFloat64(const std::vector<uint32_t> & idx) const;
557  void
558  SetPixelAsInt8(const std::vector<uint32_t> & idx, int8_t v);
559  void
560  SetPixelAsUInt8(const std::vector<uint32_t> & idx, uint8_t v);
561  void
562  SetPixelAsInt16(const std::vector<uint32_t> & idx, int16_t v);
563  void
564  SetPixelAsUInt16(const std::vector<uint32_t> & idx, uint16_t v);
565  void
566  SetPixelAsInt32(const std::vector<uint32_t> & idx, int32_t v);
567  void
568  SetPixelAsUInt32(const std::vector<uint32_t> & idx, uint32_t v);
569  void
570  SetPixelAsInt64(const std::vector<uint32_t> & idx, int64_t v);
571  void
572  SetPixelAsUInt64(const std::vector<uint32_t> & idx, uint64_t v);
573  void
574  SetPixelAsFloat(const std::vector<uint32_t> & idx, float v);
575  void
576  SetPixelAsDouble(const std::vector<uint32_t> & idx, double v);
577 
578  void
579  SetPixelAsVectorInt8(const std::vector<uint32_t> & idx, const std::vector<int8_t> & v);
580  void
581  SetPixelAsVectorUInt8(const std::vector<uint32_t> & idx, const std::vector<uint8_t> & v);
582  void
583  SetPixelAsVectorInt16(const std::vector<uint32_t> & idx, const std::vector<int16_t> & v);
584  void
585  SetPixelAsVectorUInt16(const std::vector<uint32_t> & idx, const std::vector<uint16_t> & v);
586  void
587  SetPixelAsVectorInt32(const std::vector<uint32_t> & idx, const std::vector<int32_t> & v);
588  void
589  SetPixelAsVectorUInt32(const std::vector<uint32_t> & idx, const std::vector<uint32_t> & v);
590  void
591  SetPixelAsVectorInt64(const std::vector<uint32_t> & idx, const std::vector<int64_t> & v);
592  void
593  SetPixelAsVectorUInt64(const std::vector<uint32_t> & idx, const std::vector<uint64_t> & v);
594  void
595  SetPixelAsVectorFloat32(const std::vector<uint32_t> & idx, const std::vector<float> & v);
596  void
597  SetPixelAsVectorFloat64(const std::vector<uint32_t> & idx, const std::vector<double> & v);
598 
599  void
600  SetPixelAsComplexFloat32(const std::vector<uint32_t> & idx, const std::complex<float> v);
601  void
602  SetPixelAsComplexFloat64(const std::vector<uint32_t> & idx, const std::complex<double> v);
603 
632  int8_t *
633  GetBufferAsInt8();
634  uint8_t *
635  GetBufferAsUInt8();
636  int16_t *
637  GetBufferAsInt16();
638  uint16_t *
639  GetBufferAsUInt16();
640  int32_t *
641  GetBufferAsInt32();
642  uint32_t *
643  GetBufferAsUInt32();
644  int64_t *
645  GetBufferAsInt64();
646  uint64_t *
647  GetBufferAsUInt64();
648  float *
649  GetBufferAsFloat();
650  double *
651  GetBufferAsDouble();
652  void *
653  GetBufferAsVoid();
654 
655  const int8_t *
656  GetBufferAsInt8() const;
657  const uint8_t *
658  GetBufferAsUInt8() const;
659  const int16_t *
660  GetBufferAsInt16() const;
661  const uint16_t *
662  GetBufferAsUInt16() const;
663  const int32_t *
664  GetBufferAsInt32() const;
665  const uint32_t *
666  GetBufferAsUInt32() const;
667  const int64_t *
668  GetBufferAsInt64() const;
669  const uint64_t *
670  GetBufferAsUInt64() const;
671  const float *
672  GetBufferAsFloat() const;
673  const double *
674  GetBufferAsDouble() const;
675  const void *
676  GetBufferAsVoid() const;
686  void
687  MakeUnique();
688 
691  bool
692  IsUnique() const;
693 
694  static constexpr double DefaultImageCoordinateTolerance = 1e-6;
695  static constexpr double DefaultImageDirectionTolerance = 1e-6;
696 
697 protected:
704  void
705  Allocate(const std::vector<unsigned int> & size, PixelIDValueEnum valueEnum, unsigned int numberOfComponents);
706 
711  template <class TImageType>
712  void
713  AllocateInternal(const std::vector<unsigned int> & size, unsigned int numberOfComponents);
719  template <class TImageType>
720  Image
721  ToVectorInternal(bool inPlace);
722 
723 
724  template <class TImageType>
725  Image
726  ToScalarInternal(bool inPlace);
727 
730 private:
739  void
740  InternalInitialization(PixelIDValueType type, unsigned int dimension, itk::DataObject * image);
741 
742 
743  Image(std::unique_ptr<PimpleImageBase> pimpleImage);
744 
745 
746  template <typename TImageType>
748  DispatchedInternalInitialization(itk::DataObject * image);
749 
750 
751  friend struct DispatchedInternalInitialiationAddressor;
752  friend struct AllocateMemberFunctionAddressor;
753  friend struct ToVectorAddressor;
754  friend struct ToScalarAddressor;
755 
756 
757  std::unique_ptr<PimpleImageBase> m_PimpleImage;
758 };
759 
760 } // namespace simple
761 } // namespace itk
762 
763 #endif
itk::simple::Image
The Image class for SimpleITK.
Definition: sitkImage.h:76
itk::simple::sitkLinear
@ sitkLinear
N-D linear interpolation.
Definition: sitkInterpolator.h:38
itk::simple::InterpolatorEnum
InterpolatorEnum
Definition: sitkInterpolator.h:28
sitkTemplateFunctions.h
itk::simple::PimpleImageBase
Private implementation idiom image base class.
Definition: sitkPimpleImageBase.h:38
sitkCommon.h
itk::SmartPointer
Definition: sitkImage.h:38
itk::simple::PixelIDValueEnum
PixelIDValueEnum
Enumerated values of pixelIDs.
Definition: sitkPixelIDValues.h:100
itk::simple::Image::m_PimpleImage
std::unique_ptr< PimpleImageBase > m_PimpleImage
Definition: sitkImage.h:757
itk::simple::ImageTypeToPixelIDValue
Definition: sitkPixelIDValues.h:43
sitkDetail.h
SITK_MAX_DIMENSION
#define SITK_MAX_DIMENSION
Definition: sitkConfigure.h:30
itk::simple::Image::Image
Image(itk::SmartPointer< TImageType > image)
Construct an SimpleITK Image from an pointer to an ITK image.
Definition: sitkImage.h:165
SITKCommon_EXPORT
#define SITKCommon_EXPORT
Definition: sitkCommon.h:41
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::DefaultImageDirectionTolerance
constexpr double DefaultImageDirectionTolerance
itk::DataObject
class ITK_FORWARD_EXPORT DataObject
sitkInterpolator.h
itk::simple::Image::Image
Image(TImageType *image)
Construct an SimpleITK Image from an pointer to an ITK image.
Definition: sitkImage.h:170
sitkPixelIDTokens.h
itk
itk::DefaultImageCoordinateTolerance
constexpr double DefaultImageCoordinateTolerance
itk::simple::sitkUnknown
@ sitkUnknown
Definition: sitkPixelIDValues.h:102
itk::simple::PixelIDValueType
int PixelIDValueType
Definition: sitkPixelIDValues.h:30
itk::DataObject