RTK  2.0.1
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:
161  ITK_DISALLOW_COPY_AND_ASSIGN(AdditiveGaussianNoiseImageFilter);
162 
168 
170  itkNewMacro(Self);
171 
174 
176  using OutputImageRegionType = typename Superclass::OutputImageRegionType;
177  using OutputImagePointer = typename Superclass::OutputImagePointer;
178 
180  using InputImageType = TInputImage;
181  using InputImagePointer = typename InputImageType::Pointer;
182  using InputImageConstPointer = typename InputImageType::ConstPointer;
183  using InputImageRegionType = typename InputImageType::RegionType;
184  using InputImagePixelType = typename InputImageType::PixelType;
185  using InputPixelType = typename InputImageType::PixelType;
186 
188  itkStaticConstMacro(InputImageDimension, unsigned int,
189  TInputImage::ImageDimension);
190 
191  // virtual void GenerateOutputInformation();
192 
193  void GenerateData() override;
194 
195  // Accessor & Mutator methods
196 
201  void SetMean( float mean )
202  {
203  m_NoiseFilter->GetFunctor().SetMean( mean );
204  this->Modified();
205  }
207 
212  float GetMean() const
213  {
214  return m_NoiseFilter->GetFunctor().GetMean();
215  }
216 
221  void SetStandardDeviation( float stddev )
222  {
223  m_NoiseFilter->GetFunctor().SetStandardDeviation( stddev );
224  this->Modified();
225  }
227 
232  float GetStandardDeviation() const
233  {
234  return m_NoiseFilter->GetFunctor().GetStandardDeviation();
235  }
236 
243  void SetSeed( unsigned long seed )
244  {
245  m_NoiseFilter->GetFunctor().SetSeed( seed );
246  this->Modified();
247  }
249 
252  {
253  if( min == m_NoiseFilter->GetFunctor().GetOutputMinimum() )
254  {
255  return;
256  }
257  m_NoiseFilter->GetFunctor().SetOutputMinimum( min );
258  this->Modified();
259  }
261 
264  {
265  return m_NoiseFilter->GetFunctor().GetOutputMinimum();
266  }
267 
270  {
271  if( max == m_NoiseFilter->GetFunctor().GetOutputMaximum() )
272  {
273  return;
274  }
275  m_NoiseFilter->GetFunctor().SetOutputMaximum( max );
276  this->Modified();
277  }
279 
282  {
283  return m_NoiseFilter->GetFunctor().GetOutputMaximum();
284  }
285 
286 protected:
287 
289 
290  void PrintSelf(std::ostream& os, itk::Indent indent) const override;
291 
292 public:
293 
296 
297 private:
298 
300 };
301 
302 } /* end namespace rtk */
303 
304 #ifndef ITK_MANUAL_INSTANTIATION
305 #include "rtkAdditiveGaussianNoiseImageFilter.hxx"
306 #endif
307 
308 #endif /* rtkAdditiveGaussianNoiseImageFilter_h */
typename Superclass::OutputImagePointer OutputImagePointer
typename InputImageType::ConstPointer InputImageConstPointer
typename InputImageType::RegionType InputImageRegionType
Pixel functor that adds Gaussian noise.
typename Superclass::OutputImageRegionType OutputImageRegionType
itk::Statistics::NormalVariateGenerator::Pointer m_Generator
typename InputImageType::PixelType InputImagePixelType