tools/clitkVFMerge.cxx

   1 /*=========================================================================
   2   Program:   vv                     http://www.creatis.insa-lyon.fr/rio/vv
   3
   4   Authors belong to:
   5   - University of LYON              http://www.universite-lyon.fr/
   6   - Léon Bérard cancer center       http://oncora1.lyon.fnclcc.fr
   7   - CREATIS CNRS laboratory         http://www.creatis.insa-lyon.fr
   8
   9   This software is distributed WITHOUT ANY WARRANTY; without even
  10   the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
  11   PURPOSE.  See the copyright notices for more information.
  12
  13   It is distributed under dual licence
  14
  15   - BSD        See included LICENSE.txt file
  16   - CeCILL-B   http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
  17 ======================================================================-====*/
  18 #ifndef CLITKVFMERGE_CXX
  19 #define CLITKVFMERGE_CXX
  20
  21 /**
  22  * @file   clitkVFMerge.cxx
  23  * @author Jef Vandemeulebroucke <jefvdmb@gmail.com>
  24  * @date   June 15  10:14:53 2007
  25  *
  26  * @brief  Read in one VF (ex mhd, vf)  and write to another.  Transforming from mm to voxels needed for the vf format is implemented in clitkVfImageIO.cxx .
  27  *
  28  *
  29  */
  30
  31 // clitk include
  32 #include "clitkVFMerge_ggo.h"
  33 #include "clitkIO.h"
  34 #include "clitkImageCommon.h"
  35
  36 // itk include
  37 //#include "itkReadTransform.h"
  38 #include "itkImageFileWriter.h"
  39 #include <iostream>
  40 #include "itkImageFileReader.h"
  41 //#include "itkRawImageIO.h"
  42 //#include "macro.h"
  43
  44 int main( int argc, char *argv[] )
  45 {
  46
  47   // Init command line
  48   GGO(clitkVFMerge, args_info);
  49   CLITK_INIT;
  50
  51 const unsigned int SpaceDimension = 3;
  52 const unsigned int ModelDimension = 4;
  53 typedef itk::Vector< float, SpaceDimension > Displacement;
  54 typedef itk::Image< Displacement, SpaceDimension > DeformationFieldType;
  55 typedef itk::Image< Displacement, ModelDimension > DynamicDeformationFieldType;
  56 typedef itk::ImageFileReader< DeformationFieldType > DeformationFieldReaderType;
  57 typedef itk::ImageFileWriter< DynamicDeformationFieldType > DynamicDeformationFieldWriterType;
  58
  59
  60 //declare the dynamic
  61  DynamicDeformationFieldType::Pointer  dynamicDeformationField=DynamicDeformationFieldType::New();
  62
  63
  64  //declare their iterators
  65  typedef itk::ImageRegionIterator< DynamicDeformationFieldType> DynamicDeformationFieldIteratorType;
  66  DynamicDeformationFieldIteratorType *dynamicIteratorPointer= new DynamicDeformationFieldIteratorType;
  67
  68  for (unsigned int i=0 ; i< args_info.inputs_num ; i ++)
  69    {
  70      //read in the deformation field i
  71      DeformationFieldReaderType::Pointer deformationFieldReader = DeformationFieldReaderType::New();
  72      deformationFieldReader->SetFileName( args_info.inputs[i]);
  73      if (args_info.verbose_flag) std::cout<<"Reading VF number "<< i+1 << std::endl;
  74      deformationFieldReader->Update();
  75      DeformationFieldType::Pointer currentDeformationField = deformationFieldReader->GetOutput();
  76
  77      //create an iterator for the current deformation field
  78      typedef itk::ImageRegionIterator<DeformationFieldType> FieldIteratorType;
  79      FieldIteratorType fieldIterator(currentDeformationField, currentDeformationField->GetLargestPossibleRegion());
  80
  81      //Allocate memory for the dynamic components
  82      if (i==0)
  83        {
  84          DynamicDeformationFieldType::RegionType dynamicDeformationFieldRegion;
  85          DynamicDeformationFieldType::RegionType::SizeType dynamicDeformationFieldSize;
  86          DeformationFieldType::RegionType::SizeType deformationFieldSize;
  87          deformationFieldSize= currentDeformationField->GetLargestPossibleRegion().GetSize();
  88          dynamicDeformationFieldSize[0]=deformationFieldSize[0];
  89          dynamicDeformationFieldSize[1]=deformationFieldSize[1];
  90          dynamicDeformationFieldSize[2]=deformationFieldSize[2];
  91          dynamicDeformationFieldSize[3]=args_info.inputs_num;
  92          dynamicDeformationFieldRegion.SetSize(dynamicDeformationFieldSize);
  93          DynamicDeformationFieldType::IndexType start;
  94          start.Fill(0);
  95          dynamicDeformationFieldRegion.SetIndex(start);
  96          dynamicDeformationField->SetRegions(dynamicDeformationFieldRegion);
  97          dynamicDeformationField->Allocate();
  98
  99
 100             //Set the spacing
 101              DeformationFieldType::SpacingType spacing= currentDeformationField->GetSpacing();
 102              DynamicDeformationFieldType::SpacingType dynamicSpacing;
 103             dynamicSpacing[0]=spacing[0];
 104             dynamicSpacing[1]=spacing[1];
 105             dynamicSpacing[2]=spacing[2];
 106             dynamicSpacing[3]=args_info.spacing_arg; //JV par exemple
 107             dynamicDeformationField->SetSpacing(dynamicSpacing);
 108         DynamicDeformationFieldType::PointType origin;
 109         origin[0]=args_info.xorigin_arg;
 110         origin[1]=args_info.yorigin_arg;
 111         origin[2]=args_info.zorigin_arg;
 112         origin[3]=0; //temporal origin is always 0
 113         dynamicDeformationField->SetOrigin(origin);
 114
 115             // Creating iterators for the currentDeformationField and the DynamicDeformationField
 116             DynamicDeformationFieldIteratorType *dynamicIterator= new DynamicDeformationFieldIteratorType(dynamicDeformationField, dynamicDeformationField->GetLargestPossibleRegion());
 117             dynamicIteratorPointer=dynamicIterator;
 118             dynamicIteratorPointer->GoToBegin();
 119        }
 120      if (args_info.verbose_flag) std::cout<<"Merging VF number "<< i+1 << std::endl;
 121      //Copy the current component of the input into dynamicDeformationFieldComponent
 122      fieldIterator.GoToBegin();
 123      while(!fieldIterator.IsAtEnd())
 124        {
 125          dynamicIteratorPointer->Set(fieldIterator.Get());
 126          ++fieldIterator;
 127          ++(*dynamicIteratorPointer);
 128        }
 129    }
 130
 131
 132 //Write the vector field
 133  DynamicDeformationFieldWriterType::Pointer writer = DynamicDeformationFieldWriterType::New();
 134  writer->SetInput( dynamicDeformationField );
 135  writer->SetFileName( args_info.output_arg );
 136  if (args_info.verbose_flag) std::cout<<"Writing the dynamic VF"<< std::endl;
 137
 138
 139  try {
 140
 141    writer->Update( );
 142  }
 143  catch( itk::ExceptionObject & excp ) {
 144    std::cerr << "Problem writing the output file" << std::endl;
 145    std::cerr << excp << std::endl;
 146    return EXIT_FAILURE;
 147  }
 148  return EXIT_SUCCESS;
 149 }
 150 #endif
 151
 152