前々回、一次 LPF を FIR フィルタで作ってみました。100 段ぐらいのディレイを使ってやっとサンプリング周波数の 1 / 100 のカットオフ周波数の 一次 LPF が実現できるのですから IIR に比べて効率が悪いのは事実です。ただ、時々訪問してコメントをくれる「のざわ」さんによるとポリフェーズ構成にすると大幅にハードウェアが削減できるようです。
のざわさんのブログはこちら
温泉気分で♪
http://blogs.yahoo.co.jp/nozawa_onsen_potta
ポリフェーズ構成についてはこちら
ポリフェーズフィルタの基本を知る
http://ednjapan.com/edn/articles/1009/01/news112.html
以下の記事は連載になっていますので、一通り目を通しておくことをお勧めします。
デジタルフィルターとアナログフィルター（1）――デジタルフィルター超入門編
http://ednjapan.com/edn/articles/1403/17/news023.html
さて前々回ではこんなまとめをしました。
FIR フィルタは位相直線なので一般的な伝達関数を作れない。
LPF を作ってもカットオフ周波数が低ければタップ段数が増えてしまう。
LPF はインパルス応答に対して時刻ゼロの時に応答がないので、一ディレイ分遅れた特性になる。
なんとなく当たらずといえども遠からずなのですが、一番最初の位相直線云々は違っていました。
今回は一次 HPF を作ってみてその辺りを解説してみようと思います。
実は前回の一次 LPF の場合でも通過帯域においては位相直線になりますが、通過帯域外だと位相直線になりません。連続系の LPF と同じように位相が -90°どまりになります。たださらに上になると応答遅れの影響が支配的になります。
一次 HPF は微分回路 + LPF の組み合わせの特性になっています。通常の微分回路も実際には同じ構成になっていて LPF のカットオフ周波数を限りなく高くしたもの（実用上微分と呼べる範囲）になっています。伝達関数もそのまんまで、sT / (1 + sT) と分子は微分を表し、分母は一次 LPF です。
サンプリング周波数 10KHz、カットオフ周波数 1Hz、100 段ディレイ の一次 HPF を scilab を使ってやってみます。
前回と同様インパルス応答を調べてそのパラメータを FIR の各係数に入れる方法を使ってみます。Hamming 窓関数は省略します。
スクリプトはこうです。
s=poly(0,'s')
z=poly(0,'z')
dt=1d-4
dires=[]
cutoffrate=1/100
fc=1/dt*cutoffrate
wc=2*%pi*fc
HPF=s/2/%pi/(fc+s/2/%pi)
//HPF の伝達関数
HPFsys=syslin('c',HPF);
clf(1);scf(1);clf;bode(HPFsys,0.01,1e4);
t=0:dt:dt/cutoffrate
scf(2);clf
plot(t,csim('impuls',t,HPFsys),t,csim('step',t,HPFsys),'-.')
dires=[csim('impuls',t,HPFsys)]
DHPF=0
param=dires'
n=size(param,'r')
//param(1)=-sum(dires)
//↑パラメータ修正用
for i=0:n-1
DHPF=DHPF+z^(-i)*param(i+1)*dt
end
DHPFsys=syslin(dt,DHPF);
scf(3);clf;bode(DHPFsys,1e-2,4990)
s 領域でのボーデ線図です。
インパルス応答とステップ応答です。
FIR で作ってみた結果です。
あらあら FIR の結果は変ですね。低域は位相が反転してカットオフ周波数がずれた LPF みたいです。なぜでしょうか。
もう一度インパルス応答を見てみます。
直感的な考察になってしまいますが、HPF は高い周波数成分を通して低い周波数成分はカットします。ですのでインパルスのような立ち上がりのある信号はまずはスルーするはずなのですが、ここではいきなりゼロ以下の数値で始まってしまいます。
最初の方の数値をみてみます。
0.
- 590.05479
- 554.12126
- 520.37603
- 488.68585
- 458.92555
- 430.97761
以下マイナスがずっと続きます。HPF というのはある時間の範囲の平均がゼロになるように働く（DC カット）ので、この応答では平均がゼロになっていません。よって DC レベルを持ってしまうので正しいインパルス応答を示していないのです。
実は scilab のインパルス応答は時刻 0 の数値は誤差を持つ、となっていますので、HPF のような関数のインパルス応答を調べると余計誤差が目立ちます。本来時刻 0 では無限大です。ということは平均でゼロにならなくてはいけないので、時刻 1 以降のすべての値の和もマイナス無限大になるということです。
ここでは FIR というインパルス応答が有限であるフィルタを作ろうとしていますので、応答をサンプリングした値の総和がゼロになるように時刻 0 の値を入れてあげれば、上手くいくかも知れません。
前述のスクリプトの中程にある //param(1)=-sum(dires) の // を外して処理を有効にします。時刻 0 にそれ以降の数値の合算の符号を反転させて入れます。
結果はこうなりました。
今度は上手くいったようです。
FIR フィルタは位相進みが作れないのかな、と漠然と思っていたのですが勘違いですね。ちゃんと作ることが出来ます。
scilab のインパルス応答の誤差に混乱させられるところでした。
scilab：ちょっとまって！インパルスってそもそも理論上でしか存在しないんじゃない！
画蔵 ：そういえばそうかな。
scilab：具体的に計算可能なように信号を定義してくれるならともかく、ありえない信号を定義しておいて計算ソフトになんとかしろったって、それ無理ゲー。
画蔵 ：なるほど一理あるな。
ということで、scilab を責めるのはチト酷なようです。
次回は別の方法を考えてみます。
のざわさんのブログはこちら
温泉気分で♪
http://blogs.yahoo.co.jp/nozawa_onsen_potta
ポリフェーズ構成についてはこちら
ポリフェーズフィルタの基本を知る
http://ednjapan.com/edn/articles/1009/01/news112.html
以下の記事は連載になっていますので、一通り目を通しておくことをお勧めします。
デジタルフィルターとアナログフィルター（1）――デジタルフィルター超入門編
http://ednjapan.com/edn/articles/1403/17/news023.html
さて前々回ではこんなまとめをしました。
FIR フィルタは位相直線なので一般的な伝達関数を作れない。
LPF を作ってもカットオフ周波数が低ければタップ段数が増えてしまう。
LPF はインパルス応答に対して時刻ゼロの時に応答がないので、一ディレイ分遅れた特性になる。
なんとなく当たらずといえども遠からずなのですが、一番最初の位相直線云々は違っていました。
今回は一次 HPF を作ってみてその辺りを解説してみようと思います。
実は前回の一次 LPF の場合でも通過帯域においては位相直線になりますが、通過帯域外だと位相直線になりません。連続系の LPF と同じように位相が -90°どまりになります。たださらに上になると応答遅れの影響が支配的になります。
一次 HPF は微分回路 + LPF の組み合わせの特性になっています。通常の微分回路も実際には同じ構成になっていて LPF のカットオフ周波数を限りなく高くしたもの（実用上微分と呼べる範囲）になっています。伝達関数もそのまんまで、sT / (1 + sT) と分子は微分を表し、分母は一次 LPF です。
サンプリング周波数 10KHz、カットオフ周波数 1Hz、100 段ディレイ の一次 HPF を scilab を使ってやってみます。
前回と同様インパルス応答を調べてそのパラメータを FIR の各係数に入れる方法を使ってみます。Hamming 窓関数は省略します。
スクリプトはこうです。
s=poly(0,'s')
z=poly(0,'z')
dt=1d-4
dires=[]
cutoffrate=1/100
fc=1/dt*cutoffrate
wc=2*%pi*fc
HPF=s/2/%pi/(fc+s/2/%pi)
//HPF の伝達関数
HPFsys=syslin('c',HPF);
clf(1);scf(1);clf;bode(HPFsys,0.01,1e4);
t=0:dt:dt/cutoffrate
scf(2);clf
plot(t,csim('impuls',t,HPFsys),t,csim('step',t,HPFsys),'-.')
dires=[csim('impuls',t,HPFsys)]
DHPF=0
param=dires'
n=size(param,'r')
//param(1)=-sum(dires)
//↑パラメータ修正用
for i=0:n-1
DHPF=DHPF+z^(-i)*param(i+1)*dt
end
DHPFsys=syslin(dt,DHPF);
scf(3);clf;bode(DHPFsys,1e-2,4990)
s 領域でのボーデ線図です。
インパルス応答とステップ応答です。
FIR で作ってみた結果です。
あらあら FIR の結果は変ですね。低域は位相が反転してカットオフ周波数がずれた LPF みたいです。なぜでしょうか。
もう一度インパルス応答を見てみます。
直感的な考察になってしまいますが、HPF は高い周波数成分を通して低い周波数成分はカットします。ですのでインパルスのような立ち上がりのある信号はまずはスルーするはずなのですが、ここではいきなりゼロ以下の数値で始まってしまいます。
最初の方の数値をみてみます。
0.
- 590.05479
- 554.12126
- 520.37603
- 488.68585
- 458.92555
- 430.97761
以下マイナスがずっと続きます。HPF というのはある時間の範囲の平均がゼロになるように働く（DC カット）ので、この応答では平均がゼロになっていません。よって DC レベルを持ってしまうので正しいインパルス応答を示していないのです。
実は scilab のインパルス応答は時刻 0 の数値は誤差を持つ、となっていますので、HPF のような関数のインパルス応答を調べると余計誤差が目立ちます。本来時刻 0 では無限大です。ということは平均でゼロにならなくてはいけないので、時刻 1 以降のすべての値の和もマイナス無限大になるということです。
ここでは FIR というインパルス応答が有限であるフィルタを作ろうとしていますので、応答をサンプリングした値の総和がゼロになるように時刻 0 の値を入れてあげれば、上手くいくかも知れません。
前述のスクリプトの中程にある //param(1)=-sum(dires) の // を外して処理を有効にします。時刻 0 にそれ以降の数値の合算の符号を反転させて入れます。
結果はこうなりました。
今度は上手くいったようです。
FIR フィルタは位相進みが作れないのかな、と漠然と思っていたのですが勘違いですね。ちゃんと作ることが出来ます。
scilab のインパルス応答の誤差に混乱させられるところでした。
scilab：ちょっとまって！インパルスってそもそも理論上でしか存在しないんじゃない！
画蔵 ：そういえばそうかな。
scilab：具体的に計算可能なように信号を定義してくれるならともかく、ありえない信号を定義しておいて計算ソフトになんとかしろったって、それ無理ゲー。
画蔵 ：なるほど一理あるな。
ということで、scilab を責めるのはチト酷なようです。
次回は別の方法を考えてみます。