Medical Imaging Interaction Toolkit  2021.10.00
Medical Imaging Interaction Toolkit
itkStitchImageFilter.h
Go to the documentation of this file.
1 /*============================================================================
2 
3 The Medical Imaging Interaction Toolkit (MITK)
4 
5 Copyright (c) German Cancer Research Center (DKFZ)
6 All rights reserved.
7 
8 Use of this source code is governed by a 3-clause BSD license that can be
9 found in the LICENSE file.
10 
11 ============================================================================*/
12 
13 #ifndef itkStitchImageFilter_h
14 #define itkStitchImageFilter_h
15 
16 #include "itkFixedArray.h"
17 #include "itkTransform.h"
18 #include "itkImageRegionIterator.h"
19 #include "itkImageToImageFilter.h"
20 #include "itkLinearInterpolateImageFunction.h"
21 #include "itkSize.h"
22 #include "itkDefaultConvertPixelTraits.h"
23 #include "itkDataObjectDecorator.h"
24 
25 
26 namespace itk
27 {
28  enum class StitchStrategy
29  {
30  Mean = 0, //use the mean value of all inputs that can provide a pixel vaule
31  BorderDistance = 1 //use the value that is largest minimal distance to its image borders (use e.g. if vaules tend to be not reliable at borders)
32  };
33 
34  std::ostream& operator<< (std::ostream& os, const itk::StitchStrategy& strategy)
35  {
36  if (itk::StitchStrategy::Mean == strategy)
37  os << "Mean";
38  else if (itk::StitchStrategy::BorderDistance == strategy)
39  os << "BorderDistance";
40  else
41  os << "unkown";
42 
43  return os;
44  };
45 
59 template< typename TInputImage,
60  typename TOutputImage,
61  typename TInterpolatorPrecisionType = double,
62  typename TTransformPrecisionType = TInterpolatorPrecisionType>
64  public ImageToImageFilter< TInputImage, TOutputImage >
65 {
66 public:
69  typedef ImageToImageFilter< TInputImage, TOutputImage > Superclass;
72 
73  typedef TInputImage InputImageType;
74  typedef TOutputImage OutputImageType;
75  typedef typename InputImageType::Pointer InputImagePointer;
76  typedef typename InputImageType::ConstPointer InputImageConstPointer;
77  typedef typename OutputImageType::Pointer OutputImagePointer;
78  typedef typename InputImageType::RegionType InputImageRegionType;
79 
81  itkNewMacro(Self);
82 
84  itkTypeMacro(StitchImageFilter, ImageToImageFilter);
85 
87  itkStaticConstMacro(ImageDimension, unsigned int,
88  TOutputImage::ImageDimension);
89  itkStaticConstMacro(InputImageDimension, unsigned int,
90  TInputImage::ImageDimension);
91 
93  typedef ImageBase< itkGetStaticConstMacro(ImageDimension) > ImageBaseType;
94 
98  typedef Transform< TTransformPrecisionType,
99  itkGetStaticConstMacro(ImageDimension),
100  itkGetStaticConstMacro(ImageDimension) > TransformType;
101  typedef typename TransformType::ConstPointer TransformPointerType;
102  typedef DataObjectDecorator<TransformType> DecoratedTransformType;
103  typedef typename DecoratedTransformType::Pointer DecoratedTransformPointer;
104 
105 
107  typedef InterpolateImageFunction< InputImageType,
108  TInterpolatorPrecisionType > InterpolatorType;
109  typedef typename InterpolatorType::Pointer InterpolatorPointerType;
110 
112 
113  typedef DefaultConvertPixelTraits< InterpolatorOutputType > InterpolatorConvertType;
114 
115  typedef typename InterpolatorConvertType::ComponentType ComponentType;
116 
117  typedef LinearInterpolateImageFunction< InputImageType,
118  TInterpolatorPrecisionType > LinearInterpolatorType;
119  typedef typename LinearInterpolatorType::Pointer
121 
123  typedef Size< itkGetStaticConstMacro(ImageDimension) > SizeType;
124 
126  typedef typename TOutputImage::IndexType IndexType;
127 
129  typedef typename InterpolatorType::PointType PointType;
130  //typedef typename TOutputImage::PointType PointType;
131 
133  typedef typename TOutputImage::PixelType PixelType;
134  typedef typename TInputImage::PixelType InputPixelType;
135 
136  typedef DefaultConvertPixelTraits<PixelType> PixelConvertType;
137 
138  typedef typename PixelConvertType::ComponentType PixelComponentType;
139 
141  typedef ContinuousIndex< TTransformPrecisionType, ImageDimension >
143 
145  typedef typename TOutputImage::RegionType OutputImageRegionType;
146 
148  typedef typename TOutputImage::SpacingType SpacingType;
149  typedef typename TOutputImage::PointType OriginPointType;
150  typedef typename TOutputImage::DirectionType DirectionType;
151 
152  using Superclass::GetInput;
153 
155  typedef ImageBase<ImageDimension> ReferenceImageBaseType;
156 
157  using Superclass::SetInput;
158  void SetInput(const InputImageType* image) override;
159  void SetInput(unsigned int index, const InputImageType* image) override;
162  virtual void SetInput(unsigned int index, const InputImageType* image, const TransformType* transform);
163  virtual void SetInput(unsigned int index, const InputImageType* image, const TransformType* transform, InterpolatorType* interpolator);
164 
165  const TransformType* GetTransform(unsigned int index) const;
166 
167  const InterpolatorType* GetInterpolator(unsigned int index) const;
168 
170  itkSetMacro(Size, SizeType);
171  itkGetConstReferenceMacro(Size, SizeType);
172 
175  itkSetMacro(DefaultPixelValue, PixelType);
176  itkGetConstReferenceMacro(DefaultPixelValue, PixelType);
177 
179  itkSetMacro(OutputSpacing, SpacingType);
180  virtual void SetOutputSpacing(const double *values);
181 
183  itkGetConstReferenceMacro(OutputSpacing, SpacingType);
184 
186  itkSetMacro(OutputOrigin, OriginPointType);
187  virtual void SetOutputOrigin(const double *values);
188 
190  itkGetConstReferenceMacro(OutputOrigin, OriginPointType);
191 
193  itkSetMacro(OutputDirection, DirectionType);
194  itkGetConstReferenceMacro(OutputDirection, DirectionType);
195 
197  void SetOutputParametersFromImage(const ImageBaseType *image);
198 
201  itkSetMacro(OutputStartIndex, IndexType);
202 
204  itkGetConstReferenceMacro(OutputStartIndex, IndexType);
205 
212  itkSetInputMacro(ReferenceImage, ReferenceImageBaseType);
213 
215  itkGetInputMacro(ReferenceImage, ReferenceImageBaseType);
216 
219  itkSetMacro(UseReferenceImage, bool);
220  itkBooleanMacro(UseReferenceImage);
221  itkGetConstMacro(UseReferenceImage, bool);
222 
223  itkSetMacro(StitchStrategy, StitchStrategy);
224  itkGetConstMacro(StitchStrategy, StitchStrategy);
225 
231  virtual void GenerateOutputInformation() ITK_OVERRIDE;
232 
238  virtual void GenerateInputRequestedRegion() ITK_OVERRIDE;
239 
243  virtual void BeforeThreadedGenerateData() ITK_OVERRIDE;
244 
246  virtual void AfterThreadedGenerateData() ITK_OVERRIDE;
247 
249  ModifiedTimeType GetMTime(void) const ITK_OVERRIDE;
250 
251 #ifdef ITK_USE_CONCEPT_CHECKING
252  // Begin concept checking
253  itkConceptMacro( OutputHasNumericTraitsCheck,
254  ( Concept::HasNumericTraits< PixelComponentType > ) );
255  // End concept checking
256 #endif
257 
258 protected:
260  ~StitchImageFilter() ITK_OVERRIDE {}
261  void PrintSelf(std::ostream & os, Indent indent) const ITK_OVERRIDE;
262 
268  virtual void VerifyInputInformation() ITK_OVERRIDE { }
269 
279  virtual void ThreadedGenerateData(const OutputImageRegionType & outputRegionForThread,
280  ThreadIdType threadId) ITK_OVERRIDE;
281 
283  virtual PixelType CastPixelWithBoundsChecking( const InterpolatorOutputType value,
284  const ComponentType minComponent,
285  const ComponentType maxComponent) const;
286 
287  void SetTransform(unsigned int index, const TransformType* transform);
288 
292  void EnsureTransforms();
293 
297  void EnsureInterpolators();
298 
299  static std::string GetTransformInputName(unsigned int index);
300 
301 private:
302  ITK_DISALLOW_COPY_AND_ASSIGN(StitchImageFilter);
303 
304  typedef std::vector<const InputImageType*> InputImageVectorType;
305  typedef std::map<const InputImageType*, typename TransformType::ConstPointer> TransformMapType;
306  typedef std::map<const InputImageType*, InterpolatorPointerType> InterpolatorMapType;
307 
308  InputImageVectorType GetInputs();
309  TransformMapType GetTransforms();
310 
311  InterpolatorMapType m_Interpolators; // Image function for
312  // interpolation
313  PixelType m_DefaultPixelValue; // default pixel value
314  // if the point is
315  // outside the image
316  SizeType m_Size; // Size of the output image
317  SpacingType m_OutputSpacing; // output image spacing
318  OriginPointType m_OutputOrigin; // output image origin
319  DirectionType m_OutputDirection; // output image direction cosines
320  IndexType m_OutputStartIndex; // output image start index
321  bool m_UseReferenceImage;
322  StitchStrategy m_StitchStrategy;
323 };
324 } // end namespace itk
325 
326 #ifndef ITK_MANUAL_INSTANTIATION
327 #include "itkStitchImageFilter.tpp"
328 #endif
329 
330 #endif
TOutputImage::DirectionType DirectionType
TOutputImage::PixelType PixelType
InterpolatorType::OutputType InterpolatorOutputType
InterpolatorType::PointType PointType
std::ostream & operator<<(std::ostream &os, const itk::StitchStrategy &strategy)
DataObjectDecorator< TransformType > DecoratedTransformType
ITK filter that resamples/stitches multiple images into a given reference geometry.
InterpolatorType::Pointer InterpolatorPointerType
TOutputImage::RegionType OutputImageRegionType
DefaultConvertPixelTraits< PixelType > PixelConvertType
DecoratedTransformType::Pointer DecoratedTransformPointer
PixelConvertType::ComponentType PixelComponentType
Transform< TTransformPrecisionType, itkGetStaticConstMacro(ImageDimension), itkGetStaticConstMacro(ImageDimension) > TransformType
virtual void VerifyInputInformation() ITK_OVERRIDE
SmartPointer< Self > Pointer
InterpolatorConvertType::ComponentType ComponentType
LinearInterpolatorType::Pointer LinearInterpolatorPointerType
DefaultConvertPixelTraits< InterpolatorOutputType > InterpolatorConvertType
TInputImage::PixelType InputPixelType
SmartPointer< const Self > ConstPointer
InputImageType::RegionType InputImageRegionType
InputImageType::Pointer InputImagePointer
TOutputImage::SpacingType SpacingType
InputImageType::ConstPointer InputImageConstPointer
TransformType::ConstPointer TransformPointerType
TOutputImage::IndexType IndexType
LinearInterpolateImageFunction< InputImageType, TInterpolatorPrecisionType > LinearInterpolatorType
ImageBase< ImageDimension > ReferenceImageBaseType
ContinuousIndex< TTransformPrecisionType, ImageDimension > ContinuousInputIndexType
ImageToImageFilter< TInputImage, TOutputImage > Superclass
Size< itkGetStaticConstMacro(ImageDimension) > SizeType
ImageBase< itkGetStaticConstMacro(ImageDimension) > ImageBaseType
TOutputImage::PointType OriginPointType
InterpolateImageFunction< InputImageType, TInterpolatorPrecisionType > InterpolatorType
OutputImageType::Pointer OutputImagePointer