在C中实现Goertzel算法

Question

我正在DSP处理器上实现BFSK跳频通信系统.一些论坛成员建议使用Goertzel算法来解调特定频率的跳频.我已经尝试在C中实现goertzel算法.代码如下:

float goertzel(int numSamples,int TARGET_FREQUENCY,int SAMPLING_RATE, float* data)
{
    int     k,i;
    float   floatnumSamples;
    float   omega,sine,cosine,coeff,q0,q1,q2,result,real,imag;

    floatnumSamples = (float) numSamples;
    k = (int) (0.5 + ((floatnumSamples * TARGET_FREQUENCY) / SAMPLING_RATE));
    omega = (2.0 * M_PI * k) / floatnumSamples;
    sine = sin(omega);
    cosine = cos(omega);
    coeff = 2.0 * cosine;
    q0=0;
    q1=0;
    q2=0;

    for(i=0; i<numSamples; i++)
    {
        q0 = coeff * q1 - q2 + data[i];
        q2 = q1;
        q1 = q0;
    }
    real = (q1 - q2 * cosine);
    imag = (q2 * sine);
    result = sqrtf(real*real + imag*imag);
    return result;
}

当我使用该函数计算给定数据集的特定频率的结果时,我得不到正确的结果.但是,如果我使用相同的数据集并使用MATLAB goertzel()函数计算goertzel结果,那么我会得到完美的结果.我使用C实现了算法,借助我在互联网上找到的一些在线教程.如果函数正确地实现了goertzel算法,我只想获得你们的观点.

Answer 1

如果你说Matlab实现是好的,因为它的结果与你的数据的DFT或FFT的频率的结果相匹配,那么可能是因为Matlab实现正在通过比例因子对结果进行归一化,就像使用FFT一样.

更改您的代码以将其考虑在内,看看它是否会改善您的结果.请注意,为了清楚起见,我还更改了函数和结果名称以反映您的goertzel正在计算幅度,而不是完整的复杂结果:

float goertzel_mag(int numSamples,int TARGET_FREQUENCY,int SAMPLING_RATE, float* data)
{
    int     k,i;
    float   floatnumSamples;
    float   omega,sine,cosine,coeff,q0,q1,q2,magnitude,real,imag;

    float   scalingFactor = numSamples / 2.0;

    floatnumSamples = (float) numSamples;
    k = (int) (0.5 + ((floatnumSamples * TARGET_FREQUENCY) / SAMPLING_RATE));
    omega = (2.0 * M_PI * k) / floatnumSamples;
    sine = sin(omega);
    cosine = cos(omega);
    coeff = 2.0 * cosine;
    q0=0;
    q1=0;
    q2=0;

    for(i=0; i<numSamples; i++)
    {
        q0 = coeff * q1 - q2 + data[i];
        q2 = q1;
        q1 = q0;
    }

    // calculate the real and imaginary results
    // scaling appropriately
    real = (q1 - q2 * cosine) / scalingFactor;
    imag = (q2 * sine) / scalingFactor;

    magnitude = sqrtf(real*real + imag*imag);
    return magnitude;
}