PCL
pcl::BicubicSplineInterpolation< T > Class Template Reference

Bicubic spline interpolation algorithm. More...

#include <BicubicInterpolation.h>

+ Inheritance diagram for pcl::BicubicSplineInterpolation< T >:

Public Member Functions

 BicubicSplineInterpolation (const BicubicSplineInterpolation &)=default
 
 BicubicSplineInterpolation (float clamp=__PCL_BICUBIC_SPLINE_CLAMPING_THRESHOLD)
 
 ~BicubicSplineInterpolation () override
 
float ClampingThreshold () const noexcept
 
double InterpolateVector (const double c[], double dx) const noexcept
 
void InterpolateX (double fx[], double x, double y) const noexcept
 
void InterpolateY (double fy[], double x, double y) const noexcept
 
double operator() (double x, double y) const override
 
void SetClampingThreshold (float c) noexcept
 
- Public Member Functions inherited from pcl::BicubicInterpolationBase< T >
 BicubicInterpolationBase ()=default
 
 BicubicInterpolationBase (const BicubicInterpolationBase &)=default
 
- Public Member Functions inherited from pcl::BidimensionalInterpolation< T >
 BidimensionalInterpolation ()=default
 
 BidimensionalInterpolation (const BidimensionalInterpolation &)=default
 
virtual ~BidimensionalInterpolation ()
 
const T * BeingInterpolated () const
 
double BorderFillValue () const
 
virtual void Clear ()
 
void DisableBorderFilling (bool disable=true)
 
void EnableBorderFilling (bool enable=true)
 
int Height () const
 
virtual void Initialize (const T *data, int width, int height)
 
bool IsBorderFillingEnabled () const
 
void SetBorderFillValue (double v)
 
int Width () const
 

Detailed Description

template<typename T>
class pcl::BicubicSplineInterpolation< T >

Bicubic interpolation algorithms interpolate from the nearest sixteen mapped source data items. Our implementation uses the following cubic spline function as a separable filter to perform a convolution interpolation:

f(x) = (a+2)*x^3 - (a+3)*x^2 + 1       for 0 <= x <= 1
f(x) = a*x^3 - 5*a*x^2 + 8*a*x - 4*a   for 1 < x <= 2
f(x) = 0                               otherwise

Reference: Keys, R. G. (1981), Cubic Convolution Interpolation for Digital Image Processing, IEEE Trans. Acoustics, Speech & Signal Proc., Vol. 29, pp. 1153-1160.

The a constant parameter of the cubic spline (-1 <= a < 0) controls the depth of the negative lobes of the interpolation function. In this implementation we have set a fixed value a=-1/2.

Our implementation includes an automatic linear clamping mechanism to avoid discontinuity problems when interpolating linear images. The cubic spline function has negative lobes that can cause undershoot (aka ringing) artifacts when they fall over bright pixels. This happens frequently with linear CCD images. See the documentation for the BicubicSplineInterpolation::SetClampingThreshold( float ) function for a more detailed description of the linear clamping mechanism.

Definition at line 377 of file BicubicInterpolation.h.

Constructor & Destructor Documentation

◆ BicubicSplineInterpolation() [1/2]

template<typename T >
pcl::BicubicSplineInterpolation< T >::BicubicSplineInterpolation ( float  clamp = __PCL_BICUBIC_SPLINE_CLAMPING_THRESHOLD)
inline

Constructs a BicubicSplineInterpolation instance.

The optional clamp parameter defines a threshold to trigger the linear clamping mechanism. See the documentation for the SetClampingThreshold( float ) function for a detailed description of the automatic linear clamping mechanism. The default value of clamp = 0.1 is appropriate for most linear images.

Definition at line 390 of file BicubicInterpolation.h.

◆ BicubicSplineInterpolation() [2/2]

template<typename T >
pcl::BicubicSplineInterpolation< T >::BicubicSplineInterpolation ( const BicubicSplineInterpolation< T > &  )
default

Copy constructor.

◆ ~BicubicSplineInterpolation()

template<typename T >
pcl::BicubicSplineInterpolation< T >::~BicubicSplineInterpolation ( )
inlineoverride

Virtual destructor.

Definition at line 404 of file BicubicInterpolation.h.

Member Function Documentation

◆ ClampingThreshold()

template<typename T >
float pcl::BicubicSplineInterpolation< T >::ClampingThreshold ( ) const
inlinenoexcept

Returns the current linear clamping threshold for this interpolation object.

See the documentation for SetClampingThreshold( float ) for a detailed description of the automatic linear clamping mechanism.

Definition at line 546 of file BicubicInterpolation.h.

◆ InterpolateVector()

template<typename T >
double pcl::BicubicSplineInterpolation< T >::InterpolateVector ( const double  c[],
double  dx 
) const
inlinenoexcept

Vector interpolation.

Parameters
cAddress of the first element of a vector of four data points.
dxInterpolation point in the range [0,+1[. A value of zero corresponds to the second element in the f vector.

This function interpolates four neighbor rows or columns by direct convolution with the cubic spline function. The four elements of the c vector can be raw data points, to use cubic spline interpolation in one dimension, or they can be interpolated rows or columns from an arbitrary 2-D matrix. For example, c can be generated by a previous call to InterpolateX() or InterpolateY(), respectively to provide source interpolated rows or columns, or by equivalent functions from a different separable cubic interpolation algorithm. This is useful to mix cubic spline interpolation with other interpolation algorithms, which allows for flexible interpolation schemes.

Definition at line 532 of file BicubicInterpolation.h.

◆ InterpolateX()

template<typename T >
void pcl::BicubicSplineInterpolation< T >::InterpolateX ( double  fx[],
double  x,
double  y 
) const
inlinenoexcept

Horizontal interpolation.

Parameters
x,yCoordinates of the interpolation point (horizontal, vertical).
[out]fxAddress of the first element of a vector where the four interpolated X-values will be stored.

This function interpolates horizontally the four neighbor rows for the specified x, y coordinates. This is useful when this interpolation algorithm is being used to interpolate an image on the horizontal axis exclusively.

Definition at line 457 of file BicubicInterpolation.h.

◆ InterpolateY()

template<typename T >
void pcl::BicubicSplineInterpolation< T >::InterpolateY ( double  fy[],
double  x,
double  y 
) const
inlinenoexcept

Vertical interpolation.

Parameters
x,yCoordinates of the interpolation point (horizontal, vertical).
[out]fyAddress of the first element of a vector where the four interpolated Y-values will be stored.

This function interpolates vertically the four neighbor columns for the specified x, y coordinates. This is useful when this interpolation algorithm is being used to interpolate an image on the vertical axis exclusively.

Definition at line 492 of file BicubicInterpolation.h.

◆ operator()()

template<typename T >
double pcl::BicubicSplineInterpolation< T >::operator() ( double  x,
double  y 
) const
inlineoverridevirtual

Returns an interpolated value at {x,y} location.

Parameters
x,yCoordinates of the interpolation point (horizontal, vertical).

This function performs a two-dimensional interpolation via a convolution with the separable cubic spline filter.

Implements pcl::BidimensionalInterpolation< T >.

Definition at line 417 of file BicubicInterpolation.h.

◆ SetClampingThreshold()

template<typename T >
void pcl::BicubicSplineInterpolation< T >::SetClampingThreshold ( float  c)
inlinenoexcept

Defines a threshold to trigger the linear clamping mechanism. The clamping mechanism automatically switches to linear interpolation when the differences between neighbor pixels are so large that the cubic interpolation function may cause ringing artifacts. Ringing occurs when the negative lobes of the cubic spline interpolation function fall over isolated bright source pixels or bright edges. For example, ringing problems happen frequently around stars in linear CCD images. For nonlinear images, linear clamping is almost never necessary, as in nonlinear images (e.g., stretched deep-sky images, or diurnal images) there are normally no such large variations between neighbor pixels.

Linear clamping works on a per-row or per-column basis within the interpolation neighborhood of 16 pixels. This means that when the clamping mechanism selects linear interpolation, it restricts its use to the affected (by too strong variations) row or column, without changing the bicubic interpolation scheme as a whole. This ensures artifact-free interpolated images without degrading the overall performance of bicubic spline interpolation.

The specified clamping threshold clamp must be in the [0,1] range. Lower values cause a more aggressive deringing effect.

Definition at line 574 of file BicubicInterpolation.h.

References pcl::Range().


The documentation for this class was generated from the following file: