ITK  4.3.0
Insight Segmentation and Registration Toolkit
itkMRIBiasFieldCorrectionFilter.h
Go to the documentation of this file.
1 /*=========================================================================
2  *
3  * Copyright Insight Software 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 #ifndef __itkMRIBiasFieldCorrectionFilter_h
19 #define __itkMRIBiasFieldCorrectionFilter_h
20 
21 #include <ctime>
22 
23 #include "itkImageToImageFilter.h"
24 #include "itkArray2D.h"
25 #include "itkMRASlabIdentifier.h"
30 #include "itkImageRegionIterator.h"
31 
32 namespace itk
33 {
46 template< class TImage, class TImageMask, class TBiasField >
48 {
49 public:
55 
58 
60  itkNewMacro(Self);
61 
63  typedef TImage ImageType;
64  typedef typename ImageType::Pointer ImagePointer;
65  typedef typename ImageType::PixelType ImageElementType;
66  typedef typename ImageType::IndexType ImageIndexType;
67  typedef typename ImageType::RegionType ImageRegionType;
68  typedef TImageMask MaskType;
69  typedef typename MaskType::Pointer MaskPointer;
70  typedef typename MaskType::PixelType MaskElementType;
71 
73  typedef TBiasField BiasFieldType;
74 
78 
81 
84 
85  itkStaticConstMacro(SpaceDimension, unsigned int, 3);
86 
88  typedef CompositeValleyFunction InternalEnergyFunction;
89 
91  typedef unsigned int SamplingFactorType[SpaceDimension];
92 
94  itkSetObjectMacro(Image, ImageType);
95 
97  itkSetObjectMacro(Mask, MaskType);
98 
100  itkSetMacro(Region, ImageRegionType);
101 
104  { m_BiasField = bias; }
105 
109  {
110  for ( unsigned int i = 0; i < SpaceDimension; i++ )
111  {
112  m_SamplingFactor[i] = factor[i];
113  }
114  }
116 
119  double GetEnergy0(double diff)
120  {
121  return ( *m_InternalEnergyFunction )( diff );
122  }
123 
126  MeasureType GetValue(const ParametersType & parameters) const;
127 
130  void GetDerivative( const ParametersType & itkNotUsed(parameters),
131  DerivativeType & itkNotUsed(derivative) ) const
132  {}
133 
139  Array< double > classSigmas);
140 
141  unsigned int GetNumberOfParameters(void) const;
142 
143 protected:
146 
148  virtual ~MRIBiasEnergyFunction();
149 
150 private:
151 
154 
157 
160 
163 
166 
169 
170  MRIBiasEnergyFunction(const Self &); //purposely not implemented
171  void operator=(const Self &); //purposely not implemented
172 }; // end of class
173 
224 template< class TInputImage, class TOutputImage, class TMaskImage >
226  public ImageToImageFilter< TInputImage, TOutputImage >
227 {
228 public:
234 
236  itkNewMacro(Self);
237 
240 
242  itkStaticConstMacro(ImageDimension, unsigned int,
243  TOutputImage::ImageDimension);
244 
246  typedef TOutputImage OutputImageType;
247  typedef typename TOutputImage::Pointer OutputImagePointer;
248  typedef typename TOutputImage::IndexType OutputImageIndexType;
249  typedef typename TOutputImage::PixelType OutputImagePixelType;
250  typedef typename TOutputImage::SizeType OutputImageSizeType;
251  typedef typename TOutputImage::RegionType OutputImageRegionType;
252 
253  typedef TInputImage InputImageType;
254  typedef typename TInputImage::Pointer InputImagePointer;
255  typedef typename TInputImage::IndexType InputImageIndexType;
256  typedef typename TInputImage::PixelType InputImagePixelType;
257  typedef typename TInputImage::SizeType InputImageSizeType;
258  typedef typename TInputImage::RegionType InputImageRegionType;
259 
261  typedef TMaskImage ImageMaskType;
262  typedef typename ImageMaskType::Pointer ImageMaskPointer;
263  typedef typename ImageMaskType::RegionType ImageMaskRegionType;
264 
271 
275  typedef typename SlabRegionVectorType::iterator SlabRegionVectorIteratorType;
276 
279 
285 
288 
291 
294 
298  void SetInputMask(ImageMaskType *inputMask);
299 
300  itkGetObjectMacro(InputMask, ImageMaskType);
301 
304  void SetOutputMask(ImageMaskType *outputMask);
305 
307  itkGetObjectMacro(OutputMask, ImageMaskType);
308 
312  void IsBiasFieldMultiplicative(bool flag)
313  { m_BiasMultiplicative = flag; }
314 
316  bool IsBiasFieldMultiplicative()
317  { return m_BiasMultiplicative; }
318 
322  itkSetMacro(UsingInterSliceIntensityCorrection, bool);
323  itkGetConstReferenceMacro(UsingInterSliceIntensityCorrection, bool);
325 
331  itkSetMacro(UsingSlabIdentification, bool);
332  itkGetConstReferenceMacro(UsingSlabIdentification, bool);
334 
335  itkSetMacro(SlabBackgroundMinimumThreshold, InputImagePixelType);
336  itkGetConstReferenceMacro(SlabBackgroundMinimumThreshold,
337  InputImagePixelType);
338 
339  itkSetMacro(SlabNumberOfSamples, unsigned int);
340  itkGetConstReferenceMacro(SlabNumberOfSamples, unsigned int);
341 
342  itkSetMacro(SlabTolerance, double);
343  itkGetConstReferenceMacro(SlabTolerance, double);
344 
350  itkSetMacro(UsingBiasFieldCorrection, bool);
351  itkGetConstReferenceMacro(UsingBiasFieldCorrection, bool);
353 
356  itkSetMacro(GeneratingOutput, bool);
357  itkGetConstReferenceMacro(GeneratingOutput, bool);
359 
362  itkSetMacro(SlicingDirection, int);
363 
365  itkSetMacro(BiasFieldDegree, int);
366  itkGetConstMacro(BiasFieldDegree, int);
368 
371  void SetInitialBiasFieldCoefficients(const
373  & coefficients)
374  { this->Modified(); m_BiasFieldCoefficients = coefficients; }
375 
379  BiasFieldType::CoefficientArrayType GetEstimatedBiasFieldCoefficients()
380  { return m_EstimatedBiasFieldCoefficients; }
381 
385  void SetTissueClassStatistics(const Array< double > & means,
386  const Array< double > & sigmas)
387  throw ( ExceptionObject );
388 
390  itkSetMacro(VolumeCorrectionMaximumIteration, int);
391  itkGetConstMacro(VolumeCorrectionMaximumIteration, int);
392  itkSetMacro(InterSliceCorrectionMaximumIteration, int);
393  itkGetConstMacro(InterSliceCorrectionMaximumIteration, int);
395 
397  void SetOptimizerInitialRadius(double initRadius)
398  { m_OptimizerInitialRadius = initRadius; }
399  double GetOptimizerInitialRadius()
400  { return m_OptimizerInitialRadius; }
402 
404  itkSetMacro(OptimizerGrowthFactor, double);
405  itkGetConstMacro(OptimizerGrowthFactor, double);
407 
410  itkSetMacro(OptimizerShrinkFactor, double);
411  itkGetConstMacro(OptimizerShrinkFactor, double);
412 
419  void SetNumberOfLevels(unsigned int num);
420 
422  itkGetConstMacro(NumberOfLevels, unsigned int);
423 
430  void SetSchedule(const ScheduleType & schedule);
431 
433  itkGetConstReferenceMacro(Schedule, ScheduleType);
434 
439  void SetStartingShrinkFactors(unsigned int factor);
440 
441  void SetStartingShrinkFactors(unsigned int *factors);
442 
444  const unsigned int * GetStartingShrinkFactors() const;
445 
449  static bool IsScheduleDownwardDivisible(const ScheduleType & schedule);
450 
457  void Initialize()
458  throw ( ExceptionObject );
459 
462  BiasFieldType EstimateBiasField(InputImageRegionType region,
463  unsigned int degree,
464  int maximumIteration);
465 
469  void CorrectImage(BiasFieldType & bias,
470  InputImageRegionType region);
471 
474  void CorrectInterSliceIntensityInhomogeneity(InputImageRegionType region);
475 
476 protected:
478  virtual ~MRIBiasFieldCorrectionFilter();
479  void PrintSelf(std::ostream & os, Indent indent) const;
480 
483  bool CheckMaskImage(ImageMaskType *mask);
484 
485 protected:
489  void Log1PImage(InternalImageType *source,
490  InternalImageType *target);
491 
494  void ExpImage(InternalImageType *source,
495  InternalImageType *target);
496 
499  template< class TSource, class TTarget >
500  void CopyAndConvertImage(const TSource *source,
501  TTarget *target,
502  typename TTarget::RegionType requestedRegion)
503  {
504  typedef ImageRegionConstIterator< TSource > SourceIterator;
505  typedef ImageRegionIterator< TTarget > TargetIterator;
506  typedef typename TTarget::PixelType TargetPixelType;
507 
508  SourceIterator s_iter(source, requestedRegion);
509  TargetIterator t_iter(target, requestedRegion);
510 
511  s_iter.GoToBegin();
512  t_iter.GoToBegin();
513  while ( !s_iter.IsAtEnd() )
514  {
515  t_iter.Set( static_cast< TargetPixelType >( s_iter.Get() ) );
516  ++s_iter;
517  ++t_iter;
518  }
519  }
520 
525  void GetBiasFieldSize(InputImageRegionType region,
526  BiasFieldType::DomainSizeType & domainSize);
527 
531  void AdjustSlabRegions(SlabRegionVectorType & slabs,
532  OutputImageRegionType requestedRegion);
533 
534  void GenerateData();
535 
536 private:
537  MRIBiasFieldCorrectionFilter(const Self &); //purposely not implemented
538  void operator=(const Self &); //purposely not implemented
539 
542 
545 
548 
551 
554 
557 
560 
563 
569 
570  unsigned int m_SlabNumberOfSamples;
575 
577  unsigned int m_NumberOfLevels;
578 
581 
585 
589 
592 
595 
598 
601 
604 
607 
610 };
611 
612 // ==================================
613 } // end namespace itk
614 
615 #ifndef ITK_MANUAL_INSTANTIATION
616 #include "itkMRIBiasFieldCorrectionFilter.hxx"
617 #endif
618 
619 #endif
620