前回までポンコツなベースバンドフィルタ例のアナログ・フィルタの場合を計算するためにIIRフィルタにまで立ち入って大分遠回りをしました。今回からは「ポンコツでない」例を計算するためにFIRフィルタを使ってみたいと思います。しかし信号処理素人の老人はFIRフィルタの設計などできる気がしません。どうしたらよいの？
※「手習ひデジタル信号処理」投稿順 Indexはこちら
※Windows11上で、Scilab6.1.1およびScilab上のツールボックス Scilab Communication Toolbox 0.3.1（以下comm_tbx)を使用させていただいとります。
ScilabのFIRフィルタ設計用関数
HELPを調べてみると、以下の４つのFIRフィルタ設計用の関数を見つけることができました。
-
- fsfirlin
- wfir、およびwfir_gui
- ffilt
- eqfir
４つもあるなら、きっといろいろ出来るに違いない、と思いましたが、信号処理素人の悲しさ、どれをどう使ったらよいのか皆目見当がつきませぬ。HELPファイルを手掛かりにそれぞれの特徴をみてみると以下のようでした。
-
- fsfirlin
周波数標本法（Frequency sampling method）というものを使い、線形位相のフィルタを設計する関数みたいです。目標とする周波数応答の「標本（ベクトル）」を貰ってフィルタ係数を計算する仕組みです。オプションフラグとして１型、２型などというフラグあるのですが、この関数の文脈でのその意味については何も言及ありません。今回はまずこれを実際に使ってみたいと思います。
-
- wfir、およびwfir_gui
「GUIもあるじゃん、やったね」と思ったこの関数はウインドウィング法（ WINDOWING TECHNIQUES）を使って、やはり線形位相のフィルタを設計する関数みたいです。これは次回かね。
-
- ffilt
この関数は以前にも使った朧げな記憶があるのですが、忘却力の老人には定かではありませぬ。FIRのローパスフィルタの係数を計算してくれる関数なのですが、LPFとあるだけで特に算法の記載などは見当たりませぬ。
-
- eqfir
「FIRフィルタのミニマックス近似」を求める関数なのだそうな。ううむ困った。
fsfirlinを実習
訳も分からず素人老人が関数を使ってみました。「見る前に飛べ」ってか。まあHELPファイルに処理例書いてあるので真似っこです。
まず「周波数応答標本のベクトル」を定義しないとならんのです。処理例をみると長さ64点のベクトルで定義してますが、サンプルの末尾を正規化周波数0.5として計算しているみたいです。ということはベクトルの長さはテキトーで良いのか？
処理例はバントパス・フィルタみたいですが、通過域と阻止域の間が「ほぼほぼ絶壁」急峻な例と、途中に１点挟んでちょいと勾配を緩和した例の２つの標本について計算してました。そのままコピペ。標本の波形をプロットする指令はこんな感じでよいかね。
hd=[zeros(1,15) ones(1,10) zeros(1,39)];
hd2=hd;
hd2(15)=.5;hd2(26)=.5;
hd_pas=1/prod(size(hd))*.5;
fg_pas=0:hd_pas:.5;
plot2d([1 1].*.fg_pas(1:64)',[hd' hd2'], [2,5]);
legend('hd', 'hd2');
上記で行ったプロットが以下に。青のhdの方が「ほぼほぼ絶壁」で赤のhd2が緩和したやつです。
まあねえ、標本波形をこうやって定義すれば良い、といわれると素人老人にも簡単にできるのでありがたいっす。でも後の計算はブラックボックスだけれども。
上記２つの「標本」に対して、訳の分からない「１型」「２型」のオプションプラグを組み合わせ、fsfirlin関数に合計４種類のフィルタを設計してもらいましたぞ。コマケーことは指定の必要がなく、引数２個で処理してくれるので素人には易しいですが、「周波数標本法」の原理を知らないと何がなんだか分かりませぬ。その処理が以下に。
hst11=fsfirlin(hd,1);
hst12=fsfirlin(hd,2);
hst21=fsfirlin(hd2,1);
hst22=fsfirlin(hd2,2);
pas=1/prod(size(hst11))*.5;
fg=0:pas:.5;
plot2d([1 1 1 1].*.fg(1:257)',[hst11' hst12' hst21' hst22'], [2, 3, 5, 6]);
legend('hst11', 'hst12', 'hst21', 'hst22');
生成したフィルタの結果はhstxyという名の変数群に格納されとります。xには標本波形のhdを使ったら1、hd2が2が入ります。yには「１型」なら１「２型」なら２が入ります。
計算後確認すると、以下のようにフィルタのベクトル長は全て同じです。奇数ね。
肝心の応答波形が以下に。
標本hdの方は「急峻すぎる」目標波形だったので急峻ではあるけれどリンギングがはなはだしい？緩和したhd2を指定した方がリンギングが抑えられていて良い感じに見えます。知らんけど。１型と２型の差は周波数方向への微妙なずれに現れているみたいだけれども。これから算法の違いが分かるのか？素人老人には過ぎたる疑問だ。