RTK  1.4.0
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:
54 
56  typedef DecomposedProjectionsType InputImageType;
57  typedef DecomposedProjectionsType OutputImageType;
58 
62  typedef vnl_matrix<double> DetectorResponseType;
63  typedef vnl_matrix<double> MaterialAttenuationsType;
65 
67  itkNewMacro(Self)
68 
69 
71 
73  void SetInputDecomposedProjections(const DecomposedProjectionsType* DecomposedProjections);
74  typename DecomposedProjectionsType::ConstPointer GetInputDecomposedProjections();
76 
78  void SetInputMeasuredProjections(const MeasuredProjectionsType* SpectralProjections);
79  typename MeasuredProjectionsType::ConstPointer GetInputMeasuredProjections();
81 
83  void SetDetectorResponse(const DetectorResponseImageType* DetectorResponse);
84  typename DetectorResponseImageType::ConstPointer GetDetectorResponse();
86 
88  void SetMaterialAttenuations(const MaterialAttenuationsImageType* MaterialAttenuations);
89  typename MaterialAttenuationsImageType::ConstPointer GetMaterialAttenuations();
91 
93  void SetInputIncidentSpectrum(const IncidentSpectrumImageType* IncidentSpectrum);
94  void SetInputSecondIncidentSpectrum(const IncidentSpectrumImageType* IncidentSpectrum);
95  typename IncidentSpectrumImageType::ConstPointer GetInputIncidentSpectrum();
96  typename IncidentSpectrumImageType::ConstPointer GetInputSecondIncidentSpectrum();
98 
100  itkGetMacro(NumberOfIterations, unsigned int)
101  itkSetMacro(NumberOfIterations, unsigned int)
102 
103  itkSetMacro(NumberOfEnergies, unsigned int)
104  itkGetMacro(NumberOfEnergies, unsigned int)
105 
106  itkSetMacro(NumberOfMaterials, unsigned int)
107  itkGetMacro(NumberOfMaterials, unsigned int)
108 
109  itkSetMacro(OptimizeWithRestarts, bool)
110  itkGetMacro(OptimizeWithRestarts, bool)
111 
112  itkSetMacro(Thresholds, ThresholdsType)
113  itkGetMacro(Thresholds, ThresholdsType)
114 
115  itkSetMacro(NumberOfSpectralBins, unsigned int)
116  itkGetMacro(NumberOfSpectralBins, unsigned int)
117 
118  itkSetMacro(OutputInverseCramerRaoLowerBound, bool)
119  itkGetMacro(OutputInverseCramerRaoLowerBound, bool)
120 
121  itkSetMacro(OutputFischerMatrix, bool)
122  itkGetMacro(OutputFischerMatrix, bool)
123 
124  itkSetMacro(LogTransformEachBin, bool)
125  itkGetMacro(LogTransformEachBin, bool)
126 
127  itkSetMacro(GuessInitialization, bool)
128  itkGetMacro(GuessInitialization, bool)
129 
130  itkSetMacro(IsSpectralCT, bool)
131  itkGetMacro(IsSpectralCT, bool)
132 
133 protected:
134  SimplexSpectralProjectionsDecompositionImageFilter();
135  virtual ~SimplexSpectralProjectionsDecompositionImageFilter() ITK_OVERRIDE {}
136 
137  void GenerateOutputInformation() ITK_OVERRIDE;
138 
139  void GenerateInputRequestedRegion() ITK_OVERRIDE;
140 
141  void BeforeThreadedGenerateData() ITK_OVERRIDE;
142 #if ITK_VERSION_MAJOR<5
143  void ThreadedGenerateData(const typename DecomposedProjectionsType::RegionType& outputRegionForThread, itk::ThreadIdType itkNotUsed(threadId)) ITK_OVERRIDE;
144 #else
145  void DynamicThreadedGenerateData(const typename DecomposedProjectionsType::RegionType& outputRegionForThread) ITK_OVERRIDE;
146 #endif
147 
149  typedef itk::ProcessObject::DataObjectPointerArraySizeType DataObjectPointerArraySizeType;
150  using Superclass::MakeOutput;
151  itk::DataObject::Pointer MakeOutput(DataObjectPointerArraySizeType idx) ITK_OVERRIDE;
152 
155  void VerifyInputInformation() ITK_OVERRIDE {}
156 
158  MaterialAttenuationsType m_MaterialAttenuations;
159  DetectorResponseType m_DetectorResponse;
160  ThresholdsType m_Thresholds;
161  MeanAttenuationInBinType m_MeanAttenuationInBin;
166  bool m_IsSpectralCT; //If not, it is dual energy CT
168  unsigned int m_NumberOfIterations;
169  unsigned int m_NumberOfMaterials;
170  unsigned int m_NumberOfEnergies;
172 
173 private:
174  //purposely not implemented
176  void operator=(const Self&);
177 
178 }; // end of class
179 
180 } // end namespace rtk
181 
182 
183 #ifndef ITK_MANUAL_INSTANTIATION
184 #include "rtkSimplexSpectralProjectionsDecompositionImageFilter.hxx"
185 #endif
186 
187 #endif
Decomposition of spectral projection images into material projections.
itk::ImageToImageFilter< DecomposedProjectionsType, DecomposedProjectionsType > Superclass
unsigned int ThreadIdType
#define itkSetMacro(name, type)