ITK  6.0.0
Insight Toolkit
itkMRIBiasFieldCorrectionFilter.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 itkMRIBiasFieldCorrectionFilter_h
19#define itkMRIBiasFieldCorrectionFilter_h
20
21#include <ctime>
22
24#include "itkArray2D.h"
31
32#include <memory> // For unique_ptr.
33
34
35namespace itk
36{
50template <typename TImage, typename TImageMask, typename TBiasField>
51class ITK_TEMPLATE_EXPORT MRIBiasEnergyFunction : public SingleValuedCostFunction
52{
53public:
54 ITK_DISALLOW_COPY_AND_MOVE(MRIBiasEnergyFunction);
55
61
63 itkOverrideGetNameOfClassMacro(MRIBiasEnergyFunction);
64
66 itkNewMacro(Self);
67
69 using ImageType = TImage;
71 using ImageElementType = typename ImageType::PixelType;
74 using MaskType = TImageMask;
76 using MaskElementType = typename MaskType::PixelType;
77
79 using BiasFieldType = TBiasField;
80
83 using typename Superclass::ParametersType;
84
86 using DerivativeType = Superclass::DerivativeType;
87
89 using MeasureType = Superclass::MeasureType;
90
91 static constexpr unsigned int SpaceDimension = 3;
92
94 using InternalEnergyFunction = CompositeValleyFunction;
95
97 using SamplingFactorType = unsigned int[SpaceDimension];
98
100 itkSetObjectMacro(Image, ImageType);
101
103 itkSetObjectMacro(Mask, MaskType);
104
106 itkSetMacro(Region, ImageRegionType);
107
109 void
111 {
112 m_BiasField = bias;
113 }
114
117 void
119 {
120 for (unsigned int i = 0; i < SpaceDimension; ++i)
121 {
122 m_SamplingFactor[i] = factor[i];
123 }
124 }
129 double
130 GetEnergy0(double diff)
131 {
132 return (*m_InternalEnergyFunction)(diff);
133 }
134
137 MeasureType
138 GetValue(const ParametersType & parameters) const override;
139
142 void
143 GetDerivative(const ParametersType & itkNotUsed(parameters), DerivativeType & itkNotUsed(derivative)) const override
144 {}
145
150 void
152
153 unsigned int
154 GetNumberOfParameters() const override;
155
156protected:
159
161 ~MRIBiasEnergyFunction() override = default;
162
163private:
165 BiasFieldType * m_BiasField{};
166
168 ImagePointer m_Image{};
169
171 MaskPointer m_Mask{};
172
174 ImageRegionType m_Region{};
175
177 std::unique_ptr<InternalEnergyFunction> m_InternalEnergyFunction;
178
180 SamplingFactorType m_SamplingFactor{};
181}; // end of class
182
234template <typename TInputImage, typename TOutputImage, typename TMaskImage>
235class ITK_TEMPLATE_EXPORT MRIBiasFieldCorrectionFilter : public ImageToImageFilter<TInputImage, TOutputImage>
236{
237public:
238 ITK_DISALLOW_COPY_AND_MOVE(MRIBiasFieldCorrectionFilter);
246
248 itkNewMacro(Self);
249
251 itkOverrideGetNameOfClassMacro(MRIBiasFieldCorrectionFilter);
252
254 static constexpr unsigned int ImageDimension = TOutputImage::ImageDimension;
255
257 using OutputImageType = TOutputImage;
260 using OutputImagePixelType = typename TOutputImage::PixelType;
263
264 using InputImageType = TInputImage;
267 using InputImagePixelType = typename TInputImage::PixelType;
270
272 using ImageMaskType = TMaskImage;
275
281
285 using SlabRegionVectorIteratorType = typename SlabRegionVectorType::iterator;
286
289
293
296
299
302
306 void
308 itkGetModifiableObjectMacro(InputMask, ImageMaskType);
313 void
315 itkGetModifiableObjectMacro(OutputMask, ImageMaskType);
318#if defined(ITK_LEGACY_REMOVE)
322 void
324 {
325 m_BiasFieldMultiplicative = flag;
326 }
327#endif // defined ( ITK_LEGACY_REMOVE )
328
332 itkSetMacro(BiasFieldMultiplicative, bool);
333 itkGetConstMacro(BiasFieldMultiplicative, bool);
334 itkBooleanMacro(BiasFieldMultiplicative);
337#if defined(ITK_LEGACY_REMOVE)
339 bool
341 {
342 return m_BiasFieldMultiplicative;
343 }
344#endif // defined ( ITK_LEGACY_REMOVE )
345
350 itkSetMacro(UsingInterSliceIntensityCorrection, bool);
351 itkGetConstMacro(UsingInterSliceIntensityCorrection, bool);
352 itkBooleanMacro(UsingInterSliceIntensityCorrection);
360 itkSetMacro(UsingSlabIdentification, bool);
361 itkGetConstMacro(UsingSlabIdentification, bool);
362 itkBooleanMacro(UsingSlabIdentification);
365 itkSetMacro(SlabBackgroundMinimumThreshold, InputImagePixelType);
366 itkGetConstReferenceMacro(SlabBackgroundMinimumThreshold, InputImagePixelType);
367
368 itkSetMacro(SlabNumberOfSamples, unsigned int);
369 itkGetConstReferenceMacro(SlabNumberOfSamples, unsigned int);
370
371 itkSetMacro(SlabTolerance, double);
372 itkGetConstReferenceMacro(SlabTolerance, double);
373
379 itkSetMacro(UsingBiasFieldCorrection, bool);
380 itkGetConstMacro(UsingBiasFieldCorrection, bool);
381 itkBooleanMacro(UsingBiasFieldCorrection);
386 itkSetMacro(GeneratingOutput, bool);
387 itkGetConstMacro(GeneratingOutput, bool);
388 itkBooleanMacro(GeneratingOutput);
393 itkSetMacro(SlicingDirection, int);
394 itkGetConstMacro(SlicingDirection, int);
398 itkSetMacro(BiasFieldDegree, int);
399 itkGetConstMacro(BiasFieldDegree, int);
404 void
406 {
407 this->Modified();
408 m_BiasFieldCoefficients = coefficients;
409 }
415 BiasFieldType::CoefficientArrayType
417 {
418 return m_EstimatedBiasFieldCoefficients;
419 }
420
424 void
426
429 itkSetMacro(VolumeCorrectionMaximumIteration, int);
430 itkGetConstMacro(VolumeCorrectionMaximumIteration, int);
435 itkSetMacro(InterSliceCorrectionMaximumIteration, int);
436 itkGetConstMacro(InterSliceCorrectionMaximumIteration, int);
440 itkSetMacro(OptimizerInitialRadius, double);
441 itkGetConstMacro(OptimizerInitialRadius, double);
445 itkSetMacro(OptimizerGrowthFactor, double);
446 itkGetConstMacro(OptimizerGrowthFactor, double);
451 itkSetMacro(OptimizerShrinkFactor, double);
452 itkGetConstMacro(OptimizerShrinkFactor, double);
453
460 void
461 SetNumberOfLevels(unsigned int num);
462
464 itkGetConstMacro(NumberOfLevels, unsigned int);
465
472 void
473 SetSchedule(const ScheduleType & schedule);
474
476 itkGetConstReferenceMacro(Schedule, ScheduleType);
477
482 void
483 SetStartingShrinkFactors(unsigned int factor);
484
485 void
486 SetStartingShrinkFactors(const unsigned int * factors);
487
489 const unsigned int *
491
495 static bool
497
504 void
506
510 EstimateBiasField(InputImageRegionType region, unsigned int degree, int maximumIteration);
511
515 void
517
520 void
522
523protected:
525 ~MRIBiasFieldCorrectionFilter() override = default;
526 void
527 PrintSelf(std::ostream & os, Indent indent) const override;
528
531 bool
533
534protected:
538 void
540
543 void
545
548 template <typename TSource, typename TTarget>
549 void
550 CopyAndConvertImage(const TSource * source, TTarget * target, typename TTarget::RegionType requestedRegion)
551 {
552 using SourceIterator = ImageRegionConstIterator<TSource>;
553 using TargetIterator = ImageRegionIterator<TTarget>;
554 using TargetPixelType = typename TTarget::PixelType;
555
556 SourceIterator s_iter(source, requestedRegion);
557 TargetIterator t_iter(target, requestedRegion);
558
559 s_iter.GoToBegin();
560 t_iter.GoToBegin();
561 while (!s_iter.IsAtEnd())
562 {
563 t_iter.Set(static_cast<TargetPixelType>(s_iter.Get()));
564 ++s_iter;
565 ++t_iter;
566 }
567 }
568
573 void
575
579 void
581
582 void
583 GenerateData() override;
584
585private:
587 EnergyFunctionPointer m_EnergyFunction{};
588
590 NormalVariateGeneratorType::Pointer m_NormalVariateGenerator{};
591
593 ImageMaskPointer m_InputMask{};
594
596 ImageMaskPointer m_OutputMask{};
597
599 InternalImagePointer m_InternalInput{};
600
603
605 int m_SlicingDirection{};
606
608 bool m_BiasFieldMultiplicative{ true };
609
611 bool m_UsingInterSliceIntensityCorrection{};
612 bool m_UsingSlabIdentification{ false };
613 bool m_UsingBiasFieldCorrection{ true };
614 bool m_GeneratingOutput{ true };
615
616 unsigned int m_SlabNumberOfSamples{ 200 };
617 InputImagePixelType m_SlabBackgroundMinimumThreshold{};
618 double m_SlabTolerance{ 0.0 };
619
621 int m_BiasFieldDegree{ 3 };
622
624 unsigned int m_NumberOfLevels{ 0 };
625
627 ScheduleType m_Schedule{};
628
631 BiasFieldType::CoefficientArrayType m_BiasFieldCoefficients{};
632
635 BiasFieldType::CoefficientArrayType m_EstimatedBiasFieldCoefficients{};
636
637 int m_VolumeCorrectionMaximumIteration{ 2000 };
638
639 int m_InterSliceCorrectionMaximumIteration{ 4000 };
640
642 double m_OptimizerInitialRadius{ 1.01 };
643
645 double m_OptimizerGrowthFactor{ 1.05 };
646
648 double m_OptimizerShrinkFactor{};
649
651 Array<double> m_TissueClassMeans{};
652
654 Array<double> m_TissueClassSigmas{};
655};
656
657} // end namespace itk
658
659#ifndef ITK_MANUAL_INSTANTIATION
660# include "itkMRIBiasFieldCorrectionFilter.hxx"
661#endif
662
663#endif
Array class with size defined at construction time.
Definition: itkArray.h:48
A multi-dimensional iterator templated over image type that walks a region of pixels.
A multi-dimensional iterator templated over image type that walks a region of pixels.
Base class for all process objects that output image data.
typename OutputImageType::PixelType OutputImagePixelType
typename OutputImageType::RegionType OutputImageRegionType
TOutputImage OutputImageType
typename OutputImageType::Pointer OutputImagePointer
Base class for filters that take an image as input and produce an image as output.
typename InputImageType::Pointer InputImagePointer
typename InputImageType::PixelType InputImagePixelType
typename InputImageType::RegionType InputImageRegionType
Templated n-dimensional image class.
Definition: itkImage.h:89
TPixel PixelType
Definition: itkImage.h:108
Control indentation during Print() invocation.
Definition: itkIndent.h:50
Light weight base class for most itk classes.
Identifies slabs in MR images comparing minimum intensity averages.
std::vector< ImageRegionType > SlabRegionVectorType
Represents a cost function for MRI bias field correction optimization.
typename ImageType::IndexType ImageIndexType
CompositeValleyFunction InternalEnergyFunction
void SetSamplingFactors(const SamplingFactorType factor)
unsigned int[SpaceDimension] SamplingFactorType
std::unique_ptr< InternalEnergyFunction > m_InternalEnergyFunction
~MRIBiasEnergyFunction() override=default
unsigned int GetNumberOfParameters() const override
void InitializeDistributions(Array< double > classMeans, Array< double > classSigmas)
typename ImageType::PixelType ImageElementType
typename ImageType::RegionType ImageRegionType
void GetDerivative(const ParametersType &, DerivativeType &) const override
typename MaskType::PixelType MaskElementType
MeasureType GetValue(const ParametersType &parameters) const override
typename ImageType::Pointer ImagePointer
typename MaskType::Pointer MaskPointer
void SetNumberOfLevels(unsigned int num)
typename TOutputImage::IndexType OutputImageIndexType
void SetStartingShrinkFactors(unsigned int factor)
typename TInputImage::IndexType InputImageIndexType
void SetInitialBiasFieldCoefficients(const BiasFieldType::CoefficientArrayType &coefficients)
void CorrectImage(BiasFieldType &bias, InputImageRegionType region)
void GetBiasFieldSize(InputImageRegionType region, BiasFieldType::DomainSizeType &biasSize)
static bool IsScheduleDownwardDivisible(const ScheduleType &schedule)
const unsigned int * GetStartingShrinkFactors() const
bool CheckMaskImage(ImageMaskType *mask)
void SetOutputMask(ImageMaskType *outputMask)
typename InternalImageType::PixelType InternalImagePixelType
BiasFieldType EstimateBiasField(InputImageRegionType region, unsigned int degree, int maximumIteration)
typename ImageMaskType::Pointer ImageMaskPointer
void PrintSelf(std::ostream &os, Indent indent) const override
void SetSchedule(const ScheduleType &schedule)
BiasFieldType::CoefficientArrayType GetEstimatedBiasFieldCoefficients()
typename MRASlabIdentifierType::SlabRegionVectorType SlabRegionVectorType
typename ImageMaskType::RegionType ImageMaskRegionType
typename TInputImage::SizeType InputImageSizeType
void CorrectInterSliceIntensityInhomogeneity(InputImageRegionType region)
typename TOutputImage::SizeType OutputImageSizeType
void ExpImage(InternalImageType *source, InternalImageType *target)
typename InternalImageType::Pointer InternalImagePointer
typename InternalImageType::RegionType InternalImageRegionType
void AdjustSlabRegions(SlabRegionVectorType &slabs, OutputImageRegionType requestedRegion)
typename SlabRegionVectorType::iterator SlabRegionVectorIteratorType
void SetStartingShrinkFactors(const unsigned int *factors)
void CopyAndConvertImage(const TSource *source, TTarget *target, typename TTarget::RegionType requestedRegion)
~MRIBiasFieldCorrectionFilter() override=default
typename EnergyFunctionType::Pointer EnergyFunctionPointer
void Log1PImage(InternalImageType *source, InternalImageType *target)
void SetInputMask(ImageMaskType *inputMask)
void SetTissueClassStatistics(const Array< double > &means, const Array< double > &sigmas)
2D and 3D multivariate Legendre Polynomial
The base class for all process objects (source, filters, mappers) in the Insight data processing pipe...
A region represents some portion or piece of data.
Definition: itkRegion.h:66
This class is a base for the CostFunctions returning a single value.
Superclass::ParametersType ParametersType
SmartPointer< Self > Pointer
The "itk" namespace contains all Insight Segmentation and Registration Toolkit (ITK) classes....