RTK  2.0.1
Reconstruction Toolkit
rtkSimplexSpectralProjectionsDecompositionImageFilter.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 rtkSimplexSpectralProjectionsDecompositionImageFilter_h
20 #define rtkSimplexSpectralProjectionsDecompositionImageFilter_h
21 
22 #include <itkImageToImageFilter.h>
23 #include <itkAmoebaOptimizer.h>
26 
27 namespace rtk
28 {
40 template<typename DecomposedProjectionsType,
41  typename MeasuredProjectionsType,
42  typename IncidentSpectrumImageType = itk::VectorImage<float, 2>,
43  typename DetectorResponseImageType = itk::Image<float, 2>,
44  typename MaterialAttenuationsImageType = itk::Image<float, 2> >
46  public itk::ImageToImageFilter<DecomposedProjectionsType, DecomposedProjectionsType>
47 {
48 public:
49  ITK_DISALLOW_COPY_AND_ASSIGN(SimplexSpectralProjectionsDecompositionImageFilter);
50 
56 
58  using InputImageType = DecomposedProjectionsType;
59  using OutputImageType = DecomposedProjectionsType;
60 
64  using DetectorResponseType = vnl_matrix<double>;
65  using MaterialAttenuationsType = vnl_matrix<double>;
67 
69  itkNewMacro(Self)
70 
71 
73 
75  void SetInputDecomposedProjections(const DecomposedProjectionsType* DecomposedProjections);
76  typename DecomposedProjectionsType::ConstPointer GetInputDecomposedProjections();
78 
80  void SetInputMeasuredProjections(const MeasuredProjectionsType* SpectralProjections);
81  typename MeasuredProjectionsType::ConstPointer GetInputMeasuredProjections();
83 
85  void SetDetectorResponse(const DetectorResponseImageType* DetectorResponse);
86  typename DetectorResponseImageType::ConstPointer GetDetectorResponse();
88 
90  void SetMaterialAttenuations(const MaterialAttenuationsImageType* MaterialAttenuations);
91  typename MaterialAttenuationsImageType::ConstPointer GetMaterialAttenuations();
93 
95  void SetInputIncidentSpectrum(const IncidentSpectrumImageType* IncidentSpectrum);
96  void SetInputSecondIncidentSpectrum(const IncidentSpectrumImageType* IncidentSpectrum);
97  typename IncidentSpectrumImageType::ConstPointer GetInputIncidentSpectrum();
98  typename IncidentSpectrumImageType::ConstPointer GetInputSecondIncidentSpectrum();
100 
102  itkGetMacro(NumberOfIterations, unsigned int)
103  itkSetMacro(NumberOfIterations, unsigned int)
104 
105  itkSetMacro(NumberOfEnergies, unsigned int)
106  itkGetMacro(NumberOfEnergies, unsigned int)
107 
108  itkSetMacro(NumberOfMaterials, unsigned int)
109  itkGetMacro(NumberOfMaterials, unsigned int)
110 
111  itkSetMacro(OptimizeWithRestarts, bool)
112  itkGetMacro(OptimizeWithRestarts, bool)
113 
114  itkSetMacro(Thresholds, ThresholdsType)
115  itkGetMacro(Thresholds, ThresholdsType)
116 
117  itkSetMacro(NumberOfSpectralBins, unsigned int)
118  itkGetMacro(NumberOfSpectralBins, unsigned int)
119 
120  itkSetMacro(OutputInverseCramerRaoLowerBound, bool)
121  itkGetMacro(OutputInverseCramerRaoLowerBound, bool)
122 
123  itkSetMacro(OutputFischerMatrix, bool)
124  itkGetMacro(OutputFischerMatrix, bool)
125 
126  itkSetMacro(LogTransformEachBin, bool)
127  itkGetMacro(LogTransformEachBin, bool)
128 
129  itkSetMacro(GuessInitialization, bool)
130  itkGetMacro(GuessInitialization, bool)
131 
132  itkSetMacro(IsSpectralCT, bool)
133  itkGetMacro(IsSpectralCT, bool)
134 
135 protected:
136  SimplexSpectralProjectionsDecompositionImageFilter();
137  ~SimplexSpectralProjectionsDecompositionImageFilter() override = default;
138 
139  void GenerateOutputInformation() override;
140 
141  void GenerateInputRequestedRegion() override;
142 
143  void BeforeThreadedGenerateData() override;
144 #if ITK_VERSION_MAJOR<5
145  void ThreadedGenerateData(const typename DecomposedProjectionsType::RegionType& outputRegionForThread, itk::ThreadIdType itkNotUsed(threadId)) override;
146 #else
147  void DynamicThreadedGenerateData(const typename DecomposedProjectionsType::RegionType& outputRegionForThread) override;
148 #endif
149 
151  using DataObjectPointerArraySizeType = itk::ProcessObject::DataObjectPointerArraySizeType;
152  using Superclass::MakeOutput;
153  itk::DataObject::Pointer MakeOutput(DataObjectPointerArraySizeType idx) override;
154 
157 #if ITK_VERSION_MAJOR<5
158  void VerifyInputInformation() override {}
159 #else
160  void VerifyInputInformation() const override {}
161 #endif
162 
163 
173  bool m_IsSpectralCT; //If not, it is dual energy CT
175  unsigned int m_NumberOfIterations;
176  unsigned int m_NumberOfMaterials;
177  unsigned int m_NumberOfEnergies;
179 
180 }; // end of class
181 
182 } // end namespace rtk
183 
184 
185 #ifndef ITK_MANUAL_INSTANTIATION
186 #include "rtkSimplexSpectralProjectionsDecompositionImageFilter.hxx"
187 #endif
188 
189 #endif
Decomposition of spectral projection images into material projections.
unsigned int ThreadIdType
#define itkSetMacro(name, type)