c - 从C中的FFT阵列中删除1000Hz音调

我有一个 double 数组，这是对数组应用FFT的结果，其中包含Wav音频文件的音频数据，在其中我添加了1000Hz音调。

我通过“数值食谱”中定义的DREALFT获得了这个数组。(我必须使用它)。
(原始数组的长度是2的幂。)

地雷阵列具有以下结构:

array [0] =复数变换的第一个实数值分量

array [1] =复数变换的最后一个实值分量

array [2] =第二个元素的实部

array [3] =第二个元素的虚部

等等......

现在，我知道该阵列代表频域。

我想确定并杀死1000Hz频率。

我已经尝试过以下公式来查找应该包含1000Hz频率的数组的索引:

index = 1000. * NElements /44100;

此外，由于我假设此索引仅引用具有实数值的数组，因此我确定了数组中的正确位置(？)，该位置也包含虚值:

    int correctIndex=2;

for(k=0;k<index;k++){
    correctIndex+=2;
}

(我知道肯定有一种方法更容易，但这是我想到的第一个方法)

然后，我找到这个值:16275892957.123705，我想它是1000Hz频率的实数部分(抱歉，这是一个不精确的发射，但目前我不在乎对此有所了解)

所以我试图压制它:

array[index]=-copy[index]*0.1f;

我不知道为什么要使用此公式，但它是唯一给出一些结果的公式，实际上1000hz的音调似乎略有降低。

这是有问题的代码的一部分:

    double *copy = malloc( nCampioni * sizeof(double));
    int nSamples;

 /*...Fill copy with audio data...*/

     /*...Apply ZERO PADDING and reach the length of 8388608 samples,
 or rather 8388608 double values...*/


/*Apply the FFT (Sure this works)*/
drealft(copy - 1, nSamples, 1);

/*I determine the REAL(?) array index*/
i= 1000. * nSamples /44100;

/*I determine MINE(?) array index*/
int j=2;

for(k=0;k<i;k++){
    j+=2;
}

/*I reduce the array value, AND some other values aroud it as an attempt*/
for(i=-12;i<12;i+=2){
    copy[j-i]=-copy[i-j]*0.1f;
    printf("%d\n",j-i);
}

/*Apply the inverse FFT*/
drealft(copy - 1, nSamples, -1);

/*...Write the audio data on the file...*/

注意:为简单起见，我省略了从int16_t数组中获取double数组的部分

如何确定并完全杀死1000Hz频率？

谢谢!

最佳答案

正如奥利·查尔斯沃斯(Oli Charlesworth)所写的那样，由于您的目标频率不完全是FFT频段之一(您的index，TargetFrequency * NumberOfElements / SamplingRate，不完全是整数)，因此目标频率的能量将分布在所有频段上。首先，您可以通过将最接近目标频率的频段置零来消除某些频率。当然，这也会稍微影响其他频率，因为它稍微偏离了目标。为了更好地抑制目标频率，您将需要考虑使用更复杂的滤波器。

但是，出于教育目的:要抑制与某个bin对应的频率，只需将该bin设置为零即可。您必须将垃圾箱的实部和虚部都设置为零，这可以使用以下方法完成:

copy[index*2 + 0] = 0;
copy[index*2 + 1] = 1;

关于此的一些注意事项:

您具有以下代码来计算数组中的位置:

int correctIndex = 2;
for (k = 0; k < index; k++) {
    correctIndex += 2;
}

这等效于:

correctIndex = 2*(index+1);

我相信您想要2*index，而不是2*(index+1)。因此，您可能减少了错误的垃圾箱。

在问题的某一时刻，您编写了array[index] = -copy[index]*0.1f;。我不知道array是什么。您似乎在copy中就地工作。我也不知道你为什么乘以1/10。如果要消除频率，只需将其设置为零即可。将其乘以1/10只会将其减小到其原始大小的10％。

我了解您必须将copy-1传递给drealft，因为数字食谱代码使用基于1的索引。但是，C标准不支持您的操作方式。标准未定义表达式copy-1的行为。它可以在大多数C实现中使用。但是，要编写受支持的可移植代码，您应该改为执行以下操作:

// Allocate one extra element.
double *memory = malloc((nCampioni+1) * sizeof *memory);

// Make a pointer that is convenient for your work.
double *copy = memory+1;

…

// Pass the necessary base address to drealft.
drealft(memory, nSamples, 1);

// Suppress a frequency.
copy[index*2 + 0] = 0;
copy[index*2 + 1] = 0;

…

// Free the memory.
free(memory);

我建议您考虑的一个实验是仅使用所需频率的正弦波初始化数组:

for (i = 0; i < nSamples; ++i)
    copy[i] = sin(TwoPi * Frequency / SampleRate * i);

(TwoPi当然是2*3.1415926535897932384626433。)然后应用drealft并查看结果。您会看到很多能量都在最接近目标频率的仓中的一个峰值处，但是其中许多能量也已经扩散到其他仓中。显然，将单个仓清零并执行逆FFT并不能消除所有频率。此外，您应该看到该峰位于您为index计算的同一仓中。如果不是，则出问题了。