RTK  2.7.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  using DVFSequenceImageType =
120  typename std::conditional_t<std::is_same_v<VolumeSeriesType, CPUVolumeSeriesType>,
122  itk::CudaImage<VectorForDVF, VolumeSeriesType::ImageDimension>>;
123  using DVFImageType = typename std::conditional_t<std::is_same_v<VolumeSeriesType, CPUVolumeSeriesType>,
124  itk::Image<VectorForDVF, VolumeSeriesType::ImageDimension - 1>,
125  itk::CudaImage<VectorForDVF, VolumeSeriesType::ImageDimension - 1>>;
127  typename std::conditional_t<std::is_same_v<VolumeSeriesType, CPUVolumeSeriesType>,
131  typename std::conditional_t<std::is_same_v<VolumeSeriesType, CPUVolumeSeriesType>,
134 #else
136  using DVFImageType = itk::Image<VectorForDVF, VolumeSeriesType::ImageDimension - 1>;
140 #endif
142 #ifdef RTK_USE_CUDA
144  typename std::conditional_t<std::is_same_v<VolumeSeriesType, CPUVolumeSeriesType>,
147 #else
149 #endif
150 
152  itkNewMacro(Self);
153 
156 
158  void
159  SetForwardProjectionFilter(typename Superclass::ForwardProjectionFilterType * itkNotUsed(_arg))
160  {
161  itkExceptionMacro(<< "ForwardProjection cannot be changed");
162  }
163  void
164  SetBackProjectionFilter(typename Superclass::BackProjectionFilterType * itkNotUsed(_arg))
165  {
166  itkExceptionMacro(<< "BackProjection cannot be changed");
167  }
169 
171  void
172  SetDisplacementField(const DVFSequenceImageType * DisplacementField);
173  void
174  SetInverseDisplacementField(const DVFSequenceImageType * InverseDisplacementField);
175  typename DVFSequenceImageType::ConstPointer
176  GetInverseDisplacementField();
177  typename DVFSequenceImageType::ConstPointer
178  GetDisplacementField();
180 
182  void
183  SetSignal(const std::vector<double> signal) override;
184 
186  itkSetMacro(UseCudaCyclicDeformation, bool);
187  itkGetMacro(UseCudaCyclicDeformation, bool);
189 
190 protected:
193 
195  void
196  GenerateOutputInformation() override;
197 
199  void
200  VerifyInputInformation() const override
201  {}
202 
204  void
205  GenerateData() override;
206 
210  std::vector<double> m_Signal;
212 };
213 } // namespace rtk
214 
215 
216 #ifndef ITK_MANUAL_INSTANTIATION
217 # include "rtkMotionCompensatedFourDReconstructionConjugateGradientOperator.hxx"
218 #endif
219 
220 #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
Voxel-based backprojection into warped volume implemented in CUDA.
#define itkSetMacro(name, type)
Implements part of the 4D reconstruction by conjugate gradient.
typename std::conditional_t< std::is_same_v< VolumeSeriesType, CPUVolumeSeriesType >, itk::Image< VectorForDVF, VolumeSeriesType::ImageDimension >, itk::CudaImage< VectorForDVF, VolumeSeriesType::ImageDimension > > DVFSequenceImageType
GPU version of the temporal DVF interpolator.
typename std::conditional_t< std::is_same_v< VolumeSeriesType, CPUVolumeSeriesType >, JosephForwardProjectionImageFilter< ProjectionStackType, ProjectionStackType >, CudaWarpForwardProjectionImageFilter > WarpForwardProjectionImageFilterType
typename std::conditional_t< std::is_same_v< VolumeSeriesType, CPUVolumeSeriesType >, BackProjectionImageFilter< VolumeType, VolumeType >, CudaWarpBackProjectionImageFilter > WarpBackProjectionImageFilterType
typename std::conditional_t< std::is_same_v< VolumeSeriesType, CPUVolumeSeriesType >, itk::Image< VectorForDVF, VolumeSeriesType::ImageDimension - 1 >, itk::CudaImage< VectorForDVF, VolumeSeriesType::ImageDimension - 1 > > DVFImageType
typename std::conditional_t< std::is_same_v< VolumeSeriesType, CPUVolumeSeriesType >, CPUDVFInterpolatorType, CudaCyclicDeformationImageFilter > CudaCyclicDeformationImageFilterType