ITK  4.6.0
Insight Segmentation and Registration Toolkit
itkPatchBasedDenoisingImageFilter.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 __itkPatchBasedDenoisingImageFilter_h
19 #define __itkPatchBasedDenoisingImageFilter_h
20 
22 #include "itkImageRegionIterator.h"
24 #include "itkVector.h"
25 #include "itkVectorImage.h"
26 #include "itkRGBPixel.h"
27 #include "itkRGBAPixel.h"
28 #include "itkDiffusionTensor3D.h"
29 #include "itkFixedArray.h"
30 #include "itkMatrix.h"
32 #include "itkEnableIf.h"
33 #include "itkIsSame.h"
34 
35 #include <vector>
36 
37 namespace itk
38 {
60 template <typename TInputImage, typename TOutputImage>
62  public PatchBasedDenoisingBaseImageFilter<TInputImage, TOutputImage>
63 {
64 public:
71 
73  itkNewMacro(Self);
74 
77 
81 
83  itkStaticConstMacro(ImageDimension, unsigned int,
85 
87  typedef typename InputImageType::RegionType InputImageRegionType;
88 
92 
95  typedef typename Superclass::PixelType PixelType;
97 
103 
110 
116 
124  typedef std::vector<EigenValuesArrayType*> EigenValuesCacheType;
125  typedef std::vector<EigenVectorsMatrixType*> EigenVectorsCacheType;
126 
128  {
138  };
139 
144  itkSetMacro(UseSmoothDiscPatchWeights, bool);
145  itkBooleanMacro(UseSmoothDiscPatchWeights);
146  itkGetConstMacro(UseSmoothDiscPatchWeights, bool);
148 
153  void SetKernelBandwidthSigma(const RealArrayType& kernelSigma);
154  itkGetConstMacro(KernelBandwidthSigma, RealArrayType);
156 
161  itkSetClampMacro(KernelBandwidthFractionPixelsForEstimation, double, 0.01, 1.0);
162  itkGetConstReferenceMacro(KernelBandwidthFractionPixelsForEstimation, double);
164 
167  itkSetMacro(ComputeConditionalDerivatives, bool);
168  itkBooleanMacro(ComputeConditionalDerivatives);
169  itkGetConstMacro(ComputeConditionalDerivatives, bool);
171 
185  itkSetMacro(UseFastTensorComputations, bool);
186  itkBooleanMacro(UseFastTensorComputations);
187  itkGetConstMacro(UseFastTensorComputations, bool);
189 
191  itkStaticConstMacro(MaxSigmaUpdateIterations, unsigned int,
192  20);
193 
200  itkSetClampMacro(KernelBandwidthMultiplicationFactor, double, 0.01, 100);
201  itkGetConstReferenceMacro(KernelBandwidthMultiplicationFactor, double);
203 
207  void SetNoiseSigma(const RealType& sigma);
208 
209  itkGetConstMacro(NoiseSigma, RealType);
210 
212  itkSetObjectMacro(Sampler, BaseSamplerType);
213  itkGetModifiableObjectMacro(Sampler, BaseSamplerType);
215 
216 protected:
219  virtual void PrintSelf(std::ostream& os, Indent indent) const;
220 
222  virtual void EmptyCaches();
223 
225  virtual void AllocateUpdateBuffer();
226 
227  virtual void CopyInputToOutput();
228 
229  virtual void GenerateInputRequestedRegion();
230 
237  template< typename T>
238  typename EnableIfC<
239  IsSame<T, typename NumericTraits<T>::ValueType>::Value,
240  T >::Type
241  GetComponent(const T pix,
242  unsigned int itkNotUsed( idx ) ) const
243  {
244  // The enable if idiom is used to overload this method for both
245  // scalars and multi-component types. By exploiting that
246  // NumericTraits' ValueType typedef (defines the per-element type
247  // for multi-component types ) is different then the parameterize
248  // type, the bracket operator is used only for multi-component
249  // types.
250  return pix;
251  }
252 
253  template< typename T>
254  typename DisableIfC<
255  IsSame<T, typename NumericTraits<T>::ValueType>::Value,
256  typename NumericTraits<T>::ValueType>::Type
257  GetComponent(const T& pix,
258  unsigned int idx ) const
259  {
260  return pix[idx];
261  }
262 
264  template< typename T >
265  void
266  SetComponent( T &pix,
267  unsigned int itkNotUsed( idx ),
268  typename EnableIfC< IsSame<T,
269  typename NumericTraits<T>::ValueType>::Value, RealValueType>::Type val) const
270  {
271  pix = val;
272  }
273 
274  template< typename T >
275  void
276  SetComponent( T &pix,
277  unsigned int idx,
278  typename DisableIfC< IsSame<T,
279  typename NumericTraits<T>::ValueType>::Value,
280  RealValueType>::Type val) const
281  {
282  pix[idx] = val;
283  }
284 
289  {
291  {
293  }
294  else
295  {
297  }
298  }
299 
300  template< typename TImageType>
301  typename EnableIfC<
302  IsSame<typename TImageType::PixelType, typename NumericTraits<typename TImageType::PixelType>::ValueType>::Value
303  >::Type
304  ComputeMinMax(const TImageType* img)
305  {
306  DispatchedMinMax(img);
307  }
308 
309  template< typename TImageType>
310  typename DisableIfC<
311  IsSame<typename TImageType::PixelType, typename NumericTraits<typename TImageType::PixelType>::ValueType>::Value
312  >::Type
313  ComputeMinMax(const TImageType* img)
314  {
316  }
317 
328  const RealArrayType& weight,
329  bool useCachedComputations,
330  SizeValueType cacheIndex,
331  EigenValuesCacheType& eigenValsCache,
332  EigenVectorsCacheType& eigenVecsCache,
333  RealType& diff, RealArrayType& norm)
334  {
336  {
338  useCachedComputations, cacheIndex,
339  eigenValsCache, eigenVecsCache,
340  diff, norm);
341  }
342  else
343  {
345  useCachedComputations, cacheIndex,
346  eigenValsCache, eigenVecsCache,
347  diff, norm);
348  }
349  }
351 
352  template <typename PixelT>
354  const PixelT& b,
355  const RealArrayType& weight,
356  bool useCachedComputations,
357  SizeValueType cacheIndex,
358  EigenValuesCacheType& eigenValsCache,
359  EigenVectorsCacheType& eigenVecsCache,
360  RealType& diff, RealArrayType& norm)
361  {
363  useCachedComputations, cacheIndex,
364  eigenValsCache, eigenVecsCache,
365  diff, norm);
366  }
367 
372  const RealType& b)
373  {
375  {
376  return this->AddExponentialMapUpdate(a, b);
377  }
378  else
379  {
380  return this->AddEuclideanUpdate(a, b);
381  }
382  }
384 
385  template <typename RealT>
386  RealType AddUpdate(const RealT& a,
387  const RealType& b)
388  {
389  return this->AddEuclideanUpdate(a, b);
390  }
391 
392  virtual void EnforceConstraints();
393 
394  virtual void Initialize();
395 
396  virtual void InitializeKernelSigma();
397 
398  virtual void InitializePatchWeights();
399 
400  virtual void InitializePatchWeightsSmoothDisc();
401 
402  virtual void InitializeIteration();
403 
404  virtual void ComputeKernelBandwidthUpdate(); // derived from base class;
405 
406  // define here
407 
408  virtual ThreadDataStruct ThreadedComputeSigmaUpdate(const InputImageRegionType& regionToProcess,
409  const int itkNotUsed(threadId),
410  ThreadDataStruct threadData);
411 
413 
414  virtual void ComputeImageUpdate();
415 
416  virtual ThreadDataStruct ThreadedComputeImageUpdate(const InputImageRegionType& regionToProcess,
417  const int threadId,
418  ThreadDataStruct threadData);
419 
421  typename ListAdaptorType::Pointer& inList,
422  BaseSamplerPointer& sampler,
423  ThreadDataStruct& threadData);
424 
425  virtual void ApplyUpdate();
426 
427  virtual void ThreadedApplyUpdate(const InputImageRegionType& regionToProcess,
428  const int itkNotUsed(threadId) );
429 
430  virtual void PostProcessOutput();
431 
432  virtual void SetThreadData(int threadId, const ThreadDataStruct& data);
433 
434  virtual ThreadDataStruct GetThreadData(int threadId);
435 
436  std::vector<ThreadDataStruct> m_ThreadData;
437 
439  typename OutputImageType::Pointer m_UpdateBuffer;
440 
441  unsigned int m_NumPixelComponents;
443  unsigned int m_TotalNumberPixels;
444  //
446  //
448  //
457  double m_MinSigma;
463  //
467  //
469  typename ListAdaptorType::Pointer m_SearchSpaceList;
470 
471 private:
472  PatchBasedDenoisingImageFilter(const Self&); // purposely not implemented
473  void operator=(const Self&); // purposely not implemented
474 
478 
482 
486 
487  template <typename TInputImageType>
488  void DispatchedMinMax(const TInputImageType* img);
489 
490  template <typename TInputImageType>
491  void DispatchedArrayMinMax(const TInputImageType* img);
492 
493  template <typename TInputImageType>
494  void DispatchedVectorMinMax(const TInputImageType* img);
495 
496  template <typename TInputImageType>
497  void DispatchedRiemannianMinMax(const TInputImageType* img);
498 
502 
504  const int itkNotUsed(threadId),
505  const InputImageType* img,
506  ThreadDataStruct threadData);
507 
508  virtual void ResolveRiemannianMinMax();
509 
511  const RealArrayType& weight,
512  bool useCachedComputations,
513  SizeValueType cacheIndex,
514  EigenValuesCacheType& eigenValsCache,
515  EigenVectorsCacheType& eigenVecsCache,
516  RealType& diff, RealArrayType& norm);
517 
520  const DiffusionTensor3D<PixelValueType>& spdMatrixB,
521  const RealArrayType& weight,
522  bool useCachedComputations,
523  SizeValueType cacheIndex,
524  EigenValuesCacheType& eigenValsCache,
525  EigenVectorsCacheType& eigenVecsCache,
526  RealType& symMatrixLogMap, RealArrayType& geodesicDist);
527 
528  template <typename TensorValueT>
531  Matrix< TensorValueT, 3, 3 >& eigenVecs);
532 
533  RealType AddEuclideanUpdate(const RealType& a, const RealType& b);
534 
537  const DiffusionTensor3D<RealValueType>& symMatrix);
538 
540  {
543  };
544 
545 }; // end class PatchBasedDenoisingImageFilter
546 
547 } // end namespace itk
548 
549 #ifndef ITK_MANUAL_INSTANTIATION
550 # include "itkPatchBasedDenoisingImageFilter.hxx"
551 #endif
552 
553 #endif
A templated class holding a M x N size Matrix.
Definition: itkMatrix.h:46
Derived class implementing a specific patch-based denoising algorithm, as detailed below...
virtual ThreadDataStruct ThreadedComputeSigmaUpdate(const InputImageRegionType &regionToProcess, const int, ThreadDataStruct threadData)
BaseSamplerType::InstanceIdentifier InstanceIdentifier
void ComputeDifferenceAndWeightedSquaredNorm(const PixelT &a, const PixelT &b, const RealArrayType &weight, bool useCachedComputations, SizeValueType cacheIndex, EigenValuesCacheType &eigenValsCache, EigenVectorsCacheType &eigenVecsCache, RealType &diff, RealArrayType &norm)
ConstNeighborhoodIterator< InputImageType, BoundaryConditionType > InputImagePatchIterator
Superclass::InputImagePatchIterator InputImagePatchIterator
Base class for patch-based denoisng algorithms.
DisableIfC< IsSame< T, typename NumericTraits< T >::ValueType >::Value, typename NumericTraits< T >::ValueType >::Type GetComponent(const T &pix, unsigned int idx) const
void DispatchedRiemannianMinMax(const TInputImageType *img)
static ITK_THREAD_RETURN_TYPE RiemannianMinMaxThreaderCallback(void *arg)
RealType AddEuclideanUpdate(const RealType &a, const RealType &b)
ThreadDataStruct ThreadedRiemannianMinMax(const InputImageRegionType &regionToProcess, const int, const InputImageType *img, ThreadDataStruct threadData)
Base class for all process objects that output image data.
virtual ThreadDataStruct ThreadedComputeImageUpdate(const InputImageRegionType &regionToProcess, const int threadId, ThreadDataStruct threadData)
unsigned long SizeValueType
Definition: itkIntTypes.h:143
#define ITK_THREAD_RETURN_TYPE
static ITK_THREAD_RETURN_TYPE ApplyUpdateThreaderCallback(void *arg)
RealType AddExponentialMapUpdate(const DiffusionTensor3D< RealValueType > &spdMatrix, const DiffusionTensor3D< RealValueType > &symMatrix)
Simulate a standard C array with copy semnatics.
Definition: itkFixedArray.h:50
std::vector< EigenValuesArrayType * > EigenValuesCacheType
void ComputeLogMapAndWeightedSquaredGeodesicDifference(const DiffusionTensor3D< PixelValueType > &spdMatrixA, const DiffusionTensor3D< PixelValueType > &spdMatrixB, const RealArrayType &weight, bool useCachedComputations, SizeValueType cacheIndex, EigenValuesCacheType &eigenValsCache, EigenVectorsCacheType &eigenVecsCache, RealType &symMatrixLogMap, RealArrayType &geodesicDist)
void SetKernelBandwidthSigma(const RealArrayType &kernelSigma)
PatchBasedDenoisingBaseImageFilter< TInputImage, TOutputImage > Superclass
RealType AddUpdate(const DiffusionTensor3D< RealValueType > &a, const RealType &b)
std::vector< EigenVectorsMatrixType * > EigenVectorsCacheType
typename::itk::Statistics::ImageToNeighborhoodSampleAdaptor< OutputImageType, BoundaryConditionType > ListAdaptorType
virtual RealArrayType ResolveSigmaUpdate()
OutputImageType::Pointer OutputImagePointer
void DispatchedArrayMinMax(const TInputImageType *img)
void ComputeMinMax(const Image< DiffusionTensor3D< PixelValueType >, ImageDimension > *img)
InputImageType::RegionType InputImageRegionType
This an abstract subsampler that constrains subsamples to be contained within a given image region...
A multi-dimensional iterator templated over image type that walks a region of pixels.
virtual void SetThreadData(int threadId, const ThreadDataStruct &data)
void SetNoiseSigma(const RealType &sigma)
static ITK_THREAD_RETURN_TYPE ComputeSigmaUpdateThreaderCallback(void *arg)
void SetComponent(T &pix, unsigned int, typename EnableIfC< IsSame< T, typename NumericTraits< T >::ValueType >::Value, RealValueType >::Type val) const
A method to generically set a component.
virtual void ThreadedApplyUpdate(const InputImageRegionType &regionToProcess, const int)
virtual ThreadDataStruct GetThreadData(int threadId)
void DispatchedMinMax(const TInputImageType *img)
FixedArray< PixelValueType, 3 > EigenValuesArrayType
void ComputeSignedEuclideanDifferenceAndWeightedSquaredNorm(const PixelType &a, const PixelType &b, const RealArrayType &weight, bool useCachedComputations, SizeValueType cacheIndex, EigenValuesCacheType &eigenValsCache, EigenVectorsCacheType &eigenVecsCache, RealType &diff, RealArrayType &norm)
void SetComponent(T &pix, unsigned int idx, typename DisableIfC< IsSame< T, typename NumericTraits< T >::ValueType >::Value, RealValueType >::Type val) const
virtual RealType ComputeGradientJointEntropy(InstanceIdentifier id, typename ListAdaptorType::Pointer &inList, BaseSamplerPointer &sampler, ThreadDataStruct &threadData)
NumericTraits< PixelType >::RealType RealType
EnableIfC< IsSame< T, typename NumericTraits< T >::ValueType >::Value, T >::Type GetComponent(const T pix, unsigned int) const
A method to generically get a component.
void DispatchedVectorMinMax(const TInputImageType *img)
void Compute3x3EigenAnalysis(const DiffusionTensor3D< TensorValueT > &spdMatrix, FixedArray< TensorValueT, 3 > &eigenVals, Matrix< TensorValueT, 3, 3 > &eigenVecs)
Base class for filters that take an image as input and produce an image as output.
ImageRegionIterator< OutputImageType > OutputImageRegionIteratorType
virtual void InitializePatchWeightsSmoothDisc()
NumericTraits< PixelValueType >::RealType RealValueType
static ITK_THREAD_RETURN_TYPE ComputeImageUpdateThreaderCallback(void *arg)
ImageRegionConstIterator< InputImageType > InputImageRegionConstIteratorType
Control indentation during Print() invocation.
Definition: itkIndent.h:49
virtual void PrintSelf(std::ostream &os, Indent indent) const
Define additional traits for native types such as int or float.
EnableIfC< IsSame< typename TImageType::PixelType, typename NumericTraits< typename TImageType::PixelType >::ValueType >::Value >::Type ComputeMinMax(const TImageType *img)
DisableIfC< IsSame< typename TImageType::PixelType, typename NumericTraits< typename TImageType::PixelType >::ValueType >::Value >::Type ComputeMinMax(const TImageType *img)
void ComputeDifferenceAndWeightedSquaredNorm(const DiffusionTensor3D< PixelValueType > &a, const DiffusionTensor3D< PixelValueType > &b, const RealArrayType &weight, bool useCachedComputations, SizeValueType cacheIndex, EigenValuesCacheType &eigenValsCache, EigenVectorsCacheType &eigenVecsCache, RealType &diff, RealArrayType &norm)
Represent a diffusion tensor as used in DTI images.
ListAdaptorType::NeighborhoodRadiusType PatchRadiusType
Templated n-dimensional image class.
Definition: itkImage.h:75
A multi-dimensional iterator templated over image type that walks a region of pixels.
itk::Statistics::RegionConstrainedSubsampler< PatchSampleType, InputImageRegionType > BaseSamplerType
RealType AddUpdate(const RealT &a, const RealType &b)