C言語でのFFT(構造体とポインタの使用)

質問

今、VC++2008　Express　Editionを使用して先日、作成したFFTのプログラムを構造体とポインタを使用して書き換えました。
ビルド＆コンパイルは通ったのですが、実行するとエラーが起きます。何がおかしいのでしょうか。
ソースコードを掲載しますので、教えてください。
//main.c
#include <stdio.h>
#include "stdlib.h"
#include "fft.h"
#include <math.h>
#include "data.h"
#define PI 3.1415926
FILE *ifp, *ofp;
struct DATA f;
double dT, dF;
double *amp;
int main (void)
{
//波形の生成と配列入力
for(i=0;i<1024;i++){
dumx[i]=(double)i/1024;
dumy[i]=sin(2.0*PI*50*(double)i/1024);
}
DN=1024;
//FFT用データ配列に配列のコピー
time = (double*)calloc(DN,sizeof(double));
f.Re =(double*)calloc(DN,sizeof(double));
f.Im =(double*)calloc(DN,sizeof(double));
for(i=0;i<DN;i++){
time[i]=dumx[i];
f.Re[i]=dumy[i];
f.Im[i]=0.00;
}
//サンプリング周期
dT = time[1];
//サンプリング周波数
dF = 1.0/dT;
//FFT
FFT(f);
//大きさを出力するための配列の確保
amp=(double*)calloc(DN,sizeof(double));
for(i=0;i<DN;i++){
amp[i]=sqrt(Sqr(f.Re[i])*Sqr(f.Im[i]));
}
//FFT用データ配列の初期化
free(f.Re),f.Re=NULL;
free(f.Im),f.Im=NULL;
//出力データの書き込み
ofp=fopen("c:\\make\\data\\output.csv","w");
for(i=0;i<DN/2;i++){
fprintf(ofp,"%lf,%lf\n",i*dF, amp[i]);
}
fclose(ofp);
return 0;
}
//fft.h
extern void FFT(struct DATA);
#define Sqr(x) ((x)*(x))
data.h
int g, h, i, j, k, l, m, n, p, q, DN;
double a, b, tmp;
double dumx[1024];
double dumy[1024];
double *time;
struct DATA {
double *Re;
double *Im;
};
struct ROT {
double *c;
double *s;
};
//FFT解析用ソースコード
//analysis.c
#include "data.h"
#include <math.h>
#include <stdlib.h>
#define PI 3.1415926
void Table( struct ROT);
void swap( double*, double*);
void FFT(struct DATA in)
{
struct ROT r;
r.c=NULL;
r.s=NULL;
in.Re=(double*)calloc(DN,sizeof(double));
in.Im=(double*)calloc(DN,sizeof(double));
n=DN;
m = (int)(log10(n)/log10(2));
//回転因子テーブルの作成
Table(r);
l=n,h=1;
//バタフライ演算
for(g=1;g<=m;g++){
l/=2,k=0;
for(q=1;q<=h;q++){
p=0;
for(i=k;i<=l+k-1;i++){
j=i+l;
a=in.Re[i]-in.Re[j], b=in.Im[i] - in.Im[j];
in.Re[i] = in.Re[i] + in.Re[j], in.Im[i] = in.Im[i] + in.Im[j];
if(p==0){
in.Re[j]=a,in.Im[j]=b;
}else{
in.Re[j] = a * (r.c[p]) + b * (r.s[p]), in.Im[j] = b * (r.c[p]) - a * (r.s[p]);
}
p+=h;
}
k+=l+l;
}
h+=h;
}
j=n/2;
for(i=1;i<=n-1;i++){
k=n;
if(j<i){
//ビットリバース
swap(&(in.Re[i]),&(in.Re[j]));
swap(&(in.Im[i]),&(in.Im[j]));
}
k/=2;
while(j>=k){
j=j-k,k /=2;
//無限ループ回避(ここを消すと無限ループに陥ります。)ここから
if((j==0)&&(k==0)){
break;
}
//ここまで
}
j += k;
}
}
void Table( struct ROT w)
{
w.c = (double*)calloc(DN/2,sizeof(double));
w.s = (double*)calloc(DN/2,sizeof(double));
a=0.00, b = 2.0 * PI /DN;
for(i=0;i<DN/2;i++)
{
(w.c[i])=cos(a+i*b), (w.s[i])=sin(a+i*b);
}
}
void swap(double *var1, double *var2)
{
tmp = *var1, *var1=*var2, *var2=tmp;
}

回答No.4

redfox63 · ベストアンサー

nda23氏へ

yf491224氏は『C』でやりたいのかも
メインやサブのファイル名も『.c』ですし ・・・

『Data.h』での変数宣言でこけそうな気がします
『analysis.c』『main.c』両方に同名の変数宣言を行うことになるので
おかしなことになるかと

実体を main.c(またはanalysis.c)で宣言して data.hでは 各変数を externを使って外部変数ですとした方がいいように思います

関数 FFTの中で　in.Reやin.Imを再確保していますがこれでは
引数で元データを与える意味がありません

お礼

＞yf491224氏は『C』でやりたいのかも
＞メインやサブのファイル名も『.c』ですし ・・・
あはは。C言語のほうが私には、わかりやすいのでCを扱っていました。
C＋＋に関しては、基礎的なところは本を持っていますが、それ以外のところはさっぱりです。
＞『Data.h』での変数宣言でこけそうな気がします
＞『analysis.c』『main.c』両方に同名の変数宣言を行うことになるの＞実体を main.c(またはanalysis.c)で宣言して data.hでは 各変数を ＞externを使って外部変数ですとした方がいいように思います
＞関数 FFTの中で　in.Reやin.Imを再確保していますがこれでは
＞引数で元データを与える意味がありません
わかりました。ご指摘のとおりに変更して再度やってみます。
回答ありがとうございます

回答No.3

nda23

せっかくなんで、Ｃ++っぽく
void Table( struct ROT* );
　　　　↓
void Table( struct ROT &);
これで、呼び出し側がTable(r)と書いても参照渡しになります。
関数側もw.c = ～ のように記述します。つまり、現在のままで
良いということです。

あとはredfox63さんの意見を参考になさってください。

回答No.2

redfox63

すみません ・・・ これは誤りですね
Im側の計算もなされていました m(__)m
> in.Im側の計算がされていないようですので出力はすべて0になってしまうようです

関数FFTも ポインタで参照渡しにしないと FFTでの計算結果が反映されませんよ

回答No.1

redfox63

関数Tableに対して struct ROTで引数を渡していますが
これだと『値渡し』ですので ここで確保した構造体ROTのメンバーc,sは呼び出し元に反映されません
したがって r.c[p]などとアクセスすると不正メモリーのアクセスで停止するでしょう

void Table( struct ROT* );
と宣言
void Table( struct ROT* w)
{
　　// ポインター用に『.』を『->』に変更
　　w->c = (double*)calloc(DN/2,sizeof(double));
　　w->s = (double*)calloc(DN/2,sizeof(double));
　　a=0.00, b = 2.0 * PI /DN;
　　for(i=0;i<DN/2;i++)
　　{
　　　　(w->c[i])=cos(a+i*b), (w->s[i])=sin(a+i*b);
　　}
}

また、 ROT用に確保したメモリーの解放が見当たりません
in.Im側の計算がされていないようですので出力はすべて0になってしまうようです
log10の引数に intを使うオーバーロードがないので(double)などでキャストしましょう

お礼

なるほどです。ご指摘のとおりに修正して実行してみます。

回答ありがとうございます。