手習ひデジタル信号処理(65) Scilab、FIRフィルタの係数求めるついでに周波数特性

Joseph Halfmoon

前回はFIRローパスフィルタの係数を求め、求めた係数を使って「単位」インパルス入力信号をフィルタして時間応答を求め、さらにそれをFFTして周波数特性のプロットをいたしました。まどろっこしいデス。フイルタ係数を求めるときに「ついでに」周波数特性のプロットも見て、それからOKしたいです。知らんけど。今回は微妙な改善?

※「手習ひデジタル信号処理」投稿順 Indexはこちら

参照させていただいているScilabの日本語版HELPページ

以下のページに「FIRおよびIIRフィルタの振幅」を求めるための関数 frmag についての説明があります。

frmag  FIRおよびIIRフィルタの振幅

この関数を使えば、前回の結果から周波数振幅特性が得られるんでないかい。

微妙な改善後の関数

Scilabのffilt関数に「1枚かぶせて」サンプリング周波数Fsと実際の(正規化でない)カットオフ周波数で指定できるようにしたフィルタ係数算出関数を前回作成してみました。今回、振幅周波数特性も同時に得られるように改善してみたものが以下に。

// Calculate and Plot FIR Lowpass Filter
// Fs: sampling frequency [Hz]
// nFlt: number of filter order
// Fc: cutoff frequency [Hz]
// Plt: 1 for plot, 0 for no plot
function coeff=genFIRLP2(Fs, nFlt, Fc, Plt)
    Fcnorm = Fc / Fs;
    coeff = ffilt("lp", nFlt, Fcnorm);
    if Plt == 1 then
        [hzm,fr]=frmag(coeff,256)
        clf();
        plot2d(fr*Fs, 20*log10(hzm))
        xtitle('Magnitude plot','frequency [Hz]','magnitude [dB]');
    end
endfunction
上記の関数を実行してみた結果

前回同様のサンプリング周波数10kHz、FIRフィルタの次数100、カットオフ周波数500Hzで振幅周波数プロットあり(1)のときの関数呼び出し例が以下に。

Fs=10000;
coeff=genFIRLP2(Fs, 100, 500, 1);

そしてそのとき出力されるプロットはこんな感じ。

filtMagPlot

前回のは縦軸が素のままだったですが、今回はdBに改まっております。いい感じでないかい。。。ホントか?

手習ひデジタル信号処理(64) Scilab、FIRフィルタの係数求めてフィルタしてみるの回 へ戻る

手習ひデジタル信号処理(66) Scilab、IIRフィルタ設計用関数、周波数特性表示付き へ進む