フーリエ変換のC言語プログラムについて

正弦波（およびガウス性雑音)をフーリエ変換(離散)→逆フーリエ変換するというプログラムを組みました。正弦波をフーリエ変換すると実部は２回ピークがくるはずなのですが、すべて「0.000000」または「-0.000000」と表示されてしまいます。虚部は正常なようで実装の仕方もさほど違わないので、何が問題なのかわからずにいます。念のためコードはすべて載せますが、該当箇所は関数Fourierの
fp = fopen("reX.txt", "w");//書き込む
あたりです。問題点を教えていただけないでしょうか。お願いします。


//gauss.txt, sin.txt:発生させたガウス性雑音＆正弦波
//reX, imX:フーリエ変換の実部と虚部
//re-1, im-1:逆フーリエ変換の実部と虚部

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

#define PI3.14159265358979323846
#define N256

//n:長さ, w:角周波数, p:位相(phase), a:振幅
void SinCurve(int n, double w, double p, double a)
{
FILE *fp;
double x;
int t;

fp = fopen("sin.txt", "w");//書き込むので"w"
if(fp == NULL)
{
printf("file open error\n");
} else
{
for(t = 0; t < n; t++)
{
x = a * sin( w*(double)t + p );
fprintf(fp, "%f\n", x);
}
}

fclose(fp);
}

//n:長さ, s:分散, m:平均
void Gauss(int n, double s, double m)
{
FILE *fp;
double x, x1, x2, y1;
int t;

srand((unsigned) time(NULL));

fp = fopen("gauss.txt", "w");//書き込むので"w"
if(fp == NULL)
{
printf("file open error\n");
} else
{
for(t = 0; t < n; t++)
{
x1 = ( (double)rand() + 1.0 ) / ( (double)RAND_MAX + 2.0);
x2 = ( (double)rand() + 1.0 ) / ( (double)RAND_MAX + 2.0);
y1 = pow(-2.0*log(x1), 0.5) * cos(2.0*PI*x2);
y1 = s * y1 + m;
fprintf(fp, "%f\n", y1);
}
}

fclose(fp);
}

//ファイル名sのデータをフーリエ変換し、実部と虚部をreX, imXに保存
void Fourier(int num, char *s)
{

FILE *fp;
int k, n;
double largeX, x[N+100], t;

fp = fopen(s, "r"); //読み込み
if(fp == NULL)
{
printf("file open error\n");
} else
{
//printf("%s\n", s);
for(k = 0; k < num; k++)
{
fscanf(fp, "%lf", &x[k]);
printf("x[%d]=%f\n", k, x[k]);
}
}

fp = fopen("reX.txt", "w");//書き込む
if(fp == NULL)
{
printf("file open error\n");
} else
{
for(k = 0; k < num; k++)
{
largeX = 0.0;
t = 2.0*PI*(double)k / (double)N;
for(n = 0; n < num; n++)
{
largeX += x[n] * cos((double)n*t);
//printf("%f\n", largeX);
}
fprintf(fp, "%f\n", largeX);
printf("reX[%d]=%f\n", k, largeX);
}
}

fp = fopen("imX.txt", "w");//書き込む
if(fp == NULL)
{
printf("file open error\n");
} else
{
for(k = 0; k < num; k++)
{
largeX = 0.0;
t = 2.0*PI*k / (double)N;
for(n = 0; n < num; n++)
{
largeX -= x[n] * sin(n*t);
}
fprintf(fp, "%f\n", largeX);
}
}

fclose(fp);
}

void InverseFourier(int num)
{
FILE *fp;
int k, n;
double a[N+100], b[N+100], x, t;
//a:reX, b:imX

fp = fopen("reX.txt", "r"); //読み込み
if(fp == NULL)
{
printf("file open error\n");
} else
{
for(k = 0; k < num; k++)
{
fscanf(fp, "%lf", &a[k]);
//printf("a[%d]=%f\n", k, a[k]);
}
}

fp = fopen("imX.txt", "r"); //読み込み
if(fp == NULL)
{
printf("file open error\n");
} else
{
for(k = 0; k < num; k++)
{
fscanf(fp, "%lf", &b[k]);
//printf("b[%d]=%f\n", k, b[k]);
}
}

fp = fopen("re-1.txt", "w"); //読み込み
if(fp == NULL)
{
printf("file open error\n");
} else
{
for(n = 0; n < num; n++)
{
x = 0.0;
t = 2.0*PI*(double)n / (double)N;
for(k = 0; k < num; k++)
{
x +=a[k] *cos(k*t) - b[k] *sin(k*t);
}
x /= (double)N;
fprintf(fp, "%f\n", x);
//printf("x[%d]=%f\n", n, x);
}
}
/*
fp = fopen("im-1.txt", "w"); //読み込み
if(fp == NULL)
{
printf("file open error\n");
} else
{
for(n = 0; n < num; n++)
{
x = 0.0;
for(k = 0; k < num; k++)
{
t = 2.0*PI*(double)k / (double)N;
x = x + a[k] *sin(n*t) + b[k] *cos(n*t);
}
x /= (double)N;
fprintf(fp, "%f\n", x);
}
}
*/

fclose(fp);
}

int main(void)
{
SinCurve(N, PI/8.0, 0.0, 1.0);

//Gauss(N, 1.0, 0.0);

Fourier(N, "sin.txt");

//Fourier(N, "gauss.txt");

InverseFourier(N);

return 0;
}

No.2 ベストアンサー 回答者： asdf1234151 回答日時： 2014/01/06 00:56

元のデータがsin関数なら、虚部だけしか出ないので正常では?

fopen()に対応したfclose()が必要なのはその通りなので直すべきですが。

No.1 回答者： honor 回答日時： 2014/01/05 23:19

fclose()が足りないんじゃないですか

この回答への補足

fcloseを逐一つけてみましたが、特に変わりませんでした。関数の最後につけさえすれば大丈夫なようです。。。

補足日時：2014/01/06 00:08