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 clitkLineProfileGenericFilter_cxx
19 #define clitkLineProfileGenericFilter_cxx
21 /* =================================================
22 * @file clitkLineProfileGenericFilter.cxx
28 ===================================================*/
30 #include "itkBresenhamLine.h"
35 //-----------------------------------------------------------
37 //-----------------------------------------------------------
38 template<class InputImageType>
39 void LineProfileGenericFilter::UpdateWithInputImageType()
41 typedef InputImageType ImageType;
44 std::cout << "LineProfileGenericFilter::UpdateWithInputImageType" << std::endl;
47 if (m_ArgsInfo.p0_given != ImageType::ImageDimension ||
48 m_ArgsInfo.p1_given != ImageType::ImageDimension)
49 throw std::logic_error("Dimension of input points and image do not match.");
51 typename ImageType::Pointer image = this->template GetInput<InputImageType>(0);
52 typename ImageType::RegionType region = image->GetLargestPossibleRegion();
54 typedef typename ImageType::PointType PointType;
56 for (unsigned int i = 0; i < ImageType::ImageDimension; i++) {
57 p0[i] = m_ArgsInfo.p0_arg[i];
58 p1[i] = m_ArgsInfo.p1_arg[i];
61 // compute params of line between p0 and p1
62 // length (magnitude) must be an integer value, so it's
63 // the max distance along one the axes plus one to account
64 // for the last point.
65 typedef itk::BresenhamLine<ImageType::ImageDimension> LineType;
66 typename LineType::LType direction = p1 - p0;
67 typename LineType::LType::RealValueType mag = 0;
68 for (unsigned int i = 0; i < ImageType::ImageDimension; i++) {
69 if (direction[i] > mag)
75 std::cout << "Building line with direction = " << direction << " and length = " << mag << "..." << std::endl;
77 // build the line itself
79 typename LineType::OffsetArray offsets;
80 offsets = line.BuildLine(direction, mag);
83 std::cout << "Line has " << offsets.size() << " points" << std::endl;
85 // fill the output vectors
86 typedef typename ImageType::PixelType PixelType;
87 typedef typename ImageType::OffsetType OffsetType;
88 typedef typename ImageType::IndexType IndexType;
89 typedef std::vector<IndexType> IndexArrayType;
90 typedef std::vector<PixelType> ValueArrayType;
92 IndexType index0, index1;
93 for (unsigned int i = 0; i < ImageType::ImageDimension; i++) {
94 index0[i] = m_ArgsInfo.p0_arg[i];
95 index1[i] = m_ArgsInfo.p1_arg[i];
99 std::cout << "Getting profile along line..." << std::endl;
102 IndexArrayType line_indices;
103 ValueArrayType line_values;
104 for (size_t i = 0; i < offsets.size(); i++)
106 index = index0 + offsets[i];
108 std::cout << "index " << i << " = " << index << std::endl;
111 if (region.IsInside(index)) {
113 std::cout << "value " << i << " = " << image->GetPixel(index) << std::endl;
115 line_indices.push_back(index);
116 line_values.push_back(image->GetPixel(index));
119 std::cout << "index outside image" << std::endl;
123 std::cout << "I bring to you The Computed Points!" << std::endl;
126 for (size_t i = 0; i < line_indices.size(); i++) {
127 std::cout << line_indices[i] << " " << line_values[i] << std::endl;
134 #endif //#define clitkLineProfileGenericFilter_cxx