基于FFT的频率移位器

7

我一直在使用 Rosetta Code 提供的原始 FFT 算法开发频移器。理解到,为了对样本信号进行频移,需要将原始音频应用 FFT,将每个结果正弦波的频率乘以频移因子(由用户定义),然后将正弦波加在一起。当我运行算法时,输出的质量非常低,就好像算法中收集到的正弦波数量不足以正确地再现信号。该算法在头文件中实现为类,并在其他地方(正确地)调用。

#include <complex>
#include <valarray>

typedef std::complex<double> Complex;
typedef std::valarray<Complex> CArray;

class FrequencyShifter {
float sampleRate;
public:
    FrequencyShifter() {

    }
    void setSampleRate(float inSampleRate) {
        sampleRate = inSampleRate;
    }
    double abs(double in0) {
        if (in0>=0) return in0;
        else return -in0;
    }
    void fft(CArray& x)
    {
        const size_t N = x.size();
        if (N <= 1) return;

        // divide
        CArray even = x[std::slice(0, N/2, 2)];
        CArray  odd = x[std::slice(1, N/2, 2)];

        // conquer
        fft(even);
        fft(odd);

        // combine
        for (size_t k = 0; k < N/2; ++k)
        {
            Complex t = std::polar(1.0, -2 * PI * k / N) * odd[k];
            x[k    ] = even[k] + t;
            x[k+N/2] = even[k] - t;
        }
    }
    double convertToReal(double im, double re) {
        return sqrt(abs(im*im - re*re));
    }
    void processBlock(float *inBlock, const int inFramesToProcess, float scale) {
        //inFramesToProcess is the amount of samples in inBlock
        Complex *copy = new Complex[inFramesToProcess];
        for (int frame = 0; frame<inFramesToProcess; frame++) {
            copy[frame] = Complex((double)inBlock[frame], 0.0);
        }
        CArray data(copy, inFramesToProcess);
        fft(data);
        const float freqoffsets = sampleRate/inFramesToProcess;
        for (float x = 0; x<data.size()/2; x++) {
            for (float frame = 0; frame<inFramesToProcess; frame++) {
                inBlock[(int)frame] = (float)(convertToReal(data[(int)x].imag(), data[(int)x].real())*sin(freqoffsets*x*frame*scale));
            }
        }
    }
};

我猜问题的一部分在于我只包括sampleRate/inFramesToProcess频率来覆盖正弦波。发送更大的音频文件(因此更大的*inBlockinFramesToProcess)会使音频变得不那么颗粒状吗?我该如何做到这一点而不仅仅改变参数的值或长度?


“there is no output” 是什么意思? - 1201ProgramAlarm
@1201ProgramAlarm 当我测试 *inBlock 的输出时,没有电平(音频电平为0或遇到其他错误)。基本上,算法中存在一些错误,我无法检测和修复。 - Linus Rastegar
1
“convertToReal” 是正确的吗?如果 “inFramesToProcess” 为1,则“data”将具有没有虚部的复数。fft 不会对其进行任何操作,因此在转换回来时,您将尝试对负数取平方根。如果“x.size()”为奇数,则 fft 不会对“x”的最后一个元素执行任何操作。 - 1201ProgramAlarm
@1201ProgramAlarm 啊,谢谢你提醒我!我之前也没有意识到这两个问题。我会尽快修复并更新帖子。这可能是我正在寻找的解决方案。 - Linus Rastegar
@1201ProgramAlarm 我在尝试解决这两个问题时遇到了瓶颈。如果您能详细回复并给出答案,我将不胜感激。 - Linus Rastegar
1个回答

4

这是一个更新版本的 processBlock,其中包含了许多必要的调整来实现频率偏移,下面我将对此进行描述:

void processBlock(float *inBlock, const int inFramesToProcess, float scale) {
    //inFramesToProcess is the amount of samples in inBlock
    Complex *copy = new Complex[inFramesToProcess];
    for (int frame = 0; frame<inFramesToProcess; frame++) {
        copy[frame] = Complex((double)inBlock[frame], 0.0);
    }
    CArray data(copy, inFramesToProcess);
    fft(data);
    const float freqoffsets = 2.0*PI/inFramesToProcess;
    const float normfactor  = 2.0/inFramesToProcess;
    for (int frame = 0; frame<inFramesToProcess; frame++) {
        inBlock[frame] = 0.5*data[0].real();
        for (int x = 1; x<data.size()/2; x++) {
            float arg = freqoffsets*x*frame*scale;
            inBlock[frame] += data[x].real()*cos(arg) - data[x].imag()*sin(arg);
        }
        inBlock[frame] *= normfactor;
    }
}

导出

从FFT得到的频谱是复值的,可以看作是用正弦和余弦波表示信号。使用反变换可以重建时域波形,其关系如下: enter image description here

利用频谱对称性,可以将其表示为:

enter image description here

或等价地表示为:

enter image description here

你可能已经注意到,索引为0N/2的项是频域中具有纯实系数的特殊情况。为简单起见,假设频谱没有一直延伸到N/2,则可以去掉N/2项并仍然得到合理的近似值。对于其他项,你将得到一个可实现的贡献,如下:

normfactor = 2.0/inFramesToProcess;
normfactor*(data[x].real()*cos(arg) - data[x].imag()*sin(arg))

当然,你需要将所有这些贡献添加到最终缓冲区 inBlock[frame] 中,而不是仅仅覆盖之前的结果:

inBlock[frame] += normfactor*(data[x].real()*cos(arg) - data[x].imag()*sin(arg));
//             ^^ 

请注意,归一化可以在循环后对最终结果进行,从而减少乘法的数量。这样做时,我们必须特别注意索引为0的直流项(其系数为1/N而不是2/N):
inBlock[frame] = 0.5*data[0].real();
for (int x = 1; x<data.size()/2; x++) {
    float arg = freqoffsets*x*frame*scale;
    inBlock[frame] += data[x].real()*cos(arg) - data[x].imag()*sin(arg);
}
inBlock[frame] *= normfactor;

最后,在生成音调时,给sincos函数的相位参数arg应该为形式为2*pi*k*n/inFramesToProcess(在应用scale因子之前),其中n是时域样本索引,k是频率域索引。最终结果是计算出的频率增量freqoffsets应该真正为2.0*PI/inFramesToProcess注意:
  • FFT算法基于你的基础时域信号是块长度周期性的假设。因此,块与块之间可能会出现听得见的不连续性。
  • 未来的读者应该知道,这不会将频谱移动固定的量,而是通过乘法因子压缩或扩展频谱。例如,包括100-200Hz分量的信号可能会被0.75的因子压缩到75-150Hz。请注意,下限降低了25Hz,而上限降低了50Hz。

非常感谢你的回答!在我理解之前,我需要读几遍。我已经实现了你的 processBlock 版本,确实使音频听起来好多了。我也非常感激你花时间描述推导过程。 - Linus Rastegar

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接