1 /*=========================================================================
2 Program: vv http://www.creatis.insa-lyon.fr/rio/vv
5 - University of LYON http://www.universite-lyon.fr/
6 - Léon Bérard cancer center http://www.centreleonberard.fr
7 - CREATIS CNRS laboratory http://www.creatis.insa-lyon.fr
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.
13 It is distributed under dual licence
15 - BSD See included LICENSE.txt file
16 - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
17 ===========================================================================**/
18 #ifndef clitkProfileImageGenericFilter_cxx
19 #define clitkProfileImageGenericFilter_cxx
21 /* =================================================
22 * @file clitkProfileImageGenericFilter.cxx
23 * @author Thomas Baudier <thomas.baudier@creatis.insa-lyon.fr>
28 ===================================================*/
30 #include "clitkProfileImageGenericFilter.h"
33 #include <itkLineIterator.h>
36 #include <clitkCommon.h>
43 //--------------------------------------------------------------------
44 ProfileImageGenericFilter::ProfileImageGenericFilter():
45 ImageToImageGenericFilter<Self>("ProfileImage")
47 InitializeImageType<2>();
48 InitializeImageType<3>();
49 InitializeImageType<4>();
51 //--------------------------------------------------------------------
54 //--------------------------------------------------------------------
55 template<unsigned int Dim>
56 void ProfileImageGenericFilter::InitializeImageType()
58 ADD_DEFAULT_IMAGE_TYPES(Dim);
60 //--------------------------------------------------------------------
63 //--------------------------------------------------------------------
64 vtkFloatArray* ProfileImageGenericFilter::GetArrayX()
68 //--------------------------------------------------------------------
71 //--------------------------------------------------------------------
72 vtkFloatArray* ProfileImageGenericFilter::GetArrayY()
76 //--------------------------------------------------------------------
79 //--------------------------------------------------------------------
80 vtkFloatArray* ProfileImageGenericFilter::GetCoord()
84 //--------------------------------------------------------------------
87 //--------------------------------------------------------------------
88 void ProfileImageGenericFilter::SetArgsInfo(const args_info_type & a)
91 if (mArgsInfo.verbose_given)
92 SetIOVerbose(mArgsInfo.verbose_flag);
93 if (mArgsInfo.imagetypes_given && mArgsInfo.imagetypes_flag)
94 this->PrintAvailableImageTypes();
96 if (mArgsInfo.input_given) {
97 SetInputFilename(mArgsInfo.input_arg);
99 if (mArgsInfo.output_given) {
100 SetOutputFilename(mArgsInfo.output_arg);
103 //--------------------------------------------------------------------
106 //--------------------------------------------------------------------
107 // Update with the number of dimensions and the pixeltype
108 //--------------------------------------------------------------------
109 template<class InputImageType>
111 ProfileImageGenericFilter::UpdateWithInputImageType()
115 typename InputImageType::Pointer input = this->template GetInput<InputImageType>(0);
116 typedef typename InputImageType::PixelType PixelType;
117 typedef typename InputImageType::IndexType IndexType;
119 mArrayX = vtkSmartPointer<vtkFloatArray>::New();
120 mArrayY = vtkSmartPointer<vtkFloatArray>::New();
121 mCoord = vtkSmartPointer<vtkFloatArray>::New();
122 mCoord->SetNumberOfComponents(InputImageType::ImageDimension);
123 mCoordmm = vtkSmartPointer<vtkFloatArray>::New();
124 mCoordmm->SetNumberOfComponents(InputImageType::ImageDimension);
125 mDimension = InputImageType::ImageDimension;
127 /*typename InputImageType::Pointer outputImage;
128 outputImage = InputImageType::New();
130 outputImage->SetRegions(input->GetLargestPossibleRegion());
131 outputImage->Allocate();
132 outputImage->FillBuffer(0); */
135 IndexType pointBegin, pointEnd;
137 for (int i = 0; i < mArgsInfo.point1_given; ++i) {
138 pointBegin[i] = mArgsInfo.point1_arg[i];
139 pointEnd[i] = mArgsInfo.point2_arg[i];
142 itk::LineConstIterator<InputImageType> itProfile(input, pointBegin, pointEnd);
143 itProfile.GoToBegin();
147 tuple = new double[InputImageType::ImageDimension];
148 itk::Point<double, InputImageType::ImageDimension> transformedFirstPoint;
149 itk::Point<double, InputImageType::ImageDimension> transformedCurrentPoint;
151 input->TransformIndexToPhysicalPoint(itProfile.GetIndex(), transformedFirstPoint);
153 while (!itProfile.IsAtEnd())
155 // Fill in the table the intensity value
156 mArrayY->InsertNextTuple1(itProfile.Get());
158 for (int i=0; i<InputImageType::ImageDimension; ++i) {
159 tuple[i] = itProfile.GetIndex()[i];
162 input->TransformIndexToPhysicalPoint(itProfile.GetIndex(), transformedCurrentPoint);
163 distance = transformedFirstPoint.EuclideanDistanceTo(transformedCurrentPoint);
165 // Fill in the table the distance value
166 mArrayX->InsertNextTuple1(distance);
168 // Fill in the table the voxel coordinate value
169 mCoord->InsertNextTuple(tuple); //index
170 for (int i=0; i<InputImageType::ImageDimension; ++i) {
171 tuple[i] = transformedCurrentPoint[i];
173 mCoordmm->InsertNextTuple(tuple); //mm
178 if (mArgsInfo.output_given) {
179 std::string str(mArgsInfo.output_arg);
180 this->WriteOutput(str);
184 itk::LineIterator<InputImageType> otProfile(outputImage, pointBegin, pointEnd);
185 otProfile.GoToBegin();
186 while (!otProfile.IsAtEnd())
192 this->template SetNextOutput<InputImageType>(outputImage): */
196 //--------------------------------------------------------------------
199 //--------------------------------------------------------------------
200 void ProfileImageGenericFilter::WriteOutput(std::string outputFilename)
202 ofstream fileOpen(outputFilename.c_str(), std::ofstream::trunc);
205 cerr << "Error during saving" << endl;
210 tuple = new double[mDimension];
212 fileOpen << "The Bresenham algorithm is used to travel along the line. Values represent the center of each crossed voxel (in voxel and mm)" << endl;
213 fileOpen << "Id" << "\t" << "Value" << "\t" ;
214 fileOpen << "x(vox)" << "\t" << "y(vox)" << "\t";
216 fileOpen << "z(vox)" << "\t";
218 fileOpen << "t" << "\t";
219 fileOpen << "x(mm)" << "\t" << "y(mm)" << "\t";
221 fileOpen << "z(mm)" << "\t";
223 fileOpen << "t" << "\t";
226 while (i<mArrayX->GetNumberOfTuples()) {
227 fileOpen << i << "\t" << mArrayY->GetTuple(i)[0] << "\t" ;
229 mCoord->GetTuple(i, tuple);
230 for (int j=0; j<mDimension ; ++j) {
231 fileOpen << tuple[j] << "\t" ;
233 mCoordmm->GetTuple(i, tuple);
234 for (int j=0; j<mDimension ; ++j) {
235 fileOpen << tuple[j] << "\t" ;
237 if (mDimension == 4) {
238 fileOpen << tuple[3] << "\t" ;
248 //--------------------------------------------------------------------
253 #endif //#define clitkProfileImageGenericFilter_cxx