Skip to content
Snippets Groups Projects
Commit 955d9be7 authored by Emmanuel Christophe's avatar Emmanuel Christophe
Browse files

ENH: clean up BCO

parent 9fc6206e
No related branches found
No related tags found
No related merge requests found
......@@ -28,19 +28,19 @@ namespace otb
{
/** \class BCOInterpolateImageFunction
* \brief Interpolate an image at specified positions using bicubic interpolation.
*
*
* BCOInterpolateImageFunction interpolates image intensity at
* a non-integer pixel position. This class is templated
* over the input image type and the coordinate representation type
* over the input image type and the coordinate representation type
* (e.g. float or double).
*
* This function works for 2-dimensional images.
*
*
* This function works with both Images and VectorImages.
*
*
* Parameters are the interpolation window radius and the bicubic
* optimisation coefficient alpha.
* Alpha is usually set to -0.5, -0.75 or -1 (-0.5 by default).
* Alpha is usually set to -0.5, -0.75 or -1 (-0.5 by default).
* The case alpha = -0.5 (which corresponds to the cubic Hermite
* spline) is known to produce the best approximation of the original
* function.
......@@ -48,36 +48,36 @@ namespace otb
* \ingroup ImageFunctions ImageInterpolators
*/
template< class TInputImage, class TCoordRep = double >
class ITK_EXPORT BCOInterpolateImageFunctionBase :
class ITK_EXPORT BCOInterpolateImageFunctionBase :
public itk::InterpolateImageFunction<TInputImage,TCoordRep>
{
public:
/** Standard class typedefs. */
typedef BCOInterpolateImageFunctionBase Self;
typedef itk::InterpolateImageFunction<TInputImage,TCoordRep> Superclass;
/** Run-time type information (and related methods). */
itkTypeMacro(BCOInterpolateImageFunctionBase, InterpolateImageFunction);
/** OutputType typedef support. */
typedef typename Superclass::OutputType OutputType;
/** InputImageType typedef support. */
typedef typename Superclass::InputImageType InputImageType;
/** InputPixelType typedef support. */
typedef typename Superclass::InputPixelType InputPixelType;
/** RealType typedef support. */
typedef typename Superclass::RealType RealType;
/** Dimension underlying input image. */
itkStaticConstMacro(ImageDimension, unsigned int,Superclass::ImageDimension);
/** Index typedef support. */
typedef typename Superclass::IndexType IndexType;
typedef typename Superclass::IndexValueType IndexValueType;
/** Point typedef support. */
typedef typename Superclass::PointType PointType;
......@@ -90,29 +90,29 @@ public:
/** Set/Get the window radius */
virtual void SetRadius(unsigned int radius);
virtual unsigned int GetRadius() const;
/** Set/Get the optimisation coefficient (Common values are -0.5, -0.75 or -1.0) */
virtual void SetAlpha(double alpha);
virtual double GetAlpha() const;
/** Evaluate the function at a ContinuousIndex position
*
* Returns the linearly interpolated image intensity at a
* Returns the linearly interpolated image intensity at a
* specified point position. No bounds checking is done.
* The point is assume to lie within the image buffer.
*
* ImageFunction::IsInsideBuffer() can be used to check bounds before
* calling the method. */
virtual OutputType EvaluateAtContinuousIndex( const ContinuousIndexType & index ) const = 0;
/** Compute the BCO coefficients. */
virtual void EvaluateCoef( const ContinuousIndexType & index ) const;
protected:
BCOInterpolateImageFunctionBase();
~BCOInterpolateImageFunctionBase();
BCOInterpolateImageFunctionBase() : m_Radius(2), m_Alpha(-0.5) {};
virtual ~BCOInterpolateImageFunctionBase() {};
void PrintSelf(std::ostream& os, itk::Indent indent) const;
private:
BCOInterpolateImageFunctionBase( const Self& ); //purposely not implemented
void operator=( const Self& ); //purposely not implemented
......@@ -120,13 +120,13 @@ private:
/** Used radius for the BCO */
double m_Radius;
/** Optimisation Coefficient */
double m_Alpha;
double m_Alpha;
};
template < class TInputImage, class TCoordRep = double >
class ITK_EXPORT BCOInterpolateImageFunction :
public otb::BCOInterpolateImageFunctionBase< TInputImage, TCoordRep >
class ITK_EXPORT BCOInterpolateImageFunction :
public otb::BCOInterpolateImageFunctionBase< TInputImage, TCoordRep >
{
public:
/** Standard class typedefs. */
......@@ -134,9 +134,9 @@ public:
typedef BCOInterpolateImageFunctionBase<TInputImage,TCoordRep> Superclass;
typedef itk::SmartPointer<Self> Pointer;
typedef itk::SmartPointer<const Self> ConstPointer;
itkTypeMacro(BCOInterpolateImageFunction, BCOInterpolateImageFunctionBase);
itkNewMacro(Self);
itkNewMacro(Self);
itkStaticConstMacro(ImageDimension, unsigned int,Superclass::ImageDimension);
typedef typename Superclass::OutputType OutputType;
......@@ -152,8 +152,8 @@ public:
virtual OutputType EvaluateAtContinuousIndex( const ContinuousIndexType & index ) const;
protected:
BCOInterpolateImageFunction();
~BCOInterpolateImageFunction();
BCOInterpolateImageFunction() {};
virtual ~BCOInterpolateImageFunction() {};
void PrintSelf(std::ostream& os, itk::Indent indent) const;
private:
......@@ -163,19 +163,19 @@ private:
template < typename TPixel, unsigned int VImageDimension, class TCoordRep >
class ITK_EXPORT BCOInterpolateImageFunction< otb::VectorImage<TPixel,VImageDimension> , TCoordRep > :
public otb::BCOInterpolateImageFunctionBase< otb::VectorImage<TPixel,VImageDimension> , TCoordRep >
class ITK_EXPORT BCOInterpolateImageFunction< otb::VectorImage<TPixel,VImageDimension> , TCoordRep > :
public otb::BCOInterpolateImageFunctionBase< otb::VectorImage<TPixel,VImageDimension> , TCoordRep >
{
public:
/** Standard class typedefs.*/
/** Standard class typedefs.*/
typedef BCOInterpolateImageFunction Self;
typedef BCOInterpolateImageFunctionBase
< otb::VectorImage<TPixel,VImageDimension>, TCoordRep > Superclass;
typedef itk::SmartPointer<Self> Pointer;
typedef itk::SmartPointer<const Self> ConstPointer;
itkTypeMacro(BCOInterpolateImageFunction, BCOInterpolateImageFunctionBase);
itkNewMacro(Self);
itkNewMacro(Self);
itkStaticConstMacro(ImageDimension, unsigned int,Superclass::ImageDimension);
typedef typename Superclass::OutputType OutputType;
......@@ -191,8 +191,8 @@ public:
virtual OutputType EvaluateAtContinuousIndex( const ContinuousIndexType & index ) const;
protected:
BCOInterpolateImageFunction();
~BCOInterpolateImageFunction();
BCOInterpolateImageFunction() {};
virtual ~BCOInterpolateImageFunction() {};
void PrintSelf(std::ostream& os, itk::Indent indent) const;
private:
......
......@@ -25,27 +25,13 @@
namespace otb
{
/** Constructor */
template <class TInputImage, class TCoordRep>
BCOInterpolateImageFunctionBase<TInputImage, TCoordRep>
::BCOInterpolateImageFunctionBase()
{
m_Radius = 2;
m_Alpha = -0.5;
}
/** Destructor */
template <class TInputImage, class TCoordRep>
BCOInterpolateImageFunctionBase<TInputImage, TCoordRep>
::~BCOInterpolateImageFunctionBase()
{
}
template <class TInputImage, class TCoordRep>
void BCOInterpolateImageFunctionBase<TInputImage, TCoordRep>
::PrintSelf(std::ostream& os, itk::Indent indent) const
{
Superclass::PrintSelf(os, indent);
os << indent << "Radius: " << m_Radius << std::endl;
os << indent << "Alpha: " << m_Alpha << std::endl;
}
template <class TInputImage, class TCoordRep>
......@@ -54,9 +40,9 @@ void BCOInterpolateImageFunctionBase<TInputImage, TCoordRep>
{
if (radius < 2)
{
itkExceptionMacro(<< "Radius Must Be Strictly Superior to 1");
itkExceptionMacro(<< "Radius must be strictly greater than 1");
}
else
else
{
m_Radius = radius;
}
......@@ -96,22 +82,22 @@ BCOInterpolateImageFunctionBase<TInputImage, TCoordRep>
offsetX = index[0] - itk::Math::Floor<IndexValueType>(index[0]+0.5);
offsetY = index[1] - itk::Math::Floor<IndexValueType>(index[1]+0.5);
// Compute BCO coefficients
step = 4./static_cast<double>(2*m_Radius);
position = - double(m_Radius) * step;
for ( int i = -m_Radius; i <= m_Radius; i++)
{
// Compute the BCO coefficients according to alpha.
distX = vcl_abs(position - offsetX*step);
distY = vcl_abs(position - offsetY*step);
if( distX <= 2. )
{
if (distX <= 1.)
{
BCOCoefX[m_Radius+i] = (m_Alpha + 2.)*vcl_abs(vcl_pow(distX, 3))
BCOCoefX[m_Radius+i] = (m_Alpha + 2.)*vcl_abs(vcl_pow(distX, 3))
- (m_Alpha + 3.)*vcl_pow(distX, 2) + 1;
}
else
......@@ -129,7 +115,7 @@ BCOInterpolateImageFunctionBase<TInputImage, TCoordRep>
{
if (distY <= 1.)
{
BCOCoefY[m_Radius+i] = (m_Alpha + 2.)*vcl_abs(vcl_pow(distY, 3))
BCOCoefY[m_Radius+i] = (m_Alpha + 2.)*vcl_abs(vcl_pow(distY, 3))
- (m_Alpha + 3.)*vcl_pow(distY, 2) + 1;
}
else
......@@ -146,21 +132,6 @@ BCOInterpolateImageFunctionBase<TInputImage, TCoordRep>
}
}
/** Constructor */
template <class TInputImage, class TCoordRep>
BCOInterpolateImageFunction<TInputImage, TCoordRep>
::BCOInterpolateImageFunction()
{
}
/** Destructor */
template <class TInputImage, class TCoordRep>
BCOInterpolateImageFunction<TInputImage, TCoordRep>
::~BCOInterpolateImageFunction()
{
}
template <class TInputImage, class TCoordRep>
void BCOInterpolateImageFunction<TInputImage, TCoordRep>
::PrintSelf(std::ostream& os, itk::Indent indent) const
......@@ -180,7 +151,7 @@ BCOInterpolateImageFunction<TInputImage, TCoordRep>
CoefContainerType BCOCoefY = CoefContainerType(winSize, 0.);
double offsetX, offsetY, distX, distY, position, step, alpha, norma;
unsigned int dim;
IndexType baseIndex;
IndexType neighIndex;
......@@ -195,22 +166,22 @@ BCOInterpolateImageFunction<TInputImage, TCoordRep>
position = - double(radius) * step;
alpha = this->GetAlpha();
for ( int i = -radius; i <= radius; i++)
{
distX = vcl_abs(position - offsetX*step);
distY = vcl_abs(position - offsetY*step);
if( distX <= 2. )
{
if (distX <= 1.)
{
BCOCoefX[radius+i] = (alpha + 2.)*vcl_abs(vcl_pow(distX, 3))
BCOCoefX[radius+i] = (alpha + 2.)*vcl_abs(vcl_pow(distX, 3))
- (alpha + 3.)*vcl_pow(distX, 2) + 1;
}
else
{
BCOCoefX[radius+i] = alpha*vcl_abs(vcl_pow(distX, 3))
BCOCoefX[radius+i] = alpha*vcl_abs(vcl_pow(distX, 3))
- 5*alpha*vcl_pow(distX, 2) + 8*alpha*vcl_abs(distX) - 4*alpha;
}
}
......@@ -223,7 +194,7 @@ BCOInterpolateImageFunction<TInputImage, TCoordRep>
{
if (distY <= 1.)
{
BCOCoefY[radius+i] = (alpha + 2.)*vcl_abs(vcl_pow(distY, 3))
BCOCoefY[radius+i] = (alpha + 2.)*vcl_abs(vcl_pow(distY, 3))
- (alpha + 3.)*vcl_pow(distY, 2) + 1;
}
else
......@@ -238,18 +209,17 @@ BCOInterpolateImageFunction<TInputImage, TCoordRep>
}
position += step;
}
// Compute base index = closet index
for( dim = 0; dim < ImageDimension; dim++ )
{
baseIndex[dim] = itk::Math::Floor< IndexValueType >( index[dim]+0.5 );
}
for( int i = -radius; i <= radius; i++ )
{
for( int j = -radius; j <= radius; j++ )
{
{
// get neighbor index
neighIndex[0] = baseIndex[0] + i;
neighIndex[1] = baseIndex[1] + j;
......@@ -262,7 +232,6 @@ BCOInterpolateImageFunction<TInputImage, TCoordRep>
{
neighIndex[0] = this->m_StartIndex[0];
}
if( neighIndex[1] > this->m_EndIndex[1] )
{
neighIndex[1] = this->m_EndIndex[1];
......@@ -278,28 +247,13 @@ BCOInterpolateImageFunction<TInputImage, TCoordRep>
}
value += lineRes[i+radius]*BCOCoefX[i+radius];
}
norma = (vcl_log(radius)/vcl_log(2.0));
norma = norma * norma;
return ( static_cast<OutputType>( value/norma ) );
}
/** Constructor */
template < typename TPixel, unsigned int VImageDimension, class TCoordRep >
BCOInterpolateImageFunction< otb::VectorImage<TPixel,VImageDimension> , TCoordRep >
::BCOInterpolateImageFunction()
{
}
/** Destructor */
template < typename TPixel, unsigned int VImageDimension, class TCoordRep >
BCOInterpolateImageFunction< otb::VectorImage<TPixel,VImageDimension> , TCoordRep >
::~BCOInterpolateImageFunction()
{
}
template < typename TPixel, unsigned int VImageDimension, class TCoordRep >
void BCOInterpolateImageFunction< otb::VectorImage<TPixel,VImageDimension> , TCoordRep >
::PrintSelf(std::ostream& os, itk::Indent indent) const
......@@ -322,7 +276,7 @@ BCOInterpolateImageFunction< otb::VectorImage<TPixel,VImageDimension> , TCoordRe
double offsetX, offsetY, distX, distY, position, step, alpha, norma;
unsigned int dim;
unsigned int componentNumber = this->GetInputImage()->GetNumberOfComponentsPerPixel();
IndexType baseIndex;
IndexType neighIndex;
......@@ -349,7 +303,7 @@ BCOInterpolateImageFunction< otb::VectorImage<TPixel,VImageDimension> , TCoordRe
OutputType output;
output.SetSize(componentNumber);
offsetX = index[0] - itk::Math::Floor<IndexValueType>(index[0]+0.5);
offsetY = index[1] - itk::Math::Floor<IndexValueType>(index[1]+0.5);
......@@ -357,22 +311,22 @@ BCOInterpolateImageFunction< otb::VectorImage<TPixel,VImageDimension> , TCoordRe
position = - double(radius) * step;
alpha = this->GetAlpha();
for ( int i = -radius; i <= radius; i++)
{
distX = vcl_abs(position - offsetX*step);
distY = vcl_abs(position - offsetY*step);
if( distX <= 2. )
{
if (distX <= 1.)
{
BCOCoefX[radius+i] = (alpha + 2.)*vcl_abs(vcl_pow(distX, 3))
BCOCoefX[radius+i] = (alpha + 2.)*vcl_abs(vcl_pow(distX, 3))
- (alpha + 3.)*vcl_pow(distX, 2) + 1;
}
else
{
BCOCoefX[radius+i] = alpha*vcl_abs(vcl_pow(distX, 3))
BCOCoefX[radius+i] = alpha*vcl_abs(vcl_pow(distX, 3))
- 5*alpha*vcl_pow(distX, 2) + 8*alpha*vcl_abs(distX) - 4*alpha;
}
}
......@@ -385,7 +339,7 @@ BCOInterpolateImageFunction< otb::VectorImage<TPixel,VImageDimension> , TCoordRe
{
if (distY <= 1.)
{
BCOCoefY[radius+i] = (alpha + 2.)*vcl_abs(vcl_pow(distY, 3))
BCOCoefY[radius+i] = (alpha + 2.)*vcl_abs(vcl_pow(distY, 3))
- (alpha + 3.)*vcl_pow(distY, 2) + 1;
}
else
......@@ -405,12 +359,12 @@ BCOInterpolateImageFunction< otb::VectorImage<TPixel,VImageDimension> , TCoordRe
for( dim = 0; dim < ImageDimension; dim++ )
{
baseIndex[dim] = itk::Math::Floor< IndexValueType >( index[dim]+0.5 );
}
}
for( int i = -radius; i <= radius; i++ )
{
for( int j = -radius; j <= radius; j++ )
{
{
// get neighbor index
neighIndex[0] = baseIndex[0] + i;
neighIndex[1] = baseIndex[1] + j;
......@@ -423,7 +377,6 @@ BCOInterpolateImageFunction< otb::VectorImage<TPixel,VImageDimension> , TCoordRe
{
neighIndex[0] = this->m_StartIndex[0];
}
if( neighIndex[1] > this->m_EndIndex[1] )
{
neighIndex[1] = this->m_EndIndex[1];
......@@ -445,7 +398,7 @@ BCOInterpolateImageFunction< otb::VectorImage<TPixel,VImageDimension> , TCoordRe
value.at(k) += lineRes[i+radius].at(k)*BCOCoefX[i+radius];
}
}
norma = (vcl_log(radius)/vcl_log(2.0));
norma = norma * norma;
......@@ -457,7 +410,6 @@ BCOInterpolateImageFunction< otb::VectorImage<TPixel,VImageDimension> , TCoordRe
return ( output );
}
} //namespace otb
#endif
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