离散傅里叶变换给出不正确的结果

Question

我试图在C中进行离散傅立叶变换.

最初只是蛮力方法.首先,我让程序打开一个数据文件(幅度)并将数据放入一个数组(只有一个,因为我限制自己使用实值输入).

但是变换看起来不对,所以我尝试生成一个简单的波函数并检查它是否正确转换.

这是我的代码,剥去了钟声和口哨声:

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

#define M_PI 3.14159265358979323846

//the test wavefunction
double theoretical(double t)
{
    double a = sin(M_PI * t) + 2 * sin(2 * M_PI * t) + 4 * sin(4 * M_PI * t); 
    return a;
}
//-------------------------------------------------------------------------
void dftreal(double inreal[], double outreal[], double outimag[], int linecount) 
{
    int n, k;
    for (k = 0; k < linecount; k++)
    {
        double sumreal = 0;
        double sumimag = 0;

        for (n = 0; n < linecount; n++) 
        {
            double angle = 2 * M_PI * n * ( k / (double) linecount);
            sumreal +=  inreal[n] * cos(angle);
            sumimag +=  inreal[n] * sin(angle);
        }
        outreal[k] = sumreal;
        outimag[k] = sumimag;
    }
}
//=========================================================================
int main(void)
{
     int linecount = 44100;
     //creates all necessary arrays
     double inreal[linecount], outreal[linecount], outimag[linecount], p[linecount]; 

     FILE *fout = fopen("Output.txt", "w");

     for (int i = 0 ; i < linecount ; ++i)
     {
         inreal[i] = theoretical( i / (double) linecount);
     }

     //actually computes the transform
     dftreal(inreal, outreal, outimag, linecount); 

     for (int i = 0 ; i < linecount ; ++i)
     {
          p[i] = 2*(outreal[i] * outreal[i] + outimag[i] * outimag[i]);
          fprintf(fout, "%f %f \n", (i / (double) linecount), p[i]);
     }
     fclose(fout);
     printf("\nEnd of program");
     getchar();
     return 0;
}

程序编译,完成,但在功率(频率)图上没有几个尖峰,我得到这个: .

单个频率或不同频率给出完全相同的倒浴盆曲线.

我检查了几个关于DFT的消息来源,我仍然不知道出了什么问题,这个函数似乎没有任何明显的错误:

dftreal

本身.我想就可能导致问题的原因寻求帮助.我在Windows 7上使用MinGW编译器.谢谢!

Answer 1

好消息是您的实施没有任何问题dftreal。

问题（如果可以这样称呼的话）是您使用的测试波形包含相对于您的采样率而言频率非常低的频率分量linecount。相应地，DFT 显示能量集中在前几个箱中，一些频谱泄漏到更高的箱中，因为波形频率分量不是采样率的精确倍数。

如果通过使频率相对于采样频率来增加测试波形频率，例如：

//the test wavefunction
double theoretical(double t)
{
  double f = 0.1*44100;
  double a = sin(2 * M_PI * f * t) + 2 * sin(4 * M_PI * f * t) + 4 * sin(8 * M_PI * f * t); 
  return a;
}

您应该得到如下图：