dsPIC33FJ64GP802 (15) --- 無償で使えるフィルタ設計アプリケーション (2)

dsPIC

　前回示した

ディジタル・フィルタ設計プログラム (IIR_Design、FIR_Remez、FIR_Kaiser) — (三上直樹著・CQ 出版「はじめて学ぶディジタル・フィルタと高速フーリエ変換」のサポート・プログラム)

http://www.cqpub.co.jp/hanbai/books/30/30881.htm

は、上記のリンクから飛べる CQ 出版の web サイトからダウンロードすることができます。
「readme.txt」中の記述を要約すると、

  • (このプログラムは、書籍中で) 説明しているディジタル･フィルタの設計法の理解を助けるために作成したもの
  • 営利を目的とする場合を除いて、プログラムは、
  • 利用すること、および
  • (著作権表記を消去しない限り) 改変することが可能

となっています。
　書籍を購入していない場合にも使用する権利があるかどうかについては言及がありませんが、書籍未購入の場合には使用を禁止するという表明もありません。
　したがって、書籍未購入者にも使用できるものとして解釈することにします。
　IIR_Design.exe、FIR_Remez.exe、FIR_Kaiser.exe の 3 本のアプリケーション・プログラムから構成されていて、Embarcadero C++ Builder で記述されたソース・プログラムも付属しています。
　IIR_Design の実行画面を下に示します。

  • フィルタ次数
  • サンプリング周波数
  • 遮断域端周波数

などは右側のパネル部分のテキスト・ボックスへ数値を入力し、

  • フィルタ形状 (バタワース/チェビシェフ/逆チェビシェフ/楕円)
  • フィルタ種類 (LPF/HPF/BPF/BRF)

は右側パネル部分のラジオ・ボタンで選択します。
　このフィルタ仕様入力部分は、ファイルにセーブしたりファイルをロードする機能はなく、すべて GUI で指定します。　前回の指定を覚えておく機能もなく、プログラム開始時には常にデフォルト値に戻ります。
　下部のパネルの「設計」ボタンを押すと設計が始まり、フィルタ係数値が計算されます。
　「スケーリング」ボタンを押すと、フィルタ内部のピーク・ゲインを 0 dB に抑えるようなスケーリングが施された係数値が求まります。
　「結果の保存」ボタンを押すと、右下のリスト・ボックスに表示されている計算結果を (名前を付けて) ファイルに保存します。
　「バッチ処理」的に多数のフィルタを連続して設計したい場合には、「キーボード/マウス・マクロレコーダ/プレーヤ」等を利用して、自動的なGUI 操作を実現する必要があります。
　サンプリング周波数 12 kHz、中心周波数 880 Hz の 1/12 オクターブ・バンドフィルタを 6 次 IIR フィルタとして設計した例の出力を下に示します。

帯域通過フィルタ
バタワース特性
次数　　　　：6 次
標本化周波数：12.000000 kHz
遮断周波数１：0.854948 kHz
遮断周波数２：0.905786 kHz

直接形の係数
 m    分母の係数 (am)    分子の係数 (bm)
 0                      2.2960824645E-06
 1   5.3265501415E+00   0.0000000000E+00
 2  -1.2404514338E+01  -6.8882473934E-06
 3   1.6062788090E+01   0.0000000000E+00
 4  -1.2186317001E+01   6.8882473934E-06
 5   5.1408093405E+00   0.0000000000E+00
 6  -9.4815337568E-01  -2.2960824645E-06

縦続形の係数

a11 =  1.7678805237E+00
a21 = -9.7372937932E-01
b01 =  1.0000000000E+00
b11 = -0.0000000000E+00
b21 = -1.0000000000E+00

a12 =  1.7688754785E+00
a22 = -9.8647476262E-01
b02 =  1.0000000000E+00
b12 = -0.0000000000E+00
b22 = -1.0000000000E+00

a13 =  1.7897941394E+00
a23 = -9.8708452398E-01
b03 =  1.0000000000E+00
b13 = -0.0000000000E+00
b23 = -1.0000000000E+00

利得定数 = 2.2960824645E-06

スケーリングされた縦続形の係数

s0  =  1.1677556141E-02

a11 =  1.7678805237E+00
a21 = -9.7372937932E-01
b01 =  8.8513880920E-03
b11 = -0.0000000000E+00
b21 = -8.8513880920E-03

a12 =  1.7688754785E+00
a22 = -9.8647476262E-01
b02 =  1.9583452462E-02
b12 = -0.0000000000E+00
b22 = -1.9583452462E-02

a13 =  1.7897941394E+00
a23 = -9.8708452398E-01
b03 =  1.1343181387E+00
b13 = -0.0000000000E+00
b23 = -1.1343181387E+00

　これを見ると、2 次セクションに分割された伝達関数の分母のフィルタ係数の符号は、IIRCanonic() や IIRTransposed() で必要とされている極性と一致しており、符号を付け替える必要はありません。
　係数値の並び方は一致していないので、並び替える必要があります。
　IIRCanonic() で使う場合には、

  • a2x
  • a1x
  • b2x
  • b1x
  • b0x

の順に並び替えます。　(「x」は 2 次セクション番号)
　フィルタ係数並べ替えのためのプログラムを作成して処理した結果を下に示します。　(C 言語の配列の初期値としての指定のみ。　配列の定義自体は含まない。)

// converted from file "ob12C_n81.txt"

// 帯域通過フィルタ
// バタワース特性
// 次数　　　　：6 次
// 標本化周波数：12.000000 kHz
// 遮断周波数１：0.854948 kHz
// 遮断周波数２：0.905786 kHz
// 

// biquad section 1
  Q15(-0.97372938/2),      // den12
  Q15( 1.7678805/2),       // den11
  Q15(-0.0088513881/2),    // num12
  Q15( 0.0000000/2),       // num11
  Q15( 0.0088513881/2),    // num10

// biquad section 2
  Q15(-0.98647476/2),      // den22
  Q15( 1.7688755/2),       // den21
  Q15(-0.019583452/2),     // num22
  Q15( 0.0000000/2),       // num21
  Q15( 0.019583452/2),     // num20

// biquad section 3
  Q15(-0.98708452/2),      // den32
  Q15( 1.7897941/2),       // den31
  Q15(-1.1343181/2),       // num32
  Q15( 0.0000000/2),       // num31
  Q15( 1.1343181/2),       // num30

// scaling factors for IIRCanonicStruct
  Q15( 0.011677556/2),    // initialGain
  0,                      // finalShift

　この係数で、周波数スイープするフルスケール正弦波に対する応答を観測した結果を下に示します。

　スケーリングが、あまりにびったりと 0 dB に合わせてあるので、フルスケール正弦波の周波数スイープによる過渡的なレベル・オーバーが一部生じて、クリップを引き起こし、基本波に対して - 30 dB 程度の高調波が発生しています。
　入力信号レベルが低ければクリップは生じませんが、フルスケール正弦波に対してクリップしないように、さらにスケーリングすることにしました。
　並び替えプログラム中で、1 段目と 2 段目の分子の係数をそれぞれ 0.9 倍し、-1 dB 程度レベルを下げ、3 段目の分子の係数を (1/0.9)2 倍してレベルを回復するようにした結果を下に示します。

// converted from file "ob12C_n81.txt"

// 帯域通過フィルタ
// バタワース特性
// 次数　　　　：6 次
// 標本化周波数：12.000000 kHz
// 遮断周波数１：0.854948 kHz
// 遮断周波数２：0.905786 kHz
// 

// biquad section 1
  Q15(-0.97372938/2),      // den12
  Q15( 1.7678805/2),       // den11
  Q15(-0.0084088187/2),    // num12
  Q15( 0.0000000/2),       // num11
  Q15( 0.0084088187/2),    // num10

// biquad section 2
  Q15(-0.98647476/2),      // den22
  Q15( 1.7688755/2),       // den21
  Q15(-0.018604280/2),     // num22
  Q15( 0.0000000/2),       // num21
  Q15( 0.018604280/2),     // num20

// biquad section 3
  Q15(-0.98708452/2),      // den32
  Q15( 1.7897941/2),       // den31
  Q15(-1.2568622/2),       // num32
  Q15( 0.0000000/2),       // num31
  Q15( 1.2568622/2),       // num30

// scaling factors for IIRCanonicStruct
  Q15( 0.011677556/2),    // initialGain
  0,                      // finalShift

　この係数で測定した結果を下に示します。

　まだ高調波は生じていますが、レベルは -60 dB 程度に抑えられています。