Синтез импульсных поездов в режиме реального времени с использованием SDL-микшера

Question

Я пытаюсь реализовать синтезатор звука, используя эту технику:

https://ccrma.stanford.edu/~stilti/papers/blit.pdf

Я делаю это в стандартном C, используя библиотеку SDL2_Mixer.

Это моя реализация функции BLIT:

double blit(double angle, double M, double P) {
    double x = M * angle / P;
    double denom = (M * sin(M_PI * angle / P));
    if (denom < 1)
        return (M / P) * cos(M_PI * x) / cos(M_PI * x / M);
    else {
        double numerator = sin(M_PI * x);
        return (M / P) * numerator / denom;
    }
}

Идея состоит в том, чтобы объединить его для создания прямоугольной волны, следуя инструкциям на бумаге. Я настроил SDL2_mixer с этой конфигурацией:

SDL_AudioSpec *desired, *obtained;
SDL_AudioSpec *hardware_spec;

desired = (SDL_AudioSpec*)malloc(sizeof(SDL_AudioSpec));
obtained = (SDL_AudioSpec*)malloc(sizeof(SDL_AudioSpec));

desired->freq=44100;
desired->format=AUDIO_U8;
desired->channels=1;
desired->samples=2048;
desired->callback=create_rect;
desired->userdata=NULL;

А вот мой create_rect функция. Он создает биполярный импульсный ряд, а затем интегрирует его значение для генерации прямоугольной функции с ограничением по полосе.

void create_rect(void *userdata, Uint8 *stream, int len) {
    static double angle = 0;
    static double integral = 0;
    int i = 0;
    // This is the freq of my tone
    double f1 = tone_table[current_wave.note];
    // Sample rate
    double fs = 44100;
    // Pulse
    double P = fs / f1;
    int M = 2 * floor(P / 2) + 1;

    double oldbipolar = 0;
    double bipolar = 0;
    for(i = 0; i < len; i++) {
        if (++angle > P)
            angle -= P;
        double angle2 = angle + floor(P/2);
        if (angle2 > P)
            angle2 -= P;

        bipolar = blit(angle2, M, P) - blit(angle, M, P);

        integral += (bipolar + old bipolar) * 0.5;
        oldbipolar = bipolar;
        *stream++ = (integral + 0.5) * 127;
    }
}

Моя проблема: результирующая волна вполне нормальная, но через несколько секунд она начинает издавать звуки. Я попытался построить результат, и вот оно:

Любая идея?

РЕДАКТИРОВАТЬ: Вот график биполярного BLIT до его интеграции: