1 




#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;
}

谁能说，有什么问题？？？

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

另外，我需要计算卷积。

但我找不到测试示例。

更新

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



vector<Complex> tmp( N );

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

return iF( tmp );




}

4

3 回答 3

2

我真的不知道你在问什么，但你的 DFT 和 IDFT 算法对我来说看起来是正确的。可以使用循环卷积定理使用 DFT 和 IDFT 执行卷积，该定理基本上说明了f**g = IDFT(DFT(f) * DFT(g))where**是循环卷积并且*是简单的乘法。

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

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

于 2010-11-22T15:44:42.017 回答
1

检查FFTW库“用于计算离散傅里叶变换 (DFT)”及其 C#包装器；）也许这也是；）

祝你好运！

于 2010-11-22T13:12:43.300 回答
0

变换看起来不错，但程序中没有任何东西在做卷积。

更新：卷积代码需要在逐元素乘法之前先对输入进行前向变换。

于 2010-11-22T13:11:10.077 回答