RTK  2.0.0
Reconstruction Toolkit
rtkWarpSequenceImageFilter.h
Go to the documentation of this file.
1 /*=========================================================================
2  *
3  * Copyright RTK Consortium
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 
19 #ifndef rtkWarpSequenceImageFilter_h
20 #define rtkWarpSequenceImageFilter_h
21 
22 #include "rtkConstantImageSource.h"
23 
24 #include <itkExtractImageFilter.h>
25 #include <itkPasteImageFilter.h>
26 #include <itkCastImageFilter.h>
28 
29 #ifdef RTK_USE_CUDA
30  #include "rtkCudaWarpImageFilter.h"
33 #else
34  #include <itkWarpImageFilter.h>
37 #endif
38 
39 namespace rtk
40 {
88 template< typename TImageSequence,
89  typename TDVFImageSequence = itk::Image< itk::CovariantVector < typename TImageSequence::ValueType,
90  TImageSequence::ImageDimension-1 >,
91  TImageSequence::ImageDimension >,
92  typename TImage = itk::Image< typename TImageSequence::ValueType,
93  TImageSequence::ImageDimension-1 >,
94  typename TDVFImage = itk::Image<itk::CovariantVector < typename TImageSequence::ValueType,
95  TImageSequence::ImageDimension - 1 >,
96  TImageSequence::ImageDimension - 1> >
97 class WarpSequenceImageFilter : public itk::ImageToImageFilter<TImageSequence, TImageSequence>
98 {
99 public:
104 
106  itkNewMacro(Self)
107 
108 
109  typedef typename itk::Image< typename TImage::PixelType,
110  TImage::ImageDimension> CPUImageType;
111  typedef typename itk::WarpImageFilter<TImage, TImage, TDVFImage> CPUWarpFilterType;
112 #ifdef RTK_USE_CUDA
113  typedef typename std::conditional< std::is_same< TImage, CPUImageType >::value,
115  CudaWarpImageFilter >::type WarpFilterType;
116  typedef typename std::conditional< std::is_same< TImage, CPUImageType >::value,
118  CudaForwardWarpImageFilter >::type
120  typedef typename std::conditional< std::is_same< TImage, CPUImageType >::value,
122  CudaCyclicDeformationImageFilter >::type
124 #else
125  typedef CPUWarpFilterType WarpFilterType;
126  typedef ForwardWarpImageFilter<TImage, TImage, TDVFImage> ForwardWarpFilterType;
127  typedef CyclicDeformationImageFilter<TDVFImageSequence,
129 #endif
130 
133 
134 
135  void SetDisplacementField(const TDVFImageSequence* DVFs);
136 
138  typename TDVFImageSequence::Pointer GetDisplacementField();
139 
141  itkGetMacro(ForwardWarp, bool)
142  itkSetMacro(ForwardWarp, bool)
143 
145  itkSetMacro(PhaseShift, float)
146  itkGetMacro(PhaseShift, float)
147 
149  itkSetMacro(UseNearestNeighborInterpolationInWarping, bool)
150  itkGetMacro(UseNearestNeighborInterpolationInWarping, bool)
151 
153  itkSetMacro(UseCudaCyclicDeformation, bool)
154  itkGetMacro(UseCudaCyclicDeformation, bool)
155 
159  typedef itk::ExtractImageFilter<TImageSequence, TImage> ExtractFilterType;
160  typedef rtk::CyclicDeformationImageFilter<TDVFImageSequence, TDVFImage> DVFInterpolatorType;
161  typedef itk::PasteImageFilter<TImageSequence,TImageSequence> PasteFilterType;
162  typedef itk::CastImageFilter<TImage, TImageSequence> CastFilterType;
164 
165 protected:
167  virtual ~WarpSequenceImageFilter() ITK_OVERRIDE {}
168 
170  void GenerateData() ITK_OVERRIDE;
171 
173  typename CPUWarpFilterType::Pointer m_WarpFilter;
174  typename ExtractFilterType::Pointer m_ExtractFilter;
179 
181  typename TImageSequence::RegionType m_ExtractAndPasteRegion;
182 
186 
190 #if ITK_VERSION_MAJOR<5
191  void VerifyInputInformation() ITK_OVERRIDE {}
192 #else
193  void VerifyInputInformation() const ITK_OVERRIDE {}
194 #endif
195 
196 
199  void GenerateOutputInformation() ITK_OVERRIDE;
200  void GenerateInputRequestedRegion() ITK_OVERRIDE;
202 
203  bool m_UseNearestNeighborInterpolationInWarping; //Default is false, linear interpolation is used instead
205 
206 private:
207  WarpSequenceImageFilter(const Self &); //purposely not implemented
208  void operator=(const Self &); //purposely not implemented
209 };
210 } //namespace ITK
211 
212 
213 #ifndef ITK_MANUAL_INSTANTIATION
214 #include "rtkWarpSequenceImageFilter.hxx"
215 #endif
216 
217 #endif
void GenerateData() override
Generate an n-dimensional image with constant pixel values.
Applies an N-D + time Motion Vector Field to an N-D + time sequence of images.
CPUWarpFilterType::Pointer m_WarpFilter
DVFInterpolatorType::Pointer m_DVFInterpolatorFilter
PasteFilterType::Pointer m_PasteFilter
ForwardWarpImageFilter< TImage, TImage, TDVFImage > ForwardWarpFilterType
Return 3D deformation vector field according to input 4D vector field, phase signal and frame number...
ExtractFilterType::Pointer m_ExtractFilter
TImageSequence::RegionType m_ExtractAndPasteRegion
CyclicDeformationImageFilter< TDVFImageSequence, TDVFImage > CudaCyclicDeformationImageFilterType
Warps an image using splat instead of interpolation.
itk::ImageToImageFilter< TImageSequence, TImageSequence > Superclass
void operator=(const Self &)
void GenerateOutputInformation() override
ConstantImageSourceType::Pointer m_ConstantSource
Mother class for cone beam reconstruction filters which need runtime selection of their forward and b...
TDVFImageSequence::Pointer GetDisplacementField()
void SetDisplacementField(const TDVFImageSequence *DVFs)
itk::WarpImageFilter< TImage, TImage, TDVFImage > CPUWarpFilterType
void GenerateInputRequestedRegion() override
#define itkSetMacro(name, type)