RTK  1.4.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  * 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>
23 #include <itkMutexLock.h>
25 #include "rtkConfiguration.h"
26 
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:
55 
57  itkNewMacro(Self);
58 
61 
63  typedef TInputImage InputImageType;
64  typedef typename InputImageType::Pointer ImagePointer;
65  typedef typename InputImageType::ConstPointer ImageConstPointer;
66  typedef typename Superclass::OutputImageRegionType OutputImageRegionType;
67  typedef typename InputImageType::PixelType InputImagePixelType;
68 
69  itkConceptMacro(InputImagePixelTypeIsInteger, (itk::Concept::IsInteger<InputImagePixelType>) );
70 
72  itkGetMacro(I0, InputImagePixelType)
73  itkGetMacro(I0fwhm, InputImagePixelType)
74  itkGetMacro(I0rls, InputImagePixelType)
75 
79  itkSetMacro(MaxPixelValue, InputImagePixelType)
80  itkGetMacro(MaxPixelValue, InputImagePixelType)
81 
83  itkSetMacro(ExpectedI0, InputImagePixelType)
84  itkGetMacro(ExpectedI0, InputImagePixelType)
85 
87  itkSetMacro(Lambda, float)
88  itkGetMacro(Lambda, float)
89 
92  itkSetMacro(Reset, bool);
93  itkGetConstMacro(Reset, bool);
94  itkBooleanMacro(Reset);
96 
99  itkSetMacro(SaveHistograms, bool);
100  itkGetConstMacro(SaveHistograms, bool);
101  itkBooleanMacro(SaveHistograms);
103 
104 protected:
106  virtual ~I0EstimationProjectionFilter() ITK_OVERRIDE {}
107 
108  void BeforeThreadedGenerateData() ITK_OVERRIDE;
109 
110  void ThreadedGenerateData(const OutputImageRegionType & outputRegionForThread, ThreadIdType threadId) ITK_OVERRIDE;
111 
112  void AfterThreadedGenerateData() ITK_OVERRIDE;
113 
114 private:
115  I0EstimationProjectionFilter(const Self &); //purposely not implemented
116  void operator=(const Self &); //purposely not implemented
117 
118  // Input variables
119  InputImagePixelType m_ExpectedI0; // Expected I0 value (as a result of a
120  // detector calibration)
121  InputImagePixelType m_MaxPixelValue; // Maximum encodable detector value if
122  // different from (2^16-1)
123  float m_Lambda; // RLS coefficient
124  bool m_SaveHistograms; // Save histograms in a output file
125  bool m_Reset; // Reset counters
126 
127  // Secondary inputs
128  std::vector<unsigned int>::size_type m_NBins; // Histogram size, computed
129  // from m_MaxPixelValue and bitshift
130 
131  // Main variables
132  std::vector< unsigned int > m_Histogram; // compressed (bitshifted) histogram
133  InputImagePixelType m_I0; // I0 estimate with no a priori for
134  // each new image
135  InputImagePixelType m_I0rls; // Updated RLS estimate
136  InputImagePixelType m_I0fwhm; // FWHM of the I0 mode
137 
138  // Secondary variables
139  unsigned int m_Np; // Number of previously analyzed
140  // images
141  InputImagePixelType m_Imin, m_Imax; // Define the range of consistent
142  // pixels in histogram
143  unsigned int m_DynThreshold; // Detector values with a frequency of
144  // less than dynThreshold outside
145  // min/max are discarded
146  unsigned int m_LowBound, m_HighBound; // Lower/Upper bounds of the I0 mode
147  // at half width
148 
150  int m_Nsync;
152 };
153 } // end namespace rtk
154 
155 #ifndef ITK_MANUAL_INSTANTIATION
156 #include "rtkI0EstimationProjectionFilter.hxx"
157 #endif
158 
159 #endif
I0EstimationProjectionFilter< TInputImage, TOutputImage, bitShift > Self
Estimate the I0 value from the projection histograms.
itk::InPlaceImageFilter< TInputImage, TOutputImage > Superclass
itk::SmartPointer< const Self > ConstPointer
std::vector< unsigned int >::size_type m_NBins
unsigned int ThreadIdType
Superclass::OutputImageRegionType OutputImageRegionType
#define itkConceptMacro(name, concept)
#define itkSetMacro(name, type)