2 #include <itkTimeProbe.h> 6 #include <itkCommand.h> 8 #include <tclap/CmdLine.h> 13 void eventCallback (itk::Object* caller,
const itk::EventObject& event,
void* clientData)
15 itk::ProcessObject * processObject = (itk::ProcessObject*) caller;
16 std::cout<<
"\033[K\rProgression: "<<(int)(processObject->GetProgress() * 100)<<
"%"<<std::flush;
19 int main(
int argc,
char* argv[])
21 TCLAP::CmdLine cmd(
"INRIA / IRISA - VisAGeS/Empenn Team",
' ',ANIMA_VERSION);
23 TCLAP::ValueArg<std::string> dtiArg(
"i",
"dti",
"Input DT image",
true,
"",
"DT image",cmd);
24 TCLAP::ValueArg<std::string> seedMaskArg(
"s",
"seed-mask",
"Seed mask",
true,
"",
"seed",cmd);
25 TCLAP::ValueArg<std::string> fibersArg(
"o",
"fibers",
"Output fibers",
true,
"",
"fibers",cmd);
27 TCLAP::ValueArg<std::string> filterMaskArg(
"",
"filter-mask",
"Mask for filtering fibers (default: none)",
false,
"",
"filter mask",cmd);
28 TCLAP::ValueArg<std::string> cutMaskArg(
"c",
"cut-mask",
"Mask for cutting fibers (default: none)",
false,
"",
"cut mask",cmd);
29 TCLAP::ValueArg<std::string> forbiddenMaskArg(
"f",
"forbidden-mask",
"Mask for removing fibers (default: none)",
false,
"",
"remove mask",cmd);
31 TCLAP::ValueArg<double> adcThrArg(
"A",
"adc-thr",
"ADC threshold (default: 2.0e-3)",
false,2.0e-3,
"ADC threshold",cmd);
32 TCLAP::ValueArg<double> faThrArg(
"",
"fa-thr",
"FA threshold (default: 0.1)",
false,0.1,
"FA threshold",cmd);
33 TCLAP::ValueArg<double> punctureWeightArg(
"p",
"punc-weight",
"Puncture weight to go through planar tensors (default: 0.2)",
false,0.2,
"puncture weight",cmd);
35 TCLAP::ValueArg<double> stopAngleArg(
"a",
"angle-max",
"Maximum angle for tracking (default: 60)",
false,60.0,
"maximum track angle",cmd);
36 TCLAP::ValueArg<double> stepLengthArg(
"",
"step-length",
"Length of each step (default: 1)",
false,1.0,
"step length",cmd);
37 TCLAP::ValueArg<int> nbFibersArg(
"",
"nb-fibers",
"Number of starting filters (n*n*n) per voxel (default: 1)",
false,1,
"number of seeds per voxel",cmd);
39 TCLAP::ValueArg<double> minLengthArg(
"",
"min-length",
"Minimum length for a fiber to be considered for computation (default: 10mm)",
false,10.0,
"minimum length",cmd);
40 TCLAP::ValueArg<double> maxLengthArg(
"",
"max-length",
"Maximum length of a tract (default: 150mm)",
false,150.0,
"maximum length",cmd);
42 TCLAP::SwitchArg addLocalDataArg(
"L",
"local-data",
"Add local data information to output tracks",cmd);
44 TCLAP::ValueArg<unsigned int> nbThreadsArg(
"T",
"nb-threads",
"Number of threads to run on (default: all available)",
false,itk::MultiThreaderBase::GetGlobalDefaultNumberOfThreads(),
"number of threads",cmd);
50 catch (TCLAP::ArgException& e)
52 std::cerr <<
"Error: " << e.error() <<
"for argument " << e.argId() << std::endl;
57 typedef MainFilterType::ModelImageType DTIImageType;
58 typedef MainFilterType::MaskImageType MaskImageType;
60 MainFilterType::Pointer dtiTracker = MainFilterType::New();
63 LogTensorFilterType::Pointer tensorLogger = LogTensorFilterType::New();
65 tensorLogger->SetInput(anima::readImage <DTIImageType> (dtiArg.getValue()));
66 tensorLogger->SetScaleNonDiagonal(
false);
67 tensorLogger->SetNumberOfWorkUnits(nbThreadsArg.getValue());
69 tensorLogger->Update();
71 dtiTracker->SetInputImage(tensorLogger->GetOutput());
73 dtiTracker->SetSeedingMask(anima::readImage <MaskImageType> (seedMaskArg.getValue()));
75 if (filterMaskArg.getValue() !=
"")
76 dtiTracker->SetFilteringMask(anima::readImage <MaskImageType> (filterMaskArg.getValue()));
78 if (cutMaskArg.getValue() !=
"")
79 dtiTracker->SetCutMask(anima::readImage <MaskImageType> (cutMaskArg.getValue()));
81 if (forbiddenMaskArg.getValue() !=
"")
82 dtiTracker->SetForbiddenMask(anima::readImage <MaskImageType> (forbiddenMaskArg.getValue()));
84 dtiTracker->SetNumberOfWorkUnits(nbThreadsArg.getValue());
85 dtiTracker->SetNumberOfFibersPerPixel(nbFibersArg.getValue());
86 dtiTracker->SetStepProgression(stepLengthArg.getValue());
87 dtiTracker->SetStopFAThreshold(faThrArg.getValue());
88 dtiTracker->SetStopADCThreshold(adcThrArg.getValue());
89 dtiTracker->SetPunctureWeight(punctureWeightArg.getValue());
90 dtiTracker->SetMaxFiberAngle(stopAngleArg.getValue());
91 dtiTracker->SetMinLengthFiber(minLengthArg.getValue());
92 dtiTracker->SetMaxLengthFiber(maxLengthArg.getValue());
94 bool computeColors = (fibersArg.getValue().find(
".fds") != std::string::npos) && (addLocalDataArg.isSet());
95 dtiTracker->SetComputeLocalColors(computeColors);
97 itk::CStyleCommand::Pointer callback = itk::CStyleCommand::New();
99 dtiTracker->AddObserver(itk::ProgressEvent(), callback);
101 itk::TimeProbe tmpTime;
104 dtiTracker->Update();
105 dtiTracker->RemoveAllObservers();
108 std::cout <<
"Tracking time: " << tmpTime.GetTotal() <<
"s" << std::endl;
void SetInputData(vtkPolyData *data)
int main(int argc, char *argv[])
DTI tractography image filter. Simple step by step tratpography, using advection-diffusion tricks fro...
void SetFileName(std::string &name)
void eventCallback(itk::Object *caller, const itk::EventObject &event, void *clientData)