4 #include <tclap/CmdLine.h> 5 #include <itkTimeProbe.h> 12 void eventCallback (itk::Object* caller,
const itk::EventObject& event,
void* clientData)
14 itk::ProcessObject * processObject =
static_cast<itk::ProcessObject *
> (caller);
15 std::cout<<
"\033[K\rProgression: " <<
static_cast<int>(processObject->GetProgress() * 100) <<
"%" << std::flush;
18 int main(
int argc,
char **argv)
20 TCLAP::CmdLine cmd(
"INRIA / IRISA - VisAGeS/Empenn Team",
' ',ANIMA_VERSION);
22 TCLAP::ValueArg<std::string> inArg(
"i",
"inputdwi",
"dwi_volume",
true,
"",
"DWI volume",cmd);
23 TCLAP::ValueArg<std::string> resArg(
"o",
"output",
"dit_volume",
true,
"",
"result DTI image",cmd);
24 TCLAP::ValueArg<std::string> b0OutArg(
"O",
"output-b0",
"output_b0",
false,
"",
"result B0 image",cmd);
25 TCLAP::ValueArg<std::string> varOutArg(
"N",
"output-variance",
"output_variance",
false,
"",
"result noise variance image",cmd);
27 TCLAP::ValueArg<std::string> gradsArg(
"g",
"grad",
"input_gradients",
true,
"",
"Input gradients",cmd);
28 TCLAP::ValueArg<std::string> bvalArg(
"b",
"bval",
"input_b-values",
true,
"",
"Input b-values",cmd);
29 TCLAP::SwitchArg bvalueScaleArg(
"B",
"b-no-scale",
"Do not scale b-values according to gradient norm",cmd);
30 TCLAP::ValueArg<std::string> computationMaskArg(
"m",
"mask",
"Computation mask",
false,
"",
"computation mask",cmd);
32 TCLAP::ValueArg<unsigned int> b0ThrArg(
"t",
"b0thr",
"bot_treshold",
false,0,
"B0 threshold (default : 0)",cmd);
33 TCLAP::ValueArg<unsigned int> nbpArg(
"p",
"numberofthreads",
"nb_thread",
false,itk::MultiThreaderBase::GetGlobalDefaultNumberOfThreads(),
"Number of threads to run on (default: all cores)",cmd);
34 TCLAP::ValueArg<std::string> reorientArg(
"r",
"reorient",
"dwi_reoriented",
false,
"",
"Reorient DWI given as input",cmd);
35 TCLAP::ValueArg<std::string> reorientGradArg(
"R",
"reorient-G",
"gradient reoriented output",
false,
"",
"Reorient gradients so that they are in MrTrix format (in image coordinates)",cmd);
41 catch (TCLAP::ArgException& e)
43 std::cerr <<
"Error: " << e.error() <<
"for argument " << e.argId() << std::endl;
48 typedef FilterType::MaskImageType MaskImageType;
49 typedef FilterType::InputImageType InputImageType;
50 typedef FilterType::Image4DType Image4DType;
51 typedef FilterType::OutputImageType VectorImageType;
53 itk::CStyleCommand::Pointer callback = itk::CStyleCommand::New();
56 FilterType::Pointer mainFilter = FilterType::New();
59 GFReaderType gfReader;
61 gfReader.SetBValueBaseString(bvalArg.getValue());
62 gfReader.SetGradientIndependentNormalization(bvalueScaleArg.isSet());
66 GFReaderType::GradientVectorType directions = gfReader.GetGradients();
67 GFReaderType::BValueVectorType mb = gfReader.GetBValues();
69 std::string inputFile = inArg.getValue();
71 Image4DType::Pointer input = anima::readImage<Image4DType>(inArg.getValue());
73 if (reorientArg.getValue() !=
"")
75 input = anima::reorientImage<Image4DType>(input, itk::SpatialOrientation::ITK_COORDINATE_ORIENTATION_RAI);
76 anima::writeImage<Image4DType>(reorientArg.getValue(), input);
78 inputFile = reorientArg.getValue();
81 if(reorientGradArg.getValue() !=
"")
83 anima::reorientGradients <Image4DType, vnl_vector_fixed<double,3> > (input, directions);
85 std::ofstream outGrads(reorientGradArg.getValue().c_str());
86 outGrads.precision(15);
87 for (
unsigned int i = 0;i < 3;++i)
89 for (
unsigned int j = 0;j < directions.size();++j)
90 outGrads << directions[j][i] <<
" ";
91 outGrads << std::endl;
97 mainFilter->SetBValuesList(mb);
98 for(
unsigned int i = 0;i < directions.size();++i)
99 mainFilter->AddGradientDirection(i, directions[i]);
101 std::vector <InputImageType::Pointer> inputData;
102 inputData = anima::getImagesFromHigherDimensionImage<Image4DType,InputImageType>(input);
104 for (
unsigned int i = 0;i < inputData.size();++i)
105 mainFilter->SetInput(i,inputData[i]);
107 if (computationMaskArg.getValue() !=
"")
108 mainFilter->SetComputationMask(anima::readImage<MaskImageType>(computationMaskArg.getValue()));
110 mainFilter->SetB0Threshold(b0ThrArg.getValue());
111 mainFilter->SetNumberOfWorkUnits(nbpArg.getValue());
112 mainFilter->AddObserver(itk::ProgressEvent(), callback);
114 itk::TimeProbe tmpTimer;
120 mainFilter->Update();
121 }
catch (itk::ExceptionObject &e)
123 std::cerr << e << std::endl;
129 std::cout <<
"\nEstimation done in " << tmpTimer.GetTotal() <<
" s" << std::endl;
130 std::cout <<
"Writing result to : " << resArg.getValue() << std::endl;
132 VectorImageType::Pointer output = mainFilter->GetOutput();
133 output->DisconnectPipeline();
135 anima::writeImage <VectorImageType> (resArg.getValue(),output);
137 if (b0OutArg.getValue() !=
"")
138 anima::writeImage <InputImageType> (b0OutArg.getValue(),mainFilter->GetEstimatedB0Image());
140 if (varOutArg.getValue() !=
"")
141 anima::writeImage <InputImageType> (varOutArg.getValue(),mainFilter->GetEstimatedVarianceImage());
int main(int argc, char **argv)
void SetGradientFileName(std::string fName)
void eventCallback(itk::Object *caller, const itk::EventObject &event, void *clientData)