]> Creatis software - clitk.git/blob - tools/clitkVFMerge.cxx
Romulo:
[clitk.git] / 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 #include "clitkCommon.h"
36
37 // itk include
38 #include "itkImageFileWriter.h"
39 #include <iostream>
40 #include "itkImageFileReader.h"
41
42 int main( int argc, char *argv[] )
43 {
44
45   // Init command line
46   GGO(clitkVFMerge, args_info);
47   CLITK_INIT;
48
49   const unsigned int SpaceDimension = 3;
50   const unsigned int ModelDimension = 4;
51   typedef itk::Vector< float, SpaceDimension > Displacement;
52   typedef itk::Image< Displacement, SpaceDimension > DeformationFieldType;
53   typedef itk::Image< Displacement, ModelDimension > DynamicDeformationFieldType;
54   typedef itk::ImageFileReader< DeformationFieldType > DeformationFieldReaderType;
55   typedef itk::ImageFileWriter< DynamicDeformationFieldType > DynamicDeformationFieldWriterType;
56
57
58 //declare the dynamic
59   DynamicDeformationFieldType::Pointer  dynamicDeformationField=DynamicDeformationFieldType::New();
60
61
62 //declare their iterators
63   typedef itk::ImageRegionIterator< DynamicDeformationFieldType> DynamicDeformationFieldIteratorType;
64   DynamicDeformationFieldIteratorType *dynamicIteratorPointer= new DynamicDeformationFieldIteratorType;
65
66   for (unsigned int i=0 ; i< args_info.inputs_num ; i ++) {
67     //read in the deformation field i
68     DeformationFieldReaderType::Pointer deformationFieldReader = DeformationFieldReaderType::New();
69     deformationFieldReader->SetFileName( args_info.inputs[i]);
70     if (args_info.verbose_flag) std::cout<<"Reading VF number "<< i+1 << std::endl;
71     deformationFieldReader->Update();
72     DeformationFieldType::Pointer currentDeformationField = deformationFieldReader->GetOutput();
73
74     //create an iterator for the current deformation field
75     typedef itk::ImageRegionIterator<DeformationFieldType> FieldIteratorType;
76     FieldIteratorType fieldIterator(currentDeformationField, currentDeformationField->GetLargestPossibleRegion());
77
78     //Allocate memory for the dynamic components
79     if (i==0) {
80       DynamicDeformationFieldType::RegionType dynamicDeformationFieldRegion;
81       DynamicDeformationFieldType::RegionType::SizeType dynamicDeformationFieldSize;
82       DeformationFieldType::RegionType::SizeType deformationFieldSize;
83       deformationFieldSize= currentDeformationField->GetLargestPossibleRegion().GetSize();
84       dynamicDeformationFieldSize[0]=deformationFieldSize[0];
85       dynamicDeformationFieldSize[1]=deformationFieldSize[1];
86       dynamicDeformationFieldSize[2]=deformationFieldSize[2];
87       dynamicDeformationFieldSize[3]=args_info.inputs_num;
88       dynamicDeformationFieldRegion.SetSize(dynamicDeformationFieldSize);
89       DynamicDeformationFieldType::IndexType start;
90       start.Fill(0);
91       dynamicDeformationFieldRegion.SetIndex(start);
92       dynamicDeformationField->SetRegions(dynamicDeformationFieldRegion);
93       dynamicDeformationField->Allocate();
94
95
96       //Set the spacing
97       DeformationFieldType::SpacingType spacing= currentDeformationField->GetSpacing();
98       DynamicDeformationFieldType::SpacingType dynamicSpacing;
99       dynamicSpacing[0]=spacing[0];
100       dynamicSpacing[1]=spacing[1];
101       dynamicSpacing[2]=spacing[2];
102       dynamicSpacing[3]=args_info.spacing_arg; //JV par exemple
103       dynamicDeformationField->SetSpacing(dynamicSpacing);
104       DynamicDeformationFieldType::PointType origin;
105       origin[0]=args_info.xorigin_arg;
106       origin[1]=args_info.yorigin_arg;
107       origin[2]=args_info.zorigin_arg;
108       origin[3]=0; //temporal origin is always 0
109       dynamicDeformationField->SetOrigin(origin);
110
111       // Creating iterators for the currentDeformationField and the DynamicDeformationField
112       DynamicDeformationFieldIteratorType *dynamicIterator= new DynamicDeformationFieldIteratorType(dynamicDeformationField, dynamicDeformationField->GetLargestPossibleRegion());
113       dynamicIteratorPointer=dynamicIterator;
114       dynamicIteratorPointer->GoToBegin();
115     }
116     if (args_info.verbose_flag) std::cout<<"Merging VF number "<< i+1 << std::endl;
117     //Copy the current component of the input into dynamicDeformationFieldComponent
118     fieldIterator.GoToBegin();
119     while(!fieldIterator.IsAtEnd()) {
120       dynamicIteratorPointer->Set(fieldIterator.Get());
121       ++fieldIterator;
122       ++(*dynamicIteratorPointer);
123     }
124   }
125
126
127 //Write the vector field
128   DynamicDeformationFieldWriterType::Pointer writer = DynamicDeformationFieldWriterType::New();
129   writer->SetInput( dynamicDeformationField );
130   writer->SetFileName( args_info.output_arg );
131   if (args_info.verbose_flag) std::cout<<"Writing the dynamic VF"<< std::endl;
132
133
134   try {
135
136     writer->Update( );
137   } catch( itk::ExceptionObject & excp ) {
138     std::cerr << "Problem writing the output file" << std::endl;
139     std::cerr << excp << std::endl;
140     return EXIT_FAILURE;
141   }
142   return EXIT_SUCCESS;
143 }
144 #endif
145
146