高速フーリエ変換からパワースペクトルを求める (ID:52779)

名前

ホームページ(ブログ、Twitterなど)のURL (省略可)

本文

このプログラムを Visual C++ で動くようにしたいのですか。
それとも理論として log( v^2 ) でパワーになるかどうかを
質問しているだけなのですか？

後者ならば、
W = I*V = I^2*R; 
power(dB単位) = log(W)    
なんで、
時間軸であろうと、周波数軸であろうと、パワーは
それでいいのでは。
というか、すでに log10(...) が入っているからいいのでは。

前者ならば、コンパイラオプションをざっと眺めたところ、
クラシックなスタイルのＣ言語を通すオプションが見当たりません。
（見落としていたらすまない）
まず、関数宣言のところから直す必要がありますね。
というか、gcc --traditional でコンパイルしたほうがいいのでは。

見たところ結構問題の多いプログラムです。
fopen() に対する fclose() がありません。
main() からの脱出が
exit() だったり return; だったりで統一していませんね。
make_bitrev()が見当たりません。
N の正負によって振り分ける inverse もやめたほうがよいですね。
ftt() と invftt() にしたほうがよいですよ。
今の時代この程度のルーチンでプログラムサイズをケチる必要ありません。
全体の流れのなかで、
まったく内容の同じ sintbl,bitrev を２度作成しています。
テーブルの作成部分は fft() から呼ぶのではなく、main() で
作成すべきでないでしょうか。
せっかく計算した x3,y3 出力されていません。
最後の fprintf で書式の %6.3f の数と引数の数も合っていません。
使用する変数の精度を変更できるように
#define FTYPE fload
//#define FTYPE double
で FTYPE x1[N] などとしておいた方が便利です。
...というかこのくらいの精度はどうでもいいのかな。

make_sintbl() の計算で
sin() の計算が１回しかない。
私の手元にあるソースでは、
まじめに

double dTh = - atan(1.0) * 8 / (double)n;
  for( i = 0; i < n; i ++ )
  {
     double dArg = dTh * i;
     g_dTblSin[i] = sin( dArg );
     g_dTblCos[i] = cos( dArg );
  }

やってます。
初期化時に１回実行されるだけだから、
ここは時間をかけてもさして問題にならないです。

これのオリジナルソースっておそらく
http://www.edu.ics.saitama-u.ac.jp/~j9802ah/prog2/4/fft.c
でしょ？
なので、補足。(make_bitrev())
google のキャッシュから拾ってきました。

static void make_bitrev(int n, int bitrev[])
{
    int i, j, k, n2;

n2 = n / 2;  i = j = 0;
    for ( ; ; ) {
	bitrev[i] = j;
	if (++i >= n) break;
	k = n2;
	while (k <= j) {  
	    j -= k;  k /= 2;  
	}
	j += k;
    }
}

ハミング窓とか不要なのかな？

せっかく Visual C++ なんだから 
complex<doble> 型使おうね...って処理が遅くなるかな。

最後にFFT関連でお勧めのサイトを紹介しておきます。
http://momonga.t.u-tokyo.ac.jp/~ooura/fftman/

←解決時は質問者本人がここをチェックしてください。

戻る