2016-03-03 29 views
1

我想在Armadillo包的矩陣中使用FFTW。我需要在2D矩陣上執行N個獨立的FFT。按照FFTW和其他資源在線手冊,我有以下代碼:Armadillo與FFTW的接口

#include <armadillo> 
#include <iostream> 
#include "fftw3.h" 

using namespace std; 
using namespace arma; 

fftw_plan fplan; 

int main(int argc, char *argv[]) 
{ 
    cx_mat psi; 
    int NCOL=2, NROW=4; 
    int frank=1; 
    int howmany=NCOL; 
    int n[]={NROW};   
    int idist=NROW, odist=NROW; 
    int istride=1, ostride=1; 
    int *inembed=n, *onembed=n; 

    psi.resize(NROW, NCOL); 
    psi(0,0)=cx_double(1,1); psi(0,1)=cx_double(1,1); 
    psi(1,0)=cx_double(2,1); psi(1,1)=cx_double(1,2); 
    psi(2,0)=cx_double(3,1); psi(2,1)=cx_double(1,3); 
    psi(3,0)=cx_double(4,1); psi(3,1)=cx_double(1,4); 
    cout << psi << endl; // the output is correct at here 

    fftw_complex* in = reinterpret_cast<fftw_complex*>(psi.memptr()); 

    fplan = fftw_plan_many_dft(frank, n, howmany, 
          in, inembed, istride, idist, 
          in, onembed, ostride, odist, 
          FFTW_FORWARD, FFTW_MEASURE); 

    cout << endl << "after fftw plan: " << endl << psi << endl; // all zeros 

    // fftw_execute(fplan); // output zeros again if execute 
    cx_double* A_mem=psi.memptr(); 
    cout << endl << "by reference: " << endl << *A_mem << " " << *(A_mem+1) << endl; 

    return 0; 
} 

該代碼編譯沒有任何錯誤。但是,如果我在fftw_plan_many_dft之後運行代碼,它將輸出全零。如果我刪除fftw_plan_many_dft,則輸出是正確的。我甚至不執行fftw計劃,爲什麼只通過制定計劃來清理我的數據?

+1

使用'FFTW_MEASURE'標誌,緩衝區實際上用於調整FFT算法。您應該先創建計劃*,然後初始化輸入數據。 –

+1

哦。解決了這個問題。謝謝。 – user1285419

+0

沒問題 - 評論已轉換爲現在的答案,爲將來有類似問題的讀者帶來益處。 –

回答

2

對於FFTW_MEASURE標誌,緩衝區實際上用於調整FFT算法。您應該先創建計劃,然後初始化輸入數據。

注意FFTW_ESTIMATE出現此行爲,所以你也可以只改變這個標誌在當前的代碼來解決這個問題,但你的代碼可能會運行速度較慢。