RTK  1.4.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 
110 
112  void SetDisplacementField(const TDVFImageSequence* DVFs);
113 
115  typename TDVFImageSequence::Pointer GetDisplacementField();
116 
118  itkGetMacro(ForwardWarp, bool)
119  itkSetMacro(ForwardWarp, bool)
120 
122  itkSetMacro(PhaseShift, float)
123  itkGetMacro(PhaseShift, float)
124 
126  itkSetMacro(UseNearestNeighborInterpolationInWarping, bool)
127  itkGetMacro(UseNearestNeighborInterpolationInWarping, bool)
128 
130  itkSetMacro(UseCudaCyclicDeformation, bool)
131  itkGetMacro(UseCudaCyclicDeformation, bool)
132 
134 #ifdef RTK_USE_CUDA
135  typedef rtk::CudaWarpImageFilter CudaWarpFilterType;
136  typedef rtk::CudaForwardWarpImageFilter CudaForwardWarpFilterType;
137 #endif
140 
148 
149 protected:
151  virtual ~WarpSequenceImageFilter() ITK_OVERRIDE {}
152 
154  void GenerateData() ITK_OVERRIDE;
155 
157  typename WarpFilterType::Pointer m_WarpFilter;
158  typename ExtractFilterType::Pointer m_ExtractFilter;
163 
165  typename TImageSequence::RegionType m_ExtractAndPasteRegion;
166 
170 
174 #if ITK_VERSION_MAJOR<5
175  void VerifyInputInformation() ITK_OVERRIDE {}
176 #else
177  void VerifyInputInformation() const ITK_OVERRIDE {}
178 #endif
179 
180 
183  void GenerateOutputInformation() ITK_OVERRIDE;
184  void GenerateInputRequestedRegion() ITK_OVERRIDE;
186 
187  bool m_UseNearestNeighborInterpolationInWarping; //Default is false, linear interpolation is used instead
189 
190 private:
191  WarpSequenceImageFilter(const Self &); //purposely not implemented
192  void operator=(const Self &); //purposely not implemented
193 };
194 } //namespace ITK
195 
196 
197 #ifndef ITK_MANUAL_INSTANTIATION
198 #include "rtkWarpSequenceImageFilter.hxx"
199 #endif
200 
201 #endif
itk::NearestNeighborInterpolateImageFunction< TImage, double > NearestNeighborInterpolatorType
rtk::ForwardWarpImageFilter< TImage, TImage, TDVFImage > ForwardWarpFilterType
void GenerateData() override
Generate an n-dimensional image with constant pixel values.
itk::LinearInterpolateImageFunction< TImage, double > LinearInterpolatorType
Applies an N-D + time Motion Vector Field to an N-D + time sequence of images.
DVFInterpolatorType::Pointer m_DVFInterpolatorFilter
itk::PasteImageFilter< TImageSequence, TImageSequence > PasteFilterType
PasteFilterType::Pointer m_PasteFilter
Return 3D deformation vector field according to input 4D vector field, phase signal and frame number...
itk::WarpImageFilter< TImage, TImage, TDVFImage > WarpFilterType
ExtractFilterType::Pointer m_ExtractFilter
TImageSequence::RegionType m_ExtractAndPasteRegion
itk::CastImageFilter< TImage, TImageSequence > CastFilterType
itk::ExtractImageFilter< TImageSequence, TImage > ExtractFilterType
Warps an image using splat instead of interpolation.
itk::ImageToImageFilter< TImageSequence, TImageSequence > Superclass
void operator=(const Self &)
rtk::ConstantImageSource< TImageSequence > ConstantImageSourceType
void GenerateOutputInformation() override
ConstantImageSourceType::Pointer m_ConstantSource
Mother class for cone beam reconstruction filters which need runtime selection of their forward and b...
rtk::CyclicDeformationImageFilter< TDVFImageSequence, TDVFImage > DVFInterpolatorType
TDVFImageSequence::Pointer GetDisplacementField()
void SetDisplacementField(const TDVFImageSequence *DVFs)
void GenerateInputRequestedRegion() override
#define itkSetMacro(name, type)