SimpleITK
|
Implementation of the N4 bias field correction algorithm. More...
#include <sitkN4BiasFieldCorrectionImageFilter.h>
Public Types | |
using | PixelIDTypeList = RealPixelIDTypeList |
using | Self = N4BiasFieldCorrectionImageFilter |
![]() | |
using | Self = ImageFilter |
![]() | |
using | Self = ProcessObject |
Public Member Functions | |
Image | Execute (const Image &image) |
Image | Execute (const Image &image, const Image &maskImage) |
double | GetBiasFieldFullWidthAtHalfMaximum () const |
double | GetConvergenceThreshold () const |
double | GetCurrentConvergenceMeasurement () const |
uint32_t | GetCurrentLevel () const |
uint32_t | GetElapsedIterations () const |
Image | GetLogBiasFieldAsImage (Image referenceImage) const |
The computed log bias field correction. Typically, a reduced size image is used as input to the N4 filter using something like itkShrinkImageFilter. Since the output is a corrected version of the input, the user will probably want to apply the bias field correction to the full resolution image. Returns the b-spline log bias field reconstructioned onto the space of the referenceImage parameter. An input image can be corrected by: input/exp(bias_field). More... | |
uint8_t | GetMaskLabel () const |
std::vector< uint32_t > | GetMaximumNumberOfIterations () const |
std::string | GetName () const |
std::vector< uint32_t > | GetNumberOfControlPoints () const |
uint32_t | GetNumberOfHistogramBins () const |
uint32_t | GetSplineOrder () const |
bool | GetUseMaskLabel () const |
double | GetWienerFilterNoise () const |
N4BiasFieldCorrectionImageFilter () | |
Self & | SetBiasFieldFullWidthAtHalfMaximum (double BiasFieldFullWidthAtHalfMaximum) |
Self & | SetConvergenceThreshold (double ConvergenceThreshold) |
Self & | SetMaskLabel (uint8_t MaskLabel) |
Self & | SetMaximumNumberOfIterations (std::vector< uint32_t > MaximumNumberOfIterations) |
Self & | SetNumberOfControlPoints (std::vector< uint32_t > NumberOfControlPoints) |
Self & | SetNumberOfControlPoints (uint32_t value) |
Self & | SetNumberOfHistogramBins (uint32_t NumberOfHistogramBins) |
Self & | SetSplineOrder (uint32_t SplineOrder) |
Self & | SetUseMaskLabel (bool UseMaskLabel) |
Self & | SetWienerFilterNoise (double WienerFilterNoise) |
std::string | ToString () const |
Self & | UseMaskLabelOff () |
Self & | UseMaskLabelOn () |
virtual | ~N4BiasFieldCorrectionImageFilter () |
![]() | |
ImageFilter () | |
virtual | ~ImageFilter ()=0 |
![]() | |
virtual void | Abort () |
virtual int | AddCommand (itk::simple::EventEnum event, const std::function< void()> &func) |
Directly add a callback to observe an event. More... | |
virtual int | AddCommand (itk::simple::EventEnum event, itk::simple::Command &cmd) |
Add a Command Object to observer the event. More... | |
virtual float | GetProgress () const |
An Active Measurement of the progress of execution. More... | |
virtual bool | HasCommand (itk::simple::EventEnum event) const |
Query of this object has any registered commands for event. More... | |
ProcessObject () | |
virtual void | RemoveAllCommands () |
Remove all registered commands. More... | |
virtual | ~ProcessObject () |
virtual void | DebugOn () |
virtual void | DebugOff () |
virtual bool | GetDebug () const |
virtual void | SetDebug (bool debugFlag) |
virtual void | SetNumberOfThreads (unsigned int n) |
virtual unsigned int | GetNumberOfThreads () const |
virtual void | SetNumberOfWorkUnits (unsigned int n) |
virtual unsigned int | GetNumberOfWorkUnits () const |
Private Types | |
using | MemberFunctionType = Image(Self::*)(const Image *image, const Image *maskImage) |
Private Member Functions | |
template<class TImageType > | |
Image | ExecuteInternal (const Image *image, const Image *maskImage) |
Private Attributes | |
double | m_BiasFieldFullWidthAtHalfMaximum {0.15} |
double | m_ConvergenceThreshold {0.001} |
itk::ProcessObject * | m_Filter {nullptr} |
uint8_t | m_MaskLabel {1} |
std::vector< uint32_t > | m_MaximumNumberOfIterations {std::vector<uint32_t>(4,50)} |
std::unique_ptr< detail::MemberFunctionFactory< MemberFunctionType > > | m_MemberFactory |
std::vector< uint32_t > | m_NumberOfControlPoints {std::vector<uint32_t>(3, 4)} |
uint32_t | m_NumberOfHistogramBins {200u} |
std::function< double()> | m_pfGetCurrentConvergenceMeasurement |
std::function< uint32_t()> | m_pfGetCurrentLevel |
std::function< uint32_t()> | m_pfGetElapsedIterations |
std::function< Image(Image)> | m_pfGetLogBiasFieldAsImage |
uint32_t | m_SplineOrder {3u} |
bool | m_UseMaskLabel {true} |
double | m_WienerFilterNoise {0.01} |
Friends | |
struct | detail::MemberFunctionAddressor< MemberFunctionType > |
Additional Inherited Members | |
![]() | |
static bool | GetGlobalDefaultDebug () |
static void | GlobalDefaultDebugOff () |
static void | GlobalDefaultDebugOn () |
static void | SetGlobalDefaultDebug (bool debugFlag) |
static void | GlobalWarningDisplayOn () |
static void | GlobalWarningDisplayOff () |
static void | SetGlobalWarningDisplay (bool flag) |
static bool | GetGlobalWarningDisplay () |
static double | GetGlobalDefaultCoordinateTolerance () |
Access the global tolerance to determine congruent spaces. More... | |
static void | SetGlobalDefaultCoordinateTolerance (double) |
Access the global tolerance to determine congruent spaces. More... | |
static double | GetGlobalDefaultDirectionTolerance () |
Access the global tolerance to determine congruent spaces. More... | |
static void | SetGlobalDefaultDirectionTolerance (double) |
Access the global tolerance to determine congruent spaces. More... | |
static bool | SetGlobalDefaultThreader (const std::string &threader) |
Set/Get the default threader used for process objects. More... | |
static std::string | GetGlobalDefaultThreader () |
Set/Get the default threader used for process objects. More... | |
static void | SetGlobalDefaultNumberOfThreads (unsigned int n) |
static unsigned int | GetGlobalDefaultNumberOfThreads () |
Set/Get the default threader used for process objects. More... | |
![]() | |
void | CheckImageMatchingDimension (const Image &image1, const Image &image2, const std::string &image2Name) |
void | CheckImageMatchingPixelType (const Image &image1, const Image &image2, const std::string &image2Name) |
void | CheckImageMatchingSize (const Image &image1, const Image &image2, const std::string &image2Name) |
![]() | |
virtual unsigned long | AddITKObserver (const itk::EventObject &, itk::Command *) |
virtual itk::ProcessObject * | GetActiveProcess () |
virtual void | OnActiveProcessDelete () |
virtual void | onCommandDelete (const itk::simple::Command *cmd) noexcept |
virtual void | PreUpdate (itk::ProcessObject *p) |
virtual void | RemoveITKObserver (EventCommand &e) |
![]() | |
NonCopyable ()=default | |
NonCopyable (const NonCopyable &)=delete | |
NonCopyable & | operator= (const NonCopyable &)=delete |
![]() | |
template<class TImageType > | |
static void | FixNonZeroIndex (TImageType *img) |
![]() | |
template<class TImageType > | |
static TImageType::ConstPointer | CastImageToITK (const Image &img) |
template<class TPixelType , unsigned int VImageDimension, unsigned int VLength, template< typename, unsigned int > class TVector> | |
static Image | CastITKToImage (itk::Image< TVector< TPixelType, VLength >, VImageDimension > *img) |
template<unsigned int VImageDimension, unsigned int VLength, template< unsigned int > class TVector> | |
static Image | CastITKToImage (itk::Image< TVector< VLength >, VImageDimension > *img) |
template<class TImageType > | |
static Image | CastITKToImage (TImageType *img) |
static const itk::EventObject & | GetITKEventObject (EventEnum e) |
template<typename T > | |
static std::ostream & | ToStringHelper (std::ostream &os, const T &v) |
static std::ostream & | ToStringHelper (std::ostream &os, const char &v) |
static std::ostream & | ToStringHelper (std::ostream &os, const signed char &v) |
static std::ostream & | ToStringHelper (std::ostream &os, const unsigned char &v) |
Implementation of the N4 bias field correction algorithm.
The nonparametric nonuniform intensity normalization (N3) algorithm, as introduced by Sled et al. in 1998 is a method for correcting nonuniformity associated with MR images. The algorithm assumes a simple parametric model (Gaussian) for the bias field and does not require tissue class segmentation. In addition, there are only a couple of parameters to tune with the default values performing quite well. N3 has been publicly available as a set of perl scripts (https://www.bic.mni.mcgill.ca/ServicesSoftwareAdvancedImageProcessingTools/HomePage )
The N4 algorithm, encapsulated with this class, is a variation of the original N3 algorithm with the additional benefits of an improved B-spline fitting routine which allows for multiple resolutions to be used during the correction process. We also modify the iterative update component of algorithm such that the residual bias field is continually updated
Notes for the user:
The basic algorithm iterates between sharpening the intensity histogram of the corrected input image and spatially smoothing those results with a B-spline scalar field estimate of the bias field.
Contributed by Nicholas J. Tustison, James C. Gee in the Insight Journal paper: https://www.insight-journal.org/browse/publication/640
J.G. Sled, A.P. Zijdenbos and A.C. Evans. "A Nonparametric Method for Automatic Correction of Intensity Nonuniformity in Data" IEEE Transactions on Medical Imaging, Vol 17, No 1. Feb 1998.
N.J. Tustison, B.B. Avants, P.A. Cook, Y. Zheng, A. Egan, P.A. Yushkevich, and J.C. Gee. "N4ITK: Improved N3 Bias Correction" IEEE Transactions on Medical Imaging, 29(6):1310-1320, June 2010.
Definition at line 73 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
private |
Setup for member function dispatching
Definition at line 238 of file sitkN4BiasFieldCorrectionImageFilter.h.
Define the pixels types supported by this filter
Definition at line 85 of file sitkN4BiasFieldCorrectionImageFilter.h.
Definition at line 75 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
virtual |
Destructor
itk::simple::N4BiasFieldCorrectionImageFilter::N4BiasFieldCorrectionImageFilter | ( | ) |
Default Constructor that takes no arguments and initializes default parameters
Image itk::simple::N4BiasFieldCorrectionImageFilter::Execute | ( | const Image & | image, |
const Image & | maskImage | ||
) |
Execute the filter on the input image
|
private |
|
inline |
Get the full width at half maximum parameter characterizing the width of the Gaussian deconvolution. Default = 0.15.
Definition at line 116 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
inline |
Get the convergence threshold. Convergence is determined by the coefficient of variation of the difference image between the current bias field estimate and the previous estimate. If this value is less than the specified threshold, the algorithm proceeds to the next fitting level or terminates if it is at the last level.
Definition at line 96 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
inline |
Get the current convergence measurement. This is a helper function for reporting observations.
This is an active measurement. It may be accessed while the filter is being executing in command call-backs and can be accessed after execution.
Definition at line 209 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
inline |
Get the current fitting level. This is a helper function for reporting observations.
This is an active measurement. It may be accessed while the filter is being executing in command call-backs and can be accessed after execution.
Definition at line 191 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
inline |
Get the number of elapsed iterations. This is a helper function for reporting observations.
This is an active measurement. It may be accessed while the filter is being executing in command call-backs and can be accessed after execution.
Definition at line 200 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
inline |
The computed log bias field correction. Typically, a reduced size image is used as input to the N4 filter using something like itkShrinkImageFilter. Since the output is a corrected version of the input, the user will probably want to apply the bias field correction to the full resolution image. Returns the b-spline log bias field reconstructioned onto the space of the referenceImage parameter. An input image can be corrected by: input/exp(bias_field).
This is an active measurement. It may be accessed while the filter is being executing in command call-backs and can be accessed after execution.
Definition at line 219 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
inline |
Set/Get mask label value. If a binary mask image is specified and if UseMaskValue is true, only those input image voxels corresponding with mask image values equal to MaskLabel are used in estimating the bias field. If a MaskImage is specified and UseMaskLabel is false, all input image voxels corresponding to non-zero voxels in the MaskImage are used in estimating the bias field. Default = 1.
Definition at line 183 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
inline |
Get the maximum number of iterations specified at each fitting level. Default = 50.
Definition at line 106 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
inlinevirtual |
Name of this class
Implements itk::simple::ProcessObject.
Definition at line 223 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
inline |
Get the control point grid size defining the B-spline estimate of the scalar bias field. In each dimension, the B-spline mesh size is equal to the number of control points in that dimension minus the spline order. Default = 4 control points in each dimension for a mesh size of 1 in each dimension.
Definition at line 149 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
inline |
Get number of bins defining the log input intensity histogram. Default = 200.
Definition at line 136 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
inline |
Get the spline order defining the bias field estimate. Default = 3.
Definition at line 159 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
inline |
Use a mask label for identifying mask functionality. See SetMaskLabel. Defaults to true.
Definition at line 173 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
inline |
Get the noise estimate defining the Wiener filter. Default = 0.01.
Definition at line 126 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
inline |
Set the full width at half maximum parameter characterizing the width of the Gaussian deconvolution. Default = 0.15.
Definition at line 111 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
inline |
Set the convergence threshold. Convergence is determined by the coefficient of variation of the difference image between the current bias field estimate and the previous estimate. If this value is less than the specified threshold, the algorithm proceeds to the next fitting level or terminates if it is at the last level.
Definition at line 91 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
inline |
Set/Get mask label value. If a binary mask image is specified and if UseMaskValue is true, only those input image voxels corresponding with mask image values equal to MaskLabel are used in estimating the bias field. If a MaskImage is specified and UseMaskLabel is false, all input image voxels corresponding to non-zero voxels in the MaskImage are used in estimating the bias field. Default = 1.
Definition at line 178 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
inline |
Set the maximum number of iterations specified at each fitting level. Default = 50.
Definition at line 101 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
inline |
Set the control point grid size defining the B-spline estimate of the scalar bias field. In each dimension, the B-spline mesh size is equal to the number of control points in that dimension minus the spline order. Default = 4 control points in each dimension for a mesh size of 1 in each dimension.
Definition at line 141 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
inline |
Set the values of the NumberOfControlPoints vector all to value
Definition at line 144 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
inline |
Set number of bins defining the log input intensity histogram. Default = 200.
Definition at line 131 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
inline |
Set the spline order defining the bias field estimate. Default = 3.
Definition at line 154 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
inline |
Use a mask label for identifying mask functionality. See SetMaskLabel. Defaults to true.
Definition at line 164 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
inline |
Set the noise estimate defining the Wiener filter. Default = 0.01.
Definition at line 121 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
virtual |
Print ourselves out
Reimplemented from itk::simple::ProcessObject.
|
inline |
Definition at line 168 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
inline |
Set the value of UseMaskLabel to true or false respectfully.
Definition at line 167 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
friend |
Definition at line 242 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
private |
Definition at line 254 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
private |
Definition at line 248 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
private |
Definition at line 282 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
private |
Definition at line 270 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
private |
Definition at line 251 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
private |
Definition at line 244 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
private |
Definition at line 263 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
private |
Definition at line 260 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
private |
Definition at line 277 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
private |
Definition at line 273 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
private |
Definition at line 275 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
private |
Definition at line 279 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
private |
Definition at line 266 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
private |
Definition at line 268 of file sitkN4BiasFieldCorrectionImageFilter.h.
|
private |
Definition at line 257 of file sitkN4BiasFieldCorrectionImageFilter.h.