ITK
4.1.0
Insight Segmentation and Registration Toolkit
|
00001 /*========================================================================= 00002 * 00003 * Copyright Insight Software Consortium 00004 * 00005 * Licensed under the Apache License, Version 2.0 (the "License"); 00006 * you may not use this file except in compliance with the License. 00007 * You may obtain a copy of the License at 00008 * 00009 * http://www.apache.org/licenses/LICENSE-2.0.txt 00010 * 00011 * Unless required by applicable law or agreed to in writing, software 00012 * distributed under the License is distributed on an "AS IS" BASIS, 00013 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 00014 * See the License for the specific language governing permissions and 00015 * limitations under the License. 00016 * 00017 *=========================================================================*/ 00018 #ifndef __itkN4BiasFieldCorrectionImageFilter_h 00019 #define __itkN4BiasFieldCorrectionImageFilter_h 00020 00021 #include "itkImageToImageFilter.h" 00022 00023 #include "itkArray.h" 00024 #include "itkBSplineScatteredDataPointSetToImageFilter.h" 00025 #include "itkPointSet.h" 00026 #include "itkVector.h" 00027 00028 #include "vnl/vnl_vector.h" 00029 00030 namespace itk { 00031 00090 template<class TInputImage, class TMaskImage = 00091 Image<unsigned char, ::itk::GetImageDimension<TInputImage>::ImageDimension>, 00092 class TOutputImage = TInputImage> 00093 class ITK_EXPORT N4BiasFieldCorrectionImageFilter : 00094 public ImageToImageFilter<TInputImage, TOutputImage> 00095 { 00096 public: 00098 typedef N4BiasFieldCorrectionImageFilter Self; 00099 typedef ImageToImageFilter<TInputImage, TOutputImage> Superclass; 00100 typedef SmartPointer<Self> Pointer; 00101 typedef SmartPointer<const Self> ConstPointer; 00102 00104 itkTypeMacro( N4BiasFieldCorrectionImageFilter, ImageToImageFilter ); 00105 00107 itkNewMacro( Self ); 00108 00110 itkStaticConstMacro( ImageDimension, unsigned int, 00111 TInputImage::ImageDimension ); 00112 00114 typedef TInputImage InputImageType; 00115 typedef TOutputImage OutputImageType; 00116 typedef TMaskImage MaskImageType; 00117 typedef typename MaskImageType::PixelType MaskPixelType; 00118 00119 typedef float RealType; 00120 typedef Image<RealType, ImageDimension> RealImageType; 00121 typedef typename RealImageType::Pointer RealImagePointer; 00122 typedef Array<unsigned int> VariableSizeArrayType; 00123 00125 typedef Vector<RealType, 1> ScalarType; 00126 typedef PointSet<ScalarType, itkGetStaticConstMacro( ImageDimension )> PointSetType; 00127 typedef Image<ScalarType, itkGetStaticConstMacro( ImageDimension )> ScalarImageType; 00128 typedef typename PointSetType::Pointer PointSetPointer; 00129 typedef typename PointSetType::PointType PointType; 00131 00133 typedef BSplineScatteredDataPointSetToImageFilter< 00134 PointSetType, ScalarImageType> BSplineFilterType; 00135 typedef typename BSplineFilterType::PointDataImageType BiasFieldControlPointLatticeType; 00136 typedef typename BSplineFilterType::ArrayType ArrayType; 00137 00141 void SetInput1( const InputImageType *image ) { this->SetInput( image ); } 00142 00143 00149 void SetMaskImage( const MaskImageType *mask ) 00150 { 00151 this->SetNthInput( 1, const_cast<MaskImageType *>( mask ) ); 00152 } 00153 void SetInput2( const MaskImageType *mask ) { this->SetMaskImage( mask ); } 00155 00161 const MaskImageType* GetMaskImage() const 00162 { 00163 return static_cast<const MaskImageType*>( this->ProcessObject::GetInput( 1 ) ); 00164 } 00165 00176 void SetConfidenceImage( const RealImageType *image ) 00177 { 00178 this->SetNthInput( 2, const_cast<RealImageType *>( image ) ); 00179 } 00180 void SetInput3( const RealImageType *image ) { this->SetConfidenceImage( image ); } 00182 00193 const RealImageType* GetConfidenceImage() const 00194 { 00195 return static_cast<const RealImageType*>( this->ProcessObject::GetInput( 2 ) ); 00196 } 00197 00203 itkSetMacro( MaskLabel, MaskPixelType ); 00204 00210 itkGetConstMacro( MaskLabel, MaskPixelType ); 00211 00212 // Sharpen histogram parameters: in estimating the bias field, the 00213 // first step is to sharpen the intensity histogram by Wiener deconvolution 00214 // with a 1-D Gaussian. The following parameters define this operation. 00215 // These default values in N4 match the default values in N3. 00216 00221 itkSetMacro( NumberOfHistogramBins, unsigned int ); 00222 00227 itkGetConstMacro( NumberOfHistogramBins, unsigned int ); 00228 00232 itkSetMacro( WienerFilterNoise, RealType ); 00233 00237 itkGetConstMacro( WienerFilterNoise, RealType ); 00238 00243 itkSetMacro( BiasFieldFullWidthAtHalfMaximum, RealType ); 00244 00249 itkGetConstMacro( BiasFieldFullWidthAtHalfMaximum, RealType ); 00250 00251 // B-spline parameters governing the fitting routine 00252 00256 itkSetMacro( SplineOrder, unsigned int ); 00257 00261 itkGetConstMacro( SplineOrder, unsigned int ); 00262 00270 itkSetMacro( NumberOfControlPoints, ArrayType ); 00271 00279 itkGetConstMacro( NumberOfControlPoints, ArrayType ); 00280 00287 itkSetMacro( NumberOfFittingLevels, ArrayType ); 00288 00295 void SetNumberOfFittingLevels( unsigned int n ) 00296 { 00297 ArrayType nlevels; 00298 00299 nlevels.Fill( n ); 00300 this->SetNumberOfFittingLevels( nlevels ); 00301 } 00302 00309 itkGetConstMacro( NumberOfFittingLevels, ArrayType ); 00310 00315 itkSetMacro( MaximumNumberOfIterations, VariableSizeArrayType ); 00316 00321 itkGetConstMacro( MaximumNumberOfIterations, VariableSizeArrayType ); 00322 00330 itkSetMacro( ConvergenceThreshold, RealType ); 00331 00339 itkGetConstMacro( ConvergenceThreshold, RealType ); 00340 00350 itkGetConstMacro( LogBiasFieldControlPointLattice, 00351 typename BiasFieldControlPointLatticeType::Pointer ); 00352 00357 itkGetConstMacro( ElapsedIterations, unsigned int ); 00358 00363 itkGetConstMacro( CurrentConvergenceMeasurement, RealType ); 00364 00369 itkGetConstMacro( CurrentLevel, unsigned int ); 00370 00371 protected: 00372 N4BiasFieldCorrectionImageFilter(); 00373 ~N4BiasFieldCorrectionImageFilter() {} 00374 void PrintSelf( std::ostream& os, Indent indent ) const; 00375 00376 void GenerateData(); 00377 00378 private: 00379 N4BiasFieldCorrectionImageFilter( const Self& ); //purposely not 00380 // implemented 00381 void operator=( const Self& ); //purposely not 00382 // implemented 00383 00384 // N4 algorithm functions: The basic algorithm iterates between sharpening 00385 // the intensity histogram of the corrected input image and spatially 00386 // smoothing those results with a B-spline scalar field estimate of the 00387 // bias field. The former is handled by the function SharpenImage() 00388 // whereas the latter is handled by the function UpdateBiasFieldEstimate(). 00389 // Convergence is determined by the coefficient of variation of the difference 00390 // image between the current bias field estimate and the previous estimate. 00391 00397 RealImagePointer SharpenImage( const RealImageType * ) const; 00398 00404 RealImagePointer UpdateBiasFieldEstimate( RealImageType * ); 00405 00410 RealType CalculateConvergenceMeasurement( const RealImageType *, const RealImageType * ) const; 00411 00412 00413 MaskPixelType m_MaskLabel; 00414 00415 // Parameters for deconvolution with Wiener filter 00416 00417 unsigned int m_NumberOfHistogramBins; 00418 RealType m_WienerFilterNoise; 00419 RealType m_BiasFieldFullWidthAtHalfMaximum; 00420 00421 // Convergence parameters 00422 00423 VariableSizeArrayType m_MaximumNumberOfIterations; 00424 unsigned int m_ElapsedIterations; 00425 RealType m_ConvergenceThreshold; 00426 RealType m_CurrentConvergenceMeasurement; 00427 unsigned int m_CurrentLevel; 00428 00429 // B-spline fitting parameters 00430 00431 typename 00432 BiasFieldControlPointLatticeType::Pointer m_LogBiasFieldControlPointLattice; 00433 00434 unsigned int m_SplineOrder; 00435 ArrayType m_NumberOfControlPoints; 00436 ArrayType m_NumberOfFittingLevels; 00437 00438 }; 00439 00440 } // end namespace itk 00441 00442 #ifndef ITK_MANUAL_INSTANTIATION 00443 #include "itkN4BiasFieldCorrectionImageFilter.hxx" 00444 #endif 00445 00446 #endif 00447