RTK  1.4.0
Reconstruction Toolkit
rtkLookupTableImageFilter.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 rtkLookupTableImageFilter_h
20 #define rtkLookupTableImageFilter_h
21 
24 
25 namespace rtk
26 {
27 
43 namespace Functor
44 {
45 
46 template< class TInput, class TOutput >
47 class LUT
48 {
49 public:
51  typedef typename LookupTableType::Pointer LookupTablePointer;
52  typedef typename LookupTableType::PixelType* LookupTableDataPointer;
55 
56  LUT():
57  m_LookupTableDataPointer(ITK_NULLPTR),
58  m_Interpolator(InterpolatorType::New())
59  {};
60  ~LUT() {};
61 
63  LookupTablePointer GetLookupTable() {
64  return m_LookupTablePointer;
65  }
66  void SetLookupTable(LookupTablePointer lut) {
68  m_LookupTableDataPointer = lut->GetBufferPointer();
69  m_InverseLUTSpacing = 1. / m_LookupTablePointer->GetSpacing()[0];
70  m_Interpolator->SetInputImage(lut);
71  }
73 
74  LookupTableDataPointer GetLookupTable() const {
76  }
77 
78  bool operator!=( const LUT & lut ) const {
80  }
81  bool operator==( const LUT & lut ) const {
83  }
84 
85  inline TOutput operator()( const TInput & val ) const;
86 
87 private:
88  LookupTablePointer m_LookupTablePointer;
89  LookupTableDataPointer m_LookupTableDataPointer;
90  InterpolatorPointer m_Interpolator;
92 };
93 
94 template< class TInput, class TOutput >
95 TOutput
97 ::operator()( const TInput & val ) const {
98  return m_LookupTableDataPointer[val];
99 }
100 
101 template<>
102 inline
103 float
105 ::operator()( const float & val ) const {
107  index[0] = m_InverseLUTSpacing * (val - m_LookupTablePointer->GetOrigin()[0]);
108  return float(m_Interpolator->EvaluateAtContinuousIndex(index));
109 }
110 
111 template<>
112 inline
113 double
115 ::operator()( const double & val ) const {
117  index[0] = m_InverseLUTSpacing * (val - m_LookupTablePointer->GetOrigin()[0]);
118  return double(m_Interpolator->EvaluateAtContinuousIndex(index));
119 }
120 
121 } // end namespace Functor
122 
136 template <class TInputImage, class TOutputImage>
137 class ITK_EXPORT LookupTableImageFilter : public
138  itk::UnaryFunctorImageFilter< TInputImage,
139  TOutputImage,
140  Functor::LUT< typename TInputImage::PixelType,
141  typename TOutputImage::PixelType> >
142 {
143 
144 public:
148 
154 
156  itkNewMacro(Self);
157 
160 
162  virtual void SetLookupTable(LookupTableType* _arg) {
163  //Idem as itkSetObjectMacro + call to functor SetLookupTableDataPointer
164  itkDebugMacro("setting " << "LookupTable" " to " << _arg );
165  if (this->m_LookupTable != _arg || ( _arg && _arg->GetTimeStamp()>this->GetTimeStamp() ) ) {
166  this->m_LookupTable = _arg;
167  this->Modified();
168  this->GetFunctor().SetLookupTable(_arg);
169  }
170  }
172 
174  itkGetModifiableObjectMacro(LookupTable, LookupTableType);
175 
178  void BeforeThreadedGenerateData() ITK_OVERRIDE;
179 
180 protected:
182  virtual ~LookupTableImageFilter() ITK_OVERRIDE {}
183  typename LookupTableType::Pointer m_LookupTable;
184 
185 private:
186  LookupTableImageFilter(const Self&); //purposely not implemented
187  void operator=(const Self&); //purposely not implemented
188 
189 };
190 
191 } // end namespace rtk
192 
193 #ifndef ITK_MANUAL_INSTANTIATION
194 #include "rtkLookupTableImageFilter.hxx"
195 #endif
196 
197 #endif
bool operator==(const LUT &lut) const
LookupTablePointer GetLookupTable()
bool operator!=(const LUT &lut) const
LookupTableDataPointer GetLookupTable() const
void SetLookupTable(LookupTablePointer lut)
virtual void SetLookupTable(LookupTableType *_arg)
Function to do the lookup operation.
itk::SmartPointer< Self > Pointer
LookupTableDataPointer m_LookupTableDataPointer
FunctorType::LookupTableType LookupTableType
LookupTableType::Pointer m_LookupTable
itk::SmartPointer< const Self > ConstPointer
TOutput operator()(const TInput &val) const
LookupTablePointer m_LookupTablePointer
itk::UnaryFunctorImageFilter< TInputImage, TOutputImage, FunctorType > Superclass
InterpolatorPointer m_Interpolator
LookupTableType::PixelType * LookupTableDataPointer
itk::Image< TOutput, 1 > LookupTableType
typename Superclass::ContinuousIndexType ContinuousIndexType
itk::LinearInterpolateImageFunction< LookupTableType, double > InterpolatorType
Functor::LUT< typename TInputImage::PixelType, typename TOutputImage::PixelType > FunctorType
InterpolatorType::Pointer InterpolatorPointer
LookupTableType::Pointer LookupTablePointer