RTK  2.0.1
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:
100  ITK_DISALLOW_COPY_AND_ASSIGN(WarpSequenceImageFilter);
101 
106 
108  itkNewMacro(Self)
109 
110 
111  using CPUImageType = typename itk::Image< typename TImage::PixelType,
112  TImage::ImageDimension>;
113  using CPUWarpFilterType = typename itk::WarpImageFilter<TImage, TImage, TDVFImage>;
114 #ifdef RTK_USE_CUDA
115  typedef typename std::conditional< std::is_same< TImage, CPUImageType >::value,
117  CudaWarpImageFilter >::type WarpFilterType;
118  typedef typename std::conditional< std::is_same< TImage, CPUImageType >::value,
120  CudaForwardWarpImageFilter >::type
122  typedef typename std::conditional< std::is_same< TImage, CPUImageType >::value,
124  CudaCyclicDeformationImageFilter >::type
126 #else
127  using WarpFilterType = CPUWarpFilterType;
128  using ForwardWarpFilterType = ForwardWarpImageFilter<TImage, TImage, TDVFImage>;
129  using CudaCyclicDeformationImageFilterType = CyclicDeformationImageFilter<TDVFImageSequence,
130  TDVFImage>;
131 #endif
132 
135 
136 
137  void SetDisplacementField(const TDVFImageSequence* DVFs);
138 
140  typename TDVFImageSequence::Pointer GetDisplacementField();
141 
143  itkGetMacro(ForwardWarp, bool)
144  itkSetMacro(ForwardWarp, bool)
145 
147  itkSetMacro(PhaseShift, float)
148  itkGetMacro(PhaseShift, float)
149 
151  itkSetMacro(UseNearestNeighborInterpolationInWarping, bool)
152  itkGetMacro(UseNearestNeighborInterpolationInWarping, bool)
153 
155  itkSetMacro(UseCudaCyclicDeformation, bool)
156  itkGetMacro(UseCudaCyclicDeformation, bool)
157 
161  using ExtractFilterType = itk::ExtractImageFilter<TImageSequence, TImage>;
162  using DVFInterpolatorType = rtk::CyclicDeformationImageFilter<TDVFImageSequence, TDVFImage>;
163  using PasteFilterType = itk::PasteImageFilter<TImageSequence,TImageSequence>;
164  using CastFilterType = itk::CastImageFilter<TImage, TImageSequence>;
166 
167 protected:
169  ~WarpSequenceImageFilter() override = default;
170 
172  void GenerateData() override;
173 
175  typename CPUWarpFilterType::Pointer m_WarpFilter;
181 
183  typename TImageSequence::RegionType m_ExtractAndPasteRegion;
184 
188 
192 #if ITK_VERSION_MAJOR<5
193  void VerifyInputInformation() override {}
194 #else
195  void VerifyInputInformation() const override {}
196 #endif
197 
198 
201  void GenerateOutputInformation() override;
202  void GenerateInputRequestedRegion() override;
204 
205  bool m_UseNearestNeighborInterpolationInWarping; //Default is false, linear interpolation is used instead
207 
208 };
209 } //namespace ITK
210 
211 
212 #ifndef ITK_MANUAL_INSTANTIATION
213 #include "rtkWarpSequenceImageFilter.hxx"
214 #endif
215 
216 #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
Return 3D deformation vector field according to input 4D vector field, phase signal and frame number...
ExtractFilterType::Pointer m_ExtractFilter
CyclicDeformationImageFilter< TDVFImageSequence, TDVFImage > CudaCyclicDeformationImageFilterType
TImageSequence::RegionType m_ExtractAndPasteRegion
typename itk::Image< typename TImage::PixelType, TImage::ImageDimension > CPUImageType
Warps an image using splat instead of interpolation.
typename itk::WarpImageFilter< TImage, TImage, TDVFImage > CPUWarpFilterType
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)
ForwardWarpImageFilter< TImage, TImage, TDVFImage > ForwardWarpFilterType
ImageBaseType::RegionType RegionType
void GenerateInputRequestedRegion() override
#define itkSetMacro(name, type)