remez法を使ったデジタルフィルタをfreeのpython scipyを使って設計しました。
備忘録に残します。
１．remez法
よく比較される窓関数法と比べて、所望のフィルタ特性を、より少ないタップ数
で設計できます。
その分計算が複雑であるため、有償ツールであるmatlabやMathematicaを使って
設計していました。
今回フリーであるpythonとライブラリのSciPy,NumPyを使用して設計できること
がわかり、やってみました。
２．環境構築
下記のサイトを参考に、numpy, matplotlib, scipy, pandasをインストールします。
下記サイトがわかりやすかったです。
https://qiita.com/koKekkoh/items/e5c8ff85379fbac8bea4
フィルタ特性のグラフを表示させるため「matplotlib」も必要になります。
３．サンプルコード
SciPyの公式サイトに使い方が載っています。
https://docs.scipy.org/doc/scipy/reference/generated/scipy.signal.remez.html
引数は以下のように記載されています
■numtaps・・・フィルタ長(タップ数)
■band・・・帯域で配列で指定します。全て周波数Hzで指定します。
[通過帯域の開始、通過帯域の終了、遮断帯域の開始、遮断帯域の修了]
■desired・・・理想ゲインを指定します。配列で指定します。
[通過帯域の理想ゲイン、遮断帯域の理想ゲイン]
■weight・・・各帯域の重みづけを配列で指定します。
■Hz・・・・・サンプリング周波数 int型。
（注 不使用が推奨され、代わりにfsの使用が推奨されています）
■type・・・フィルタタイプを選択します。デフォルトはbandpass
■maxiter・・・アルゴリズムの最大反復回数。デフォルトは25回
■grid_density・・アルゴリズムで使用されるグリッド密度。デフォルトは16
■fs・・・・・サンプリング周波数。Hzよりこちらを推奨。float型。
実際の使用例は、以下のような形で使用します。
fs = 500 # サンプリングレート
cutoff = 2 # カットオフ周波数
trans_width = 20 # 通過領域から遮断領域までの周波数
numtaps = 64
desired = [1.0, 0.00] # 通過領域ゲイン１、遮断領域ゲイン0
remez_impres = sp.signal.remez(
numtaps, [0, cutoff, cutoff + trans_width, 0.5 * fs], desired, fs=fs
)
Hz、weight、maxiter、grid_density に関しては特に指定しません。
４．フィルタ設計ツール
上記の使用方法で、簡単にフィルタ設計できるツールを作成しました。
GitHubにコードをアップします。
4つの値を入れるだけでローパスフィルタを作成するツールになります。
UIはpython標準のUIアプリであるtkinterを使っています。
4つの値を入れるだけでローパスフィルタを作成するツールになります。
UIはpython標準のUIアプリであるtkinterを使っています。
各値を入力し、calcボタンを押すとフィルタの周波数特性グラフが出力されます。 正規表現を比較するガードを付けていますが、英数字以外を入れたり、空白のまま だとエラーになります。
ターミナル上には、フィルタ係数を表示するようにしています。
得られたフィルタ特性です。グラフはmatplotlibで描いています。
５．設計したフィルタが正しいのか
実際に設計したフィルタですが、他のツールを使用しても同等のフィルタ特性
となるかを確認しました。
digitalfilter.com社様のdsplinksというすばらしいツールと比較させていただきました。
http://digitalfilter.com/products/dsplinks/jpdsplinks.html
若干の違いはあるのですが、同等のフィルタ係数と周波数特性を得ることが出来ました。
最後にエクセルにて、フィルタの多項式に係数を代入して作成したグラフを載せます。
当たり前ですが同様の周波数特性が得られました。
PythonやSciPyといった無料のツールで複雑なデジタルフィルタ設計が簡単に出来ました。
今後もSciPyやNumpyを使用したツールを作成していきたいです。