FFTによるフィルタの実験
FIRフィルタの演算量は、普通はタップ数Nに比例すると考えられています。ところが、視点を変えると演算量はサンプリング周波数Fsの2乗に比例します。つまり、48kHzサンプルから96kHzサンプルに変更すると、演算量は4倍になります。
なぜでしょうか。
この理由はフィルタの特性がFsに対する相対値でなく、物理的な周波数を元に決められるからです。ある周波数をカットオフとして、ある特性を実現すると、その特性はHzとdBを軸とするグラフで表現されます。サンプリング周波数は出てきません。そして、このようなフィルタの特性の良さ*1は、そのインパルスの長さで決まります。この長さは秒であらわされます。つまり、一定の特性のフィルタを作ろうとすると、どのようなサンプル周波数でも、同じ秒数で示される長さのインパルス応答が必要になります。このことから、必然的にサンプル周波数が倍になるときには、インパルス応答の長さを倍にしなければならないことが分かります。サンプル周波数が倍になることでフィルタの適用回数が倍になり、且つ一回のフィルタの演算長が倍になることから、フィルタの演算量は倍の倍、つまり4倍になります。
これが、フィルタの演算量がサンプリング周波数の2乗に比例する理由です。
この問題に対抗するためにFFTを使ったフィルタが用いられることがあります。あらかじめインパルス応答をFFTにかけて周波数応答に変換しておきます。入力信号列を一定ブロックごとにFFTにかけ、先のFFTしたインパルス応答と要素間の積をとります。その積を逆FFTにかけるとあら不思議、フィルタのかかった出力が得られます。
ブロックの継ぎ目を嫌って、信号に三角形の窓関数を適用しつつ、半周期ごとのブロック・オーバーラップをかける方法などもあります。
さて、この方法は周波数フィルタをそのまま実行するので、それほど不思議ではないと思うかもしれません。しかし、私はどうしてもこの方法が時間のシフトに対して不変であると信じられなかったのです。
つまり、
F(a) = w f(a)
とおいたときに*2
F(a+b) = w f(a+b)
という線形性を保ちうるか、疑っていました。
さて、今日はハードを触らないと決め付けましたので、この点をScilabで実験してみたわけですが…どうやら線形性が保たれているようです。とても不思議です。
直感的には保たれない気がするんだけどねぇ。
