DMRITool  v0.1.1-139-g860d86b4
Diffusion MRI Tool
itkMeshFromSHCoefficientsImageFilter.h
Go to the documentation of this file.
1 /*=========================================================================
2 
3  Program: Mesh From SH Coefficients Image Filter
4 
5  Copyright (c) Pew-Thian Yap, Jian Cheng. All rights reserved.
6  See http://www.unc.edu/~ptyap/ for details.
7 
8  This software is distributed WITHOUT ANY WARRANTY; without even
9  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
10  PURPOSE. See the above copyright notices for more information.
11 
12  =========================================================================*/
13 #ifndef __itkMeshFromSHCoefficientsImageFilter_h
14 #define __itkMeshFromSHCoefficientsImageFilter_h
15 
16 
18 
19 namespace itk
20 {
21 
29 template < class TInputImage, class TOutputMesh=vtkPolyData>
31  : public MeshFromContinuousSphericalFunctionImageFilter<TInputImage, TOutputMesh>
32 {
33 public:
37  typedef SmartPointer<Self> Pointer;
38  typedef SmartPointer<const Self> ConstPointer;
39 
41  itkNewMacro(Self);
42 
45 
46  //** Convenient typedefs */
47  typedef TInputImage InputImageType;
48  typedef typename InputImageType::Pointer InputImagePointer;
49  typedef typename InputImageType::ConstPointer InputImageConstPointer;
50  typedef typename InputImageType::IndexType InputImageIndexType;
51  typedef typename InputImageType::SizeType InputImageSizeType;
52  typedef typename InputImageType::SizeValueType InputImageSizeValueType;
53  typedef typename InputImageType::SpacingType InputImageSpacingType;
54  typedef typename InputImageType::PixelType InputImagePixelType;
55  typedef typename InputImageType::RegionType InputImageRegionType;
56  typedef typename InputImageType::PointType InputImagePointType;
57 
58  typedef TOutputMesh OutputMeshPolyDataType;
59 
65 
66 
71 
72 
80 
81 
83  itkSetMacro(MaxOrder, int);
84  itkGetMacro(MaxOrder, int);
85 
87  {
88  utlShowPosition(this->GetDebug());
89 
90  InputImageConstPointer inputPtr = this->GetInput();
91  unsigned int maxNumberOfBasis = inputPtr->GetNumberOfComponentsPerPixel();
92  if (m_MaxOrder<=0)
93  m_MaxOrder = utl::DimToRankSH(maxNumberOfBasis);
94 
95  utlGlobalException(this->m_Orientations->Rows()==0, "need to set m_Orientations");
96  auto orientationsFliped = utl::FlipOrientations(*this->m_Orientations, this->m_Flip);
97  this->m_BasisMatrix = utl::ComputeSHMatrix(m_MaxOrder, orientationsFliped, CARTESIAN_TO_SPHERICAL);
98  // utl::PrintVnlMatrix(*this->m_BasisMatrix, "m_BasisMatrix 0");
99  return this->m_BasisMatrix;
100  }
101 
102  VectorType NormalizeUnitIntegral(const VectorType& x) const ITK_OVERRIDE
103  {
104  VectorType result(x);
105  if (std::fabs(x[0]>1e-8))
106  result /= (std::sqrt(4.0 * M_PI) * x[0]);
107  return result;
108  }
109 
110 protected:
112  {
113  m_MaxOrder = -1;
114  }
116  {};
117 
118  void PrintSelf(std::ostream& os, Indent indent) const ITK_OVERRIDE
119  {
120  Superclass::PrintSelf(os, indent);
121  PrintVar1(true, m_MaxOrder, os<<indent);
122  }
123 
124  typename LightObject::Pointer InternalClone() const ITK_OVERRIDE
125  {
126  typename LightObject::Pointer loPtr = Superclass::InternalClone();
127 
128  typename Self::Pointer rval = dynamic_cast<Self *>(loPtr.GetPointer());
129  if(rval.IsNull())
130  {
131  itkExceptionMacro(<< "downcast to type " << this->GetNameOfClass()<< " failed.");
132  }
133 
134  rval->m_MaxOrder = m_MaxOrder;
135 
136  return loPtr;
137  }
138 
139 
140  // virtual void GenerateData();
141  // void ThreadedGenerateData(const typename TInputImage::RegionType& regionForThread,ThreadIdType threadId );
142  // void BeforeThreadedGenerateData();
143 
145 
146 
147 private:
148  MeshFromSHCoefficientsImageFilter(const Self&); //purposely not implemented
149  void operator=(const Self&);//purposely not implemented
150 
151 
152 };
153 
154 } // end namespace itk
155 
156 
157 #endif
#define PrintVar1(cond, var, os)
Definition: utlCoreMacro.h:447
MeshFromContinuousSphericalFunctionImageFilter< TInputImage, TOutputMesh > Superclass
int DimToRankSH(const int dimm)
Definition: utlDMRI.h:182
Superclass::OutputMeshCellArrayType OutputMeshCellArrayType
NDArray< T, 2 > FlipOrientations(const NDArray< T, 2 > &in, const std::vector< int > &flip)
#define ITK_OVERRIDE
Definition: utlITKMacro.h:46
std::shared_ptr< NDArray< T, 2 > > ComputeSHMatrix(const unsigned int rank, const NDArray< T, 2 > &grad, const int mode)
Definition: utl.h:171
#define M_PI
Definition: utlCoreMacro.h:57
#define utlGlobalException(cond, expout)
Definition: utlCoreMacro.h:372
SphereTessellatorType::BasicShapeType BasicShapeType
Superclass::OutputMeshPolyDataPointer OutputMeshPolyDataPointer
VectorType NormalizeUnitIntegral(const VectorType &x) const ITK_OVERRIDE
Superclass::OutputMeshPolyDataPointer OutputMeshPolyDataPointer
Compute mesh from continuous spherical function linearly represented by a set of bases.
#define utlShowPosition(cond)
Definition: utlCoreMacro.h:554
void PrintSelf(std::ostream &os, Indent indent) const ITK_OVERRIDE
SmartPointer< Self > Pointer
Tesselates a sphere via subdivision of tetrahedron, octahedron, or icosahedron.
LightObject::Pointer InternalClone() const ITK_OVERRIDE