RTK  1.4.0
Reconstruction Toolkit
rtkAdditiveGaussianNoiseImageFilter.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 /*Insight Segmentation & Registration Toolkit
20  Module: $RCSfile: $
21  Language: C++
22  Date: $Date: 2007-08-31 22:17:25 +0200 (Fri, 31 Aug 2007) $
23  Version: $Revision: 2 $
24  Author: Gavin Baker <gavinb@cs.mu.oz.au>
25 
26  Copyright (c) 2004 Insight Software Consortium. All rights reserved.
27  See ITKCopyright.txt or http://www.itk.org/HTML/Copyright.htm for details.
28 
29  This software is distributed WITHOUT ANY WARRANTY; without even
30  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
31  PURPOSE. See the above copyright notices for more information.
32 
33 =========================================================================*/
34 #ifndef rtkAdditiveGaussianNoiseImageFilter_h
35 #define rtkAdditiveGaussianNoiseImageFilter_h
36 
37 #include <itkImageToImageFilter.h>
40 
41 #include "rtkMacro.h"
42 
43 namespace rtk
44 {
45 
54 template < class TPixel >
56 {
57 public:
58 
60  {
61  m_Mean = 0.0;
62  m_StandardDeviation = 1.0;
63 
65 
66  this->SetSeed( 42 );
67  }
68 
69  float GetMean() const
70  {
71  return m_Mean;
72  }
73 
74  void SetMean( float mean )
75  {
76  m_Mean = mean;
77  }
78 
79  float GetStandardDeviation() const
80  {
81  return m_StandardDeviation;
82  }
83 
84  void SetStandardDeviation( float stddev )
85  {
86  m_StandardDeviation = stddev;
87  }
88 
89  void SetSeed( unsigned long seed )
90  {
91  m_Generator->Initialize( seed );
92  }
93 
94  void SetOutputMinimum( TPixel min ) {
95  m_OutputMinimum = min;
96  }
97 
98  void SetOutputMaximum( TPixel max ) {
99  m_OutputMaximum = max;
100  }
101 
102  TPixel GetOutputMinimum() const {
103  return m_OutputMinimum;
104  }
105 
106  TPixel GetOutputMaximum() const {
107  return m_OutputMaximum;
108  }
109 
110  TPixel operator()( TPixel input )
111  {
112  // Get the minimum and maximum output values
113  static const float min = static_cast<float>( m_OutputMinimum );
114  static const float max = static_cast<float>( m_OutputMaximum );
115 
116  // Compute the output
117  float output = static_cast<float>( input ) +
118  m_Mean + m_StandardDeviation * m_Generator->GetVariate();
119 
120  // Clamp the output value in valid range
121  output = ( output < min ? min : output );
122  output = ( output > max ? max : output );
123 
124  return static_cast< TPixel > ( output );
125  }
126 
127 private:
130  float m_Mean;
133 };
134 
135 
156 template <class TInputImage >
158  public itk::ImageToImageFilter< TInputImage, TInputImage >
159 {
160 public:
166 
168  itkNewMacro(Self);
169 
172 
174  typedef typename Superclass::OutputImageRegionType OutputImageRegionType;
175  typedef typename Superclass::OutputImagePointer OutputImagePointer;
176 
178  typedef TInputImage InputImageType;
179  typedef typename InputImageType::Pointer InputImagePointer;
180  typedef typename InputImageType::ConstPointer InputImageConstPointer;
181  typedef typename InputImageType::RegionType InputImageRegionType;
182  typedef typename InputImageType::PixelType InputImagePixelType;
183  typedef typename InputImageType::PixelType InputPixelType;
184 
186  itkStaticConstMacro(InputImageDimension, unsigned int,
187  TInputImage::ImageDimension);
188 
189  // virtual void GenerateOutputInformation();
190 
191  void GenerateData() ITK_OVERRIDE;
192 
193  // Accessor & Mutator methods
194 
199  void SetMean( float mean )
200  {
201  m_NoiseFilter->GetFunctor().SetMean( mean );
202  this->Modified();
203  }
205 
210  float GetMean() const
211  {
212  return m_NoiseFilter->GetFunctor().GetMean();
213  }
214 
219  void SetStandardDeviation( float stddev )
220  {
221  m_NoiseFilter->GetFunctor().SetStandardDeviation( stddev );
222  this->Modified();
223  }
225 
230  float GetStandardDeviation() const
231  {
232  return m_NoiseFilter->GetFunctor().GetStandardDeviation();
233  }
234 
241  void SetSeed( unsigned long seed )
242  {
243  m_NoiseFilter->GetFunctor().SetSeed( seed );
244  this->Modified();
245  }
247 
249  void SetOutputMinimum( InputImagePixelType min )
250  {
251  if( min == m_NoiseFilter->GetFunctor().GetOutputMinimum() )
252  {
253  return;
254  }
255  m_NoiseFilter->GetFunctor().SetOutputMinimum( min );
256  this->Modified();
257  }
259 
261  InputImagePixelType GetOutputMinimum( )
262  {
263  return m_NoiseFilter->GetFunctor().GetOutputMinimum();
264  }
265 
267  void SetOutputMaximum( InputImagePixelType max )
268  {
269  if( max == m_NoiseFilter->GetFunctor().GetOutputMaximum() )
270  {
271  return;
272  }
273  m_NoiseFilter->GetFunctor().SetOutputMaximum( max );
274  this->Modified();
275  }
277 
279  InputImagePixelType GetOutputMaximum( )
280  {
281  return m_NoiseFilter->GetFunctor().GetOutputMaximum();
282  }
283 
284 protected:
285 
287 
288  void PrintSelf(std::ostream& os, itk::Indent indent) const ITK_OVERRIDE;
289 
290 private:
291 
292  AdditiveGaussianNoiseImageFilter(const Self&); // intentionally not implemented
293  void operator=(const Self&); // intentionally not implemented
294 
295 public:
296 
297  typedef itk::UnaryFunctorImageFilter< InputImageType, InputImageType,
300 
301 private:
302 
304 };
305 
306 } /* end namespace rtk */
307 
308 #ifndef ITK_MANUAL_INSTANTIATION
309 #include "rtkAdditiveGaussianNoiseImageFilter.hxx"
310 #endif
311 
312 #endif /* rtkAdditiveGaussianNoiseImageFilter_h */
itk::ImageToImageFilter< TInputImage, TInputImage > Superclass
Pixel functor that adds Gaussian noise.
itk::Statistics::NormalVariateGenerator::Pointer m_Generator
itk::UnaryFunctorImageFilter< InputImageType, InputImageType, NormalVariateNoiseFunctor< typename InputImageType::PixelType > > NoiseFilterType