RTK  2.0.1
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 
24 
25 #include <itkInPlaceImageFilter.h>
26 
27 namespace rtk
28 {
42 template<typename DecomposedProjectionsType,
43  typename MeasuredProjectionsType,
44  typename IncidentSpectrumImageType = itk::VectorImage<float, 2>,
45  typename DetectorResponseImageType = itk::Image<float, 2>,
46  typename MaterialAttenuationsImageType = itk::Image<float, 2> >
48  public itk::InPlaceImageFilter<MeasuredProjectionsType, MeasuredProjectionsType>
49 {
50 public:
51  ITK_DISALLOW_COPY_AND_ASSIGN(SpectralForwardModelImageFilter);
52 
58 
60  using InputImageType = MeasuredProjectionsType;
61  using OutputImageType = MeasuredProjectionsType;
62 
65  using DetectorResponseType = vnl_matrix<double>;
66  using MaterialAttenuationsType = vnl_matrix<double>;
67 
69  itkNewMacro(Self)
70 
71 
73 
75  void SetInputIncidentSpectrum(const IncidentSpectrumImageType* IncidentSpectrum);
76  void SetInputSecondIncidentSpectrum(const IncidentSpectrumImageType* IncidentSpectrum);
77  typename IncidentSpectrumImageType::ConstPointer GetInputIncidentSpectrum();
78  typename IncidentSpectrumImageType::ConstPointer GetInputSecondIncidentSpectrum();
80 
82  void SetInputDecomposedProjections(const DecomposedProjectionsType* DecomposedProjections);
83  typename DecomposedProjectionsType::ConstPointer GetInputDecomposedProjections();
85 
87  void SetInputMeasuredProjections(const MeasuredProjectionsType* SpectralProjections);
88  typename MeasuredProjectionsType::ConstPointer GetInputMeasuredProjections();
90 
92  void SetDetectorResponse(const DetectorResponseImageType* DetectorResponse);
93  typename DetectorResponseImageType::ConstPointer GetDetectorResponse();
95 
97  void SetMaterialAttenuations(const MaterialAttenuationsImageType* MaterialAttenuations);
98  typename MaterialAttenuationsImageType::ConstPointer GetMaterialAttenuations();
100 
101  itkSetMacro(Thresholds, ThresholdsType)
102  itkGetMacro(Thresholds, ThresholdsType)
103 
104  itkSetMacro(NumberOfSpectralBins, unsigned int)
105  itkGetMacro(NumberOfSpectralBins, unsigned int)
106 
107  itkSetMacro(NumberOfMaterials, unsigned int)
108  itkGetMacro(NumberOfMaterials, unsigned int)
109 
110  itkSetMacro(NumberOfEnergies, unsigned int)
111  itkGetMacro(NumberOfEnergies, unsigned int)
112 
113  itkSetMacro(IsSpectralCT, bool)
114  itkGetMacro(IsSpectralCT, bool)
115 
116  itkSetMacro(ComputeVariances, bool)
117  itkGetMacro(ComputeVariances, bool)
118 
119 protected:
120  SpectralForwardModelImageFilter();
121  ~SpectralForwardModelImageFilter() override = default;
122 
125  using Superclass::MakeOutput;
126  itk::DataObject::Pointer MakeOutput(DataObjectPointerArraySizeType idx) override;
127 
128  void GenerateOutputInformation() override;
129 
130  void GenerateInputRequestedRegion() override;
131 
132  void BeforeThreadedGenerateData() override;
133 #if ITK_VERSION_MAJOR<5
134  void ThreadedGenerateData(const typename OutputImageType::RegionType& outputRegionForThread, itk::ThreadIdType itkNotUsed(threadId)) override;
135 #else
136  void DynamicThreadedGenerateData(const typename OutputImageType::RegionType& outputRegionForThread) override;
137 #endif
138 
141 #if ITK_VERSION_MAJOR<5
142  void VerifyInputInformation() override {}
143 #else
144  void VerifyInputInformation() const override {}
145 #endif
146 
147 
150 
153  unsigned int m_NumberOfEnergies;
154 
156  unsigned int m_NumberOfIterations;
157  unsigned int m_NumberOfMaterials;
159  bool m_IsSpectralCT; // If not, it is dual energy CT
160  bool m_ComputeVariances; // Only implemented for dual energy CT
161 
162 }; // end of class
163 
164 // Function to bin a detector response matrix according to given energy thresholds
165 template<typename OutputElementType, typename DetectorResponseImageType, typename ThresholdsType>
166 vnl_matrix<OutputElementType>
167 SpectralBinDetectorResponse(const DetectorResponseImageType *drm,
168  const ThresholdsType &thresholds,
169  const unsigned int numberOfEnergies);
170 
171 } // end namespace rtk
172 
173 
174 #ifndef ITK_MANUAL_INSTANTIATION
175 #include "rtkSpectralForwardModelImageFilter.hxx"
176 #endif
177 
178 #endif
Forward model for the decomposition of spectral projection images into material projections.
itk::ProcessObject::DataObjectPointerArraySizeType DataObjectPointerArraySizeType
unsigned int ThreadIdType
vnl_matrix< OutputElementType > SpectralBinDetectorResponse(const DetectorResponseImageType *drm, const ThresholdsType &thresholds, const unsigned int numberOfEnergies)
#define itkSetMacro(name, type)