calc_fourier.c

/*-----------------------------------------------------------------------------

離散フーリエ変換と離散フーリエ逆変換

discrete Fourier transform and discrete Fourier inverse transform

calc_fourier.c

-----------------------------------------------------------------------------*/

#include "FT_surrogate.h"

void fourier( int num, complex *x, complex *F )

{

int k, n;

double real_sum;

double imag_sum;

for(k = 0; k < num; k++){

real_sum = 0.0;

imag_sum = 0.0;

for( n = 0; n < num; n++){

/* 実部の計算 */

real_sum += ( x[n].real * cos(2 * M_PI * n * k/(double)num)

+ x[n].imag * sin(2 * M_PI * n * k/(double)num) );

/* 虚部の計算 */

imag_sum += ( -x[n].real * sin(2 * M_PI * n * k/(double)num)

+ x[n].imag * cos(2 * M_PI * n * k/(double)num) );

}

F[k].real = real_sum;

F[k].imag = imag_sum;

}

}

void inv_fourier( int num, complex *x, complex *F )

{

int k, n;

double real_sum;

double imag_sum;

for(n = 0; n < num; n++){

real_sum = 0.0;

imag_sum = 0.0;

for( k = 0; k < num; k++){

/* 実部の計算 */

real_sum += F[k].real * cos( 2 * M_PI * n * k/(double)num )

- F[k].imag * sin( 2 * M_PI * n * k/(double)num );

/* 虚部の計算 */

imag_sum += F[k].real * sin( 2 * M_PI * n * k/(double)num )

+ F[k].imag * cos( 2 * M_PI * n * k/(double)num );

}

x[n].real = real_sum/(double)num;

x[n].imag = imag_sum/(double)num;

}

}