   // this initialization corresponds to mirror boundaries
   // accelerated loop
   float zn = IU_BSPLINE_POLE;
-  floatN Sum = *c;
+  PixelType Sum = *c;
   for (unsigned int n = 1; n < Horizon; n++) {
     c += step;
     Sum += zn * *c;


#define IU_BSPLINE_POLE (sqrtf(3.0f)-2.0f)  // pole for cubic b-spline

template<class floatN> __device__ floatN InitialCausalCoefficient(
  floatN* c,        // coefficients
  unsigned int DataLength,  // number of coefficients
  int step)
{

  return(Sum);
}

template<class PixelType> __device__ PixelType InitialAntiCausalCoefficient(
  PixelType* c,        // last coefficient
  unsigned int DataLength,  // number of samples or coefficients
  int step)
{
  // this initialization corresponds to mirror boundaries
  return((IU_BSPLINE_POLE / (IU_BSPLINE_POLE * IU_BSPLINE_POLE - 1.0f)) *
         (IU_BSPLINE_POLE * c[-step] + *c));
}

//-----------------------------------------------------------------------------
template<class PixelType> __device__ void ConvertToInterpolationCoefficients(
  PixelType* coeffs,   // input samples --> output coefficients
  unsigned int DataLength,  // number of samples or coefficients
  int step)
{
  // compute the overall gain
  const float lambda = (1.0f - IU_BSPLINE_POLE) *
      (1.0f - 1.0f / IU_BSPLINE_POLE);

  // causal initialization
  PixelType* c = coeffs;
  PixelType previous_c;  //cache the previously calculated c rather than look it up again (faster!)
  *c = previous_c = lambda * InitialCausalCoefficient(c, DataLength, step);
  // causal recursion
  for (unsigned int n = 1; n < DataLength; ++n) {
    c += step;
    *c = previous_c = lambda * (*c) + IU_BSPLINE_POLE * previous_c;
  }

  // anticausal initialization
