如何理解DFT结果

4

我使用这个DFT实现:

/*
Direct fourier transform
*/
int DFT(int dir,int m,double *x1,double *y1)
{
    long i,k;
    double arg;
    double cosarg,sinarg;
    double *x2=NULL,*y2=NULL;

    x2 = malloc(m*sizeof(double));
    y2 = malloc(m*sizeof(double));
    if (x2 == NULL || y2 == NULL)
       return(FALSE);

    for (i=0;i<m;i++) {
       x2[i] = 0;
       y2[i] = 0;
       arg = - dir * 2.0 * 3.141592654 * (double)i / (double)m;
       for (k=0;k<m;k++) {
          cosarg = cos(k * arg);
          sinarg = sin(k * arg);
          x2[i] += (x1[k] * cosarg - y1[k] * sinarg);
          y2[i] += (x1[k] * sinarg + y1[k] * cosarg);
       }
    }

    /* Copy the data back */
    if (dir == 1) {
      for (i=0;i<m;i++) {
         x1[i] = x2[i] / (double)m;
         y1[i] = y2[i] / (double)m;
      }
   } else {
      for (i=0;i<m;i++) {
         x1[i] = x2[i];
         y1[i] = y2[i];
      }
   }

   free(x2);
   free(y2);
   return(TRUE);
}

这里放置了一个链接:http://paulbourke.net/miscellaneous/dft/

第一个问题是为什么在应用直接变换(dir=1)后,我们需要对值进行缩放?我读了一些关于DFT实现的想法,但没有找到有关此问题的任何信息。

我使用的输入是1024采样频率的cos。

#define SAMPLES 2048
#define ZEROES_NUMBER 512

double step = PI_2/(SAMPLES-2*ZEROES_NUMBER);
for(int i=0; i<SAMPLES; i++)
{
    /*
     * Fill in the beginning and end with zeroes 
     */
    if(i<ZEROES_NUMBER || i > SAMPLES-ZEROES_NUMBER)
    {
        samplesReal[i] = 0;
        samplesImag[i] = 0;
    }
    /*
     *  Generate one period cos with 1024 samples
     */
    else
    {
        samplesReal[i] = cos(step*(double)(i-ZEROES_NUMBER));
        samplesImag[i] = 0;
    }
}

为了绘图,我删除了之前问过的缩放操作,因为输出值变得非常小,无法绘制图形。

我得到了幅度和相位的图表:

enter image description here enter image description here

如您所见,相位总是为0,幅度频谱被颠倒。为什么?

以下是更易读且不使用缩放操作的版本,其产生相同的结果:

void DFT_transform(double complex* samples, int num, double complex*   res)
{
    for(int k=0; k<num; k++)
    {
        res[k] = 0;
        for(int n=0; n<num; n++)
        {
            double complex Wkn = cos(PI_2*(double)k*(double)n/(double)num) -
            I*sin(PI_2*(double)k*(double)n/(double)num);

            res[k] += samples[n]*Wkn;
        }
    }
}

1
每个问题一个问题。缩放只是一种约定,有不同的约定-通常是正向方向为1 / N,反向方向没有缩放。要查看是否您的DFT正常工作,请从纯实正弦波开始,并查看输出箱的幅度(sqrt(re ^ 2 + im ^ 2))。 - Paul R
@LongSmith:但是你在问关于C代码,对吗? - Lightness Races in Orbit
@PaulR 谢谢。我尝试了输入 cos(i)。这是我得到的 振幅相位。振幅只有一点点不同,但相位出现了一些新的东西。这是在开头和结尾填充零时发生的。 - Long Smith
1
@LongSmith: 你所说的振幅是指大小,即 sqrt(re^2 + im^2) 吗?请注意,输入应该在 x1 中有正弦波,在 y1 中填充零点 - 不要像你在问题代码中那样在开头和结尾填充零点。您可以使用已知的良好DFT或FFT实现(例如Octave(免费的MATLAB克隆版))检查结果。 - Paul R
@PaulR 是的,我指的是振幅,抱歉搞混了。这里有一个没有零点的振幅相位的图表。它们显示出一些奇怪的东西。稍后会尝试使用MATLAB。 - Long Smith
显示剩余2条评论
1个回答

3

好的,大家好。我很高兴地说这个实现成功了。问题是绘制图形的方式错误以及对公式的理解不足。

它是如何工作的

DFT Formulas

正如您所看到的,有一个k变量用于改变频率。因此,频率为ν = k / T,其中T是取样时间。 T = N/S,其中S是您的采样频率。然后,您可以将您的频率找到为v = S*k/N

因此,当您获得结果时,您应该为每个点计算频率,并删除所有高于S/2的内容,然后仅绘制图形幅度 = 幅度(频率)。这就是我以前不理解的地方。希望对某人有所帮助。

我得到的一些图表。

正弦100HZ。

Sin 100Hz

Sin 100HZ + Cos 200HZ。 sin(100x)+cos(200x) Sin 100HZ + (Cos 200HZ)/2 sin(100x)+cos(200x)/2 如您所见,频率及其相关幅度被显示出来。虽然存在缩放问题,但如果我们要确定信号中呈现的频率,这并不重要。
感谢 @PaulR

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