/usr/include/ITK-4.5/itkBresenhamLine.hxx is in libinsighttoolkit4-dev 4.5.0-3.
This file is owned by root:root, with mode 0o644.
The actual contents of the file can be viewed below.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 | /*=========================================================================
*
* Copyright Insight Software Consortium
*
* Licensed under the Apache License, Version 2.0 (the "License");
* you may not use this file except in compliance with the License.
* You may obtain a copy of the License at
*
* http://www.apache.org/licenses/LICENSE-2.0.txt
*
* Unless required by applicable law or agreed to in writing, software
* distributed under the License is distributed on an "AS IS" BASIS,
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
* limitations under the License.
*
*=========================================================================*/
#ifndef __itkBresenhamLine_hxx
#define __itkBresenhamLine_hxx
#include "itkBresenhamLine.h"
#include "itkPoint.h"
#include "itkMath.h"
#include "vnl/vnl_math.h"
namespace itk
{
template< unsigned int VDimension >
typename BresenhamLine< VDimension >::OffsetArray BresenhamLine< VDimension >
::BuildLine(LType Direction, unsigned int length)
{
// copied from the line iterator
/** Variables that drive the Bresenham-Algorithm */
// The dimension with the largest difference between start and end
unsigned int m_MainDirection;
// Accumulated error for the other dimensions
IndexType m_AccumulateError;
// Increment for the error for each step. Two times the difference between
// start and end
IndexType m_IncrementError;
// If enough is accumulated for a dimension, the index has to be
// incremented. Will be the number of pixels in the line
IndexType m_MaximalError;
// Direction of increment. -1 or 1
IndexType m_OverflowIncrement;
// After an overflow, the accumulated error is reduced again. Will be
// two times the number of pixels in the line
IndexType m_ReduceErrorAfterIncrement;
OffsetArray result(length);
IndexType m_CurrentImageIndex, StartIndex, LastIndex;
Direction.Normalize();
// we are going to start at 0
m_CurrentImageIndex.Fill(0);
StartIndex.Fill(0);
for ( unsigned i = 0; i < VDimension; i++ )
{
LastIndex[i] = (IndexValueType)( length * Direction[i] );
}
// Find the dominant direction
IndexValueType maxDistance = 0;
unsigned int maxDistanceDimension = 0;
for ( unsigned i = 0; i < VDimension; i++ )
{
IndexValueType distance = static_cast<long>(vnl_math_abs(LastIndex[i]));
if ( distance > maxDistance )
{
maxDistance = distance;
maxDistanceDimension = i;
}
m_IncrementError[i] = 2 * distance;
m_OverflowIncrement[i] = ( LastIndex[i] < 0 ? -1 : 1 );
}
m_MainDirection = maxDistanceDimension;
m_MaximalError.Fill(maxDistance);
m_ReduceErrorAfterIncrement.Fill(2 * maxDistance);
m_AccumulateError.Fill(0);
unsigned int steps = 1;
result[0] = m_CurrentImageIndex - StartIndex;
while ( steps < length )
{
// This part is from ++ in LineConstIterator
// We need to modify m_AccumulateError, m_CurrentImageIndex, m_IsAtEnd
for ( unsigned int i = 0; i < VDimension; ++i )
{
if ( i == m_MainDirection )
{
m_CurrentImageIndex[i] += m_OverflowIncrement[i];
}
else
{
m_AccumulateError[i] += m_IncrementError[i];
if ( m_AccumulateError[i] >= m_MaximalError[i] )
{
m_CurrentImageIndex[i] += m_OverflowIncrement[i];
m_AccumulateError[i] -= m_ReduceErrorAfterIncrement[i];
}
}
}
result[steps] = m_CurrentImageIndex - StartIndex; // produce an offset
++steps;
}
return ( result );
}
template< unsigned int VDimension >
typename BresenhamLine< VDimension >::IndexArray BresenhamLine< VDimension >
::BuildLine(IndexType p0, IndexType p1)
{
itk::Point<float,VDimension> point0;
itk::Point<float,VDimension> point1;
for(unsigned int i = 0; i < VDimension; i++)
{
point0[i] = p0[i];
point1[i] = p1[i];
}
const unsigned int distance = itk::Math::RoundHalfIntegerToEven<unsigned int, double>( point0.EuclideanDistanceTo(point1) );
OffsetArray offsets = this->BuildLine(point1-point0, distance);
IndexArray indices(offsets.size());
for(unsigned int i = 0; i < offsets.size(); i++)
{
indices[i] = p0 + offsets[i];
}
return indices;
}
} // namespace itk
#endif
|