需要帮助了解FFT算法中的这一行

问题描述

在我的程序中，我有一个执行快速傅立叶变换的函数。我知道有很多很好的免费实现，但是这是一个学习的事情，所以我不想使用那些。我最终通过以下实现找到了this comment（它源自FFT的意大利语条目）：

void transform(complex<double>* f,int N) //
{
  ordina(f,N);    //first: reverse order
  complex<double> *W;
  W = (complex<double> *)malloc(N / 2 * sizeof(complex<double>));
  W[1] = polar(1.,-2. * M_PI / N);
  W[0] = 1;
  for(int i = 2; i < N / 2; i++)
    W[i] = pow(W[1],i);
  int n = 1;
  int a = N / 2;
  for(int j = 0; j < log2(N); j++) {
    for(int k = 0; k < N; k++) {
      if(!(k & n)) {
        complex<double> temp = f[k];
        complex<double> Temp = W[(k * a) % (n * a)] * f[k + n];
        f[k] = temp + Temp;
        f[k + n] = temp - Temp;
      }
    }
    n *= 2;
    a = a / 2;
  }
  free(W);
}

到目前为止，我已经进行了很多更改，但这是我的出发点。我所做的更改之一是不缓存旋转因素，因为我决定先查看是否需要。现在，我决定要缓存它们。此实现的执行方式似乎是使用长度为W的数组N/2，其中每个索引k的值为

。我不明白的是这个表情：

W[(k * a) % (n * a)]

请注意，n * a始终等于N/2。我知道这应该等于

，并且我可以看到它依赖于

。我也得到模数可以在这里使用，因为旋转因子是循环的。但是，我没有得到一件事：这是一个长度为N的DFT，但仅计算了N/2旋转因子。数组的长度不应该为N，取模的长度是否应为N？

解决方法

但是我没有一件事：这是一个长度为N的DFT，但只计算了N / 2个旋转因子。数组的长度不应该为N，取模的长度是否应为N？

旋转因子是单位圆上等距的点，并且由于N是2的幂，所以点数偶数。在绕过圆的一半（从1开始，沿X轴逆时针旋转）之后，下半部分是上半部分的重复，但是这次是在X轴以下（点可以通过原点反映出来） ）。这就是第二次减去Temp的原因。这种减法是对旋转因子的否定。