Skip to content
Snippets Groups Projects
Commit 311382ea authored by Christophe Palmann's avatar Christophe Palmann
Browse files

ENH: BCO interpolation (1)

parent 7a3579b0
No related branches found
No related tags found
No related merge requests found
...@@ -107,20 +107,23 @@ public: ...@@ -107,20 +107,23 @@ public:
virtual OutputType EvaluateAtContinuousIndex( const ContinuousIndexType & index ) const = 0; virtual OutputType EvaluateAtContinuousIndex( const ContinuousIndexType & index ) const = 0;
protected: protected:
BCOInterpolateImageFunctionBase() : m_Radius(2), m_Alpha(-0.5) {}; BCOInterpolateImageFunctionBase() : m_Radius(2), m_WinSize(5), m_Alpha(-0.5) {};
virtual ~BCOInterpolateImageFunctionBase() {}; virtual ~BCOInterpolateImageFunctionBase() {};
void PrintSelf(std::ostream& os, itk::Indent indent) const; void PrintSelf(std::ostream& os, itk::Indent indent) const;
/** Compute the BCO coefficients. */ /** Compute the BCO coefficients. */
virtual CoefContainerType EvaluateCoef( const ContinuousIndexValueType & indexValue ) const; virtual CoefContainerType EvaluateCoef( const ContinuousIndexValueType & indexValue ) const;
/** Used radius for the BCO */
unsigned int m_Radius;
/** Used winsize for the BCO */
unsigned int m_WinSize;
/** Optimisation Coefficient */
double m_Alpha;
private: private:
BCOInterpolateImageFunctionBase( const Self& ); //purposely not implemented BCOInterpolateImageFunctionBase( const Self& ); //purposely not implemented
void operator=( const Self& ); //purposely not implemented void operator=( const Self& ); //purposely not implemented
/** Used radius for the BCO */
unsigned int m_Radius;
/** Optimisation Coefficient */
double m_Alpha;
}; };
......
...@@ -45,6 +45,7 @@ void BCOInterpolateImageFunctionBase<TInputImage, TCoordRep> ...@@ -45,6 +45,7 @@ void BCOInterpolateImageFunctionBase<TInputImage, TCoordRep>
else else
{ {
m_Radius = radius; m_Radius = radius;
m_WinSize = 2*m_Radius+1;
} }
} }
...@@ -76,21 +77,21 @@ BCOInterpolateImageFunctionBase<TInputImage, TCoordRep> ...@@ -76,21 +77,21 @@ BCOInterpolateImageFunctionBase<TInputImage, TCoordRep>
::EvaluateCoef( const ContinuousIndexValueType & indexValue ) const ::EvaluateCoef( const ContinuousIndexValueType & indexValue ) const
{ {
// Init BCO coefficient container // Init BCO coefficient container
int radius = this->GetRadius();
unsigned int winSize = 2*radius+1; CoefContainerType BCOCoef(m_WinSize, 0.);
CoefContainerType BCOCoef = CoefContainerType(winSize, 0.);
double offset, dist, position, step; double offset, dist, position, step;
offset = indexValue - itk::Math::Floor<IndexValueType>(indexValue+0.5); offset = indexValue - itk::Math::Floor<IndexValueType>(indexValue+0.5);
// Compute BCO coefficients // Compute BCO coefficients
step = 4./static_cast<double>(2*radius); step = 4./static_cast<double>(2*m_Radius);
position = - double(radius) * step; position = - static_cast<double>(m_Radius) * step;
double sum = 0.0; double sum = 0.0;
for ( int i = -radius; i <= radius; ++i) for ( unsigned int i = 0; i < m_WinSize; ++i)
{ {
// Compute the BCO coefficients according to alpha. // Compute the BCO coefficients according to alpha.
dist = vcl_abs(position - offset*step); dist = vcl_abs(position - offset*step);
...@@ -98,25 +99,25 @@ BCOInterpolateImageFunctionBase<TInputImage, TCoordRep> ...@@ -98,25 +99,25 @@ BCOInterpolateImageFunctionBase<TInputImage, TCoordRep>
{ {
if (dist <= 1.) if (dist <= 1.)
{ {
BCOCoef[radius+i] = (m_Alpha + 2.)*vcl_abs(dist * dist * dist) BCOCoef[i] = (m_Alpha + 2.)*vcl_abs(dist * dist * dist)
- (m_Alpha + 3.)*dist*dist + 1; - (m_Alpha + 3.)*dist*dist + 1;
} }
else else
{ {
BCOCoef[radius+i] = m_Alpha*vcl_abs(dist * dist * dist) - 5 BCOCoef[i] = m_Alpha*vcl_abs(dist * dist * dist) - 5
*m_Alpha*dist*dist + 8*m_Alpha*vcl_abs(dist) - 4*m_Alpha; *m_Alpha*dist*dist + 8*m_Alpha*vcl_abs(dist) - 4*m_Alpha;
} }
} }
else else
{ {
BCOCoef[m_Radius+i] = 0; BCOCoef[i] = 0;
} }
sum += BCOCoef[m_Radius+i]; sum += BCOCoef[i];
position += step; position += step;
} }
for ( unsigned int i = 0; i < winSize; ++i) for ( unsigned int i = 0; i < m_WinSize; ++i)
BCOCoef[i] = BCOCoef[i] / sum; BCOCoef[i] = BCOCoef[i] / sum;
return BCOCoef; return BCOCoef;
...@@ -269,10 +270,11 @@ BCOInterpolateImageFunction< otb::VectorImage<TPixel, VImageDimension> , TCoordR ...@@ -269,10 +270,11 @@ BCOInterpolateImageFunction< otb::VectorImage<TPixel, VImageDimension> , TCoordR
neighIndex[1] = this->m_StartIndex[1]; neighIndex[1] = this->m_StartIndex[1];
} }
const InputPixelType & pixel = this->GetInputImage()->GetPixel( neighIndex );
for( unsigned int k = 0; k<componentNumber; ++k) for( unsigned int k = 0; k<componentNumber; ++k)
{ {
lineRes[i+radius].at(k) = lineRes[i+radius].at(k) lineRes[i+this->m_Radius].at(k) = lineRes[i+this->m_Radius].at(k)
+ this->GetInputImage()->GetPixel( neighIndex ).GetElement(k) * BCOCoefY[j+radius]; + pixel.GetElement(k) * BCOCoefY[j+this->m_Radius];
} }
} }
for( unsigned int k = 0; k<componentNumber; ++k) for( unsigned int k = 0; k<componentNumber; ++k)
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment