ANIMA  4.0
animaSmoothingRecursiveYvvGaussianImageFilter.hxx
Go to the documentation of this file.
1 #pragma once
2 
4 #include <itkImageRegionIteratorWithIndex.h>
5 #include <itkProgressAccumulator.h>
6 
7 namespace anima
8 {
9 
13 template <typename TInputImage, typename TOutputImage >
16 {
17  m_NormalizeAcrossScale = false;
18 
19  m_FirstSmoothingFilter = FirstGaussianFilterType::New();
20  m_FirstSmoothingFilter->SetDirection(ImageDimension - 1);
21  m_FirstSmoothingFilter->SetNormalizeAcrossScale( m_NormalizeAcrossScale );
22  m_FirstSmoothingFilter->ReleaseDataFlagOn();
23 
24  for( unsigned int i = 0; i<ImageDimension-1; i++ )
25  {
26  m_SmoothingFilters[i] = InternalGaussianFilterType::New();
27  m_SmoothingFilters[i]->SetNormalizeAcrossScale( m_NormalizeAcrossScale );
28  m_SmoothingFilters[i]->SetDirection(i);
29  m_SmoothingFilters[i]->ReleaseDataFlagOn();
30  m_SmoothingFilters[i]->InPlaceOn();
31  }
32 
33  m_SmoothingFilters[0]->SetInput( m_FirstSmoothingFilter->GetOutput() );
34  for( unsigned int i = 1; i<ImageDimension-1; i++ )
35  {
36  m_SmoothingFilters[i]->SetInput(m_SmoothingFilters[i-1]->GetOutput());
37  }
38 
39  m_CastingFilter = CastingFilterType::New();
40  m_CastingFilter->SetInput(m_SmoothingFilters[ImageDimension-2]->GetOutput());
41  m_CastingFilter->InPlaceOn();
42 
43  this->InPlaceOff();
44 
45  //
46  // NB: We must call SetSigma in order to initialize the smoothing
47  // filters with the default scale. However, m_Sigma must first be
48  // initialized (it is used inside SetSigma) and it must be different
49  // from 1.0 or the call will be ignored.
50  this->m_Sigma.Fill(0.0);
51  this->SetSigma( 1.0 );
52 }
53 
54 
55 template <typename TInputImage, typename TOutputImage>
56 void
58 ::SetNumberOfWorkUnits( itk::ThreadIdType nb )
59 {
60  Superclass::SetNumberOfWorkUnits( nb );
61  for( unsigned int i = 0; i<ImageDimension-1; i++ )
62  {
63  m_SmoothingFilters[ i ]->SetNumberOfWorkUnits( nb );
64  }
65  m_FirstSmoothingFilter->SetNumberOfWorkUnits( nb );
66 
67 }
68 
69 template< typename TInputImage, typename TOutputImage >
70 bool
73 {
74  // Note: There are two different ways this filter may try to run
75  // in-place:
76  // 1) Similar to the standard way, when the input and output image
77  // are of the same type, they can share the bulk data. The output
78  // will be grafted onto the last filter. In this fashion the input
79  // and output will be the same bulk data, but the intermediate
80  // mini-pipeline will use different data.
81  // 2) If the input image is the same type as the RealImage used for
82  // the mini-pipeline, then all the filters may re-use the same
83  // bulk data, stealing it from the input then moving it down the
84  // pipeline filter by filter. Additionally, if the output is also
85  // RealType then the last filter will run in-place making the entire
86  // pipeline in-place and only utilizing on copy of the bulk data.
87 
88  return m_FirstSmoothingFilter->CanRunInPlace() || this->Superclass::CanRunInPlace();
89 }
90 
91 // Set value of Sigma (isotropic)
92 
93 template <typename TInputImage, typename TOutputImage>
94 void
97 {
98  SigmaArrayType sigmas(sigma);
99  this->SetSigmaArray(sigmas);
100 }
101 
102 
103 // Set value of Sigma (an-isotropic)
104 
105 template <typename TInputImage, typename TOutputImage>
106 void
109 {
110  if (this->m_Sigma != sigma)
111  {
112  this->m_Sigma = sigma;
113  for( unsigned int i = 0; i<ImageDimension-1; i++ )
114  {
115  m_SmoothingFilters[ i ]->SetSigma( m_Sigma[i] );
116  }
117  m_FirstSmoothingFilter->SetSigma(m_Sigma[ImageDimension - 1]);
118 
119  this->Modified();
120  }
121 }
122 
123 
124 // Get the sigma array.
125 template <typename TInputImage, typename TOutputImage>
129 {
130  return m_Sigma;
131 }
132 
133 
134 // Get the sigma scalar. If the sigma is anisotropic, we will just
135 // return the sigma along the first dimension.
136 template <typename TInputImage, typename TOutputImage>
139 ::GetSigma() const
140 {
141  return m_Sigma[0];
142 }
143 
144 
148 template <typename TInputImage, typename TOutputImage>
149 void
151 ::SetNormalizeAcrossScale( bool normalize )
152 {
153  m_NormalizeAcrossScale = normalize;
154 
155  for( unsigned int i = 0; i<ImageDimension-1; i++ )
156  {
157  m_SmoothingFilters[i]->SetNormalizeAcrossScale(normalize);
158  }
159  m_FirstSmoothingFilter->SetNormalizeAcrossScale(normalize);
160 
161  this->Modified();
162 }
163 
164 template <typename TInputImage, typename TOutputImage>
165 void
168 {
169  // call the superclass' implementation of this method. this should
170  // copy the output requested region to the input requested region
171  Superclass::GenerateInputRequestedRegion();
172 
173  // This filter needs all of the input
175  if( image )
176  {
177  image->SetRequestedRegion( this->GetInput()->GetLargestPossibleRegion() );
178  }
179 }
180 
181 template <typename TInputImage, typename TOutputImage>
182 void
184 ::EnlargeOutputRequestedRegion(itk::DataObject *output)
185 {
186  TOutputImage *out = dynamic_cast<TOutputImage*>(output);
187 
188  if (out)
189  {
190  out->SetRequestedRegion( out->GetLargestPossibleRegion() );
191  }
192 }
193 
197 template <typename TInputImage, typename TOutputImage >
198 void
201 {
202 
203  itkDebugMacro(<< "SmoothingRecursiveYvvGaussianImageFilter generating data ");
204 
205  const typename TInputImage::ConstPointer inputImage( this->GetInput() );
206 
207  const typename TInputImage::RegionType region = inputImage->GetRequestedRegion();
208  const typename TInputImage::SizeType size = region.GetSize();
209 
210  for( unsigned int d=0; d < ImageDimension; d++)
211  {
212  if(size[d] < 4)
213  {
214  itkExceptionMacro("The number of pixels along dimension " << d << " is less than 4. This filter requires a minimum of four pixels along the dimension to be processed.");
215  }
216  }
217 
218  // If this filter is running in-place, then set the first smoothing
219  // filter to steal the bulk data, by running in-place.
220  if (this->CanRunInPlace() && this->GetInPlace())
221  {
222  m_FirstSmoothingFilter->InPlaceOn();
223 
224  // To make this filter's input and out share the same data, call
225  // the InPlace's/Superclass's allocate methods, which takes care
226  // of the needed bulk data sharing.
227  this->AllocateOutputs();
228  }
229  else
230  {
231  m_FirstSmoothingFilter->InPlaceOff();
232  }
233 
234  // If the last filter is running in-place then this bulk data is not
235  // needed, release it to save memory.
236  if ( m_CastingFilter->CanRunInPlace() )
237  {
238  this->GetOutput()->ReleaseData();
239  }
240 
241  // Create a process accumulator for tracking the progress of this minipipeline
242  itk::ProgressAccumulator::Pointer progress = itk::ProgressAccumulator::New();
243  progress->SetMiniPipelineFilter(this);
244 
245  // Register the filter with the with progress accumulator using
246  // equal weight proportion
247  for(unsigned int i = 0;i < ImageDimension - 1;++i)
248  {
249  progress->RegisterInternalFilter(m_SmoothingFilters[i],1.0 / (ImageDimension));
250  }
251 
252  progress->RegisterInternalFilter(m_FirstSmoothingFilter,1.0 / (ImageDimension));
253  m_FirstSmoothingFilter->SetInput( inputImage );
254  // graft our output to the internal filter to force the proper regions
255  // to be generated
256  m_CastingFilter->GraftOutput(this->GetOutput());
257  m_CastingFilter->Update();
258  this->GraftOutput(m_CastingFilter->GetOutput());
259 
260 }
261 
262 
263 template <typename TInputImage, typename TOutputImage>
264 void
266 ::PrintSelf(std::ostream& os, itk::Indent indent) const
267 {
268  Superclass::PrintSelf(os,indent);
269 
270  os << "NormalizeAcrossScale: " << m_NormalizeAcrossScale << std::endl;
271  os << "Sigma: " << m_Sigma << std::endl;
272 }
273 
274 
275 } // end namespace anima
itk::FixedArray< ScalarRealType, itkGetStaticConstMacro(ImageDimension) > SigmaArrayType
void PrintSelf(std::ostream &os, itk::Indent indent) const ITK_OVERRIDE