ITK 6.0.0
Insight Toolkit
 
Loading...
Searching...
No Matches
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
22
23#include "itkArray.h"
25#include "itkPointSet.h"
26#include "itkVector.h"
27
28#include "vnl/vnl_vector.h"
29
30namespace itk
31{
32
84
85template <typename TInputImage,
87 class TOutputImage = TInputImage>
88class ITK_TEMPLATE_EXPORT N4BiasFieldCorrectionImageFilter : public ImageToImageFilter<TInputImage, TOutputImage>
89{
90public:
91 ITK_DISALLOW_COPY_AND_MOVE(N4BiasFieldCorrectionImageFilter);
92
98
100 itkOverrideGetNameOfClassMacro(N4BiasFieldCorrectionImageFilter);
101
103 itkNewMacro(Self);
104
106 static constexpr unsigned int ImageDimension = TInputImage::ImageDimension;
107
109 using InputImageType = TInputImage;
110 using OutputImageType = TOutputImage;
111 using MaskImageType = TMaskImage;
112 using MaskPixelType = typename MaskImageType::PixelType;
113
114 using RealType = float;
118
125
130
132 void
134
138 void
140 {
141 this->SetInput(image);
142 }
143
151 void
153 {
154 this->SetMaskImage(mask);
155 }
156
157
164
174 itkSetMacro(MaskLabel, MaskPixelType);
175 itkGetConstMacro(MaskLabel, MaskPixelType);
177
182 itkSetMacro(UseMaskLabel, bool);
183 itkGetConstMacro(UseMaskLabel, bool);
184 itkBooleanMacro(UseMaskLabel);
186
199 void
201 {
202 this->SetConfidenceImage(image);
203 }
204
205
217
218 // Sharpen histogram parameters: in estimating the bias field, the
219 // first step is to sharpen the intensity histogram by Wiener deconvolution
220 // with a 1-D Gaussian. The following parameters define this operation.
221 // These default values in N4 match the default values in N3.
222
227 itkSetMacro(NumberOfHistogramBins, unsigned int);
228
233 itkGetConstMacro(NumberOfHistogramBins, unsigned int);
234
238 itkSetMacro(WienerFilterNoise, RealType);
239
243 itkGetConstMacro(WienerFilterNoise, RealType);
244
249 itkSetMacro(BiasFieldFullWidthAtHalfMaximum, RealType);
250
255 itkGetConstMacro(BiasFieldFullWidthAtHalfMaximum, RealType);
256
257 // B-spline parameters governing the fitting routine
258
262 itkSetMacro(SplineOrder, unsigned int);
263
267 itkGetConstMacro(SplineOrder, unsigned int);
268
276 itkSetMacro(NumberOfControlPoints, ArrayType);
277
285 itkGetConstMacro(NumberOfControlPoints, ArrayType);
286
293 itkSetMacro(NumberOfFittingLevels, ArrayType);
294
301 void
303 {
304 auto nlevels = MakeFilled<ArrayType>(n);
305 this->SetNumberOfFittingLevels(nlevels);
306 }
307
314 itkGetConstMacro(NumberOfFittingLevels, ArrayType);
315
320 itkSetMacro(MaximumNumberOfIterations, VariableSizeArrayType);
321
326 itkGetConstMacro(MaximumNumberOfIterations, VariableSizeArrayType);
327
335 itkSetMacro(ConvergenceThreshold, RealType);
336
344 itkGetConstMacro(ConvergenceThreshold, RealType);
345
355 itkGetConstObjectMacro(LogBiasFieldControlPointLattice, BiasFieldControlPointLatticeType);
356
361 itkGetConstMacro(ElapsedIterations, unsigned int);
362
367 itkGetConstMacro(CurrentConvergenceMeasurement, RealType);
368
373 itkGetConstMacro(CurrentLevel, unsigned int);
374
380
381protected:
384 void
385 PrintSelf(std::ostream & os, Indent indent) const override;
386
387 void
388 GenerateData() override;
389
390private:
391 // N4 algorithm functions: The basic algorithm iterates between sharpening
392 // the intensity histogram of the corrected input image and spatially
393 // smoothing those results with a B-spline scalar field estimate of the
394 // bias field. The former is handled by the function SharpenImage()
395 // whereas the latter is handled by the function UpdateBiasFieldEstimate().
396 // Convergence is determined by the coefficient of variation of the difference
397 // image between the current bias field estimate and the previous estimate.
398
404 void
405 SharpenImage(const RealImageType * unsharpenedImage, RealImageType * sharpenedImage) const;
406
414
421
423 bool m_UseMaskLabel{ false };
424
425 // Parameters for deconvolution with Wiener filter
426
427 unsigned int m_NumberOfHistogramBins{ 200 };
428 RealType m_WienerFilterNoise{ static_cast<RealType>(0.01) };
430
431 // Convergence parameters
432
434 unsigned int m_ElapsedIterations{ 0 };
437 unsigned int m_CurrentLevel{ 0 };
438
439 // B-spline fitting parameters
440
441 typename BiasFieldControlPointLatticeType::Pointer m_LogBiasFieldControlPointLattice{};
442
443 unsigned int m_SplineOrder{ 3 };
446};
447
448} // end namespace itk
449
450#ifndef ITK_MANUAL_INSTANTIATION
451# include "itkN4BiasFieldCorrectionImageFilter.hxx"
452#endif
453
454#endif
Array class with size defined at construction time.
Definition itkArray.h:48
Image filter which provides a B-spline output approximation.
Base class for all data objects in ITK.
virtual void SetInput(const InputImageType *input)
Templated n-dimensional image class.
Definition itkImage.h:89
Control indentation during Print() invocation.
Definition itkIndent.h:50
itkGetInputMacro(ConfidenceImage, RealImageType)
virtual void SetNumberOfFittingLevels(ArrayType _arg)
Image< ScalarType, Self::ImageDimension > ScalarImageType
void SharpenImage(const RealImageType *unsharpenedImage, RealImageType *sharpenedImage) const
itkGetInputMacro(MaskImage, MaskImageType)
itkSetInputMacro(MaskImage, MaskImageType)
BSplineScatteredDataPointSetToImageFilter< PointSetType, ScalarImageType > BSplineFilterType
BiasFieldControlPointLatticeType::Pointer m_LogBiasFieldControlPointLattice
void PrintSelf(std::ostream &os, Indent indent) const override
itkSetInputMacro(ConfidenceImage, RealImageType)
void EnlargeOutputRequestedRegion(DataObject *) override
RealImagePointer UpdateBiasFieldEstimate(RealImageType *, vcl_size_t)
typename BSplineFilterType::PointDataImageType BiasFieldControlPointLatticeType
RealImagePointer ReconstructBiasField(const BiasFieldControlPointLatticeType *)
~N4BiasFieldCorrectionImageFilter() override=default
ImageToImageFilter< TInputImage, TOutputImage > Superclass
RealType CalculateConvergenceMeasurement(const RealImageType *, const RealImageType *) const
PointSet< ScalarType, Self::ImageDimension > PointSetType
A superclass of the N-dimensional mesh structure; supports point (geometric coordinate and attribute)...
Definition itkPointSet.h:82
Implements transparent reference counting.
A templated class holding a n-Dimensional vector.
Definition itkVector.h:63
The "itk" namespace contains all Insight Segmentation and Registration Toolkit (ITK) classes....
constexpr TContainer MakeFilled(typename TContainer::const_reference value)