RTK  2.6.0
Reconstruction Toolkit
rtkMotionCompensatedFourDReconstructionConjugateGradientOperator.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 #ifndef rtkMotionCompensatedFourDReconstructionConjugateGradientOperator_h
19 #define rtkMotionCompensatedFourDReconstructionConjugateGradientOperator_h
20 
24 
25 #ifdef RTK_USE_CUDA
29 #endif
30 
31 namespace rtk
32 {
99 template <typename VolumeSeriesType, typename ProjectionStackType>
101  : public FourDReconstructionConjugateGradientOperator<VolumeSeriesType, ProjectionStackType>
102 {
103 public:
105 
110 
112  using VolumeType = ProjectionStackType;
113  using VectorForDVF = itk::CovariantVector<typename VolumeSeriesType::ValueType, VolumeSeriesType::ImageDimension - 1>;
114 
116  using CPUVolumeSeriesType =
118 #ifdef RTK_USE_CUDA
119  typedef typename std::conditional<std::is_same<VolumeSeriesType, CPUVolumeSeriesType>::value,
121  itk::CudaImage<VectorForDVF, VolumeSeriesType::ImageDimension>>::type
123  typedef
124  typename std::conditional<std::is_same<VolumeSeriesType, CPUVolumeSeriesType>::value,
125  itk::Image<VectorForDVF, VolumeSeriesType::ImageDimension - 1>,
126  itk::CudaImage<VectorForDVF, VolumeSeriesType::ImageDimension - 1>>::type DVFImageType;
127  typedef typename std::conditional<std::is_same<VolumeSeriesType, CPUVolumeSeriesType>::value,
129  CudaWarpForwardProjectionImageFilter>::type WarpForwardProjectionImageFilterType;
130  typedef typename std::conditional<std::is_same<VolumeSeriesType, CPUVolumeSeriesType>::value,
132  CudaWarpBackProjectionImageFilter>::type WarpBackProjectionImageFilterType;
133 #else
135  using DVFImageType = itk::Image<VectorForDVF, VolumeSeriesType::ImageDimension - 1>;
139 #endif
141 #ifdef RTK_USE_CUDA
142  typedef typename std::conditional<std::is_same<VolumeSeriesType, CPUVolumeSeriesType>::value,
144  CudaCyclicDeformationImageFilter>::type CudaCyclicDeformationImageFilterType;
145 #else
147 #endif
148 
150  itkNewMacro(Self);
151 
153 #ifdef itkOverrideGetNameOfClassMacro
155 #else
158 #endif
159 
160 
162  void
164  {
165  itkExceptionMacro(<< "ForwardProjection cannot be changed");
166  }
167  void
169  {
170  itkExceptionMacro(<< "BackProjection cannot be changed");
171  }
173 
175  void
176  SetDisplacementField(const DVFSequenceImageType * DisplacementField);
177  void
178  SetInverseDisplacementField(const DVFSequenceImageType * InverseDisplacementField);
179  typename DVFSequenceImageType::ConstPointer
180  GetInverseDisplacementField();
181  typename DVFSequenceImageType::ConstPointer
182  GetDisplacementField();
184 
186  void
187  SetSignal(const std::vector<double> signal) override;
188 
190  itkSetMacro(UseCudaCyclicDeformation, bool);
191  itkGetMacro(UseCudaCyclicDeformation, bool);
193 
194 protected:
197 
199  void
200  GenerateOutputInformation() override;
201 
203  void
204  VerifyInputInformation() const override
205  {}
206 
208  void
209  GenerateData() override;
210 
214  std::vector<double> m_Signal;
216 };
217 } // namespace rtk
218 
219 
220 #ifndef ITK_MANUAL_INSTANTIATION
221 # include "rtkMotionCompensatedFourDReconstructionConjugateGradientOperator.hxx"
222 #endif
223 
224 #endif
Return 3D deformation vector field according to input 4D vector field, phase signal and frame number...
typename itk::Image< typename VolumeSeriesType::PixelType, VolumeSeriesType::ImageDimension > CPUVolumeSeriesType
#define itkSetMacro(name, type)
Implements part of the 4D reconstruction by conjugate gradient.