2016-09-05 14 views
0

FFTとバックトランスフォームを計算して、同じ動作をしているかどうかをチェックしたいと思います。私は* 4 4 * 4行列でそれをテストしようとした私のコードには大きな3Dマトリックスのアプリケーションがあり、ここに私のコードですCでの3D FFTと逆FFTの計算

`

#include <stdio.h> 
#include <stdlib.h> 
#include <complex.h> 
#include <time.h> 
#include <math.h> 
#include <fftw3.h> 


int main() 
{ 
int N = 4; //Dimension of matrix 
unsigned int seed = 1; 
double *in = (double*)malloc(sizeof(double)*N*N*N); 
fftw_complex *out = fftw_malloc(sizeof(fftw_complex)*N*N*N); 
double *out1 = (double*)malloc(sizeof(double)*N*N*N); 

fftw_plan plan_backward; 
fftw_plan plan_forward; 

srand (seed); 

for(int i = 0; i < N; i++) 
{ 
    for(int j = 0; j < N; j++) 
    { 
     for(int k = 0; k < N; k++) 
     { 
      in[i*(N*N) + j*N + k] = rand (); 
     } 
    } 
} 

printf(" Given matrix in\n"); 


for(int i = 0; i < N; i++) 
{ 
    for(int j = 0; j < N; j++) 
    { 
     for(int k = 0; k < N; k++) 
     { 
      printf("%f\n", in[i*(N*N) + j*N + k]); 
     } 
    } 
} 


printf("\n"); 

plan_backward = fftw_plan_dft_r2c_3d (N, N, N, in, out, FFTW_ESTIMATE); 

fftw_execute (plan_backward); 

fftw_destroy_plan (plan_backward); 

printf("out matrix\n"); 

for(int i = 0; i < N; i++) 
{ 
    for(int j = 0; j < N; j++) 
    { 
     for(int k = 0; k < N; k++) 
     { 
      printf("%f\t%f\n", creal(out[i*(N*N) + j*N + k]), cimag(out[i*(N*N) + j*N + k])); 
     } 
    } 
} 

printf("\n"); 

plan_forward = fftw_plan_dft_c2r_3d (N, N, N, out, out1, FFTW_ESTIMATE); 

fftw_execute (plan_forward); 

fftw_destroy_plan (plan_forward); 

printf("out1 matrix\n"); 


for(int i = 0; i < N; i++) 
{ 
    for(int j = 0; j < N; j++) 
    { 
     for(int k = 0; k < N; k++) 
     { 
      printf("%f\n", out1[i*(N*N) + j*N + k]); 
     } 
    } 
} 

fftw_free(in); 
free(out); 
fftw_free(out1); 

return 0; 

}` 

どうやら私の変換結果は同じではありません。何がうまくいかないのか分かりませんか?

答えて

0

FFTが正規化されていません。入力と出力の間には一定の要因があります。

これらの変換は、それほどC2R続いR2C変換(またはその逆)、非正規化されている外観をhere

を取るは、ある要素、実データの数でスケーリングされた元のデータをもたらします実データの(論理)次元の積。

したがって、因子はN * N * Nである必要があります。データをこの要素で分けるだけで、入力と同じデータを取得する必要があります。

+0

ありがとうございました:)ちょうど明確化私は逆変換されているデータにのみこの正規化を行う必要がありますか?複雑なデータに対して最初にc2rを適用し、r2c /(N * N * N)を返すとしますか? – verito

+0

フーリエ変換は線形なので、 'c2r - > 1/N^3 - > r2c'、' 1/N^3 - > c2r-> r2c'などどこでも正規化を行うことができます。 –

0

FFTW3のマニュアルに記載されているように:

これらの変換は、非正規化されているので、C2R続いR2C変換(またはその逆)は、実際のデータ要素の数でスケーリングされた元のデータをもたらしますすなわち、実データの(論理)次元の積である。

本当の次元数は、あなたのケースで4^3であり、その数で最初の結果115474520512を分割することは、第1の入力115474520512/(4^3) = 1804289383をバック与えます。

+0

ありがとう私は答えがあります:) – verito