DMRITool  v0.1.1-139-g860d86b4
Diffusion MRI Tool
mexReadDWIList.cxx
Go to the documentation of this file.
1 
18 #include "utl.h"
19 #include "utlMEX.h"
20 
21 #include "itkDWIReader.h"
22 
23 template <typename T>
24  inline void callFunction(mxArray* plhs[], const mxArray* prhs[],
25  const int nlhs,const int nrhs)
26 {
27  utlGlobalException(mxGetClassID(prhs[0])!=mxCHAR_CLASS, "the first input has to be a string");
28  utlGlobalException(nrhs==2 && mxGetClassID(prhs[1])!=mxSTRUCT_CLASS, "the second input has to be a struct");
29 
30  std::string dwiFile;
31  utl::GetString(prhs[0], dwiFile);
32 
33 
34  bool normalize = utl::GetScalarStructDef<bool>(nrhs==2?prhs[1]:NULL,"normalize",true);
35  bool correctDWI = utl::GetScalarStructDef<bool>(nrhs==2?prhs[1]:NULL,"correctDWI",true);
36  bool warn = utl::GetScalarStructDef<bool>(nrhs==2?prhs[1]:NULL,"warn",true);
37  double bThreshold = utl::GetScalarStructDef<bool>(nrhs==2?prhs[1]:NULL,"bThreshold",-1);
38 
39  // utlPrintVar3(true, normalize, correctDWI, bThreshold);
40  typedef itk::DWIReader<T,3> DWIReaderType;
41  typename DWIReaderType::Pointer dwiReader = DWIReaderType::New();
42 
43  dwiReader->SetConfigurationFile(dwiFile);
44  dwiReader->SetNormalizeDWI(normalize);
45  if (bThreshold>0)
46  {
47  dwiReader->GetSamplingSchemeQSpace()->SetBThresholdSingleShell(bThreshold);
48  }
49  dwiReader->SetCorrectDWIValues(correctDWI);
50  dwiReader->SetShowWarnings(warn);
51 
52  typedef itk::Image<T,3> Image3DType;
53  typename Image3DType::Pointer b0Image = Image3DType::New();
54  mxArray* b0Array = utl::GetArrayStruct(prhs[1], "b0Image" );
55  if (b0Array)
56  {
57  itk::GetITKImageFromMXArray(b0Array, b0Image);
58  // itk::PrintImage(b0Image, "b0");
59  dwiReader->SetB0Image(b0Image);
60  }
61 
62  dwiReader->Update();
63 
64  typedef itk::VectorImage<T,3> VectorImageType;
65  typename VectorImageType::Pointer dwiImage = dwiReader->GetOutput();
66  // itk::PrintVectorImage(dwiImage, "dwiImage");
67  itk::GetMXArrayFromITKVectorImage(dwiImage, plhs[0]);
68 
69  utl_shared_ptr<std::vector<double> > bVec = dwiReader->GetSamplingSchemeQSpace()->GetBVector();
70  utl::GetMXArrayFromSTDVector<double>(bVec.get(), plhs[1]);
71 
72  utl_shared_ptr<utl::NDArray<double,2> > grad = dwiReader->GetSamplingSchemeQSpace()->GetOrientationsCartesian();
73  utl::GetMXArrayFromUtlMatrix<double>(grad.get(), plhs[2]);
74 
75  if (nlhs>3)
76  {
77  typename Image3DType::Pointer b0Image = dwiReader->GetB0Image();
78  // itk::PrintImage(b0Image, "b02");
79  itk::GetMXArrayFromITKImage(b0Image, plhs[3]);
80  }
81 }
82 
83 void mexFunction(int nlhs, mxArray *plhs[],
84  int nrhs, const mxArray *prhs[])
85 {
86  utl_shared_ptr<std::vector<double> > grad(new std::vector<double>());
87  utlGlobalException(nrhs!=1 && nrhs!=2, "Bad number of inputs arguments");
88  utlGlobalException(nlhs!=3 && nlhs!=4, "Bad number of outputs arguments");
89 
90  // utlPrintVar1(true, mxGetClassID(prhs[0]));
91  // utlPrintVar2(true, mxDOUBLE_CLASS, mxSINGLE_CLASS);
92 
93  // if (mxGetClassID(prhs[0]) == mxSINGLE_CLASS)
94  // callFunction<float>(plhs,prhs,nlhs,nrhs);
95  // else
96  callFunction<double>(plhs,prhs,nlhs,nrhs);
97 
98 }
mxArray * GetArrayStruct(const mxArray *pr_struct, const char *name)
Definition: mexutils.h:333
helper functions specifically used in dmritool
void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
void GetMXArrayFromITKVectorImage(const SmartPointer< VectorImage< T, VImageDimension > > &image, mxArray *&pr)
Definition: mexITK.h:96
void callFunction(mxArray *plhs[], const mxArray *prhs[], const int nlhs, const int nrhs)
#define utlGlobalException(cond, expout)
Definition: utlCoreMacro.h:372
Load gradient file, b values and DWI files (with optional index file)
Definition: itkDWIReader.h:59
void GetString(const mxArray *pr, std::string &str)
Definition: mexutils.h:275
void GetMXArrayFromITKImage(const SmartPointer< Image< T, VImageDimension > > &image, mxArray *&pr)
Definition: mexITK.h:32
itk::VectorImage< ScalarType, 3 > VectorImageType
Definition: 4DImageMath.cxx:30
void GetITKImageFromMXArray(const mxArray *pr, SmartPointer< Image< T, VImageDimension > > &image)
Definition: mexITK.h:55