RTK  2.0.0
Reconstruction Toolkit
rtkMechlemOneStepSpectralReconstructionFilter.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 rtkMechlemOneStepSpectralReconstructionFilter_h
20 #define rtkMechlemOneStepSpectralReconstructionFilter_h
21 
26 #include "rtkConstantImageSource.h"
30 
31 #include <itkExtractImageFilter.h>
32 #include <itkAddImageFilter.h>
33 #include <itkMultiplyImageFilter.h>
34 
35 #ifdef RTK_USE_CUDA
37 #endif
38 
39 namespace rtk
40 {
125 template< typename TOutputImage, typename TPhotonCounts, typename TSpectrum >
127 {
128 public:
133 
135  itkNewMacro(Self)
136 
137 
139 
141  itkStaticConstMacro(nBins, unsigned int, TPhotonCounts::PixelType::Dimension);
142  itkStaticConstMacro(nMaterials, unsigned int, TOutputImage::PixelType::Dimension);
143  typedef typename TOutputImage::PixelType::ValueType dataType;
145 
147  typedef typename itk::Image< typename TOutputImage::PixelType,
148  TOutputImage::ImageDimension> CPUOutputImageType;
149 #ifdef RTK_USE_CUDA
150  typedef typename std::conditional< std::is_same< TOutputImage, CPUOutputImageType >::value,
151  itk::Image< itk::Vector<dataType, nMaterials * nMaterials>, TOutputImage::ImageDimension >,
152  itk::CudaImage< itk::Vector<dataType, nMaterials * nMaterials>, TOutputImage::ImageDimension > >::type
154  typedef typename std::conditional< std::is_same< TOutputImage, CPUOutputImageType >::value,
156  itk::CudaImage<dataType, TOutputImage::ImageDimension> >::type
158 #else
159  typedef typename itk::Image< itk::Vector<dataType, nMaterials * nMaterials>,
160  TOutputImage::ImageDimension > THessiansImage;
161  typedef typename itk::Image<dataType, TOutputImage::ImageDimension> SingleComponentImageType;
162 #endif
163 
164 #if !defined( ITK_WRAPPING_PARSER )
165 #ifdef RTK_USE_CUDA
166  typedef typename std::conditional< std::is_same< TOutputImage, CPUOutputImageType >::value,
168  CudaWeidingerForwardModelImageFilter<TOutputImage, TPhotonCounts, TSpectrum> >::type
170  typedef typename std::conditional< std::is_same< TOutputImage, CPUOutputImageType >::value,
172  CudaForwardProjectionImageFilter<SingleComponentImageType, SingleComponentImageType> >::type
174  typedef typename std::conditional< std::is_same< TOutputImage, CPUOutputImageType >::value,
176  CudaBackProjectionImageFilter<THessiansImage> >::type
178 #else
179  typedef WeidingerForwardModelImageFilter<TOutputImage, TPhotonCounts, TSpectrum> WeidingerForwardModelType;
182  typedef BackProjectionImageFilter<THessiansImage, THessiansImage> CudaHessiansBackProjectionImageFilterType;
183 #endif
184  typedef TOutputImage TGradientsImage;
185 #endif
186 
189 
190 #if !defined( ITK_WRAPPING_PARSER )
191 
195  SingleComponentImageType > SingleComponentForwardProjectionFilterType;
208 #endif
209 
211  void SetForwardProjectionFilter (ForwardProjectionType _arg) ITK_OVERRIDE;
212 
214  void SetBackProjectionFilter (BackProjectionType _arg) ITK_OVERRIDE;
215 
217  itkSetConstObjectMacro(Geometry, ThreeDCircularProjectionGeometry)
218 
219  itkSetMacro(NumberOfIterations, int)
220  itkGetMacro(NumberOfIterations, int)
221 
223  itkSetMacro(NumberOfSubsets, int)
224  itkGetMacro(NumberOfSubsets, int)
225 
229  itkSetMacro(ResetNesterovEvery, int)
230  itkGetMacro(ResetNesterovEvery, int)
231 
233  void SetInputMaterialVolumes(const TOutputImage* materialVolumes);
234  void SetInputPhotonCounts(const TPhotonCounts* photonCounts);
235  void SetInputSpectrum(const TSpectrum* spectrum);
236  void SetSupportMask(const SingleComponentImageType* support);
238 
240  itkSetMacro(RegularizationWeights, typename TOutputImage::PixelType)
241  itkGetMacro(RegularizationWeights, typename TOutputImage::PixelType)
242 
244  itkSetMacro(RegularizationRadius, typename TOutputImage::RegionType::SizeType)
245  itkGetMacro(RegularizationRadius, typename TOutputImage::RegionType::SizeType)
246 
247 // itkSetMacro(IterationCosts, bool)
248 // itkGetMacro(IterationCosts, bool)
249 
252  typedef vnl_matrix<dataType> BinnedDetectorResponseType;
253  typedef vnl_matrix<dataType> MaterialAttenuationsType;
254  virtual void SetBinnedDetectorResponse(const BinnedDetectorResponseType & detResp);
255  virtual void SetMaterialAttenuations(const MaterialAttenuationsType & matAtt);
257 
258 protected:
259  MechlemOneStepSpectralReconstructionFilter();
260  virtual ~MechlemOneStepSpectralReconstructionFilter() ITK_OVERRIDE {}
261 
263  void GenerateData() ITK_OVERRIDE;
264 
265 #if !defined( ITK_WRAPPING_PARSER )
266 
267  typename ExtractPhotonCountsFilterType::Pointer m_ExtractPhotonCountsFilter;
284 #endif
285 
289 #if ITK_VERSION_MAJOR<5
290  void VerifyInputInformation() ITK_OVERRIDE {}
291 #else
292  void VerifyInputInformation() const ITK_OVERRIDE {}
293 #endif
294 
295 
298  void GenerateInputRequestedRegion() ITK_OVERRIDE;
299  void GenerateOutputInformation() ITK_OVERRIDE;
301 
303  typename TOutputImage::ConstPointer GetInputMaterialVolumes();
304  typename TPhotonCounts::ConstPointer GetInputPhotonCounts();
305  typename TSpectrum::ConstPointer GetInputSpectrum();
306  typename SingleComponentImageType::ConstPointer GetSupportMask();
308 
309 #if !defined( ITK_WRAPPING_PARSER )
310 
314 #endif
315 
316 
318 
325 
326  typename TOutputImage::PixelType m_RegularizationWeights;
327  typename TOutputImage::RegionType::SizeType m_RegularizationRadius;
328 
329 private:
330  MechlemOneStepSpectralReconstructionFilter(const Self &); //purposely not implemented
331  void operator=(const Self &); //purposely not implemented
332 
333 // bool m_IterationCosts;
334 };
335 } //namespace ITK
336 
337 
338 #ifndef ITK_MANUAL_INSTANTIATION
339 #include "rtkMechlemOneStepSpectralReconstructionFilter.hxx"
340 #endif
341 
342 #endif
rtk::SeparableQuadraticSurrogateRegularizationImageFilter< TGradientsImage > SQSRegularizationType
void SetInputPhotonCounts(const TPhotonCounts *photonCounts)
Base class for forward projection, i.e. accumulation along x-ray lines.
SingleComponentForwardProjectionFilterType::Pointer InstantiateSingleComponentForwardProjectionFilter(int fwtype)
HessiansBackProjectionFilterType::Pointer InstantiateHessiansBackProjectionFilter(int bptype)
Applies Nesterov&#39;s momentum technique.
Generate an n-dimensional image with constant pixel values.
void SetSupportMask(const SingleComponentImageType *support)
JosephForwardProjectionImageFilter< SingleComponentImageType, SingleComponentImageType > CudaSingleComponentForwardProjectionImageFilterType
rtk::ForwardProjectionImageFilter< TOutputImage, TOutputImage > ForwardProjectionFilterType
void SetBackProjectionFilter(BackProjectionType _arg) override
BackProjectionImageFilter< THessiansImage, THessiansImage > CudaHessiansBackProjectionImageFilterType
virtual void SetMaterialAttenuations(const MaterialAttenuationsType &matAtt)
Projection geometry for a source and a 2-D flat panel.
ImageBaseType::SizeType SizeType
WeidingerForwardModelImageFilter< TOutputImage, TPhotonCounts, TSpectrum > WeidingerForwardModelType
Computes update from gradient and Hessian in Newton&#39;s method.
void SetInputMaterialVolumes(const TOutputImage *materialVolumes)
rtk::ForwardProjectionImageFilter< SingleComponentImageType, SingleComponentImageType > SingleComponentForwardProjectionFilterType
rtk::BackProjectionImageFilter< THessiansImage, THessiansImage > HessiansBackProjectionFilterType
rtk::ConstantImageSource< SingleComponentImageType > SingleComponentImageSourceType
itk::ExtractImageFilter< TPhotonCounts, TPhotonCounts > ExtractPhotonCountsFilterType
void SetInputSpectrum(const TSpectrum *spectrum)
Performs intermediate computations in Weidinger2016.
rtk::AddMatrixAndDiagonalImageFilter< TGradientsImage, THessiansImage > AddMatrixAndDiagonalFilterType
For each vector-valued pixel, adds a vector to the diagonal of a matrix.
For one-step inversion of spectral CT data by the method Mechlem2017, computes regularization term&#39;s ...
Mother class for cone beam reconstruction filters which need runtime selection of their forward and b...
virtual void SetBinnedDetectorResponse(const BinnedDetectorResponseType &detResp)
SingleComponentImageType::ConstPointer GetSupportMask()
rtk::BackProjectionImageFilter< TGradientsImage, TGradientsImage > GradientsBackProjectionFilterType
Implements the one-step spectral CT inversion method described by Mechlem et al.
rtk::GetNewtonUpdateImageFilter< TGradientsImage, THessiansImage > NewtonFilterType
ImageBaseType::RegionType RegionType
SingleComponentForwardProjectionFilterType::Pointer m_SingleComponentForwardProjectionFilter
TPhotonCounts::ConstPointer GetInputPhotonCounts()
TOutputImage::ConstPointer GetInputMaterialVolumes()
void SetForwardProjectionFilter(ForwardProjectionType _arg) override
itk::Image< itk::Vector< dataType, nMaterials *nMaterials >, TOutputImage::ImageDimension > THessiansImage
#define itkSetMacro(name, type)
itk::Image< dataType, TOutputImage::ImageDimension > SingleComponentImageType
IterativeConeBeamReconstructionFilter< TOutputImage, TOutputImage > Superclass
itk::MultiplyImageFilter< TOutputImage, SingleComponentImageType > MultiplyFilterType