Band-limited wavetable playback via sinc interpolation and mipmapping

1

hey, ive recreated the band-limited wavetable approach from the go book.
Its using windowed sinc interpolation and mipmapping.
The result is pretty good anti-aliased even for high index FM.

// prepare buffers
(
var getSinc = { |numRipples|
    var x = Array.interpolation(4096, -1.0, 1.0);
    sincPi(x * 0.5pi * numRipples);
};
var sinc = getSinc.(8);
var kaiser = Signal.kaiserWindow(4096);
var windowedSinc = sinc * kaiser;
~sincBuf = Buffer.loadCollection(s, windowedSinc);
)

(
~sndBuf = Buffer.alloc(s, 2048, 1, { |buf|
    buf.sine1Msg((1..512).reciprocal, asWavetable: false)
});
)

// run example!
(
var rampToSlope = { |phase|
    var history = Delay1.ar(phase);
    var delta = (phase - history);
    delta.wrap(-0.5, 0.5);
};

var sincInterpolated = { |phase, sndBuf, sincBuf, sampleSpacing, numSamples = 8|

    var sampleIndex, fracPart, intPart;
    var samples, windows;
    var halfSamples = numSamples.div(2);

    // Calculate the base sample position
    sampleIndex = phase * BufFrames.kr(sndBuf) / sampleSpacing;
    fracPart = sampleIndex.wrap(0, 1);  // fractional part for interpolation
    intPart = sampleIndex - fracPart;   // integer part for base position

    // Get sample values at integer offsets
    samples = (halfSamples.neg..(halfSamples - 1)).collect{ |offset|
        var readPos = ((intPart + offset) * sampleSpacing).wrap(0, BufFrames.kr(sndBuf));
        BufRd.ar(1, sndBuf, readPos, interpolation: 1);
    };

    // Get window values from sinc function
    windows = (1..numSamples).collect{ |i|
        var sincPos = (i / numSamples) - (fracPart / numSamples) * BufFrames.kr(sincBuf);
        BufRd.ar(1, sincBuf, sincPos, interpolation: 4);
    };

    (samples * windows).sum;
};

var mipmapInterpolate = { |phase, sndBuf, sincBuf|

    var slope, samplesPerFrame, octave, layer;
    var spacing1, spacing2, sig1, sig2;

	// Calculate mipmap parameters
    slope = rampToSlope.(phase);
    samplesPerFrame = slope.abs * BufFrames.kr(sndBuf);
    octave = max(0, log2(samplesPerFrame));
    layer = octave.ceil;

	// Calculate spacings for adjacent mipmap levels
    spacing1 = 2 ** layer;
    spacing2 = 2 ** (layer + 1);

	// Get and crossfade interpolated signals
    sig1 = sincInterpolated.(phase, sndBuf, sincBuf, spacing1);
    sig2 = sincInterpolated.(phase, sndBuf, sincBuf, spacing2);
    LinXFade2.ar(sig1, sig2, octave.wrap(0, 1) * 2 - 1);
};

{
    var rate, phase, sig;

	rate = SinOsc.ar(0.1, 1.5pi).linlin(-1, 1, 100, 8000);
    phase = Phasor.ar(DC.ar(0), rate * SampleDur.ir);
    sig = mipmapInterpolate.(phase, ~sndBuf, ~sincBuf);

    sig = LeakDC.ar(sig);
    sig!2 * 0.1;
}.play;
)

// check the freqscope!
s.freqscope;
6 Likes
2

will add the 1-D/2-D/3-D/wave terrain morphing capabilities in the upcoming days, which could be a blueprint for a C++ implementation. Will learn more C++ in the upcoming year, i promise haha :slight_smile:

2 Likes
3

Very nice! However, in the Oversampling Oscillators library there is already a 3d implementation that you may want to look at - OscOS3. This is a real kitchen sink oscillator (I don’t know what else you would want it to do) - 3d wavetable loading, mipmap building, oversampling, triggered sync, variable phasor shape.

I don’t have sinc interpolation, but you have inspired me to add it! Will try to do so in the coming weeks.

Sam

2 Likes
4

Hey, thanks! oversampling/downsampling is adding a small delay, i think this should therefore be the last instance of anti-aliasing. Additionally I dont like the wavetable format, im just using BufRd for everything buffer based.

I have not anything specific in mind, i just share what i come up with.

5

I also don’t like the wavetable format, so I don’t use it. Normal buffers are used in these objects. The wavetable format is only for interpolation, and since you are using sinc interpolation, you don’t need it anyway.

You are also welcome to not use oversampling in my objects! Set the oversampling index to 0 and the up/down sampling is bypassed.

So I think we agree here.

1 Like
6

i think what would be cool is to implement the sinc interpolation in the c++ code and keep the oversampling optional. The sinc interpolation (8 ripples + kaiser window) + mipmapping is already doing a good job of anti-aliasing. Havent done any benchmark but could be the case that its not as costly on CPU.

but then its not anti-aliased

7

It is if you make a mipmap first. Then it is just a sine wave in the top octave. You build a multichannel buffer which contains the mipmap and play back from that.

8

thanks, interesting!

9

How significant of a delay are you measuring?

I’m getting about 2 samples AFAICS, which I wouldn’t lose sleep over (though it may be critical in some applications, it wouldn’t be critical to me).

b = Buffer.alloc(s, 2048, 1, { |buf| buf.sine1Msg((1 .. 200).reciprocal, asWavetable: false) });

a = {
	var trig = Impulse.ar(0);
	var phase = Phasor.ar(0, 440 * SampleDur.ir, 0, 1);
	var sig = OscOS.ar(b,
		phase,
		bufdivs: 1, bufloc: 0, oversample: 8
	);
	var sig2 = BufRd.ar(1, b, phase * BufFrames.ir(b), interpolation: 0);
	[trig, sig, sig2]
}.plot(duration: 0.002);

hjh

10

hey, i havent measured any delay i just know that oversampling always causes a delay.
In my use case for granulation it does matter :slight_smile:

Unfortunately, you also have a delay when you disable the oversampling. I dont know if thats intented.

b = Buffer.alloc(s, 2048, 1, { |buf| buf.sine1Msg((1 .. 200).reciprocal, asWavetable: false) });

(
a = {
	var trig = Impulse.ar(0);
	var phase = Phasor.ar(0, 440 * SampleDur.ir, 0, 1);
	var sig = OscOS.ar(b,
		phase,
		buf_divs: 1, buf_loc: 0, oversample: 0
	);
	var sig2 = BufRd.ar(1, b, phase * BufFrames.ir(b), interpolation: 0);
	[trig, sig, sig2]
}.plot(duration: 0.002);
)
11

AFAICS the delay is constant, and independent of oscillator frequency.

I’ve tried with the band-limited saw buffer, and also with a Dirac impulse in the buffer (positioned halfway through, and offsetting the phase by half a cycle). In both cases, at different oscillator frequencies, OscOS generates the peak in the same place (while BufRd “smears” that corner at lower frequencies).

So (I guess) the effect for granulation would be that the phase position of the grain’s waveform would be very slightly differently aligned with the window (comparing OscOS and BufRd), but I wouldn’t expect additional phase cancellation coming from the oversampling. You might even get less phase cancellation from OscOS – here, the anti-aliasing wiggly bits are all in phase (in that their positive vs negative samples are all aligned – however the inter-sample peaks would reflect the relative oscillator frequencies).

oversampling-phase
oversampling-phase840Ă—524 30.1 KB

hjh

12

Thanks. Im not really happy with the phase offset between carrier and window, but maybe not that big of a deal. Ideally the oversampling could be disabled for no delay.

This should be correct. The Delay just depends on the Order of the Filter. Im not an expert here, but currently investigating the topic in greater Detail and picked up some Information Reading some threads on the gen Discord. Oversampling requires up and downsampling filters of typically tens of samples. Blep/blamp kernels imply latency too. There are minimum-phase/IIR variant options for anti-aliasing with lower latency, but they change the phase structure, which means the latency is not a fixed amount but varies by frequency.

The rules of thumb for AA are: 1.) if you can change your algorithm to not generate audible frequencies beyond available bandwidth, focus on that. 2.) if you can isolate aliasing to sporadic moments (e.g. corners of square/saw etc. waves) focus on those (e.g. blep, blamp etc.) 3.) only after 1 and 2 add oversampling.

The sinc interpolation + midmapping does a good job on 1.) without adding any delay from a biquad which is used here for up and downsampling.
When oversampling is disabled, there should be no delay.

13

I have compared OscOS with the sincInterpolation. The cutoff could be steeper, but Its doing an even better job on anti-aliasing then 16x oversampling.

(
var getSinc = { |numRipples|
    var x = Array.interpolation(4096, -1.0, 1.0);
    sincPi(x * 0.5pi * numRipples);
};
var sinc = getSinc.(8);
var kaiser = Signal.kaiserWindow(4096);
var windowedSinc = sinc * kaiser;
~sincBuf = Buffer.loadCollection(s, windowedSinc);
)

(
~sndBuf = Buffer.alloc(s, 2048, 1, { |buf|
    buf.sine1Msg((1..512).reciprocal, asWavetable: false)
});
)

(
var rampToSlope = { |phase|
    var history = Delay1.ar(phase);
    var delta = (phase - history);
    delta.wrap(-0.5, 0.5);
};

var sincInterpolated = { |phase, sndBuf, sincBuf, sampleSpacing, numSamples = 8|

    var sampleIndex, fracPart, intPart;
    var samples, windows;
    var halfSamples = numSamples.div(2);

    // Calculate the base sample position
    sampleIndex = phase * BufFrames.kr(sndBuf) / sampleSpacing;
    fracPart = sampleIndex.wrap(0, 1);  // fractional part for interpolation
    intPart = sampleIndex - fracPart;   // integer part for base position

    // Get sample values at integer offsets
    samples = (halfSamples.neg..(halfSamples - 1)).collect{ |offset|
        var readPos = ((intPart + offset) * sampleSpacing).wrap(0, BufFrames.kr(sndBuf));
        BufRd.ar(1, sndBuf, readPos, interpolation: 1);
    };

    // Get window values from sinc function
    windows = (1..numSamples).collect{ |i|
        var sincPos = (i / numSamples) - (fracPart / numSamples) * BufFrames.kr(sincBuf);
        BufRd.ar(1, sincBuf, sincPos, interpolation: 4);
    };

    (samples * windows).sum;
};

var mipmapInterpolate = { |phase, sndBuf, sincBuf|

    var slope, samplesPerFrame, octave, layer;
    var spacing1, spacing2, sig1, sig2;

	// Calculate mipmap parameters
    slope = rampToSlope.(phase);
    samplesPerFrame = slope.abs * BufFrames.kr(sndBuf);
    octave = max(0, log2(samplesPerFrame));
    layer = octave.ceil;

	// Calculate spacings for adjacent mipmap levels
    spacing1 = 2 ** layer;
    spacing2 = 2 ** (layer + 1);

	// Get and crossfade interpolated signals
    sig1 = sincInterpolated.(phase, sndBuf, sincBuf, spacing1);
    sig2 = sincInterpolated.(phase, sndBuf, sincBuf, spacing2);
    LinXFade2.ar(sig1, sig2, octave.wrap(0, 1) * 2 - 1);
};

{
	var freq = SinOsc.ar(0.1, 1.5pi).linlin(-1, 1, 100, 8000);
	var phase = Phasor.ar(DC.ar(0), freq * SampleDur.ir);
	var osc = OscOS.ar(~sndBuf, phase, buf_divs: 1, buf_loc: 0, oversample: 4);
	var sinc = mipmapInterpolate.(phase, ~sndBuf, ~sincBuf);

	// compare OscOS with sincInterpolation
	//var sig = osc!2 * 0.1;
	var sig = sinc!2 * 0.1;
	LeakDC.ar(sig);
}.play;
)

// check the freqscope!
s.freqscope;