How do you interpolate/avoid distortion/aliasing when changing the frequency of the oscillator?

Hi Folks!

I have a mystical and maybe a stupid question. But after a couple of years of developing my synthesizer, when I came to the functional that changes the frequency of the oscillator (pitch, fine, semi).

JFYI: Sample is already bandlimited, so there is no aliasing. All harmonics removed

I was surprised to find that there are distortion when it moves. Although I do a pretty good interpolation. I use Hermite interpolation (with 4 dots) And as you can see it’s really sound and looks clean on analyzer. But it’s clean only when the frequency is static and not moving.

This is how it looks with static frequency:

This is when frequency is moving up slowly. Which seems to be looks fine.

But here is what happening if it will be 2khz or higher. You can see there is some spike’s

I understand that there is a slight phase shift and it seems to be torn relative to the previous delta. But I don’t understand how this can be avoided? I tried to transfer to the interpolation function previous two dots from the buffer (not from the sample), but the result is worse.

What is the practice in this case? Or am I worried in vain and this is inevitable? Although it seemed to me that these spikes in other synthesizers are either completely absent or less noticeable (perhaps because they do not update each sample frequency like I’m doing). But even if I’m updating frequency each 1024 or 2048 sample, then I still see this spikes.

Another strange thing, that this spikes is absent if I’m using large buffer: 1024 or 2048. But they now appear if I’m playing from 100 hz to 1 khz :frowning:

But once I’m using buffer size: 32. It is significantly noticeable:

Second question. How often do you update the frequency? each sample, each 4-8 sample or by buffer size? What is the common practice here?

Thank you!!!

Presumably you’ve implemented band-limited wavetable? There could be an issue with how you’re dealing with the band-limiting. Are you using mip-maps for your wavetables at semitone or octave intervals or are you generating the band-limited wavetable using inverse FFT every ‘n’ samples based on oscillator frequency?

Are you still seeing the issue if you change your wavetable to a sine wavetable? The artefacts don’t look like aliasing but they do look like there could be an issue with your Hermite interpolation. Are you making sure that the points you select are handling wrap around correctly in the wavetable?

In any case, the artefacts are well below any audio range but I appreciate it’s annoying when you’re going for clean.

As to your second question. I update the frequency every sample but you need to make sure your frequency calculations are optimal. Using the standard c++ library log/exp etc won’t cut it performance wise. You could probably calculate every 4, 8 or 16 samples but even then you should interpolate in-between.

1 Like

Thanks for reply @Nitsuj70

Presumably you’ve implemented band-limited wavetable?

Yes, correct and I even use Double instead of Float. it doesn’t matter, but still )

Are you using mip-maps for your wavetables at semitone or octave intervals or are you generating the band-limited wavetable using inverse FFT every ‘n’ samples based on oscillator frequency?

Actually, I tried both: 11 octaves and FFT inverse each N sample (Basically each 4 sample)

Are you still seeing the issue if you change your wavetable to a sine wavetable?

Yes, with a Sine wave the same and even more noticeable spikes.

The artefacts don’t look like aliasing but they do look like there could be an issue with your Hermite interpolation

Correct. it’s not seems to be aliasing. I see the same spikes when I’m not interpolating or using linear even if the frequency stays the same/static. So I assume I’m using wrong Hermite interpolation when I’m changing frequency increments.

but even then you should interpolate in-between.

That’s the moment which I’m afraid I don’t understand. For now, I’m just interpolating the sample by double phasePosition with range [0.0 - 2048]. But it looks like interpolation is still needed between how I write data to the output buffer and the increment points?

In any case, the artefacts are well below any audio range but I appreciate it’s annoying when you’re going for clean.

yes, it’s a matter of principle here) I didn’t check other synths a lot, not many of them seems to be updating each sample, so it’s hard to say how they deal with it, but it seemed to me that they have less spikes

The most strange thing is why this spikes is missing when I’m using larger buffer.
I know someone will probably say that this is a bug in my code, but no. I even tested it on a new empty project with simple few lines code: std::sine(…); increment+= 0.0001;
By the way, it is likely that the analyzers like Fabfilter or Ableton Spectrum does not have enough data with buffer size: 32 or 64 to show with hi resolution. Hmm…

Doing the FFT inverse every 4 samples is going to melt your CPU. At an absolute minimum I’d go for ~44 samples which is every 1ms if you’re at 44100Hz. Even that will be a CPU hit and you may want to set it to be less frequent. Samples intervals of 64, 128 and 256 aren’t uncommon.

Here’s how I’m handling the phasor. The phasor is a float [0…1] where 0 is inclusive and 1 is exclusive. I convert it to a wavetable position int [0…2047] ‘p’, and the fractional ‘f’. Your hermite points are the wavetable values at p-1, p, p+1 and p+2.

Then I just call: sample = hermite(wt[p-1], wt[p], wt[p+1], wt[p+2], f);

Where ‘wt’ is your wavetable. You must make sure that various offsets of p don’t extend before the wavetable start and end. You can either use a mask on the indexes for example (p-1)&2047 or you can have a slightly larger wavetable and copy wrapped values to the start and end.

FWIW, here’s the hermite function I use:

template <typename T>
inline T hermite(T xm1, T x0, T x1, T x2, T frac)
{
    const T c = (x1 - xm1) * T(0.5);
    const T v = x0 - x1;
    const T w = c + v;
    const T a = w + v + (x2 - x0) * T(0.5);
    const T b_neg = w + a;

    return ((((a * frac) - b_neg) * frac + c) * frac + x0);
}

Then I just increase the phasor and make sure it wraps to stay in the [0…1] range where 0 is inclusive and 1 is exclusive (i.e. if (phasor >= 1.f) phasor -= 1.f; )

1 Like

Thank you @Nitsuj70 looks like we are doing same thing :slight_smile:

Yes, agree 4 samples is too much, but it’s just for the test, on the iPad A12X it’s not even a problem. I found that having 11 octaves pretty enough if I have a steep cut of Chebyshev type II filter (120db and 48 order) and 4x oversampling.Even too much, but for visual pleasure quite :slight_smile:

{
  DoubleType fmPhase[2]{DoubleType(0.0), DoubleType(0.0)};
  
  if (hasFM && assignedFMSource != nullptr) {
      const auto fmDelta = assignedFMSource->outputBuffer_[0][i];
      fmPhase[0] = fmDelta * state.fmAmountA;
      fmPhase[1] = fmDelta * state.fmAmountB;
  }
  
  // Sound
  const auto soundA = HermiteInterpolation::interpolate<PhaseLengthMask>(state.sampleA, state.phase + fmPhase[0]);
  
  // Cross-fade between sample A-B
  if (state.transitionFrac < 0.0) {
      outputBuffer_[0][i] = soundA;
  } else {
      const auto soundB = HermiteInterpolation::interpolate<PhaseLengthMask>(state.sampleB, state.phase + fmPhase[1]);
      outputBuffer_[0][i] = soundA * (DoubleType(1.0) - state.transitionFrac) + soundB * state.transitionFrac;
  }
  
  // Set frequency and move phase
  state.phase += state.increment * frequency;
  }

And here is interpolation implementation:

class HermiteInterpolation {
    public:

        template<unsigned long PhaseLengthMask>
        FORCE_INLINE static DoubleType interpolate(const DoubleType sample[], const DoubleType pos) noexcept {
            const int posInt = int(pos);
            const int pos0 = (posInt - 1) & PhaseLengthMask;
            const int pos1 = posInt & PhaseLengthMask;
            const int pos2 = (pos1 + 1) & PhaseLengthMask;
            const int pos3 = (pos1 + 2) & PhaseLengthMask;
            const DoubleType frac = pos - DoubleType(posInt);  // Warning! do not use here std::abs, value should be negative as well.

            return hermiteInterpolation(sample[pos0], sample[pos1], sample[pos2], sample[pos3], frac);
        }

    private:

        static FORCE_INLINE DoubleType hermiteInterpolation(DoubleType value0, DoubleType value1, DoubleType value2, DoubleType value3, DoubleType sampleOffset) {
            DoubleType slope0 = (value2 - value0) * DoubleType(0.5);
            DoubleType slope1 = (value3 - value1) * DoubleType(0.5);
            DoubleType v = value1 - value2;
            DoubleType w = slope0 + v;
            DoubleType a = w + v + slope1;
            DoubleType bNeg = w + a;
            DoubleType stage1 = a * sampleOffset - bNeg;
            DoubleType stage2 = stage1 * sampleOffset + slope0;
            return DoubleType(stage2 * sampleOffset + value1);
        }

    };

I can’t see anything wrong with the code you’ve posted. I’d recommend wrapping state.phase so that it doesn’t run away over time, but other than that it looks okay.

I’m guessing that state.sampleA and state.sampleB are two wavetables but I’m not sure what the scheme is here, for example, why you’ve got two fmPhases or why state.transitionFrac could be lower than 0.0.

In my code I have two wavetables and implement double buffering. Let’s call these wavetables ‘src’ and ‘dst’. Based on the current frequency, I generate the ‘dst’ band-limited wavetable using inverse FFT transform. Then, over the next 64 samples or so I interpolate between ‘src’ and ‘dst’ samples. At the end of 64 samples I swap the ‘src’ and ‘dst’ buffer pointers and re-trigger the band-limited wavetable generation into ‘dst’. So the code is always transitioning between ‘src’ and ‘dst’ wavetables.

This way, any transition to wavetable for changing frequency remains smooth.

Is this what your cross-fade A-B is doing in the code you provided? Are state.sampleA and state.sampleB pointers to wavetable buffers the same way that ‘src’ and ‘dst’ are in my description above? If that’s the case then state.transitionFrac should never be below 0 or above 1.

Unless…you’re using state.transitionFrac < 0 to indicate that you don’t need to transition any more? If that’s the case then make sure that state.sampleA is pointing at your transitioned wavetable when your state.transitionFrac < 0.

1 Like

I’m doing this at the end of loop, to save a cpu resources.

void process(...) {
  for (loop) {
    ...
  }
  
  state.phase = state.phase - PhaseLengthDouble * int(state.phase * PhaseLengthDoubleMult); // same as std::fmod
}

I’m guessing that state.sampleA and state.sampleB are two wavetables but I’m not sure what the scheme is here, for example, why you’ve got two fmPhases or why state.transitionFrac could be lower than 0.0.

Two fmAmount and it’s because I need smooth crossfade between FM modulation. For example if LFO/or user set FM Amount to 1.0 and it was before 0.0 or 0.5, then Sine at the low notes may crackle due to abrupt transition of amplitude.

Same actually for stateA and stateB it’s for morphing between two different waves in the wavetable. If user or LFO changing sample, then to avoid crackles I’m doing crossfade for about 0.005 or 0.007 seconds (it’s transition about 220.5 sample at 44100 sample rate). These crackles are especially audible if it is, for example, a transition from a wave where the phase significantly changes its amplitude: From Sine to Square or Square to Triangle, etc.

the state.transitionFrac < 0 in my code mean’s NoTransitionRequired, sorry for confusion. state.transitionFrac can be -2.0, -1.0 or 0.0-1.0.

-2.0: Sample need to be initialized. although I’ll move it to noteOn or can use simply sampleA == nullptr.
-1.0: No transition. i.e. I’m using just sampleA from state.

That’s just to not to create extra bool values like bool transitionRequired, sampleInitialized, since I’m also copying osc State from Heap memory to the Stack before main loop. That is, I squeeze it to the maximum and save processor resources.

Here is code. It looks a little ugly for now, but I’m still testing, in the future of course it will be the norm:

struct OscState {
        DoubleType increment{};
        DoubleType phase{};
        DoubleType morphRatioA{};
        DoubleType fmAmountA{};
        DoubleType morphRatioB{};
        DoubleType fmAmountB{};
        DoubleType transitionFrac{-2.0}; // -2 sample init -1 no transition
        DoubleType *sampleA{};
        DoubleType *sampleB{};
    } JUCE_PACKED;

    OscState state = state_; // copy on stack

    for (loop){
       ...
    }

    state.phase = state.phase - PhaseLengthDouble * int(state.phase * PhaseLengthDoubleMult); // same as std::fmod
            state_ = state; // copy back to the HEAP memory

Then, over the next 64 samples or so I interpolate between ‘src’ and ‘dst’ samples. At the end of 64 samples I swap the ‘src’ and ‘dst’ buffer pointers and re-trigger the band-limited wavetable generation into ‘dst’

When you say you are interpolating between src and dst. Are you doing the same as me with transitionFrac ? or do you also do an additional hermite interpolation somehow between these two waves? Or just linear transition between two samples like me? soundA * a + soundB * b

I’m also do swap between samples at the end of transition. Simply changing pointers. But I also check that the values in B state did not receive new values until the transition ended.
So I have something like stateA, stateB, stateReal.

auto newOctave = ...
auto fmAmount = ...
auto morphRatio = ...

// Init sample TODO: Move this to the noteOn()
if (state.transitionFrac < -1.0) {
    state.transitionFrac = -1.0;
    state.morphRatioA = morphRatio;
    state.fmAmountA = fmAmount;
    state.sampleA = oscConfiguration_.waveTable.getSample(newOctave, morphRatio);
    state.sampleB = oscConfiguration_.waveTable.getSample(newOctave, morphRatio);
    currentOctave = newOctave;
} else {

    // Has changes?
    if (state.transitionFrac < 0.0) {

        const bool hasMorphChanges = morphRatio != state.morphRatioA;
        const bool hasFMChanges = fmAmount != state.fmAmountA;

        if (hasMorphChanges || hasFMChanges) {
            state.morphRatioB = morphRatio;
            state.fmAmountB = fmAmount;
            state.transitionFrac = 0.0; // initiate transition process
            currentOctave = newOctave;

            if (hasMorphChanges) {
                state.sampleB = oscConfiguration_.waveTable.getSample(newOctave, morphRatio);
            }
        }

    }
}

// Transition in progress
if (state.transitionFrac > -1.0) {
    state.transitionFrac += transitionIncrement;
    if (state.transitionFrac >= 1.0) {
        state.transitionFrac = -1.0; // end of transition
        state.morphRatioA = state.morphRatioB;
        state.fmAmountA = state.fmAmountB;
        std::swap(state.sampleA, state.sampleB);
    }
}

if (currentOctave != newOctave) {
    currentOctave = newOctave;
    state.sampleA = oscConfiguration_.waveTable.getSample(newOctave, state.morphRatioA);
    state.sampleB = oscConfiguration_.waveTable.getSample(newOctave, state.morphRatioB);
}

Ah, I see. In principle I’m doing something similar but not quite the same.

So, as stated I have ‘src’ and ‘dst’ (your state.sampleA and state.sampleB). I’m always transitioning by interpolating between them using linear interpolation. When I’ve fully transitioned to ‘dst’, I swap ‘src’ and ‘dst’ and then for ‘dst’ I generate a band-limited wavetable using inverse FFT transform that takes into account the current frequency of the oscillator and the current required waveform. Note, I’m using the current frequency of the oscillator to generate the band-limited wavetable, not the octave. This works extremely well without any oversampling at all.

As a suggestion, maybe try:

    // Has changes?
    if (state.transitionFrac < 0.0) {

        const bool hasMorphChanges = morphRatio != state.morphRatioA;
        const bool hasFMChanges = fmAmount != state.fmAmountA;
        const bool hasOctaveChanges = currentOctave != newOctave;

        if (hasMorphChanges || hasFMChanges || hasOctaveChanges) {
            state.morphRatioB = morphRatio;
            state.fmAmountB = fmAmount;
            state.transitionFrac = 0.0; // initiate transition process
            currentOctave = newOctave;

            if (hasMorphChanges || hasOctaveChanges) {
                state.sampleB = oscConfiguration_.waveTable.getSample(currentOctave, morphRatio);
            }
        }
    }

And get rid of the if statement that begins:

if (currentOctave != newOctave) {

As a suggestion, maybe try:
And get rid of the if statement that begins:
if (currentOctave != newOctave) {

In this case sampleA and sampleB risk not being band-limited if there is no crossfade.

Note, I’m using the current frequency of the oscillator to generate the band-limited wavetable, not the octave. This works extremely well without any oversampling at all.

I’m using oversampling because I didn’t find another way out for the ladder filter that I’m using after osc on each voice. And also for the FM modulation. As soon as FM for example reaches amount 6.0 then I hear critically lacks of resolution to play without aliasing. So it’s reduce it a little bit.

In this case sampleA and sampleB risk not being band-limited if there is no crossfade.

You’re still checking for state.transitionFrac < 0.0 which it would be when there is no transition. I introduced ‘hasOctaveChanges’ so that you’d generate state.sampleB in this case and a transition would begin. So, if there’s a morph OR an FM change OR an octave change, then transition.

Doing this would ensure that you don’t change the wavetables arbitrarily whilst in transition which may well cause issues.

Another thing to try with the above would be to check for semitone rather than octave, just to see if it helps with the problem.