RTK  2.0.1
Reconstruction Toolkit
rtkUnwarpSequenceImageFilter.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 rtkUnwarpSequenceImageFilter_h
20 #define rtkUnwarpSequenceImageFilter_h
21 
25 #include "rtkConstantImageSource.h"
26 
27 #ifdef RTK_USE_CUDA
30 #endif
31 
32 namespace rtk
33 {
71  template< typename TImageSequence,
72  typename TDVFImageSequence = itk::Image< itk::CovariantVector < typename TImageSequence::ValueType,
73  TImageSequence::ImageDimension-1 >,
74  TImageSequence::ImageDimension >,
75  typename TImage = itk::Image< typename TImageSequence::ValueType,
76  TImageSequence::ImageDimension-1 >,
77  typename TDVFImage = itk::Image<itk::CovariantVector < typename TImageSequence::ValueType,
78  TImageSequence::ImageDimension - 1 >,
79  TImageSequence::ImageDimension - 1> >
80 class UnwarpSequenceImageFilter : public itk::ImageToImageFilter<TImageSequence, TImageSequence>
81 {
82 public:
83  ITK_DISALLOW_COPY_AND_ASSIGN(UnwarpSequenceImageFilter);
84 
89 
91  itkNewMacro(Self)
92 
93 
95 
97  TDVFImageSequence,
98  TImage,
99  TDVFImage>;
100  using WarpForwardFilterType = WarpSequenceImageFilter< TImageSequence,
101  TDVFImageSequence,
102  TImage,
103  TDVFImage>;
105 
107  using CPUImageSequence = typename itk::Image< typename TImageSequence::PixelType,
108  TImageSequence::ImageDimension>;
109 #ifdef RTK_USE_CUDA
110  typedef typename std::conditional< std::is_same< TImageSequence, CPUImageSequence >::value,
112  CudaConstantVolumeSeriesSource >::type
114  typedef typename std::conditional< std::is_same< TImageSequence, CPUImageSequence >::value,
116  CudaConjugateGradientImageFilter<TImageSequence> >::type
118 #else
119  using ConstantSourceType = ConstantImageSource<TImageSequence>;
120  using CudaConjugateGradientType = ConjugateGradientFilterType;
121 #endif
122 
124  void SetDisplacementField(const TDVFImageSequence* DVFs);
125 
127  typename TDVFImageSequence::Pointer GetDisplacementField();
128 
130  itkSetMacro(NumberOfIterations, float)
131  itkGetMacro(NumberOfIterations, float)
132 
134  itkSetMacro(PhaseShift, float)
135  itkGetMacro(PhaseShift, float)
136 
137  itkSetMacro(UseNearestNeighborInterpolationInWarping, bool)
138  itkGetMacro(UseNearestNeighborInterpolationInWarping, bool)
139 
140  itkSetMacro(CudaConjugateGradient, bool)
141  itkGetMacro(CudaConjugateGradient, bool)
142 
144  itkSetMacro(UseCudaCyclicDeformation, bool)
145  itkGetMacro(UseCudaCyclicDeformation, bool)
146 
147 protected:
148  UnwarpSequenceImageFilter();
149  ~UnwarpSequenceImageFilter() override = default;
150 
152  void GenerateData() override;
153 
155  typename ConjugateGradientFilterType::Pointer m_ConjugateGradientFilter;
158  typename ConstantSourceType::Pointer m_ConstantSource;
159 
162 
166 #if ITK_VERSION_MAJOR<5
167  void VerifyInputInformation() override {}
168 #else
169  void VerifyInputInformation() const override {}
170 #endif
171 
172 
175  void GenerateInputRequestedRegion() override;
176  void GenerateOutputInformation() override;
178 
179  bool m_UseNearestNeighborInterpolationInWarping; //Default is false, linear interpolation is used instead
182 
183 private:
184  unsigned int m_NumberOfIterations;
185 };
186 } //namespace ITK
187 
188 
189 #ifndef ITK_MANUAL_INSTANTIATION
190 #include "rtkUnwarpSequenceImageFilter.hxx"
191 #endif
192 
193 #endif
ConjugateGradientFilterType CudaConjugateGradientType
Implements the operator A used in the conjugate gradient unwarp sequence filter.
CGOperatorFilterType::Pointer m_CGOperator
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.
ConjugateGradientImageFilter< TImageSequence > ConjugateGradientFilterType
Finds the image sequence that, once warped, equals the input image sequence.
void SetDisplacementField(const TDVFImageSequence *DVFs)
ConstantSourceType::Pointer m_ConstantSource
void GenerateInputRequestedRegion() override
void GenerateOutputInformation() override
ConjugateGradientFilterType::Pointer m_ConjugateGradientFilter
Solves AX = B by conjugate gradient.
typename itk::Image< typename TImageSequence::PixelType, TImageSequence::ImageDimension > CPUImageSequence
WarpForwardFilterType::Pointer m_WarpForwardFilter
#define itkSetMacro(name, type)
TDVFImageSequence::Pointer GetDisplacementField()
ConstantImageSource< TImageSequence > ConstantSourceType