FFT その8 - Maximaでこうぞうりきがく

前回に引き続き，C言語のFFTプログラムのお話です
コンパイラにはBorland C++ Compiler 5.5 (BCC32) を使いました
環境依存する様なものは書いてないので何でコンパイルしても問題ない・・・筈・・・
fft.zip

fft.exe　プログラム本体
fft.dat　入力データ（テキストファイル）

入力データは気象庁主な地震の強震観測データで公開されている，ある地域のUD方向の地震応答加速度振幅を使いました
サンプル数は2のべき乗でしか受け付けないので適当な所（2^15）で切り捨ててます

【fft.dat】

0.01
0.000	0.000
-0.002	0.000
-0.004	0.000
-0.005	0.000
-0.007	0.000
-0.011	0.000
-0.013	0.000
.
.
.

1行目はサンプリング間隔[s]
2行目以降は複素数振幅[cm/s2]で，これが2^15 = 32768個並びます

fft.exeを実行するとDOS窓とgnuplotウインドが開きます
（DOS窓で何かキーを押すとどちらも閉じます）

【gnuplot】

gnuplotは横軸に周波数[Hz]，縦軸に加速度振幅[cm/s2]として対数軸でプロットします

またfft.exe実行後に以下のテキストファイルを生成します

fft.out　FFT計算結果出力データ
gnu.out　gnuplot描画データ

【fft.out】

 FFT : fast Fourier transform with Cooley-Tukey algorithm

 number of samples	: 32768
 sampling interval	: 0.01 [s]
 sampling time		: 327.68 [s]
 sampling rate		: 100 [Hz]
 nyquist frequency	: 50 [Hz]
 resolution of freq	: 16384
 minimum frequency	: 0.00305176 [Hz]

 TIME	INPUT			FFT
	 real	    image	 real	    image
    0 |  0.000e+00  0.000e+00 |  7.051e-03  0.000e+00
    1 | -2.000e-03  0.000e+00 | -6.505e-03 -5.252e-03
    2 | -4.000e-03  0.000e+00 |  5.341e-03  1.099e-02
    3 | -5.000e-03  0.000e+00 |  6.526e-03 -1.317e-02
    4 | -7.000e-03  0.000e+00 | -1.425e-02  4.924e-03
    5 | -1.100e-02  0.000e+00 |  1.706e-02  1.154e-02
    6 | -1.300e-02  0.000e+00 |  5.532e-03 -2.305e-02
    7 | -6.000e-03  0.000e+00 | -2.617e-02  4.342e-03
    8 |  0.000e+00  0.000e+00 |  2.115e-02  2.639e-02
.
.
.

4〜10行はFFT計算の主要なパラメータ
11行以降は入力データと計算結果の複素数振幅[cm/s2]で，これが2^15 = 32768個並びます
（16385個以降の計算結果は折り返し分になりますので注意して下さい）

【gnu.out】

 0.000e+00	 0.000e+00
 3.052e-03	 1.672e-02
 6.104e-03	 2.443e-02
 9.155e-03	 2.940e-02
 1.221e-02	 3.016e-02
 1.526e-02	 4.118e-02
 1.831e-02	 4.742e-02
 2.136e-02	 5.305e-02
.
.
.

gnuplot描画用の周波数[Hz]と加速度振幅[cm/s2]で，これが2^(15-1) = 16384個並びます