RTK  2.6.0
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  * 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 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>
46 class ITK_TEMPLATE_EXPORT I0EstimationProjectionFilter : public itk::InPlaceImageFilter<TInputImage, TOutputImage>
47 {
48 public:
49  ITK_DISALLOW_COPY_AND_MOVE(I0EstimationProjectionFilter);
50 
56 
58  itkNewMacro(Self);
59 
61 #ifdef itkOverrideGetNameOfClassMacro
62  itkOverrideGetNameOfClassMacro(I0EstimationProjectionFilter);
63 #else
65 #endif
66 
67 
69  using InputImageType = TInputImage;
70  using ImagePointer = typename InputImageType::Pointer;
71  using ImageConstPointer = typename InputImageType::ConstPointer;
72  using OutputImageRegionType = typename Superclass::OutputImageRegionType;
73  using InputImagePixelType = typename InputImageType::PixelType;
74 
75  itkConceptMacro(InputImagePixelTypeIsInteger, (itk::Concept::IsInteger<InputImagePixelType>));
76 
78  itkGetMacro(I0, InputImagePixelType);
79  itkGetMacro(I0fwhm, InputImagePixelType);
80  itkGetMacro(I0rls, InputImagePixelType);
82 
86  itkSetMacro(MaxPixelValue, InputImagePixelType);
87  itkGetMacro(MaxPixelValue, InputImagePixelType);
89 
91  itkSetMacro(ExpectedI0, InputImagePixelType);
92  itkGetMacro(ExpectedI0, InputImagePixelType);
94 
96  itkSetMacro(Lambda, float);
97  itkGetMacro(Lambda, float);
99 
102  itkSetMacro(Reset, bool);
103  itkGetConstMacro(Reset, bool);
104  itkBooleanMacro(Reset);
106 
109  itkSetMacro(SaveHistograms, bool);
110  itkGetConstMacro(SaveHistograms, bool);
111  itkBooleanMacro(SaveHistograms);
113 
114 protected:
116  ~I0EstimationProjectionFilter() override = default;
117 
118  void
119  BeforeThreadedGenerateData() override;
120 
121  void
122  ThreadedGenerateData(const OutputImageRegionType & outputRegionForThread, ThreadIdType threadId) override;
123 
124  void
125  AfterThreadedGenerateData() override;
126 
127 private:
128  // Input variables
129  InputImagePixelType m_ExpectedI0; // Expected I0 value (as a result of a
130  // detector calibration)
131  InputImagePixelType m_MaxPixelValue; // Maximum encodable detector value if
132  // different from (2^16-1)
133  float m_Lambda; // RLS coefficient
134  bool m_SaveHistograms; // Save histograms in a output file
135  bool m_Reset; // Reset counters
136 
137  // Secondary inputs
138  std::vector<unsigned int>::size_type m_NBins; // Histogram size, computed
139  // from m_MaxPixelValue and bitshift
140 
141  // Main variables
142  std::vector<unsigned int> m_Histogram; // compressed (bitshifted) histogram
143  InputImagePixelType m_I0; // I0 estimate with no a priori for
144  // each new image
145  InputImagePixelType m_I0rls; // Updated RLS estimate
146  InputImagePixelType m_I0fwhm; // FWHM of the I0 mode
147 
148  // Secondary variables
149  unsigned int m_Np; // Number of previously analyzed
150  // images
151  InputImagePixelType m_Imin, m_Imax; // Define the range of consistent
152  // pixels in histogram
153  unsigned int m_DynThreshold; // Detector values with a frequency of
154  // less than dynThreshold outside
155  // min/max are discarded
156  unsigned int m_LowBound, m_HighBound; // Lower/Upper bounds of the I0 mode
157  // at half width
158 
159  std::mutex m_Mutex;
160  int m_Nsync;
162 };
163 } // end namespace rtk
164 
165 #ifndef ITK_MANUAL_INSTANTIATION
166 # include "rtkI0EstimationProjectionFilter.hxx"
167 #endif
168 
169 #endif
Estimate the I0 value from the projection histograms.
typename InputImageType::ConstPointer ImageConstPointer
TInputImage InputImageType
typename InputImageType::PixelType InputImagePixelType
#define itkSetMacro(name, type)
typename OutputImageType::RegionType OutputImageRegionType
std::vector< unsigned int >::size_type m_NBins
typename InputImageType::Pointer ImagePointer
unsigned int ThreadIdType
#define itkConceptMacro(name, concept)