FFTの倍速化

例によってFFTライブラリのチューニング。

今日は、前に作った表計算ソフト用のFFTシートを
ちょこっと作り直して、128個の実数を2組同時に
FFT掛けちゃうとどうなるか…っていうシミュレート。

2つの実数データの組（128個×2組）を、片方を実数データ
として、もう片方を虚数データとして入力しFFTを掛け、
出力データを偶関数・奇関数の性質をうまく使うと
分離できるっていうあれ。
http://www.kurims.kyoto-u.ac.jp/~ooura/fftman/ftmn2_1.html#sec2_1

で、これを愚直に表計算シートに反映していきたいところ
なんだけど、最初で躓く…。

データ数が128点なら、N=128ってこと。それなのに、
データは0～127の計128個。/A(N-k）なら最初の要素は
/A(128)だから、129番目の要素ってことになっちゃう
ジャン！

まぁ、入力するデータは循環していることというのが
FFTの大前提だから、ちょっと考えればすぐに答えは
わかるんだけど、オイラ的にはどこかに何かちゃんと
書いてないかな…って思って探してみる。あった↓。
http://homepage2.nifty.com/m_kamada/math/fftmul.htm#parallelfft

ここの(m=0のとき…）という件。

添え字Nのデータと添え字0のデータは一緒ってこと。
やっぱね。

ということで、表計算シートのFFTに加味してみる。
早速いくつかの波形を入力してみると…　やっぱり
ちゃんと2つの波形の周波数成分に分離できる。ｏｋ！


さて、それが出来たところで何を確認したいのかって
いうと…

実数側と虚数側にどんな波形を入れた場合、丁度山と
山がかち合って、レンジが広がっちゃうのかってこと。

さすがに周波数があまりに違っちゃってるとかち合わ
無いのは想像できるんだけど、同じ周波数のデータを
ためしに入れてみると…　普通に計算できちゃう…。

うーん。位相を90度ずらしてみる…。出た。

どうやら位相を90度ずらしてみる（sinとcosとか）と、
計算途中の数値に実数側、虚数側それぞれに入力した
データが丁度ぶつかってしまうっぽい。

その点だけ気をつけておけば、まぁ今回の用途では
実用上問題無さそう。

で、イザぶつかっちゃった場合の対処方法（処理方法）
なんだけど、とりあえずレンジを越えた分は
　　　「一定のレンジに圧縮しちゃう」
っていう無理やりな方法にしたいんだけど、例によって
数値計算の誤差で問題になるのは「差を取る時」
なので、FFT結果から実数側・虚数側に分離するときに
誤差を拡大してしまう恐れがあるわけ。

で、そういう処理をさせてみた場合を表計算上で
シミュレートしてみると…　まぁ、それなりに
大丈夫っぽい答えが出てきた。よし、よし。


例の出力レンジを広げる（S/N比やダイナミックレンジを
もうちょっと広げたい）っていう話と絡めて、破綻しない
方法に落とし込んで行きたい所。

そこまで出来ればカスタマイズ版のFFTライブラリと
してはバッチリって感じじゃぁないかなぁと。

こうやって2つのデータ群をいっぺんにFFT掛けられちゃう
なら、128点FFT（1回あたり約62000クロック）でも4弦を
4chのまま取り込むことが出来るし、なにより2ch同時だから
2回のFFT計算でｏｋ。
計算上は1秒間に130回近く処理できるはず。まぁ、2chに
分離する処理とか、そのほかメインルーチン側の処理など
もアレコレあるので、実際は100回程度かな。

それでも遅延は1/100秒程度と考えれば、せいぜい3m先
まで手を伸ばして弾いている感覚。怪物くんウクレレ
とまでは行かないんじゃないかな。