FastCubicIpol.hpp

//
// $Id: FastCubicIpol.hpp,v 1.1 2006/07/07 14:10:53 werner Exp $
//
// $Log: FastCubicIpol.hpp,v $
// Revision 1.1  2006/07/07 14:10:53  werner
// Intermediate state: multidimensional interpolation with on-demand
// fractional loading appears to work, but still problem with vectorfield
// interpolation.
//
// Revision 1.7  2004/08/12 16:08:45  werner
// Various improvements for implementing the tangential space
//
// Revision 1.6  2004/07/09 20:15:26  werner
// Added local alignment capability when interpolating vector fields and
// option to compute the interpolated derivative.
//
// Revision 1.5  2004/07/04 17:24:03  werner
// Recursive multidimensional Interpolation interface implemented.
//
// Revision 1.4  2004/06/18 23:29:46  werner
// Intermediate version that supports the recursive interpolation scheme.
// Need to beautify the interface, fix some indexing bug, and make it to
// work with the ViewPtr.
//
// Revision 1.3  2004/06/15 23:31:00  werner
// Support for interpolation of quaternions within the interpolation interfaces.
//
// Revision 1.2  2004/06/15 22:45:31  werner
// Enhanced the fixed array member functions.
// Using better name for interpolation template parameters.
//
// Revision 1.1  2004/06/14 22:42:00  werner
// Moved interpolation files from light++ to VecAl and Multindex from Fiber to VecAl.
//
//
#ifndef __IPOL_FASTCUBICIPOL_PP
#define __IPOL_FASTCUBICIPOL_HPP "Created 11.06.2004 22:28:21 by bzfbenge"

#include "IpolDelimiter.hpp"

namespace Fiber
{

template <class T>
class   FastCubicIpol
{
public:

        template <class Storage1D, class Limiter>
static  T interpolate(const Storage1D&Data, double t, index_t size, const Limiter&L)
        {
                // TODO: Check for Samplings.size() < 2 !
                assert(size > 2);

                if (t<0     ) return Data[0];
                if (t>=size ) return Data[size-1];

        index_t i = index_t(t); 
                t = t - floor(t); 

        //      3 x^2 - 2 x^3
                t = t*t*(3-2*t);

                return Data[i] * (1-t) + Data[i+1] * (t);
        }

        template <class Storage1D>
static  T derivative(const Storage1D&Data, double t, index_t size)
        {
        index_t i = index_t(t); 
                t = t - floor(t); 

                //      3 x^2 - 2 x^3
                //      6 x   - 6 x^2
        double  x = 6*t*(1-t);

                return (Data[i+1] - Data[i] )*x;
        }
};

} /* namespace VecAl */ 

#endif /* __IPOL_FASTCUBICIPOL_HPP */