audio - 为循环 wav 样本添加抗锯齿/带宽限制(非傅立叶变换)

Question

如何使用 C++ 代码构建抗锯齿插值？我有一个简单的 4096 或 1024 缓冲区。当然，当我以高频播放时，我会遇到混叠问题。为避免这种情况，信号必须受到高频带宽的限制。粗略地说，高频的“锯齿”波应该看起来像一个正弦波。这就是我想要得到的，这样我的声音就不会像你在车里的旧 FM/AM radio 中移动旋钮一样脏。

我知道如何用傅里叶变换构建带限正方形、三角形、锯齿形。所以我的问题只是关于波表

Answer 1

在 AudioKit 源中找到解决方案。一个缓冲区将被分成 10 个缓冲区/ Octave 。因此，当您播放声音时，您不会播放原始波形，而是播放为特定 Octave 音程准备的样本。

导入您的项目 WaveStack.hpp

namespace AudioKitCore
{

    // WaveStack represents a series of progressively lower-resolution sampled versions of a
    // waveform. Client code supplies the initial waveform, at a resolution of 1024 samples,
    // equivalent to 43.6 Hz at 44.1K samples/sec (about 23.44 cents below F1, midi note 29),
    // and then calls initStack() to create the filtered higher-octave versions.
    // This provides a basis for anti-aliased oscillators; see class WaveStackOscillator.

    struct WaveStack
    {
        // Highest-resolution rep uses 2^maxBits samples
        static constexpr int maxBits = 10;  // 1024

        // maxBits also defines the number of octave levels; highest level has just 2 samples
        float *pData[maxBits];

        WaveStack();
        ~WaveStack();

        // Fill pWaveData with 1024 samples, then call this
        void initStack(float *pWaveData, int maxHarmonic=512);

        void init();
        void deinit();

        float interp(int octave, float phase);
    };

}

WaveStack.cpp

#include "WaveStack.hpp"
#include "kiss_fftr.h"

namespace AudioKitCore
{

    WaveStack::WaveStack()
    {
        int length = 1 << maxBits;                  // length of level-0 data
        pData[0] = new float[2 * length];           // 2x is enough for all levels
        for (int i=1; i<maxBits; i++)
        {
            pData[i] = pData[i - 1] + length;
            length >>= 1;
        }
    }

    WaveStack::~WaveStack()
    {
        delete[] pData[0];
    }

    void WaveStack::initStack(float *pWaveData, int maxHarmonic)
    {
        // setup
        int fftLength = 1 << maxBits;
        float *buf = new float[fftLength];
        kiss_fftr_cfg fwd = kiss_fftr_alloc(fftLength, 0, 0, 0);
        kiss_fftr_cfg inv = kiss_fftr_alloc(fftLength, 1, 0, 0);

        // copy supplied wave data for octave 0
        for (int i=0; i < fftLength; i++) pData[0][i] = pWaveData[i];

        // perform initial forward FFT to get spectrum
        kiss_fft_cpx spectrum[fftLength / 2 + 1];
        kiss_fftr(fwd, pData[0], spectrum);

        float scaleFactor = 1.0f / (fftLength / 2);

        for (int octave = (maxHarmonic==512) ? 1 : 0; octave < maxBits; octave++)
        {
            // zero all harmonic coefficients above new Nyquist limit
            int maxHarm = 1 << (maxBits - octave - 1);
            if (maxHarm > maxHarmonic) maxHarm = maxHarmonic;
            for (int h=maxHarm; h <= fftLength/2; h++)
            {
                spectrum[h].r = 0.0f;
                spectrum[h].i = 0.0f;
            }

            // perform inverse FFT to get filtered waveform
            kiss_fftri(inv, spectrum, buf);

            // resample filtered waveform
            int skip = 1 << octave;
            float *pOut = pData[octave];
            for (int i=0; i < fftLength; i += skip) *pOut++ = scaleFactor * buf[i];
        }

        // teardown
        kiss_fftr_free(inv);
        kiss_fftr_free(fwd);
        delete[] buf;
    }

    void WaveStack::init()
    {
    }

    void WaveStack::deinit()
    {
    }

    float WaveStack::interp(int octave, float phase)
    {
        while (phase < 0) phase += 1.0;
        while (phase >= 1.0) phase -= 1.0f;

        int nTableSize = 1 << (maxBits - octave);
        float readIndex = phase * nTableSize;
        int ri = int(readIndex);
        float f = readIndex - ri;
        int rj = ri + 1; if (rj >= nTableSize) rj -= nTableSize;

        float *pWaveTable = pData[octave];
        float si = pWaveTable[ri];
        float sj = pWaveTable[rj];
        return (float)((1.0 - f) * si + f * sj);
    }

}

然后这样使用:

//wave and outputWave should be float[1024];
void getSample(int octave, float* wave, float* outputWave){

    uint_fast32_t impulseCount = 1024;

    if (octave == 0){
        impulseCount = 737;
    }else if (octave == 1){
        impulseCount = 369;
    }
    else if (octave == 2){
        impulseCount = 185;
    }
    else if (octave == 3){
        impulseCount = 93;
    }
    else if (octave == 4){
        impulseCount = 47;
    }
    else if (octave == 5){
        impulseCount = 24;
    }
    else if (octave == 6){
        impulseCount = 12;
    }
    else if (octave == 7){
        impulseCount = 6;
    }
    else if (octave == 8){
        impulseCount = 3;
    }
    else if (octave == 9){
        impulseCount = 2;
    }

    //Get sample for octave
    stack->initStack(wave, impulseCount);

    for (int i = 0; i < 1024;i++){
        float phase = (1.0/float(1024))*i;

        //get interpolated wave and apply volume compensation
        outputWave[i] = stack->interp(0, phase) / 2.0;
    }

}

然后...当 10 个缓冲区准备就绪时。您可以在播放声音时使用它们。使用此代码，您可以获得缓冲区/ Octave 音阶的索引，具体取决于您的频率

uint_fast8_t getBufferIndex(const float& frequency){

    if (frequency >= 0 && frequency < 40){
        return 0;
    }
    else if (frequency >= 40 && frequency < 80){
        return 1;
    }else if (frequency >= 80 && frequency < 160){
        return 2;
    }else if (frequency >= 160 && frequency < 320){
        return 3;
    }else if (frequency >= 320 && frequency < 640){
        return 4;
    }else if (frequency >= 640 && frequency < 1280){
        return 5;
    }else if (frequency >= 1280 && frequency < 2560){
        return 6;
    }else if (frequency >= 2560 && frequency < 5120){
        return 7;
    }else if (frequency >= 5120 && frequency < 10240){
        return 8;
    }else if (frequency >= 10240){
        return 9;
    }
    return 0;
}

所以如果我知道我的音符频率是 440hz。然后对于这个笔记，我以这种方式得到波:

float notInterpolatedSound[1024];
float interpolatedSound[1024];

uint_fast8_t octaveIndex = getBufferIndex(440.0);
getSample(octaveIndex, notInterpolatedSound, interpolatedSound);

//tada! 

ps. the code above is a low pass filter. I also tried to do sinc interpolation. But sinc worked for me very expensive and not exactly. Although maybe I did it wrong.