RTK  1.4.0
Reconstruction Toolkit
rtkSpectralForwardModelImageFilter.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 rtkSpectralForwardModelImageFilter_h
20 #define rtkSpectralForwardModelImageFilter_h
21 
25 
26 #include <itkInPlaceImageFilter.h>
27 
28 namespace rtk
29 {
43 template<typename DecomposedProjectionsType,
44  typename MeasuredProjectionsType,
45  typename IncidentSpectrumImageType = itk::VectorImage<float, 2>,
46  typename DetectorResponseImageType = itk::Image<float, 2>,
47  typename MaterialAttenuationsImageType = itk::Image<float, 2> >
49  public itk::InPlaceImageFilter<MeasuredProjectionsType, MeasuredProjectionsType>
50 {
51 public:
57 
59  typedef MeasuredProjectionsType InputImageType;
60  typedef MeasuredProjectionsType OutputImageType;
61 
64  typedef vnl_matrix<double> DetectorResponseType;
65  typedef vnl_matrix<double> MaterialAttenuationsType;
66 
68  itkNewMacro(Self)
69 
70 
72 
74  void SetInputIncidentSpectrum(const IncidentSpectrumImageType* IncidentSpectrum);
75  void SetInputSecondIncidentSpectrum(const IncidentSpectrumImageType* IncidentSpectrum);
76  typename IncidentSpectrumImageType::ConstPointer GetInputIncidentSpectrum();
77  typename IncidentSpectrumImageType::ConstPointer GetInputSecondIncidentSpectrum();
79 
81  void SetInputDecomposedProjections(const DecomposedProjectionsType* DecomposedProjections);
82  typename DecomposedProjectionsType::ConstPointer GetInputDecomposedProjections();
84 
86  void SetInputMeasuredProjections(const MeasuredProjectionsType* SpectralProjections);
87  typename MeasuredProjectionsType::ConstPointer GetInputMeasuredProjections();
89 
91  void SetDetectorResponse(const DetectorResponseImageType* DetectorResponse);
92  typename DetectorResponseImageType::ConstPointer GetDetectorResponse();
94 
96  void SetMaterialAttenuations(const MaterialAttenuationsImageType* MaterialAttenuations);
97  typename MaterialAttenuationsImageType::ConstPointer GetMaterialAttenuations();
99 
100  itkSetMacro(Thresholds, ThresholdsType)
101  itkGetMacro(Thresholds, ThresholdsType)
102 
103  itkSetMacro(NumberOfSpectralBins, unsigned int)
104  itkGetMacro(NumberOfSpectralBins, unsigned int)
105 
106  itkSetMacro(NumberOfMaterials, unsigned int)
107  itkGetMacro(NumberOfMaterials, unsigned int)
108 
109  itkSetMacro(NumberOfEnergies, unsigned int)
110  itkGetMacro(NumberOfEnergies, unsigned int)
111 
112  itkSetMacro(IsSpectralCT, bool)
113  itkGetMacro(IsSpectralCT, bool)
114 
115  itkSetMacro(ComputeVariances, bool)
116  itkGetMacro(ComputeVariances, bool)
117 
118 protected:
119  SpectralForwardModelImageFilter();
120  virtual ~SpectralForwardModelImageFilter() ITK_OVERRIDE {}
121 
123  typedef itk::ProcessObject::DataObjectPointerArraySizeType DataObjectPointerArraySizeType;
124  using Superclass::MakeOutput;
125  itk::DataObject::Pointer MakeOutput(DataObjectPointerArraySizeType idx) ITK_OVERRIDE;
126 
127  void GenerateOutputInformation() ITK_OVERRIDE;
128 
129  void GenerateInputRequestedRegion() ITK_OVERRIDE;
130 
131  void BeforeThreadedGenerateData() ITK_OVERRIDE;
132 #if ITK_VERSION_MAJOR<5
133  void ThreadedGenerateData(const typename OutputImageType::RegionType& outputRegionForThread, itk::ThreadIdType itkNotUsed(threadId)) ITK_OVERRIDE;
134 #else
135  void DynamicThreadedGenerateData(const typename OutputImageType::RegionType& outputRegionForThread) ITK_OVERRIDE;
136 #endif
137 
140 #if ITK_VERSION_MAJOR<5
141  void VerifyInputInformation() ITK_OVERRIDE {}
142 #else
143  void VerifyInputInformation() const ITK_OVERRIDE {}
144 #endif
145 
146 
147  ThresholdsType m_Thresholds;
149 
150  MaterialAttenuationsType m_MaterialAttenuations;
151  DetectorResponseType m_DetectorResponse;
152  unsigned int m_NumberOfEnergies;
153 
155  unsigned int m_NumberOfIterations;
156  unsigned int m_NumberOfMaterials;
158  bool m_IsSpectralCT; // If not, it is dual energy CT
159  bool m_ComputeVariances; // Only implemented for dual energy CT
160 
161 private:
162  //purposely not implemented
163  SpectralForwardModelImageFilter(const Self&);
164  void operator=(const Self&);
165 
166 }; // end of class
167 
168 } // end namespace rtk
169 
170 
171 #ifndef ITK_MANUAL_INSTANTIATION
172 #include "rtkSpectralForwardModelImageFilter.hxx"
173 #endif
174 
175 #endif
Forward model for the decomposition of spectral projection images into material projections.
itk::ProcessObject::DataObjectPointerArraySizeType DataObjectPointerArraySizeType
unsigned int ThreadIdType
itk::ImageToImageFilter< MeasuredProjectionsType, MeasuredProjectionsType > Superclass
#define itkSetMacro(name, type)