在FFT算法

Question

我的朋友們的一個錯誤，我現在編程FFT算法，我有以下程序：

#include<iostream> 
#include<vector> 
#include<complex> 
#include<cmath> 
#define _USE_MATH_DEFINES 
#include<math.h> 
#include<stdlib.h> 

using namespace std; 

const complex<double> I(0,1); 
const int pi = M_PI; 

//This function will check if a number is a power of certain number 
bool checkpower(float base,float num){ 
    float a = ceil(log(num)/log(base))-floor(log(num)/log(base)); 
    if (a==0){ 
     return 1; 
    } 
    else{ 
    return 0; 
    } 
} 


//Fast Fourier Transform for DFT 
vector<complex<double>> FFT_DFT(vector<complex<double>> samples){ 
    int N = samples.size(); 
    cout << N << endl; 
    if(!checkpower(2,N)){ 
     cout << "Please make sure the sample size is of 2^n!" << endl; 
     exit (EXIT_FAILURE); 
    } 
    else{ 
     vector<complex<double>> F(N); 
     if(N==1){ 
      F.push_back(samples[0]); 
     } 
     else{ 
      int M = N/2; 
      cout << M << endl; 
      vector<complex<double>> O; //vector to store the odd elements 
      vector<complex<double>> E; //vector to store the even elements 
      //Reoder the samples 
      for(int l=0;l<M;l++){ 
       E.push_back(samples[2*l]); 
       O.push_back(samples[2*l+1]); 
      } 
      vector<complex<double>> ODFT; 
      cout << "Start recursive for odd" << endl; 
      ODFT = FFT_DFT(O); 
      vector<complex<double>> EDFT; 
      cout << "Start recursive for even" << endl; 
      EDFT = FFT_DFT(E); 
      for(int k=0;k<M;k++){ 
       cout << real(EDFT[k]) << " + "<< imag(EDFT[k]) << "I" << endl; 
       cout << real(ODFT[k]) << " + "<< imag(ODFT[k]) << "I" << endl; 
       F[k] = EDFT[k]+exp(-2.0*pi*k*I/(double)N)*ODFT[k]; 
       F[k+M] = EDFT[k]-exp(-2.0*pi*k*I/(double)N)*ODFT[k]; 
       cout << real(F[k]) << " + "<< imag(F[k]) << "I" << endl; 
       cout << real(F[k+M]) << " + "<< imag(F[k+M]) << "I" << endl; 
      } 
     } 
     return F; 
    } 
} 


int main(){ 
    vector<complex<double>> samples; 
    samples.push_back(8.0); 
    samples.push_back(4.0); 
    samples.push_back(8.0); 
    samples.push_back(0.0); 

    vector<complex<double>> dft = FFT_DFT(samples); 
    vector<complex<double>>::iterator item; 

    for(item=dft.begin();item!=dft.end();item++){ 
     cout << real(*item) << " + "<< imag(*item) << "I" << endl; 
    } 


    return 0; 
} 

我使用Visual Studio 2010專業的編譯器。我不知道我的遞歸算法有什麼問題，所以我在VS中使用了調試模式，並且我一行一行地檢查了過程，但它似乎總是給我所有的值爲0。我用普通的FFT算法測試它，它工作得很完美。那麼，任何人都可以幫我看看我的程序嗎？我已經調試了大約4個小時，仍然找不到錯誤。也許，我做了一件非常愚蠢的事情，我沒有注意到。

（只是爲了您的通知，在FFT功能，爲了看到底發生了什麼，我也增添了許多COUT線）

謝謝你幫我出去！

Answer 1

更換

F.push_back(samples[0]); 

它增加了一個元素F - 這已經有一個元素，零 - 與

F[0] = samples[0]; 

你與一個具有零矢量結束你的遞歸第一個元素和第二個樣本。
由於您以後只使用第一個元素，所有內容都將變爲零。

Answer 2

我不能立即說出你的問題是什麼，但你至少有一個嚴重的錯誤。這裏：不是事實

float a = ceil(log(num)/log(base))-floor(log(num)/log(base)); 
if (a==0){ 
    return 1; 
} 

其他一個有效的比較將0.0f，你不應做一些浮點運算，並期待值是一個具體的數字。如果您期望0.0，它可能實際上是0.00000000132。因爲那是how floating point numbers work。不幸的是，他們總是有這樣的工作。至於你的FFT計算，我會更仔細地看看你分配的線F[k]和F[k+M]。在那些，你乘以ODFT[k]一個值，我的懷疑是ODFT[k]可能0那裏。