<直線位相FIRフィルター>

直線位相とは耳慣れない言葉かもしれない。

アナログフィルターでは実現するのは難しいのであるが、デジタルフィルターの世界では伝達しても波形の形が崩れないフィルターを構成することができる。そ れは、波形をある有限な時間だけ遅 延させて伝達する様なイメージである。波形といっても、単一のSIN波だけの信号ならば、動作範囲が線形のフィルターならば形崩れしないで SIN波をSIN波の形にして伝達してくれるが、ここでの波形とは、単一のSIN波ではなく、複雑な形をした波形を言っている。通常のフィルターでは、波 形に含まれる周波数成分とその位相変化量がまちまちとなるため、フィルターの伝送時に波形の崩れが生じる。

デジタルフィルターの設計の理論は数式を使った高等な学問である。しかし、今では希望の周波数特性を入力するとデジタルフィルター係数を自動的に計算して 答えを教えてくれる便利なWEBサイトもある。そこで、ここでは、参考用に、直線位相FIRフィルターのインパルス応答を求めるC++のプログラムの例を 示す。


直線位相FIRフィルターのインパルス応答を求めるサンプルプログラム

LPF(ローパスフィルター)と HPF(ハイパスフィルター)の2つを兼用している例である。

#include <math.h>


void  TCFIR::FIR_Filter(int fs, int fc, int nzi, int filter_type)
//
//    fs    fs サンプリング周波数[Hz]
//    fc    カットオフ周波数[Hz]
//    nzi   フィルターの次数。 ここでは、次数は奇数である必要がある。
//    filter_type を LPF_TYPEに設定するとLPFのインパルス応答を求める。
//                            又は、HPF_TYPEに設定するとHPFのインパルス応答を求める。
{
   int i;
   double x;
                            // フィルターは直線位相FIRを想定すると、
   m_N= nzi;                // フィルターの次数は奇数である必要があり、
   m_L= (m_N -1)/2;      &   // フィルターの次数のまんなかを中心にしてインパルス応答が対称となる。
                                  // m_L x m_T 時間 遅延する。                     

   x= (double)fc/ (double) fs;  // fc カットオフ周波数[Hz], fs サンプリング周波数[Hz]
   x= x * 2.0 * M_PI;             // サンプリング周波数が2πになる。
   m_wc=x;                        // カットオフ周波数の単位をラジアンに変換する。

   m_T= 1.0/ (double) fs;         // サンプリング間隔の時間[秒]を求める。

  
   // 係数計算のための作業用メモリーを確保する。
   m_h = (double *)malloc( m_N * sizeof(double));
 
  
  switch(filter_type)
  {
   case LPF_TYPE:  // 理想的なLPFの周波数特性を逆フーリエ級数変換して、インパルス応答を計算する。
      for(i=0;i<m_L;++i)
      {
           x= m_wc * (double)(i - m_L) ;
           x= sin(x);
           x = x / ((double)(i - m_L) * M_PI);

           // m_Lを中心に対称なので、同じ値をいれる。
           *(m_h + i) =x;
           *(m_h + (m_N -1 -i))=x;

      }
     *(m_h + m_L) = m_wc / M_PI;   // 対称中心の計算する。
     break;
                 
    case HPF_TYPE:    // 理想的なHPFの周波数特性を逆フーリエ級数変換して、インパルス応答を計算する。
      for(i=0;i<m_L;++i)
      {
           x= m_wc * (double)(i - m_L) ;
           x=  -sin(x);
           x = x / ((double)(i - m_L) * M_PI);

           // m_Lを中心に対称なので、同じ値をいれる。
           *(m_h + i) =x;
           *(m_h + (m_N -1 -i))=x;

      }
     *(m_h + m_L) = ( M_PI - m_wc )/ M_PI;   // 対称中心の計算する。
     break;

 }

   // ギブス現象を抑えて、減衰量を多くとるために
   // 上記で求めたインパルス応答に
   //  更に ハミングの窓をかける。
   for(i=0;i<m_L;++i)
   {
       // ハミングの窓の係数の計算
        x= 2.0 * M_PI * (double)i;
        x= x / (double)(m_N -1);
        x= cos(x);
        x= 0.54 - 0.46*x;

        // 上記で求めたインパルス応答に、ハミングの係数を乗じる。
        x= x * *(m_h + i);

       // m_Lを中心に対称なので、同じ値をいれる。
       *(m_h + i) =x;
       *(m_h + (m_N -1 -i))=x;

   }
 
    // m_h に求めたインパルス応答の値が格納される。

  ・・・・・

   
}

int  TCFIR::return_delay_points( )
{
      // フィルターをかけたときに 入力の波形が遅延するポイント数を返す。
     return m_L;
}



class TCFIR
{
public:
        void FIR_Filter(int fs, int fc, int nzi, int filter_type);
        int  return_delay_points( );
        ・・・・・

private:
        int    m_N;
        int    m_L;
        double m_wc;
        double m_T;
        double *m_h;

};

#define LPF_TYPE   1
#define HPF_TYPE   2

このサンプルを実際に動作させるためには、おそらく、
ご使用されるプログラミング環境に合わせて修正・変更・追加などが必要です。
万一、このサンプルを動作させる場合は、あなたの責任においておこなってください。


Home page
No.4   2007年4月10日