Loading [MathJax]/extensions/tex2jax.js
ITK  6.0.0
Insight Toolkit
itkN4BiasFieldCorrectionImageFilter.h
Go to the documentation of this file.
1 /*=========================================================================
2  *
3  * Copyright NumFOCUS
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  * https://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 itkN4BiasFieldCorrectionImageFilter_h
19 #define itkN4BiasFieldCorrectionImageFilter_h
20 
21 #include "itkImageToImageFilter.h"
22 
23 #include "itkArray.h"
25 #include "itkPointSet.h"
26 #include "itkVector.h"
27 
28 #include "vnl/vnl_vector.h"
29 
30 namespace itk
31 {
32 
93 template <typename TInputImage,
94  typename TMaskImage = Image<unsigned char, TInputImage::ImageDimension>,
95  class TOutputImage = TInputImage>
96 class ITK_TEMPLATE_EXPORT N4BiasFieldCorrectionImageFilter : public ImageToImageFilter<TInputImage, TOutputImage>
97 {
98 public:
99  ITK_DISALLOW_COPY_AND_MOVE(N4BiasFieldCorrectionImageFilter);
100 
106 
108  itkOverrideGetNameOfClassMacro(N4BiasFieldCorrectionImageFilter);
109 
111  itkNewMacro(Self);
112 
114  static constexpr unsigned int ImageDimension = TInputImage::ImageDimension;
115 
117  using InputImageType = TInputImage;
118  using OutputImageType = TOutputImage;
119  using MaskImageType = TMaskImage;
120  using MaskPixelType = typename MaskImageType::PixelType;
121 
122  using RealType = float;
126 
133 
138 
140  void
141  EnlargeOutputRequestedRegion(DataObject *) override;
142 
146  void
147  SetInput1(const InputImageType * image)
148  {
149  this->SetInput(image);
150  }
151 
157  itkSetInputMacro(MaskImage, MaskImageType);
158  void
159  SetInput2(const MaskImageType * mask)
160  {
161  this->SetMaskImage(mask);
162  }
170  itkGetInputMacro(MaskImage, MaskImageType);
171 
180  itkSetMacro(MaskLabel, MaskPixelType);
181  itkGetConstMacro(MaskLabel, MaskPixelType);
187  itkSetMacro(UseMaskLabel, bool);
188  itkGetConstMacro(UseMaskLabel, bool);
189  itkBooleanMacro(UseMaskLabel);
202  itkSetInputMacro(ConfidenceImage, RealImageType);
203  void
204  SetInput3(const RealImageType * image)
205  {
206  this->SetConfidenceImage(image);
207  }
220  itkGetInputMacro(ConfidenceImage, RealImageType);
221 
222  // Sharpen histogram parameters: in estimating the bias field, the
223  // first step is to sharpen the intensity histogram by Wiener deconvolution
224  // with a 1-D Gaussian. The following parameters define this operation.
225  // These default values in N4 match the default values in N3.
226 
231  itkSetMacro(NumberOfHistogramBins, unsigned int);
232 
237  itkGetConstMacro(NumberOfHistogramBins, unsigned int);
238 
242  itkSetMacro(WienerFilterNoise, RealType);
243 
247  itkGetConstMacro(WienerFilterNoise, RealType);
248 
253  itkSetMacro(BiasFieldFullWidthAtHalfMaximum, RealType);
254 
259  itkGetConstMacro(BiasFieldFullWidthAtHalfMaximum, RealType);
260 
261  // B-spline parameters governing the fitting routine
262 
266  itkSetMacro(SplineOrder, unsigned int);
267 
271  itkGetConstMacro(SplineOrder, unsigned int);
272 
280  itkSetMacro(NumberOfControlPoints, ArrayType);
281 
289  itkGetConstMacro(NumberOfControlPoints, ArrayType);
290 
297  itkSetMacro(NumberOfFittingLevels, ArrayType);
298 
305  void
306  SetNumberOfFittingLevels(unsigned int n)
307  {
308  auto nlevels = MakeFilled<ArrayType>(n);
309  this->SetNumberOfFittingLevels(nlevels);
310  }
319  itkGetConstMacro(NumberOfFittingLevels, ArrayType);
320 
325  itkSetMacro(MaximumNumberOfIterations, VariableSizeArrayType);
326 
331  itkGetConstMacro(MaximumNumberOfIterations, VariableSizeArrayType);
332 
340  itkSetMacro(ConvergenceThreshold, RealType);
341 
349  itkGetConstMacro(ConvergenceThreshold, RealType);
350 
360  itkGetConstObjectMacro(LogBiasFieldControlPointLattice, BiasFieldControlPointLatticeType);
361 
366  itkGetConstMacro(ElapsedIterations, unsigned int);
367 
372  itkGetConstMacro(CurrentConvergenceMeasurement, RealType);
373 
378  itkGetConstMacro(CurrentLevel, unsigned int);
379 
383  RealImagePointer
384  ReconstructBiasField(const BiasFieldControlPointLatticeType *);
385 
386 protected:
388  ~N4BiasFieldCorrectionImageFilter() override = default;
389  void
390  PrintSelf(std::ostream & os, Indent indent) const override;
391 
392  void
393  GenerateData() override;
394 
395 private:
396  // N4 algorithm functions: The basic algorithm iterates between sharpening
397  // the intensity histogram of the corrected input image and spatially
398  // smoothing those results with a B-spline scalar field estimate of the
399  // bias field. The former is handled by the function SharpenImage()
400  // whereas the latter is handled by the function UpdateBiasFieldEstimate().
401  // Convergence is determined by the coefficient of variation of the difference
402  // image between the current bias field estimate and the previous estimate.
403 
409  void
410  SharpenImage(const RealImageType * unsharpenedImage, RealImageType * sharpenedImage) const;
411 
417  RealImagePointer
418  UpdateBiasFieldEstimate(RealImageType *, size_t);
419 
424  RealType
425  CalculateConvergenceMeasurement(const RealImageType *, const RealImageType *) const;
426 
427  MaskPixelType m_MaskLabel{};
428  bool m_UseMaskLabel{ false };
429 
430  // Parameters for deconvolution with Wiener filter
431 
432  unsigned int m_NumberOfHistogramBins{ 200 };
433  RealType m_WienerFilterNoise{ static_cast<RealType>(0.01) };
434  RealType m_BiasFieldFullWidthAtHalfMaximum{ static_cast<RealType>(0.15) };
435 
436  // Convergence parameters
437 
438  VariableSizeArrayType m_MaximumNumberOfIterations{};
439  unsigned int m_ElapsedIterations{ 0 };
440  RealType m_ConvergenceThreshold{ static_cast<RealType>(0.001) };
441  RealType m_CurrentConvergenceMeasurement{};
442  unsigned int m_CurrentLevel{ 0 };
443 
444  // B-spline fitting parameters
445 
446  typename BiasFieldControlPointLatticeType::Pointer m_LogBiasFieldControlPointLattice{};
447 
448  unsigned int m_SplineOrder{ 3 };
449  ArrayType m_NumberOfControlPoints{};
450  ArrayType m_NumberOfFittingLevels{};
451 };
452 
453 } // end namespace itk
454 
455 #ifndef ITK_MANUAL_INSTANTIATION
456 # include "itkN4BiasFieldCorrectionImageFilter.hxx"
457 #endif
458 
459 #endif
itk::N4BiasFieldCorrectionImageFilter
Implementation of the N4 bias field correction algorithm.
Definition: itkN4BiasFieldCorrectionImageFilter.h:96
Pointer
SmartPointer< Self > Pointer
Definition: itkAddImageFilter.h:93
itk::N4BiasFieldCorrectionImageFilter::RealImagePointer
typename RealImageType::Pointer RealImagePointer
Definition: itkN4BiasFieldCorrectionImageFilter.h:124
itk::PointSet
A superclass of the N-dimensional mesh structure; supports point (geometric coordinate and attribute)...
Definition: itkPointSet.h:81
itk::GTest::TypedefsAndConstructors::Dimension2::PointType
ImageBaseType::PointType PointType
Definition: itkGTestTypedefsAndConstructors.h:51
itk::BSplineScatteredDataPointSetToImageFilter
Image filter which provides a B-spline output approximation.
Definition: itkBSplineScatteredDataPointSetToImageFilter.h:132
itk::N4BiasFieldCorrectionImageFilter::MaskPixelType
typename MaskImageType::PixelType MaskPixelType
Definition: itkN4BiasFieldCorrectionImageFilter.h:120
itk::Vector
A templated class holding a n-Dimensional vector.
Definition: itkVector.h:62
itk::SmartPointer< Self >
itk::Indent
Control indentation during Print() invocation.
Definition: itkIndent.h:49
itk::N4BiasFieldCorrectionImageFilter::SetInput1
void SetInput1(const InputImageType *image)
Definition: itkN4BiasFieldCorrectionImageFilter.h:147
itk::ImageToImageFilter
Base class for filters that take an image as input and produce an image as output.
Definition: itkImageToImageFilter.h:108
itk::ImageSource
Base class for all process objects that output image data.
Definition: itkImageSource.h:67
itk::N4BiasFieldCorrectionImageFilter::PointType
typename PointSetType::PointType PointType
Definition: itkN4BiasFieldCorrectionImageFilter.h:132
itk::N4BiasFieldCorrectionImageFilter::SetNumberOfFittingLevels
void SetNumberOfFittingLevels(unsigned int n)
Definition: itkN4BiasFieldCorrectionImageFilter.h:306
itk::N4BiasFieldCorrectionImageFilter::BiasFieldControlPointLatticeType
typename BSplineFilterType::PointDataImageType BiasFieldControlPointLatticeType
Definition: itkN4BiasFieldCorrectionImageFilter.h:136
itkBSplineScatteredDataPointSetToImageFilter.h
itk::ImageToImageFilter::InputImageType
TInputImage InputImageType
Definition: itkImageToImageFilter.h:129
itkImageToImageFilter.h
itk::FixedArray< unsigned int, Self::ImageDimension >
itkArray.h
itk::N4BiasFieldCorrectionImageFilter::SetInput2
void SetInput2(const MaskImageType *mask)
Definition: itkN4BiasFieldCorrectionImageFilter.h:159
itk::N4BiasFieldCorrectionImageFilter::PointSetPointer
typename PointSetType::Pointer PointSetPointer
Definition: itkN4BiasFieldCorrectionImageFilter.h:131
itk
The "itk" namespace contains all Insight Segmentation and Registration Toolkit (ITK) classes....
Definition: itkAnatomicalOrientation.h:29
itk::N4BiasFieldCorrectionImageFilter::RealType
float RealType
Definition: itkN4BiasFieldCorrectionImageFilter.h:122
itk::ProcessObject
The base class for all process objects (source, filters, mappers) in the Insight data processing pipe...
Definition: itkProcessObject.h:139
itkVector.h
itk::N4BiasFieldCorrectionImageFilter::SetInput3
void SetInput3(const RealImageType *image)
Definition: itkN4BiasFieldCorrectionImageFilter.h:204
itk::N4BiasFieldCorrectionImageFilter::ArrayType
typename BSplineFilterType::ArrayType ArrayType
Definition: itkN4BiasFieldCorrectionImageFilter.h:137
itk::Array< unsigned int >
itk::Image
Templated n-dimensional image class.
Definition: itkImage.h:88
itkPointSet.h
itk::ImageSource::OutputImageType
TOutputImage OutputImageType
Definition: itkImageSource.h:90
itk::DataObject
Base class for all data objects in ITK.
Definition: itkDataObject.h:293
itk::N4BiasFieldCorrectionImageFilter::MaskImageType
TMaskImage MaskImageType
Definition: itkN4BiasFieldCorrectionImageFilter.h:119