Goertzel 算法给出无限结果

Posted 2023-02-25

【中文标题】Goertzel 算法给出无限结果

【英文标题】：Goertzel algorithm giving infinite result

【发布时间】：2015-03-09 18:17:25

【问题描述】：

我有一个使用 Audacity 软件创建的 20hz - 1 振幅的正弦波。也只有500ms。

我正在使用以下算法来检测频率。

我只想检测音调幅度是否超过阈值并在 20 hz 频率周期给我肯定的结果。

static 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;

调用函数

// If there's more packets, read them
        inCompleteAQBuffer->mAudioDataByteSize = numBytes;
        CheckError(AudioQueueEnqueueBuffer(inAQ,
                                           inCompleteAQBuffer,
                                           (sound->packetDescs?nPackets:0),
                                           sound->packetDescs),
                   "couldn't enqueue buffer");
        sound->packetPosition += nPackets;

        NSLog(@"number of packets %i",nPackets);

        float *data=(float*)inCompleteAQBuffer->mAudioData;
        int nn = sizeof(data)/sizeof(float);
        float gort = goertzel_mag(nn, 20, 44100, data);
        NSLog(@"gort:%f", gort);
        if (gort == INFINITY)
             NSLog(@"positive infinity");  

函数内部的断点

输出

number of packets 8192
gort:36029896530591744.000000
number of packets 8192
gort:inf
positive infinity
number of packets 5666
gort:inf
positive infinity 

为什么我得到 inf 结果？我不知道如何读取返回值，我知道幅度总是必须是正值，但我创建了 1 幅度的文件，我不应该得到 0 到 1 的结果吗？

编辑； 音频信息

afinfo 500ms.aiff
File:           500ms.aiff
File type ID:   AIFF
Num Tracks:     1
----
Data format:     1 ch,  44100 Hz, 'lpcm' (0x0000000E) 16-bit big-endian signed integer
                no channel layout.
estimated duration: 0.500000 sec
audio bytes: 44100
audio packets: 22050
bit rate: 705600 bits per second
packet size upper bound: 2
maximum packet size: 2
audio data file offset: 54
optimized
source bit depth: I16

【参考方案1】：

我认为这些行有一个问题：

float *data=(float*)inCompleteAQBuffer->mAudioData;
int nn = sizeof(data)/sizeof(float);

我相信这是为了告诉你样本的数量。我没有重现您的代码的信息或资源，但可以通过以下方式重现错误：

#include <stdio.h>
#include <stdlib.h>

int main(void) 
    float *data=malloc(sizeof(float) * 10);
    printf ("Sizeof 'data' = %d\n", sizeof(data));
    return 0;

程序输出

Sizeof 'data' = 4

在我的 32 位编译中是数组指针的大小，不是数组。并且使用sizeof(*data) 不会让您有任何收获，因为它只会告诉您数据类型float 的大小，而不是数组。

您无法从其指针确定数组的大小或元素的数量，所以很遗憾，我的回答是，您需要更多信息，也许是numBytes？还是numBytes/sizeof(float)？

【讨论】：

我已经编辑了问题（添加了音频信息），根据我的日志，nn 总是 2

@MordFustang nn=2 是因为你有一个 64 位系统所以 64 位指针所以 int nn = sizeof(data)/sizeof(float); 是 8/4 即 2。这被传递给 goertzel_mag(int numSamples...) 的第一个参数，它运行一个过滤器for(i=0; i<numSamples; i++)... 只用 2 个样本去哪里？