RTK  2.4.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  * 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 
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>>
47 class ITK_TEMPLATE_EXPORT SpectralForwardModelImageFilter
48  : public itk::InPlaceImageFilter<MeasuredProjectionsType, MeasuredProjectionsType>
49 {
50 public:
51  ITK_DISALLOW_COPY_AND_MOVE(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 
73 
75  void
76  SetInputIncidentSpectrum(const IncidentSpectrumImageType * IncidentSpectrum);
77  void
78  SetInputSecondIncidentSpectrum(const IncidentSpectrumImageType * SecondIncidentSpectrum);
79  typename IncidentSpectrumImageType::ConstPointer
80  GetInputIncidentSpectrum();
81  typename IncidentSpectrumImageType::ConstPointer
82  GetInputSecondIncidentSpectrum();
84 
86  void
87  SetInputDecomposedProjections(const DecomposedProjectionsType * DecomposedProjections);
88  typename DecomposedProjectionsType::ConstPointer
89  GetInputDecomposedProjections();
91 
93  void
94  SetInputMeasuredProjections(const MeasuredProjectionsType * SpectralProjections);
95  typename MeasuredProjectionsType::ConstPointer
96  GetInputMeasuredProjections();
98 
100  void
101  SetDetectorResponse(const DetectorResponseImageType * DetectorResponse);
102  typename DetectorResponseImageType::ConstPointer
103  GetDetectorResponse();
105 
107  void
108  SetMaterialAttenuations(const MaterialAttenuationsImageType * MaterialAttenuations);
109  typename MaterialAttenuationsImageType::ConstPointer
110  GetMaterialAttenuations();
112 
113  itkSetMacro(Thresholds, ThresholdsType);
114  itkGetMacro(Thresholds, ThresholdsType);
115 
116  itkSetMacro(NumberOfSpectralBins, unsigned int);
117  itkGetMacro(NumberOfSpectralBins, unsigned int);
118 
119  itkSetMacro(NumberOfMaterials, unsigned int);
120  itkGetMacro(NumberOfMaterials, unsigned int);
121 
122  itkSetMacro(NumberOfEnergies, unsigned int);
123  itkGetMacro(NumberOfEnergies, unsigned int);
124 
125  itkSetMacro(IsSpectralCT, bool);
126  itkGetMacro(IsSpectralCT, bool);
127 
128  itkSetMacro(ComputeVariances, bool);
129  itkGetMacro(ComputeVariances, bool);
130 
131 protected:
133  ~SpectralForwardModelImageFilter() override = default;
134 
137  using Superclass::MakeOutput;
139  MakeOutput(DataObjectPointerArraySizeType idx) override;
140 
141  void
142  GenerateOutputInformation() override;
143 
144  void
145  GenerateInputRequestedRegion() override;
146 
147  void
148  BeforeThreadedGenerateData() override;
149  void
150  DynamicThreadedGenerateData(const typename OutputImageType::RegionType & outputRegionForThread) override;
151 
154  void
155  VerifyInputInformation() const override
156  {}
157 
160 
163  unsigned int m_NumberOfEnergies;
164 
166  unsigned int m_NumberOfIterations;
167  unsigned int m_NumberOfMaterials;
169  bool m_IsSpectralCT; // If not, it is dual energy CT
170  bool m_ComputeVariances; // Only implemented for dual energy CT
171 
172 }; // end of class
173 
174 // Function to bin a detector response matrix according to given energy thresholds
175 template <typename OutputElementType, typename DetectorResponseImageType, typename ThresholdsType>
176 vnl_matrix<OutputElementType>
177 SpectralBinDetectorResponse(const DetectorResponseImageType * drm,
178  const ThresholdsType & thresholds,
179  const unsigned int numberOfEnergies);
180 
181 } // end namespace rtk
182 
183 
184 #ifndef ITK_MANUAL_INSTANTIATION
185 # include "rtkSpectralForwardModelImageFilter.hxx"
186 #endif
187 
188 #endif
Superclass::DataObjectPointerArraySizeType DataObjectPointerArraySizeType
DataObjectPointerArray::size_type DataObjectPointerArraySizeType
Forward model for the decomposition of spectral projection images into material projections.
TInputImage InputImageType
#define itkSetMacro(name, type)
TOutputImage OutputImageType
vnl_matrix< OutputElementType > SpectralBinDetectorResponse(const DetectorResponseImageType *drm, const ThresholdsType &thresholds, const unsigned int numberOfEnergies)