Medical Imaging Interaction Toolkit  2016.11.0
Medical Imaging Interaction Toolkit
itkOdfMaximaExtractionFilter.h
Go to the documentation of this file.
1 /*=========================================================================
2 
3  Program: Insight Segmentation & Registration Toolkit
4  Module: $RCSfile: itkDiffusionTensor3DReconstructionImageFilter.h,v $
5  Language: C++
6  Date: $Date: 2006-03-27 17:01:06 $
7  Version: $Revision: 1.12 $
8 
9  Copyright (c) Insight Software Consortium. All rights reserved.
10  See ITKCopyright.txt or http://www.itk.org/HTML/Copyright.htm for details.
11 
12  This software is distributed WITHOUT ANY WARRANTY; without even
13  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
14  PURPOSE. See the above copyright notices for more information.
15 
16 =========================================================================*/
17 #ifndef __itkOdfMaximaExtractionFilter_h_
18 #define __itkOdfMaximaExtractionFilter_h_
19 
20 #include <mitkFiberBundle.h>
21 
22 #include <itkImageToImageFilter.h>
23 #include <itkVectorContainer.h>
26 
27 namespace itk{
32 template< class TOdfPixelType >
33 class OdfMaximaExtractionFilter : public ProcessObject
34 {
35 
36 public:
37 
42  };
43 
47  typedef ProcessObject Superclass;
48 
50  itkFactorylessNewMacro(Self)
51  itkCloneMacro(Self)
52 
54  itkTypeMacro(OdfMaximaExtractionFilter, ImageToImageFilter)
55 
57  typedef QballReconstructionFilterType::CoefficientImageType CoefficientImageType;
58  typedef CoefficientImageType::PixelType CoefficientPixelType;
59  typedef itk::VectorImage< short, 3 > DiffusionImageType;
60  typedef vnl_vector_fixed< double, 3 > Vector3D;
61  typedef VectorContainer< unsigned int, Vector3D > DirectionContainerType;
62  typedef VectorContainer< unsigned int, DirectionContainerType::Pointer > ContainerType;
63  typedef itk::Image<unsigned char, 3> ItkUcharImgType;
64  typedef itk::Image< itk::Vector< float, 3 >, 3> ItkDirectionImage;
65  typedef itk::VectorContainer< unsigned int, ItkDirectionImage::Pointer > ItkDirectionImageContainer;
66 
67  // output
68  itkGetMacro( OutputFiberBundle, mitk::FiberBundle::Pointer)
69  itkGetMacro( NumDirectionsImage, ItkUcharImgType::Pointer)
70  itkGetMacro( DirectionImageContainer, ItkDirectionImageContainer::Pointer)
71 
72  // input
73  itkSetMacro( MaskImage, ItkUcharImgType::Pointer)
74  itkSetMacro( NormalizationMethod, NormalizationMethods)
75  itkSetMacro( DiffusionGradients, DirectionContainerType::Pointer)
76  itkSetMacro( DiffusionImage, DiffusionImageType::Pointer)
77  itkSetMacro( Bvalue, float)
78  itkSetMacro( ShCoeffImage, CoefficientImageType::Pointer)
79  itkSetMacro( MaxNumPeaks, unsigned int)
80  itkSetMacro( PeakThreshold, double)
81 
82  void GenerateData() override;
83 
84 protected:
85  OdfMaximaExtractionFilter();
86  ~OdfMaximaExtractionFilter(){}
87 
89  bool ReconstructQballImage();
90 
92  std::vector<double> SolveCubic(const double& a, const double& b, const double& c, const double& d);
93 
95  double ODF_dtheta2(const double& sn, const double& cs, const double& A, const double& B, const double& C, const double& D, const double& E, const double& F, const double& G, const double& H);
96 
98  double ODF_dphi2(const double& sn, const double& cs, const double& A, const double& B, const double& C, const double& D, const double& E, const double& F, const double& G, const double& H);
99 
101  double ODF_dtheta(const double& sn, const double& cs, const double& A, const double& B, const double& C, const double& D, const double& E, const double& F, const double& G, const double& H);
102 
104  void FindCandidatePeaks(const CoefficientPixelType& SHcoeff);
105 
107  std::vector< Vector3D > ClusterPeaks(const CoefficientPixelType& shCoeff);
108 
109  void Cart2Sph(const std::vector< Vector3D >& dir, vnl_matrix<double>& sphCoords);
110  vnl_matrix<double> CalcShBasis(vnl_matrix<double>& sphCoords, const int& shOrder);
111 
112  // diffusion weighted image (mandatory input)
115  float m_Bvalue;
116 
117  // binary mask image (optional input)
119 
120  // input parameters
123  unsigned int m_MaxNumPeaks;
124 
125  // intermediate results
127  std::vector< vnl_vector_fixed< double, 2 > > m_CandidatePeaks;
128 
129  // output data
133 
134 private:
135 
136 };
137 
138 }
139 
140 #ifndef ITK_MANUAL_INSTANTIATION
142 #endif
143 
144 #endif //__itkOdfMaximaExtractionFilter_h_
145 
unsigned int m_MaxNumPeaks
if more peaks are found, only the largest are kept
itk::SmartPointer< Self > Pointer
double ODF_dtheta(const double &sn, const double &cs, const double &A, const double &B, const double &C, const double &D, const double &E, const double &F, const double &G, const double &H)
QballReconstructionFilterType::CoefficientImageType CoefficientImageType
itk::VectorContainer< unsigned int, ItkDirectionImage::Pointer > ItkDirectionImageContainer
itk::Image< itk::Vector< float, 3 >, 3 > ItkDirectionImage
NormalizationMethods m_NormalizationMethod
normalization for peaks
double ODF_dphi2(const double &sn, const double &cs, const double &A, const double &B, const double &C, const double &D, const double &E, const double &F, const double &G, const double &H)
DataCollection - Class to facilitate loading/accessing structured data.
CoefficientImageType::PixelType CoefficientPixelType
void FindCandidatePeaks(const CoefficientPixelType &SHcoeff)
vnl_vector_fixed< double, 3 > Vector3D
ItkUcharImgType::Pointer m_NumDirectionsImage
number of peaks per voxel
CoefficientImageType::Pointer m_ShCoeffImage
conatins spherical harmonic coefficients
std::vector< double > SolveCubic(const double &a, const double &b, const double &c, const double &d)
vnl_matrix< double > CalcShBasis(vnl_matrix< double > &sphCoords, const int &shOrder)
float m_Bvalue
input for qball reconstruction
class ITK_EXPORT Image
This class takes as input one or more reference image (acquired in the absence of diffusion sensitizi...
double m_PeakThreshold
threshold on peak length relative to the largest peak in the current voxel
ItkDirectionImageContainer::Pointer m_DirectionImageContainer
output peaks
std::vector< vnl_vector_fixed< double, 2 > > m_CandidatePeaks
container for candidate peaks (all extrema, also minima)
DiffusionImageType::Pointer m_DiffusionImage
input for qball reconstruction
double ODF_dtheta2(const double &sn, const double &cs, const double &A, const double &B, const double &C, const double &D, const double &E, const double &F, const double &G, const double &H)
#define QBALL_ODFSIZE
mitk::FiberBundle::Pointer m_OutputFiberBundle
vector field (peak sizes rescaled for visualization purposes)
void Cart2Sph(const std::vector< Vector3D > &dir, vnl_matrix< double > &sphCoords)
itk::VectorImage< short, 3 > DiffusionImageType
std::vector< Vector3D > ClusterPeaks(const CoefficientPixelType &shCoeff)
VectorContainer< unsigned int, Vector3D > DirectionContainerType
VectorContainer< unsigned int, DirectionContainerType::Pointer > ContainerType
unsigned short PixelType
itk::Image< unsigned char, 3 > ItkUcharImgType
ItkUcharImgType::Pointer m_MaskImage
only voxels inside the binary mask are processed
DirectionContainerType::Pointer m_DiffusionGradients
input for qball reconstruction
largest peak is normalized to length 1, other peaks relative to it