RTK  2.0.1
Reconstruction Toolkit
rtkI0EstimationProjectionFilter.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 rtkI0EstimationProjectionFilter_h
20 #define rtkI0EstimationProjectionFilter_h
21 
22 #include <itkInPlaceImageFilter.h>
24 #include "rtkConfiguration.h"
25 
26 #include <mutex>
27 #include <vector>
28 #include <string>
29 
30 namespace rtk
31 {
43 template<class TInputImage = itk::Image< unsigned short, 3>,
44  class TOutputImage = TInputImage,
45  unsigned char bitShift = 2 >
47  public itk::InPlaceImageFilter<TInputImage, TOutputImage>
48 {
49 public:
50  ITK_DISALLOW_COPY_AND_ASSIGN(I0EstimationProjectionFilter);
51 
57 
59  itkNewMacro(Self);
60 
63 
65  using InputImageType = TInputImage;
66  using ImagePointer = typename InputImageType::Pointer;
67  using ImageConstPointer = typename InputImageType::ConstPointer;
68  using OutputImageRegionType = typename Superclass::OutputImageRegionType;
69  using InputImagePixelType = typename InputImageType::PixelType;
70 
71  itkConceptMacro(InputImagePixelTypeIsInteger, (itk::Concept::IsInteger<InputImagePixelType>) );
72 
74  itkGetMacro(I0, InputImagePixelType)
75  itkGetMacro(I0fwhm, InputImagePixelType)
76  itkGetMacro(I0rls, InputImagePixelType)
77 
81  itkSetMacro(MaxPixelValue, InputImagePixelType)
82  itkGetMacro(MaxPixelValue, InputImagePixelType)
83 
85  itkSetMacro(ExpectedI0, InputImagePixelType)
86  itkGetMacro(ExpectedI0, InputImagePixelType)
87 
89  itkSetMacro(Lambda, float)
90  itkGetMacro(Lambda, float)
91 
94  itkSetMacro(Reset, bool);
95  itkGetConstMacro(Reset, bool);
96  itkBooleanMacro(Reset);
98 
101  itkSetMacro(SaveHistograms, bool);
102  itkGetConstMacro(SaveHistograms, bool);
103  itkBooleanMacro(SaveHistograms);
105 
106 protected:
108  ~I0EstimationProjectionFilter() override = default;
109 
110  void BeforeThreadedGenerateData() override;
111 
112  void ThreadedGenerateData(const OutputImageRegionType & outputRegionForThread, ThreadIdType threadId) override;
113 
114  void AfterThreadedGenerateData() override;
115 
116 private:
117  // Input variables
118  InputImagePixelType m_ExpectedI0; // Expected I0 value (as a result of a
119  // detector calibration)
120  InputImagePixelType m_MaxPixelValue; // Maximum encodable detector value if
121  // different from (2^16-1)
122  float m_Lambda; // RLS coefficient
123  bool m_SaveHistograms; // Save histograms in a output file
124  bool m_Reset; // Reset counters
125 
126  // Secondary inputs
127  std::vector<unsigned int>::size_type m_NBins; // Histogram size, computed
128  // from m_MaxPixelValue and bitshift
129 
130  // Main variables
131  std::vector< unsigned int > m_Histogram; // compressed (bitshifted) histogram
132  InputImagePixelType m_I0; // I0 estimate with no a priori for
133  // each new image
134  InputImagePixelType m_I0rls; // Updated RLS estimate
135  InputImagePixelType m_I0fwhm; // FWHM of the I0 mode
136 
137  // Secondary variables
138  unsigned int m_Np; // Number of previously analyzed
139  // images
140  InputImagePixelType m_Imin, m_Imax; // Define the range of consistent
141  // pixels in histogram
142  unsigned int m_DynThreshold; // Detector values with a frequency of
143  // less than dynThreshold outside
144  // min/max are discarded
145  unsigned int m_LowBound, m_HighBound; // Lower/Upper bounds of the I0 mode
146  // at half width
147 
148  std::mutex m_Mutex;
149  int m_Nsync;
150  int m_Nthreads;
151 };
152 } // end namespace rtk
153 
154 #ifndef ITK_MANUAL_INSTANTIATION
155 #include "rtkI0EstimationProjectionFilter.hxx"
156 #endif
157 
158 #endif
STL namespace.
Estimate the I0 value from the projection histograms.
typename InputImageType::ConstPointer ImageConstPointer
typename Superclass::OutputImageRegionType OutputImageRegionType
unsigned int ThreadIdType
typename InputImageType::Pointer ImagePointer
#define itkConceptMacro(name, concept)
typename InputImageType::PixelType InputImagePixelType
#define itkSetMacro(name, type)