c++ - DFT算法和卷积。怎么了？

Question

<pre><code>#include <vector>
using std::vector;

#include <complex>
using std::complex;
using std::polar;

typedef complex<double> Complex;

#define Pi 3.14159265358979323846

// direct Fourier transform
vector<Complex> dF( const vector<Complex>& in )
{
 const int N = in.size();

 vector<Complex> out( N );

 for (int k = 0; k < N; k++)
 {
  out[k] = Complex( 0.0, 0.0 );

  for (int n = 0; n < N; n++)
  {
   out[k] += in[n] * polar<double>( 1.0, - 2 * Pi * k * n / N );
  }
 }

 return out;
}

// inverse Fourier transform
vector<Complex> iF( const vector<Complex>& in )
{
 const int N = in.size();

 vector<Complex> out( N );

 for (int k = 0; k < N; k++)
 {
  out[k] = Complex( 0.0, 0.0 );

  for (int n = 0; n < N; n++)
  {
   out[k] += in[n] * polar<double>(1, 2 * Pi * k * n / N );
  }

  out[k] *= Complex(  1.0 / N , 0.0 );
 }

 return out;
}
</code></pre>

<p></p>

谁能说，哪里错了？？？

也许我不明白这个算法的实现细节......但我找不到它)))

另外，我需要计算卷积。

但是我找不到测试例子。

更新

// convolution. I suppose that x0.size == x1.size
vector convolution( const vector& x0, const vector& x1) 
{
    const int N = x0.size();<p></p>

<pre><code>vector<Complex> tmp( N );

for ( int i = 0; i < N; i++ )
{
    tmp[i] = x0[i] * x1[i];
}

return iF( tmp );
</code></pre>

<p>}
</p>

Answer 1

我真的不知道你到底在问什么，但你的 DFT 和 IDFT 算法对我来说看起来是正确的。可以使用 circular convolution theorem 使用 DFT 和 IDFT 执行卷积它基本上表明 f**g = IDFT(DFT(f) * DFT(g)) 其中 ** 是循环卷积和 *是简单的乘法。

要使用 DFT 计算线性卷积(非循环)，您必须对每个输入进行零填充，以便循环环绕只发生在零值样本上，不会影响输出。每个输入序列需要补零到长度为 N >= L+M-1 其中 L 和 M 是序列的长度输入序列。然后执行如上所示的循环卷积，第一个 L+M-1 样本是线性卷积输出(超出此范围的样本应为零)。

注意:使用您展示的 DFT 和 IDFT 算法执行卷积比直接计算要低效得多。只有在使用 FFT 和 IFFT(O(NlogN)) 算法代替 DFT 和 IDFT (O(N^2)) 时才会有优势。