RTK  2.6.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  * https://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 {
89 template <typename TImageSequence,
90  typename TDVFImageSequence =
91  itk::Image<itk::CovariantVector<typename TImageSequence::ValueType, TImageSequence::ImageDimension - 1>,
92  TImageSequence::ImageDimension>,
93  typename TImage = itk::Image<typename TImageSequence::ValueType, TImageSequence::ImageDimension - 1>,
94  typename TDVFImage =
95  itk::Image<itk::CovariantVector<typename TImageSequence::ValueType, TImageSequence::ImageDimension - 1>,
96  TImageSequence::ImageDimension - 1>>
97 class ITK_TEMPLATE_EXPORT WarpSequenceImageFilter : public itk::ImageToImageFilter<TImageSequence, TImageSequence>
98 {
99 public:
100  ITK_DISALLOW_COPY_AND_MOVE(WarpSequenceImageFilter);
101 
106 
108  itkNewMacro(Self);
109 
113 #ifdef RTK_USE_CUDA
114  typedef
115  typename std::conditional<std::is_same<TImage, CPUImageType>::value, CPUWarpFilterType, CudaWarpImageFilter>::type
117  typedef typename std::conditional<std::is_same<TImage, CPUImageType>::value,
119  CudaForwardWarpImageFilter>::type ForwardWarpFilterType;
120  typedef typename std::conditional<std::is_same<TImage, CPUImageType>::value,
122  CudaCyclicDeformationImageFilter>::type CudaCyclicDeformationImageFilterType;
123 #else
127 #endif
128 
130 #ifdef itkOverrideGetNameOfClassMacro
131  itkOverrideGetNameOfClassMacro(WarpSequenceImageFilter);
132 #else
134 #endif
135 
136 
138  void
139  SetDisplacementField(const TDVFImageSequence * DVFs);
140 
142  typename TDVFImageSequence::Pointer
143  GetDisplacementField();
144 
146  itkGetMacro(ForwardWarp, bool);
147  itkSetMacro(ForwardWarp, bool);
149 
151  itkSetMacro(PhaseShift, float);
152  itkGetMacro(PhaseShift, float);
154 
156  itkSetMacro(UseNearestNeighborInterpolationInWarping, bool);
157  itkGetMacro(UseNearestNeighborInterpolationInWarping, bool);
159 
161  itkSetMacro(UseCudaCyclicDeformation, bool);
162  itkGetMacro(UseCudaCyclicDeformation, bool);
164 
173 
174 protected:
176  ~WarpSequenceImageFilter() override = default;
177 
179  void
180  GenerateData() override;
181 
183  typename CPUWarpFilterType::Pointer m_WarpFilter;
189 
191  typename TImageSequence::RegionType m_ExtractAndPasteRegion;
192 
196 
200  void
201  VerifyInputInformation() const override
202  {}
203 
206  void
207  GenerateOutputInformation() override;
208  void
209  GenerateInputRequestedRegion() override;
211 
212  bool m_UseNearestNeighborInterpolationInWarping; // Default is false, linear interpolation is used instead
214 };
215 } // namespace rtk
216 
217 
218 #ifndef ITK_MANUAL_INSTANTIATION
219 # include "rtkWarpSequenceImageFilter.hxx"
220 #endif
221 
222 #endif
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.
typename itk::Image< typename TImage::PixelType, TImage::ImageDimension > CPUImageType
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
TImageSequence::RegionType m_ExtractAndPasteRegion
#define itkSetMacro(name, type)
Warps an image using splat instead of interpolation.
void VerifyInputInformation() const override
typename itk::WarpImageFilter< TImage, TImage, TDVFImage > CPUWarpFilterType
ConstantImageSourceType::Pointer m_ConstantSource
Mother class for cone beam reconstruction filters which need runtime selection of their forward and b...