Writing a sinc-interpolating playback object

Hello,

After a lot of reading and head-scratching, I've begun to write a sample playback engine using windowed sinc interpolation. The idea is, it will work a bit like [groove~] (taking a signal to determine playback speed), but it will be click-triggerable and have built-in ADSRs and a variety of looping options.

The ADSR and loop logic is all done, but the interpolator is giving me grief. Maybe somebody on here is savvy enough to see where I'm going wrong...

When I pass integer values into my "interp(playhead)" function, everything is hunky-dory. But when I pass fractional values I get hideous ringing and the interpolated values are WAY off.

First off, I populate a table with the right-hand wing of a sinc function. I then multiply this table with the right-hand wing of a windowing function (Blackman-Harris in my case). I'm pretty damn sure that I'm populating and windowing my table correctly. There are (SINC_POINTS + 1) lobes in it.

Here is my naive implementation so far:

#define SINC_RES 512   // resolution of each lobe of the sinc table
#define SINC_POINTS 16 // number of samples to use for interpolation
                       // on either side of the sample in question

float interp(double playhead) {

  long playheadint = floor(playhead);
  double playheadfrac = playhead - playheadint;

  double accum = 0.0;

  // Loop through every sample we're interested in...
  int i;
  for(i = -SINC_POINTS; i < SINC_POINTS + 1; i ++) {

    long thissample = playheadint + i;

    //  HACK: To handle start and end points,
    //  we just duplicate the start and end samples.
    //  POSSIBLE FIX: mirror start/end of input?
    if(thissample < 0) thissample = 0;
    if(thissample >= frames) thissample = frames - 1;

    double moo = samplebuffer[thissample];

    // Get index into sinc table
    // This could be moved outside the loop but it's left here for clarity
    long sincindex = round(playheadfrac * SINC_RES);

    // Offset index relative to centre sample
    sincindex += abs(i) * SINC_RES;

    moo *= sinctable[sincindex];
    accum += moo;
  }

  return (float)accum;
}

If anybody can help figure out where I'm going wrong, I'd be extremely grateful. Here's what I'm trying to do:

1. Conceptually, replace each input sample with a sinc function.
2. Take the fractional part of my desired 'sub-sample' and scale this so it can be used as an offset into the first (or main) lobe of the sinc function.
3. Take the relevant value from the sinc function and multiply this by the input sample (as determined by the integer part).
4. Do the same for an arbitrary number of samples on either side of the sample we're interested in, but offset the index into the sinc function accordingly; i.e. if we're looking at one sample away from the main sample, we have to offset the sinc function index by one zero-crossing.
5. Sum the results of the above table lookups.

Is this the right approach? I think my code above reflects this, but maybe I've just gone code-blind...

P.S. Maybe I should mention that the eventual goal for this is to create a vari-speed sample player. I don't think this ought to affect the results of my simple function, but perhaps there are other issues to worry about further down the line? Some of the source I've seen for resampling libraries processes the source in chunks, rather than updating on a per-sample basis. However, I suspect this may be for performance reasons. For now I'd just like a slow-but-easy-to-follow piece of code that works...

If debugging the above is too much of an ask, does anyone have a similar interpolation function they wouldn't mind sharing? One that reports the interpolated value at an arbitrary location in an input buffer?

Thanks for looking,

Chris