Fixing bug in cubic spline interpolation routine

parent 4f5a748b
......@@ -4,11 +4,13 @@
#include <ITASampleBuffer.h>
#include <cassert>
#include <vector>
// Calculate cubic spline set (second derivatives, ypp) for equidistant data
void spline_cubic_set_equidistant( const int n, const float* y, float* ypp );
void spline_cubic_set_equidistant( const int n, const float* y, std::vector< float >& ypp );
// Evaluate cubic spline interpolation point using only base values (y) and second derivatives (ypp) for equidistant data
float spline_cubic_val_equidistant( const int n, const float fX, const float* y, const float* ypp );
float spline_cubic_val_equidistant( const float fX, const float* y, const std::vector< float >& ypp );
bool CITASampleLinearInterpolation::Interpolate( const ITASampleBuffer* pInput, const int iInputLength, const int iInputStartOffset, ITASampleBuffer* pOutput, const int iOutputLength, const int iOutputOffset /*=0*/ ) const
{
......@@ -159,8 +161,8 @@ bool CITASampleCubicSplineInterpolation::Interpolate( const ITASampleBuffer* pIn
float* pOutputData = pOutput->data();
// Return value
float* pInputSecondDerivatives = new float[ iInputLength ];
spline_cubic_set_equidistant( iInputLength, pInputData, pInputSecondDerivatives ); // Natural spline interpolation
std::vector< float > vfInputSecondDerivatives( iInputLength );
spline_cubic_set_equidistant( iInputLength, pInputData, vfInputSecondDerivatives ); // Natural spline interpolation
for( int i = 0; i < iOutputLength; i++ )
{
......@@ -171,14 +173,12 @@ bool CITASampleCubicSplineInterpolation::Interpolate( const ITASampleBuffer* pIn
assert( iOutputLength + iDeltaDelay == iInputLength - iInputStartOffset );
// Spline-Segment auswerten
float fYOutput = spline_cubic_val_equidistant( iInputLength, fXInput + iInputStartOffset, pInputData, pInputSecondDerivatives );
float fYOutput = spline_cubic_val_equidistant( fXInput + iInputStartOffset, pInputData, vfInputSecondDerivatives );
pOutputData[ i + iOutputOffset ] = fYOutput;
}
delete[] pInputSecondDerivatives;
return true;
}
......@@ -246,7 +246,7 @@ bool CITASampleWindowedSincInterpolation::Interpolate( const ITASampleBuffer* pI
return true;
}
void spline_cubic_set_equidistant( const int n, const float* y, float* ypp )
void spline_cubic_set_equidistant( const int n, const float* y, std::vector< float >& ypp )
{
assert( n > 2 );
......@@ -287,13 +287,16 @@ void spline_cubic_set_equidistant( const int n, const float* y, float* ypp )
delete[] a;
}
float spline_cubic_val_equidistant( const int n, const float fX, const float* y, const float* ypp )
float spline_cubic_val_equidistant( const float fX, const float* y, const std::vector< float >& ypp )
{
assert( n > fX );
assert( ypp.size() > fX );
int iXInt = ( int ) std::floor( fX );
float fFrac = fX - iXInt;
if( fFrac < ITAConstants::EPS_F_L )
return y[ iXInt ];
float y_interp = y[ iXInt ]
+ fFrac*( ( y[ iXInt + 1 ] - y[ iXInt ] ) - ( ypp[ iXInt + 1 ] / 6.0f + ypp[ iXInt ] / 3.0f )
+ fFrac * ( 0.5f * ypp[ iXInt ] + fFrac * ( ypp[ iXInt + 1 ] - ypp[ iXInt ] ) / 6.0f ) );
......
Markdown is supported
0%
or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment