DMRITool  v0.1.1-139-g860d86b4
Diffusion MRI Tool
itkMeshFromSphericalFunctionImageFilter.h
Go to the documentation of this file.
1 
18 #ifndef __itkMeshFromSphericalFunctionImageFilter_h
19 #define __itkMeshFromSphericalFunctionImageFilter_h
20 
21 #include "itkSphereTessellator.h"
23 
24 namespace itk
25 {
26 
34 template < class TInputImage, class TOutputMesh=vtkPolyData>
36  : public MeshFromImageImageFilter<TInputImage, TOutputMesh>
37 {
38 public:
42  typedef SmartPointer<Self> Pointer;
43  typedef SmartPointer<const Self> ConstPointer;
44 
46  itkNewMacro(Self);
47 
50 
51  //** Convenient typedefs */
52  typedef TInputImage InputImageType;
53  typedef typename InputImageType::Pointer InputImagePointer;
54  typedef typename InputImageType::ConstPointer InputImageConstPointer;
55  typedef typename InputImageType::IndexType InputImageIndexType;
56  typedef typename InputImageType::SizeType InputImageSizeType;
57  typedef typename InputImageType::SizeValueType InputImageSizeValueType;
58  typedef typename InputImageType::SpacingType InputImageSpacingType;
59  typedef typename InputImageType::PixelType InputImagePixelType;
60  typedef typename InputImageType::RegionType InputImageRegionType;
61  typedef typename InputImageType::PointType InputImagePointType;
62 
63  typedef TOutputMesh OutputMeshPolyDataType;
64 
65  typedef typename Superclass::OutputMeshPointsType OutputMeshPointsType;
66  typedef typename Superclass::OutputMeshCellArrayType OutputMeshCellArrayType;
67  typedef typename Superclass::OutputMeshPolyDataPointer OutputMeshPolyDataPointer;
68  typedef typename Superclass::OutputMeshScalarType OutputMeshScalarType;
69  typedef typename Superclass::OutputMeshRGBType OutputMeshRGBType;
70 
71 
76 
78  typedef typename Superclass::MatrixType MatrixType;
80  typedef typename Superclass::VectorType VectorType;
81  typedef typename Superclass::VectorPointer VectorPointer;
82  typedef typename Superclass::STDVectorType STDVectorType;
83  typedef typename Superclass::STDVectorPointer STDVectorPointer;
84 
86  typedef enum {NONE=0, MIN_MAX, UNIT_MAX, UNIT_INTEGRAL} NormalizationType;
88  itkSetMacro(Normalization, NormalizationType);
89  itkGetMacro(Normalization, NormalizationType);
90 
91  itkSetMacro(Orientations, MatrixPointer);
92  itkGetMacro(Orientations, MatrixPointer);
93 
95  itkSetMacro(Pow, double);
96  itkGetMacro(Pow, double);
97 
99  itkSetMacro(RemoveNegativeValues, bool);
100  itkGetMacro(RemoveNegativeValues, bool);
101  itkBooleanMacro(RemoveNegativeValues);
102 
103  unsigned int GetNumberOfOrientations() const
104  {
105  return m_Orientations->Rows();
106  }
107 
108 
109 protected:
111  m_Orientations(new MatrixType())
112  {
113  m_Pow = 1.0;
114  m_RemoveNegativeValues = false;
115  m_Normalization = NONE;
116  this->m_RemoveNegativeValues = false;
117 
118  this->m_ColorScheme = Superclass::DIRECTION;
119  }
120 
122  {};
123 
124  void PrintSelf(std::ostream& os, Indent indent) const ITK_OVERRIDE
125  {
126  PrintVar(true, os<<indent, m_Pow, m_RemoveNegativeValues);
127  PrintEnum4(true, m_Normalization, NONE, MIN_MAX, UNIT_MAX, UNIT_INTEGRAL, os<<indent);
128  utl::PrintUtlMatrix(*m_Orientations, "m_Orientations", " ", os<<indent);
129  }
130 
131  typename LightObject::Pointer InternalClone() const ITK_OVERRIDE
132  {
133  typename LightObject::Pointer loPtr = Superclass::InternalClone();
134 
135  typename Self::Pointer rval = dynamic_cast<Self *>(loPtr.GetPointer());
136  if(rval.IsNull())
137  {
138  itkExceptionMacro(<< "downcast to type " << this->GetNameOfClass()<< " failed.");
139  }
140 
141  rval->m_Pow = m_Pow;
142  rval->m_RemoveNegativeValues = m_RemoveNegativeValues;
143  rval->m_Orientations = m_Orientations;
144  rval->m_Normalization = m_Normalization;
145 
146  return loPtr;
147  }
148 
149  // virtual void GenerateData();
150  // void ThreadedGenerateData(const typename TInputImage::RegionType& regionForThread,ThreadIdType threadId );
151  // void BeforeThreadedGenerateData();
152 
154 
155  double m_Pow;
156 
157  MatrixPointer m_Orientations;
158 
159  NormalizationType m_Normalization;
160 
161 private:
162  MeshFromSphericalFunctionImageFilter(const Self&); //purposely not implemented
163  void operator=(const Self&);//purposely not implemented
164 
165 };
166 
167 } // end namespace itk
168 
169 
170 #endif
utl_shared_ptr< MatrixType > MatrixPointer
auto Pow(const TLeft &lhs, const TRight &rhs) -> decltype(utl::F< utl::Functor::Pow< Expr2ValueType< TLeft, TRight >> >(lhs, rhs))
Definition: utlFunctors.h:385
void PrintUtlMatrix(const NDArray< T, 2 > &mat, const std::string &str="", const char *separate=" ", std::ostream &os=std::cout)
#define PrintEnum4(cond, var, val1, val2, val3, val4, os)
Definition: utlCoreMacro.h:423
MeshFromImageImageFilter< TInputImage, TOutputMesh > Superclass
void PrintSelf(std::ostream &os, Indent indent) const ITK_OVERRIDE
#define ITK_OVERRIDE
Definition: utlITKMacro.h:46
Compute mesh from SH coefficients.
#define PrintVar(cond, os,...)
Definition: utlCoreMacro.h:242
Superclass::OutputMeshPolyDataPointer OutputMeshPolyDataPointer
Compute mesh from general spherical function.
SmartPointer< Self > Pointer
Tesselates a sphere via subdivision of tetrahedron, octahedron, or icosahedron.