ANIMA  4.0
animaDistortionCorrection.cxx
Go to the documentation of this file.
3 
4 #include <tclap/CmdLine.h>
5 
6 #include <itkImage.h>
7 
8 int main(int ac, const char** av)
9 {
10  std::string descriptionMessage;
11  descriptionMessage += "Compute a vector field in order to correct a distorted EPI\n";
12  descriptionMessage += "INRIA / IRISA - VisAGeS/Empenn Team";
13 
14  TCLAP::CmdLine cmd(descriptionMessage, ' ',ANIMA_VERSION);
15 
16  TCLAP::ValueArg<std::string> forwardArg("f", "forwardImage", "Forward image", true, "", "Forward image", cmd);
17  TCLAP::ValueArg<std::string> backwardArg("b", "backwardImage", "Backward image", true, "", "Backward image", cmd);
18  TCLAP::ValueArg<unsigned int> distortionDirection("d", "dir", "Direction of distortion (0,1,2)", false, 1, "number of the direction of distortion", cmd);
19  TCLAP::ValueArg<std::string> outArg("o", "outputVectorField", "Distortion vector field", true, "","Output vector field", cmd );
20 
21  TCLAP::ValueArg<double> sigmaArg("s", "sigma-smooth", "Sigma for gaussian smoothing of the transformation (in pixels, default: 2)",false,2,"gaussian smoothing sigma", cmd );
22  TCLAP::ValueArg<unsigned int> nbpArg("T","numberofthreads","Number of threads to run on (default : all cores)",false,itk::MultiThreaderBase::GetGlobalDefaultNumberOfThreads(),"number of threads",cmd);
23 
24  try
25  {
26  cmd.parse(ac,av);
27  }
28  catch (TCLAP::ArgException& e)
29  {
30  std::cerr << "Error: " << e.error() << "for argument " << e.argId() << std::endl;
31  return EXIT_FAILURE;
32  }
33 
34  const unsigned int Dimension = 3;
35 
36  typedef double PixelType;
37 
38  typedef itk::Image <PixelType, Dimension> ImageType;
39  typedef itk::VectorImage <PixelType, Dimension> VectorFieldType;
40 
42 
43  ImageType::Pointer backwardImage = anima::readImage<ImageType>(backwardArg.getValue());
44  backwardImage->DisconnectPipeline();
45 
46  double meanSpacing = 0;
47  for (unsigned int i = 0;i < Dimension;++i)
48  meanSpacing += backwardImage->GetSpacing()[i];
49 
50  meanSpacing /= Dimension;
51 
52  FilterType::Pointer filter = FilterType::New();
53  filter->SetInput(0, anima::readImage<ImageType>(forwardArg.getValue()));
54  filter->SetInput(1, backwardImage);
55  filter->SetDirection(distortionDirection.getValue());
56  filter->SetFieldSmoothingSigma(sigmaArg.getValue() * meanSpacing);
57  filter->SetNumberOfWorkUnits(nbpArg.getValue());
58  filter->Update();
59 
60  anima::writeImage<VectorFieldType>(outArg.getValue(),filter->GetOutput());
61 
62  return EXIT_SUCCESS;
63 }
int main(int ac, const char **av)