RTK  1.4.0
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:
87 
89  itkNewMacro(Self)
90 
91 
93 
94  typedef rtk::UnwarpSequenceConjugateGradientOperator<TImageSequence,
95  TDVFImageSequence,
96  TImage,
98  typedef rtk::WarpSequenceImageFilter< TImageSequence,
99  TDVFImageSequence,
100  TImage,
103 // typedef rtk::CyclicDeformationImageFilter<TDVFImage> DVFInterpolatorType;
104  typedef rtk::ConstantImageSource<TImageSequence> ConstantSourceType;
105 
107  void SetDisplacementField(const TDVFImageSequence* DVFs);
108 
110  typename TDVFImageSequence::Pointer GetDisplacementField();
111 
113  itkSetMacro(NumberOfIterations, float)
114  itkGetMacro(NumberOfIterations, float)
115 
117  itkSetMacro(PhaseShift, float)
118  itkGetMacro(PhaseShift, float)
119 
120  itkSetMacro(UseNearestNeighborInterpolationInWarping, bool)
121  itkGetMacro(UseNearestNeighborInterpolationInWarping, bool)
122 
123  itkSetMacro(CudaConjugateGradient, bool)
124  itkGetMacro(CudaConjugateGradient, bool)
125 
127  itkSetMacro(UseCudaCyclicDeformation, bool)
128  itkGetMacro(UseCudaCyclicDeformation, bool)
129 
130 protected:
131  UnwarpSequenceImageFilter();
132  virtual ~UnwarpSequenceImageFilter() ITK_OVERRIDE {}
133 
135  void GenerateData() ITK_OVERRIDE;
136 
142 
145 
149  void VerifyInputInformation() ITK_OVERRIDE {}
150 
153  void GenerateInputRequestedRegion() ITK_OVERRIDE;
154  void GenerateOutputInformation() ITK_OVERRIDE;
156 
157  bool m_UseNearestNeighborInterpolationInWarping; //Default is false, linear interpolation is used instead
160 
161 private:
162  UnwarpSequenceImageFilter(const Self &); //purposely not implemented
163  void operator=(const Self &); //purposely not implemented
164 
165  unsigned int m_NumberOfIterations;
166 };
167 } //namespace ITK
168 
169 
170 #ifndef ITK_MANUAL_INSTANTIATION
171 #include "rtkUnwarpSequenceImageFilter.hxx"
172 #endif
173 
174 #endif
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.
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
itk::ImageToImageFilter< TImageSequence, TImageSequence > Superclass
ConjugateGradientFilterType::Pointer m_ConjugateGradientFilter
void operator=(const Self &)
Solves AX = B by conjugate gradient.
WarpForwardFilterType::Pointer m_WarpForwardFilter
#define itkSetMacro(name, type)
TDVFImageSequence::Pointer GetDisplacementField()