Band Limited Table Look Up Oscillator C++

Question

I would like to create a Band Limited Table Lookup Oscillator that functions much like the sin() function in C++ works where i would inout a phase value and receive a sample value. for example Sample s = saw(x); or Sample s = square(x);.

To fill the table I would first generate an oversampled table by a factor of "X" that contains the wave in a non band limited state. Then I would run the table through a filter of some sort FIR or IIR most likely. Then i would decimate the table by a factor of "X" and hopefully then have a table that i could read from that would produce a band limited signal. I would also have to interpolate across the table to ensure that all phase angles will be accommodated for.

I would like to know if this theory would work? Meaning would this produce a signal that is Band limited.

I will mock up the basic code below.

double saw(double phase){
    return interpolate(phase);//i know this is incomplete
}

double fillTable(){
    double table[size*oversampling];

    //fill table here. skipping for length of example

    for(int i =0;i<tableSize;++i){
        table[i]=filter(table[i]);
    }
    double table2[size];
    int x=0;
    for(int i =0;i<tableSize;++i){
        table2[i]=table[x];
        x+=oversampling;
    }
}

//usage
sample x = saw(phase);

I know the code is incomplete and syntactically incorrect but i just wanted to give an example of what i would like to be able to do.

So again. Will this be Band limited?

Answer 1

Please note that you cannot implement a band-limited oscillator with a pure (stateless) function like:

Sample s = saw(x);

The reason is that the function must be aware of the rate at which x change to generate a signal with the correct bandwidth

For example, let us assume your sample rate is 48kHz. If you call your function in this order:

saw(0.00000000);
saw(0.00000001);
saw(0.00000002);

You can see that the argument is slowly changing to the point that aliasing is not an issue, so your function could return something like x / M_PI - 1 (assuming input is in [0, 2 pi[ ; and output in [-1, 1]) without any lookup table. But now, if you call it this way:

saw(0.0);
saw(3.0);
saw(6.0);

We are getting close to Nyquist frequency and the output must then be heavily band-limited - actually closer to sin.

So a band-limited sawtooth or square cannot be generated by mere waveshaping of the naive saw. You need your function to be aware of the target frequency so that it can make the right decision.

Another related misconception that can be found in your code is that you assume that you would need only one table. There is just not one single "band-limited sawtooth" waveform that would work for all frequencies. For low frequencies, there is very little different between the band-limited saw and the naive saw. For high frequencies, the band-limited saw is effectively a sine. So you have to define how much memory you want to dedicate to your task, and generate many band-limited tables corresponding to various frequencies. At run-time, you decide which one to use depending on the frequency of the signal. If this is for a musical application in which you don't want to hear the transition from one table to the next as notes rise in pitch, you'll need to add one further layer of interpolation on top of that to crossfade from one table to the next. A typical implementation would use one table by octave or per 8 semitones.

Answer 2

Note that if the function with which you fill the table has spectral content at or above the Nyquist frequency for even the longer table length, then there may be (likely will be) aliasing noise that subsequent low pass filtering will not remove. Perhaps a lesser amount of aliasing noise than with just using a shorter table, but this depends on the table filling function used.

If you plan on using this table to produce a periodic oscillation, you will have to wrap the filtering around the table circularly to account for the impulse response of the low pass filter crossing period boundaries, or else you will end up with filter truncation noise.

A better approach might be to calculate the Fourier coefficients of the periodic function with which you plan of filling the table, and use those coefficients to synthesize the table by summing just those sinusoidal components (exactly periodic in the table length) that fall below the Nyquist frequency for the table length.