Slicer 5.8
Slicer is a multi-platform, free and open source software package for visualization and medical image computing
Loading...
Searching...
No Matches
vtkDiffusionTensorMathematics.h
Go to the documentation of this file.
1/*=auto=========================================================================
2
3 Portions (c) Copyright 2005 Brigham and Women's Hospital (BWH) All Rights Reserved.
4
5 See COPYRIGHT.txt
6 or http://www.slicer.org/copyright/copyright.txt for details.
7
8 Program: 3D Slicer
9 Module: $RCSfile: vtkDiffusionTensorMathematics.h,v $
10 Date: $Date: 2006/12/19 17:14:44 $
11 Version: $Revision: 1.20 $
12
13=========================================================================auto=*/
18//
23//
24
25
26#ifndef __vtkDiffusionTensorMathematics_h
27#define __vtkDiffusionTensorMathematics_h
28
29// vtkTeem includes
30#include "vtkTeemConfigure.h"
31
32// VTK includes
33#include <vtkThreadedImageAlgorithm.h>
34
35class vtkMatrix4x4;
36class vtkImageData;
37class VTK_Teem_EXPORT vtkDiffusionTensorMathematics : public vtkThreadedImageAlgorithm
38{
39public:
41 vtkTypeMacro(vtkDiffusionTensorMathematics,vtkThreadedImageAlgorithm);
42 void PrintSelf(ostream& os, vtkIndent indent) override;
43
45 enum
46 {
47 VTK_TENS_TRACE = 0,
48 VTK_TENS_DETERMINANT = 1,
49 VTK_TENS_RELATIVE_ANISOTROPY = 2,
50 VTK_TENS_FRACTIONAL_ANISOTROPY = 3,
51 VTK_TENS_MAX_EIGENVALUE = 4,
52 VTK_TENS_MID_EIGENVALUE = 5,
53 VTK_TENS_MIN_EIGENVALUE = 6,
54 VTK_TENS_LINEAR_MEASURE = 7,
55 VTK_TENS_PLANAR_MEASURE = 8,
56 VTK_TENS_SPHERICAL_MEASURE = 9,
57 VTK_TENS_COLOR_ORIENTATION = 10,
58 VTK_TENS_D11 = 11,
59 VTK_TENS_D22 = 12,
60 VTK_TENS_D33 = 13,
61 VTK_TENS_MODE = 14,
62 VTK_TENS_COLOR_MODE =15,
63 VTK_TENS_MAX_EIGENVALUE_PROJX = 16,
64 VTK_TENS_MAX_EIGENVALUE_PROJY = 17,
65 VTK_TENS_MAX_EIGENVALUE_PROJZ = 18,
66 VTK_TENS_RAI_MAX_EIGENVEC_PROJX = 19,
67 VTK_TENS_RAI_MAX_EIGENVEC_PROJY = 20,
68 VTK_TENS_RAI_MAX_EIGENVEC_PROJZ = 21,
69 VTK_TENS_MAX_EIGENVEC_PROJX = 22,
70 VTK_TENS_MAX_EIGENVEC_PROJY = 23,
71 VTK_TENS_MAX_EIGENVEC_PROJZ = 24,
72 VTK_TENS_PARALLEL_DIFFUSIVITY = 25,
73 VTK_TENS_PERPENDICULAR_DIFFUSIVITY = 26,
74 VTK_TENS_COLOR_ORIENTATION_MIDDLE_EIGENVECTOR = 27,
75 VTK_TENS_COLOR_ORIENTATION_MIN_EIGENVECTOR = 28,
76 VTK_TENS_MEAN_DIFFUSIVITY = 29
77 };
80 vtkGetMacro(Operation,int);
81 vtkSetClampMacro(Operation,int, VTK_TENS_TRACE, VTK_TENS_MEAN_DIFFUSIVITY);
82
86 {this->SetOperation(VTK_TENS_TRACE);};
87
91 {this->SetOperation(VTK_TENS_DETERMINANT);};
92
96 {this->SetOperation(VTK_TENS_RELATIVE_ANISOTROPY);};
98 {this->SetOperation(VTK_TENS_FRACTIONAL_ANISOTROPY);};
100 {this->SetOperation(VTK_TENS_LINEAR_MEASURE);};
102 {this->SetOperation(VTK_TENS_PLANAR_MEASURE);};
104 {this->SetOperation(VTK_TENS_SPHERICAL_MEASURE);};
108 {this->SetOperation(VTK_TENS_MODE);};
110 {this->SetOperation(VTK_TENS_PARALLEL_DIFFUSIVITY);};
112 {this->SetOperation(VTK_TENS_PERPENDICULAR_DIFFUSIVITY);};
114 {this->SetOperation(VTK_TENS_MEAN_DIFFUSIVITY);};
115
119 {this->SetOperation(VTK_TENS_MAX_EIGENVALUE);};
121 {this->SetOperation(VTK_TENS_MID_EIGENVALUE);};
123 {this->SetOperation(VTK_TENS_MIN_EIGENVALUE);};
124
128 {this->SetOperation(VTK_TENS_MAX_EIGENVALUE_PROJX);};
130 {this->SetOperation(VTK_TENS_MAX_EIGENVALUE_PROJY);};
132 {this->SetOperation(VTK_TENS_MAX_EIGENVALUE_PROJZ);};
133
137 {this->SetOperation(VTK_TENS_RAI_MAX_EIGENVEC_PROJX);}
139 {this->SetOperation(VTK_TENS_RAI_MAX_EIGENVEC_PROJY);}
141 {this->SetOperation(VTK_TENS_RAI_MAX_EIGENVEC_PROJZ);}
142
146 {this->SetOperation(VTK_TENS_MAX_EIGENVEC_PROJX);}
148 {this->SetOperation(VTK_TENS_MAX_EIGENVEC_PROJY);}
150 {this->SetOperation(VTK_TENS_MAX_EIGENVEC_PROJZ);}
151
152
156 {this->SetOperation(VTK_TENS_D11);};
158 {this->SetOperation(VTK_TENS_D22);};
160 {this->SetOperation(VTK_TENS_D33);};
161
167 {this->SetOperation(VTK_TENS_COLOR_ORIENTATION);};
168
175 {this->SetOperation(VTK_TENS_COLOR_MODE);};
176
180 vtkSetMacro(ScaleFactor,double);
181 vtkGetMacro(ScaleFactor,double);
182
185 vtkSetMacro(ExtractEigenvalues,int);
186 vtkBooleanMacro(ExtractEigenvalues,int);
187 vtkGetMacro(ExtractEigenvalues,int);
188
193 //
196 //
207 //
208 virtual void SetTensorRotationMatrix(vtkMatrix4x4*);
209 vtkGetObjectMacro(TensorRotationMatrix, vtkMatrix4x4);
210
214 vtkBooleanMacro(MaskWithScalars, int);
215 vtkSetMacro(MaskWithScalars, int);
216 vtkGetMacro(MaskWithScalars, int);
217
218 vtkBooleanMacro(FixNegativeEigenvalues, int);
219 vtkSetMacro(FixNegativeEigenvalues, int);
220 vtkGetMacro(FixNegativeEigenvalues, int);
221
224 virtual void SetScalarMask(vtkImageData*);
225 vtkGetObjectMacro(ScalarMask, vtkImageData);
226
229 vtkSetMacro(MaskLabelValue, int);
230 vtkGetMacro(MaskLabelValue, int);
231
233 static void ModeToRGB(double Mode, double FA,
234 double &R, double &G, double &B);
235
236 static void RGBToIndex(double R, double G,
237 double B, double &index);
238
241 static int FixNegativeEigenvaluesMethod(double w[3]);
242 static double Determinant(double D[3][3]);
243 static double Trace(double D[3][3]);
244 static double Trace(double w[3]);
245 static double RelativeAnisotropy(double w[3]);
246 static double FractionalAnisotropy(double w[3]);
247 static double LinearMeasure(double w[3]);
248 static double PlanarMeasure(double w[3]);
249 static double SphericalMeasure(double w[3]);
250 static double MaxEigenvalue(double w[3]);
251 static double MiddleEigenvalue(double w[3]);
252 static double ParallelDiffusivity(double w[3]);
253 static double PerpendicularDiffusivity(double w[3]);
254 static double MeanDiffusivity(double w[3]);
255 static double MinEigenvalue(double w[3]);
256 static double RAIMaxEigenvecX(double **v, double w[3]);
257 static double RAIMaxEigenvecY(double **v, double w[3]);
258 static double RAIMaxEigenvecZ(double **v, double w[3]);
259 static double MaxEigenvecX(double **v, double w[3]);
260 static double MaxEigenvecY(double **v, double w[3]);
261 static double MaxEigenvecZ(double **v, double w[3]);
262 static double MaxEigenvalueProjectionX(double **v, double w[3]);
263 static double MaxEigenvalueProjectionY(double **v, double w[3]);
264 static double MaxEigenvalueProjectionZ(double **v, double w[3]);
265 static double Mode(double w[3]);
266 static void ColorByMode(double w[3], double &R,double &G, double &B);
267
268 //Description
269 //Wrap function to teem eigen solver
270 static int TeemEigenSolver(double **m, double *w, double **v);
271 void ComputeTensorIncrements(vtkImageData *imageData, vtkIdType incr[3]);
272
273protected:
276
278 double ScaleFactor;
280
282 vtkImageData *ScalarMask;
284
285 vtkMatrix4x4 *TensorRotationMatrix;
287
288 int RequestInformation (vtkInformation*,
289 vtkInformationVector**,
290 vtkInformationVector*) override;
291
292 void ThreadedRequestData(vtkInformation *request,
293 vtkInformationVector **inputVector,
294 vtkInformationVector *outputVector,
295 vtkImageData ***inData,
296 vtkImageData **outData,
297 int extent[6], int threadId) override;
298
299 int FillInputPortInformation(int port, vtkInformation* info) override;
300
301 // Reimplemented to delete the tensor array of the output.
302 int RequestData(vtkInformation* request,
303 vtkInformationVector** inputVector,
304 vtkInformationVector* outputVector) override;
305private:
307 void operator=(const vtkDiffusionTensorMathematics&) = delete;
308};
309
310#endif
int FillInputPortInformation(int port, vtkInformation *info) override
void PrintSelf(ostream &os, vtkIndent indent) override
static void ModeToRGB(double Mode, double FA, double &R, double &G, double &B)
Public for access from threads.
~vtkDiffusionTensorMathematics() override
static double ParallelDiffusivity(double w[3])
static void RGBToIndex(double R, double G, double B, double &index)
int RequestData(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector) override
static double FractionalAnisotropy(double w[3])
static double PlanarMeasure(double w[3])
static double MaxEigenvalueProjectionX(double **v, double w[3])
static double MaxEigenvalue(double w[3])
virtual void SetTensorRotationMatrix(vtkMatrix4x4 *)
int RequestInformation(vtkInformation *, vtkInformationVector **, vtkInformationVector *) override
void SetOperationToDeterminant()
Output the determinant.
static double Mode(double w[3])
static double RelativeAnisotropy(double w[3])
static double MinEigenvalue(double w[3])
static double MeanDiffusivity(double w[3])
void SetOperationToD11()
Output a matrix (tensor) component.
int MaskWithScalars
Boolean controls eigenfunction extraction.
static double RAIMaxEigenvecZ(double **v, double w[3])
static vtkDiffusionTensorMathematics * New()
static double Trace(double w[3])
static double MaxEigenvalueProjectionY(double **v, double w[3])
static double MiddleEigenvalue(double w[3])
void SetOperationToMaxEigenvalueProjectionX()
Output Maxeigenvalue*Maxeigenvec_projection also known as L1Z.
static double MaxEigenvecZ(double **v, double w[3])
static double MaxEigenvecY(double **v, double w[3])
static double RAIMaxEigenvecY(double **v, double w[3])
void SetOperationToMaxEigenvecX()
Output Relative_anisotropy*Maxeigenvec_projection also known as L1z.
static int TeemEigenSolver(double **m, double *w, double **v)
static int FixNegativeEigenvaluesMethod(double w[3])
Helper functions to perform operations pixel-wise.
void SetOperationToMaxEigenvalue()
Output a selected eigenvalue.
static double SphericalMeasure(double w[3])
double ScaleFactor
math operation to perform
int ExtractEigenvalues
Scale factor for output scalars.
void SetOperationToRelativeAnisotropy()
Output various anisotropy and shape measures.
void SetOperationToRAIMaxEigenvecX()
Output Relative_anisotropy*Maxeigenvec_projection also known as L1z.
static double Trace(double D[3][3])
static double PerpendicularDiffusivity(double w[3])
void ThreadedRequestData(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector, vtkImageData ***inData, vtkImageData **outData, int extent[6], int threadId) override
static double RAIMaxEigenvecX(double **v, double w[3])
static double LinearMeasure(double w[3])
virtual void SetScalarMask(vtkImageData *)
Scalar mask.
static void ColorByMode(double w[3], double &R, double &G, double &B)
static double MaxEigenvalueProjectionZ(double **v, double w[3])
static double MaxEigenvecX(double **v, double w[3])
void ComputeTensorIncrements(vtkImageData *imageData, vtkIdType incr[3])
void SetOperationToTrace()
Output the trace (sum of eigenvalues = sum along diagonal)
static double Determinant(double D[3][3])