RTK  1.4.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 
34 #ifdef RTK_USE_CUDA
36 #endif
37 
38 namespace rtk
39 {
119 template< typename TOutputImage, typename TPhotonCounts, typename TSpectrum >
121 {
122 public:
127 
129  itkNewMacro(Self)
130 
131 
133 
135  itkStaticConstMacro(nBins, unsigned int, TPhotonCounts::PixelType::Dimension);
136  itkStaticConstMacro(nMaterials, unsigned int, TOutputImage::PixelType::Dimension);
137  typedef typename TOutputImage::PixelType::ValueType dataType;
138 #if !defined( ITK_WRAPPING_PARSER )
139 #ifdef RTK_USE_CUDA
140  typedef itk::CudaImage< itk::Vector<dataType, nMaterials * nMaterials>, TOutputImage::ImageDimension > THessiansImage;
141  typedef itk::CudaImage<dataType, TOutputImage::ImageDimension> TSingleComponentImage;
142 #else
145 #endif
146 
147 
148  typedef TOutputImage TGradientsImage;
149 #endif
150 
153 
154 #if !defined( ITK_WRAPPING_PARSER )
155 
168 #ifdef RTK_USE_CUDA
169  typedef rtk::CudaWeidingerForwardModelImageFilter<TOutputImage, TPhotonCounts, TSpectrum> WeidingerForwardModelType;
170 #else
172 #endif
175 #endif
176 
178  void SetForwardProjectionFilter (ForwardProjectionType _arg) ITK_OVERRIDE;
179 
181  void SetBackProjectionFilter (BackProjectionType _arg) ITK_OVERRIDE;
182 
184  itkSetConstObjectMacro(Geometry, ThreeDCircularProjectionGeometry)
185 
186  itkSetMacro(NumberOfIterations, int)
187  itkGetMacro(NumberOfIterations, int)
188  itkSetMacro(NumberOfSubsets, int)
189  itkGetMacro(NumberOfSubsets, int)
190 
192  void SetInputMaterialVolumes(const TOutputImage* materialVolumes);
193  void SetInputPhotonCounts(const TPhotonCounts* photonCounts);
194  void SetInputSpectrum(const TSpectrum* spectrum);
196 
198  itkSetMacro(RegularizationWeights, typename TOutputImage::PixelType)
199  itkGetMacro(RegularizationWeights, typename TOutputImage::PixelType)
200 
202  itkSetMacro(RegularizationRadius, typename TOutputImage::RegionType::SizeType)
203  itkGetMacro(RegularizationRadius, typename TOutputImage::RegionType::SizeType)
204 
205 // itkSetMacro(IterationCosts, bool)
206 // itkGetMacro(IterationCosts, bool)
207 
210  typedef vnl_matrix<dataType> BinnedDetectorResponseType;
211  typedef vnl_matrix<dataType> MaterialAttenuationsType;
212  virtual void SetBinnedDetectorResponse(const BinnedDetectorResponseType & detResp);
213  virtual void SetMaterialAttenuations(const MaterialAttenuationsType & matAtt);
215 
216 protected:
217  MechlemOneStepSpectralReconstructionFilter();
218  virtual ~MechlemOneStepSpectralReconstructionFilter() ITK_OVERRIDE {}
219 
221  void GenerateData() ITK_OVERRIDE;
222 
223 #if !defined( ITK_WRAPPING_PARSER )
224 
225  typename ExtractPhotonCountsFilterType::Pointer m_ExtractPhotonCountsFilter;
241 #endif
242 
246 #if ITK_VERSION_MAJOR<5
247  void VerifyInputInformation() ITK_OVERRIDE {}
248 #else
249  void VerifyInputInformation() const ITK_OVERRIDE {}
250 #endif
251 
252 
255  void GenerateInputRequestedRegion() ITK_OVERRIDE;
256  void GenerateOutputInformation() ITK_OVERRIDE;
258 
260  typename TOutputImage::ConstPointer GetInputMaterialVolumes();
261  typename TPhotonCounts::ConstPointer GetInputPhotonCounts();
262  typename TSpectrum::ConstPointer GetInputSpectrum();
264 
265 #if !defined( ITK_WRAPPING_PARSER )
266 
270 #endif
271 
272 
274 
280 
281  typename TOutputImage::PixelType m_RegularizationWeights;
282  typename TOutputImage::RegionType::SizeType m_RegularizationRadius;
283 
284 private:
285  MechlemOneStepSpectralReconstructionFilter(const Self &); //purposely not implemented
286  void operator=(const Self &); //purposely not implemented
287 
288 // bool m_IterationCosts;
289 };
290 } //namespace ITK
291 
292 
293 #ifndef ITK_MANUAL_INSTANTIATION
294 #include "rtkMechlemOneStepSpectralReconstructionFilter.hxx"
295 #endif
296 
297 #endif
rtk::SeparableQuadraticSurrogateRegularizationImageFilter< TGradientsImage > SQSRegularizationType
itk::Image< dataType, TOutputImage::ImageDimension > TSingleComponentImage
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.
rtk::ForwardProjectionImageFilter< TOutputImage, TOutputImage > ForwardProjectionFilterType
void SetBackProjectionFilter(BackProjectionType _arg) override
rtk::ForwardProjectionImageFilter< TSingleComponentImage, TSingleComponentImage > SingleComponentForwardProjectionFilterType
virtual void SetMaterialAttenuations(const MaterialAttenuationsType &matAtt)
Projection geometry for a source and a 2-D flat panel.
ImageBaseType::SizeType SizeType
Computes update from gradient and Hessian in Newton&#39;s method.
void SetInputMaterialVolumes(const TOutputImage *materialVolumes)
rtk::BackProjectionImageFilter< THessiansImage, THessiansImage > HessiansBackProjectionFilterType
itk::ExtractImageFilter< TPhotonCounts, TPhotonCounts > ExtractPhotonCountsFilterType
void SetInputSpectrum(const TSpectrum *spectrum)
Performs intermediate computations in Weidinger2016.
rtk::ConstantImageSource< TSingleComponentImage > SingleComponentImageSourceType
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)
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)
rtk::WeidingerForwardModelImageFilter< TOutputImage, TPhotonCounts, TSpectrum > WeidingerForwardModelType
IterativeConeBeamReconstructionFilter< TOutputImage, TOutputImage > Superclass