ANIMA  4.0
animaConvertImage.cxx
Go to the documentation of this file.
1 #include <tclap/CmdLine.h>
2 
3 #include <itkImage.h>
4 #include <itkOrientImageFilter.h>
5 
6 #include <animaReorientation.h>
9 
10 struct arguments
11 {
12  std::string input, output, reorient, space;
13  std::string gradientFileName;
15 };
16 
17 template <class ImageType>
18 void
19 convert(const arguments &args)
20 {
21  typename ImageType::Pointer input = anima::readImage<ImageType>(args.input);
22 
23  if(args.space != "")
24  {
25  typename ImageType::Pointer spaceImage = anima::readImage<ImageType>(args.space);
26 
27  input->CopyInformation(spaceImage);
28  }
29 
30  std::vector < vnl_vector_fixed <double,3> > gradients;
31 
32  if (args.gradientFileName != "")
33  {
34  typedef anima::GradientFileReader < vnl_vector_fixed<double,3>, double > GFReaderType;
35  GFReaderType gfReader;
37  gfReader.SetGradientIndependentNormalization(false);
38 
39  gfReader.Update();
40 
41  gradients = gfReader.GetGradients();
42  }
43 
44  if (args.output != "")
45  {
46  if(args.reorient == "AXIAL")
47  input = anima::reorientImage(input,itk::SpatialOrientation::ITK_COORDINATE_ORIENTATION_RAI);
48  else if(args.reorient == "CORONAL")
49  input = anima::reorientImage(input,itk::SpatialOrientation::ITK_COORDINATE_ORIENTATION_RSA);
50  else if(args.reorient == "SAGITTAL")
51  input = anima::reorientImage(input,itk::SpatialOrientation::ITK_COORDINATE_ORIENTATION_ASL);
52  }
53 
54  if (args.gradientFileName != "")
55  {
56  anima::reorientGradients(input,gradients);
57 
58  std::string outGradFileName = args.input;
59 
60  if (args.output != "")
61  outGradFileName = args.output;
62 
63  std::size_t pointLocation = outGradFileName.find_last_of('.');
64  outGradFileName = outGradFileName.substr(0,pointLocation);
65  std::string tmpStr;
66 
67  if (args.output != "")
68  tmpStr = args.output.substr(pointLocation + 1);
69  else
70  tmpStr = args.input.substr(pointLocation + 1);
71 
72  if (tmpStr == "gz")
73  {
74  pointLocation = outGradFileName.find_last_of('.');
75  outGradFileName = outGradFileName.substr(0,pointLocation);
76  }
77 
78  outGradFileName += ".bvec";
79  std::ofstream outGrads(outGradFileName);
80  outGrads.precision(15);
81  for (unsigned int i = 0;i < 3;++i)
82  {
83  for (unsigned int j = 0;j < gradients.size();++j)
84  outGrads << gradients[j][i] << " ";
85  outGrads << std::endl;
86  }
87  outGrads.close();
88  }
89 
90  if (args.output != "")
91  {
92  input->SetMetaDataDictionary(itk::MetaDataDictionary());
93  anima::writeImage<ImageType>(args.output, input);
94  }
95 }
96 
97 template <class ComponentType, int dimension>
98 void
99 checkIfComponentsAreVectors(itk::ImageIOBase::Pointer imageIO, const arguments &args)
100 {
101  ANIMA_CHECK_IF_COMPONENTS_ARE_VECTORS(imageIO, ComponentType, dimension, convert, args)
102 }
103 
104 template <class ComponentType>
105 void
106 retrieveNbDimensions(itk::ImageIOBase::Pointer imageIO, const arguments &args)
107 {
108  ANIMA_RETRIEVE_NUMBER_OF_DIMENSIONS(imageIO, ComponentType, checkIfComponentsAreVectors, imageIO, args);
109 }
110 
111 int main(int ac, const char** av)
112 {
113  TCLAP::CmdLine cmd("animaConvertImage can be used to rewrite an image in a new compatible format ie: nii to nrrd. "
114  "It can also be used to display generic information on an image such as size, origin etc. using the --info option. "
115  "If an image file is given to the --space option it will be used to rewrite the input in the same "
116  "real coordinates repair. "
117  "The --reorient option allow to reorient the image in either the AXIAL, CORONAL or SAGITALL orientation. "
118  "Note that the reorientation is performed after the change of coordinate space if --space and --reorient "
119  "are given together.\n"
120  "INRIA / IRISA - VisAGeS/Empenn Team",' ',ANIMA_VERSION);
121 
122  TCLAP::ValueArg<std::string> inputArg("i",
123  "input",
124  "input_filename",
125  true,
126  "",
127  "Input image.",
128  cmd);
129 
130  TCLAP::ValueArg<std::string> outputArg("o",
131  "output",
132  "output_filename",
133  false,
134  "",
135  "Output image.",
136  cmd);
137 
138  TCLAP::SwitchArg infoArg("I",
139  "info",
140  "Display info on the image or not",
141  cmd);
142 
143  TCLAP::ValueArg<std::string> reorientArg("R",
144  "reorient",
145  "orientation",
146  false,
147  "",
148  "Reorient the image in 'AXIAL' or 'CORONAL' or 'SAGITTAL' direction. [defalut: No reorientation]",
149  cmd);
150 
151  TCLAP::ValueArg<std::string> gradsArg("g","grad","input gradients (apply converted image orientation matrix to get to MrTrix compatible format (iamge coordinates)",false,"","Input gradients",cmd);
152 
153  TCLAP::ValueArg<std::string> spaceArg("s",
154  "space",
155  "space_filename",
156  false,
157  "",
158  "Image used as space reference.",
159  cmd);
160 
161  try
162  {
163  cmd.parse(ac,av);
164  }
165  catch (TCLAP::ArgException& e)
166  {
167  std::cerr << "Error: " << e.error() << "for argument " << e.argId() << std::endl;
168  return EXIT_FAILURE;
169  }
170 
171  // Find out the type of the image in file
172  itk::ImageIOBase::Pointer imageIO = itk::ImageIOFactory::CreateImageIO(inputArg.getValue().c_str(),
173  itk::ImageIOFactory::ReadMode);
174 
175  if( !imageIO )
176  {
177  std::cerr << "Itk could not find suitable IO factory for the input" << std::endl;
178  return EXIT_FAILURE;
179  }
180 
181  // Now that we found the appropriate ImageIO class, ask it to read the meta data from the image file.
182  imageIO->SetFileName(inputArg.getValue());
183  imageIO->ReadImageInformation();
184 
185  arguments args;
186  args.input = inputArg.getValue(); args.output = outputArg.getValue(); args.displayInfo = infoArg.getValue();
187  args.reorient = reorientArg.getValue(); args.space = spaceArg.getValue();
188  args.gradientFileName = gradsArg.getValue();
189 
190  int numDimensions = imageIO->GetNumberOfDimensions() - 1;
191  if(args.displayInfo)
192  {
193  std::cout << "NUMBER OF DIMENSIONS: " << imageIO->GetNumberOfDimensions() << std::endl;
194 
195  std::cout << "SIZE : [";
196  for(int dim = 0; dim < numDimensions; ++dim)
197  std::cout << imageIO->GetDimensions(dim) << ", ";
198  std::cout << imageIO->GetDimensions(numDimensions) << "]" << std::endl;
199 
200  std::cout << "ORIGIN : [";
201  for(int dim = 0; dim < numDimensions; ++dim)
202  std::cout << imageIO->GetOrigin(dim) << ", ";
203  std::cout << imageIO->GetOrigin(numDimensions) << "]" << std::endl;
204 
205  std::cout << "SPACING : [";
206  for(int dim = 0; dim < numDimensions; ++dim)
207  std::cout << imageIO->GetSpacing(dim) << ", ";
208  std::cout << imageIO->GetSpacing(numDimensions) << "]" << std::endl;
209 
210  std::cout << "DIRECTIONS :" << std::endl;
211  for(unsigned int dir = 0; dir < imageIO->GetNumberOfDimensions(); ++dir)
212  {
213  std::cout << "\tDIM " << dir << " : [";
214  for (int dim = 0; dim < numDimensions; ++dim)
215  std::cout << imageIO->GetDirection(dim)[dir] << ", ";
216  std::cout << imageIO->GetDirection(numDimensions)[dir] << "]" << std::endl;
217  }
218  std::cout << "NUMBER OF COMPONENTS: "<< imageIO->GetNumberOfComponents()<< std::endl;
219  std::cout << "COMPONENT TYPE: "<< itk::ImageIOBase::GetComponentTypeAsString(imageIO->GetComponentType()) << std::endl;
220  }
221 
222  if ((args.output != "")||(args.gradientFileName != ""))
223  {
224  try
225  {
226  ANIMA_RETRIEVE_COMPONENT_TYPE(imageIO, retrieveNbDimensions, imageIO, args)
227  }
228  catch ( itk::ExceptionObject & err )
229  {
230  std::cerr << "Itk cannot convert, be sure to use valid arguments..." << std::endl;
231  std::cerr << err << std::endl;
232  return EXIT_FAILURE;
233  }
234  }
235 
236  return EXIT_SUCCESS;
237 }
std::string reorient
std::string output
std::string space
void checkIfComponentsAreVectors(itk::ImageIOBase::Pointer imageIO, const arguments &args)
itk::SmartPointer< ImageType > reorientImage(typename itk::SmartPointer< ImageType > input, itk::SpatialOrientation::ValidCoordinateOrientationFlags orientation)
void SetGradientFileName(std::string fName)
itk::ImageIOBase::Pointer imageIO
void reorientGradients(typename itk::SmartPointer< ImageType > input, std::vector< GradientType > &gradients)
#define ANIMA_RETRIEVE_COMPONENT_TYPE(imageIO, function,...)
void convert(const arguments &args)
#define ANIMA_RETRIEVE_NUMBER_OF_DIMENSIONS(imageIO, ComponentType, function,...)
std::string input
std::string gradientFileName
void retrieveNbDimensions(itk::ImageIOBase::Pointer imageIO, const arguments &args)
#define ANIMA_CHECK_IF_COMPONENTS_ARE_VECTORS(imageIO, ComponentType, Dimension, function,...)
int main(int ac, const char **av)