インターポレータの実装

さて、いよいよ実装です。入力を整数倍するインターポレーターを実装してみましょう。せっかくなので汎用性をもたせてC++のクラスに仕立てます。

簡単な実装法

いちばん簡単な実装は図1の通りにインターポレーターを実装してしまうことです。

インターポレーター

図1 ポリフェーズフィルタによるインターポレーター

これは簡単です。L倍のインターポレーターは内部にL個のFIRフィルタを持たせ、全部に同じ入力を与えて出力1サンプル毎にフィルターを実行していけばそれで足ります。

確かに簡単ですし、プログラムの再利用という点からも生産性という点からも好ましく見えます。が、面白くないですねぇ。基本的に風呂上りの楽しみである2191空挺団としてはもうちょっとひねってみたいところです。そこで少し凝った構成に挑戦してみます。その前にもう一度ポリフェーズフィルタとアップサンプラーの関係を見てみましょう。

直感的な理解

アップサンプラの後に置かれたインターポレーションフィルタについて考えてみます。このフィルタをFIR構成で作るとして、状態変数(あるいはディレイライン)の中を流れていくデータは図2のようにあらわすことができます。

データの流れ

図2 インターポレーションフィルタの中のデータの流れ

ここでhxはインパルス応答、xnは入力データです。この図はL=3のアップサンプラの場合を説明しています。アップサンプラは入力データデータの間に0をつめて出力します。アップサンプラの出力1サンプル毎に図2のようデータは状態変数の中を流れていきますが、その際0に対応するインパルス応答は無視できます。その結果、アップサンプラによる0挿入は図1のコミュテーターと同様、インパルス応答の中の1/Lを抜き出す操作としても働きます。

なかなか文章で書くとわかりづらいのですが、ポリフェーズ分解によるインターポレーターはひとつのインターポレーションフィルタをL位相に分割して、出力1サンプル毎に1位相を処理する作業と等価です。したがって、L個のフィルタを用意せずとも、ひとつのフィルタを1位相ごとに処理できるようプログラムを書けば、きれいにまとまります。

さて、図2では内部変数の中をアップサンプルされたデータが1サンプル毎に流れていく様子を説明しましたが、アップサンプルせずに、入力したデータをインパルス応答側で切り替えて(Comutateして)処理する場合は図3のようになります。これは図1とまったく同じ処理だとわかります。この場合、Nタップのインターポレーションフィルタの実装には長さN/Lの状態変数を用意し、N/LタップのフィルタをL回実行することでLサンプル分の出力を生成します。そして長さNのインパルス応答と長さN/Lの状態変数をひとつずつ用意して巧みにポインタを操ることにより、簡単に任意のLに対応するインターポレーターが作れそうだということもわかります。

コミュテーション

図3 インパル応答によるコミュテーション

アセンブリ言語での実装

アセンブリ言語で書いたポリフェーズフィルタのプログラムを紹介します。これまでにない長さのプログラムですが、ポリフェーズフィルタ処理そのものは16行しかありません。プログラムの2/3程度はコメントと引数の読み込みなどの初期化コードです。

このプログラムは汎用性を持たせるために次のようにしました。

結局全部自由なのであらゆる場面で使えます。アップサンプル比にすれば普通のFIRフィルタとしても使えるはずです。この柔軟性をささえるために引数の数は全部で11もあります。このうち最初の5つは入出力バッファと処理サンプル数を指定しています(下のプログラム)。6つ目の引数ratioはアップサンプル比Lを与えます。7番目から10番目は内部のFIRフィルタそのものに関する引数です。ただし、tapsはポリフェーズ分解前のFIRフィルタのタップ数Nです。

11番目は少し毛色が違います。引数tapsDivRatioはFIRフィルタのタップ数をアップサンプル比で割った値です。これは一番内側のFIRフィルタが実際に行うフィルター処理のタップ数N/Lです。

ratio, taps, tapsDivRatioを三つとも指定するのは冗長です。たとえばtapsはあとの二つの積なので簡単に計算できます。引数の受け渡しの手間と内部の乗算の手間の比較をしてみると実際にはほとんど差がでません。そこで可読性をとって冗長な引数を渡すことにしました。

この関数に渡す引数tapsはratioの整数倍でなければなりません。

#include <asm_sprt.h>
.section/pm program;
.global _polyphaseUpSamplerFunc;
// void polyphaseUpSamplerFunc
// ( 
//      short input[],          : first input data address
//      short inputModify,      : distance of data ( dxxdxxdxx = 3 )
//      short output[],         : first output data address
//      short outputModify,     : distance of adta
//      short count,            : input data count to be processed
//      short ratio,            : upsampling ratio
//      short delayLine[],      : delayLine for FIR filter
//      short h[],              : impuls response
//      short taps,             : TAPs of FIR filter ( mutiple of ratio ).
//      short **curPtr,         : current Pointer to delayline
//      short tapsDivRatio      : taps / ratio. 
// )

  // offset of parameter from frame pointer
#define input 1
#define inputModify 2
#define output 3
#define outputModify 4
#define count 5
#define ratio 6
#define delayLine 7
#define h 8
#define taps 9
#define curPtr 10
#define tapsDivRatio 11

プログラムの入り口ではスタックフレームを確立し、保護すべきレジスタを退避します。アセンブリ・ルーチンの呼び方を参照にしてスクラッチレジスタを優先して使うことにより、退避しなければならないレジスタを抑えています。しかし手間を考えるとsave_regマクロとrestore_regマクロを使って全部退避したほうが楽だったかもしれません。

MSTATを退避したらDIS m_modeを実行して固定小数点乗算モードに変更します。

_polyphaseUpSamplerFunc:
    //******************** Stack frame
    entry;                  // establish the frame.
    pushs(i7);              // save i7
    pushs(m0);              // save m0
    pushs(mstat);           // save mstat
    dis     m_mode;         // go into fractional mode.

プログラムの前半はポインタ群の設定に費やします。ここではディレイライン、インパルス応答に循環バッファを使い、入出力のバッファアクセスもポインタでおこないます。これらいずれも後更新アクセスを行うので最初のポインタの設定は大事です。途中で道を見失わないように、あらかじめレジスタの配分を決めてコメントに書いておきます。これは後々のデバッグでも役に立ちます。

プログラム中で使っているget_argマクロとget_arg_dregマクロはフレームポインタを使って引数にアクセスできる便利なマクロです。この二つはasm_sprt.hで定義されています。

get_argマクロはm6レジスタを破壊するので注意してください。

// Register Allocation
//     i0,l0,b0,m0 : input,   0,             ,           inputModify
//     i1,l1,b1,m1 : h,       taps,          h,          ratio
//             ,m3 :                                     1
//     i6,l6,b6,m6 : output,  0,             ,           outputModify
//     i7,l7,b7,m7 : curptr*, tapsDivRatio,  delayline,  1
//
//     mx1 : tapsDivRatio-1
//     my1 : h

// Regsiters to be saved&restored
//     i7, m0, mstat
// Registers to be restored to 0
//     l1, l7
     
    //******************** Initialization of pointer
                      // set i0,l0,b0 & m0 = (input,0,-) & inputModify
    get_arg( i0, input );                // i0 : input
    get_arg( m0, inputModify );          // m0 : inputModify

                      // set i1,l1,b1 &m1 = (h,taps,h) & ratio
    get_arg_dreg( my1, h );              // my1 : h
    i1=my1;                              // i1 : h
    get_arg( l1, taps );                 // l1 : taps
    get_arg( m1, ratio );                // m1 : ratio
    reg(b1)=my1;                         // b1 : h

                      // set i7,l7,b7 & m7 = (*curPtr,tapsDivRatio,delayLine ) & 1
    get_arg( i7, curPtr );
    ax0=dm(i7+0);
    i7=ax0;                              // i7 : curPtr;
    get_arg_dreg( ax0, delayLine );
    reg(b7)=ax0;                         // b7 : delayLine
    m7=1;                                // m7 : 1
    get_arg_dreg(mx1, tapsDivRatio);     // mx1 : tapsDivRatio
    l7=mx1;                              // l7 : tapsDivRatio

                      // set i6,l6,b6 & m6 = ( output,0,-) & outputModify
    get_arg( i6, output );               // i6 : output
    get_arg_dreg( ax0, outputModify );
    m6=ax0;                              // m6 : outPutModify

        //******************** Set global       
    ar=mx1-1;
    mx1=ar;                                // mx1 : tapsDivRatio-1
    m3=1;                                  // m3 : 1;

ループは三重の入れ子になっています。一番内側はFIRフィルタでN/Lタップを実行し、一番外側は入力サンプルを処理します。中間のループは図3のコミュテーションを行うためにL回の繰り返しを行います。

一番外側のループはデータを持ってきて長さN/Lのディレイラインに挿入し、内部のポリフェーズフィルターを実行します。

L分割ポリフェーズフィルターは二番目のループ(innerLoop)が実行します。この部分には普通のFIRフィルターと異なる部分がいくつかあります。

  1. FIR部のインパルス応答のポインタ(i1)は、1ずつではなくratio(L)ずつ進む。
  2. FIRフィルタの終わりで最後のMACに重ねてインパルス応答のポインタを1だけ進めている。これで次回のFIRは次の特性に切り替わりる。
  3. ポリフェーズ分解したフィルターの処理が完了すると、2の処理のせいでインパルス応答のポインタがratio(L)だけ進みすぎている。その分の補正をおこなう。

この三つの処理のうち、オーバーヘッドとして表れるのは3だけです。2番の処理は複合命令で実行しますのでサイクル数の増加はありません。

        // ポリフェーズフィルター本体
        //******************** Outer loop. acquire input "count" times


    get_arg_dreg(ax0, count);       
    cntr=ax0;
    do outerloop until ce;

        modify(i7+=m5);                    // (*countPtr)-- : m5 is always -1
        ax0=dm(i0+=m0);                    // get input. input++ by inputModify
        dm(i7+0)=ax0;                      // store input into delayline.

           //**************** Inner loop. Polyphase filter
        cntr=m1;                                // cntr = ratio
        do innerLoop until ce;

            //************ FIR loop. tapsDivRatio taps filter
            cntr=mx1;                           // cntr = tapsDivRatio-1;
            mr=0, mx0=dm(i1,m1), my0=pm(i7,m7); // *(h++ratio), *(delayline++)
            do firLoop until ce;                // Let's Roll!
firLoop:        mr=mr+mx0*my0(ss), mx0=dm(i1,m1), my0=pm(i7,m7);  // MAC!!
            mr=mr+mx0*my0(rnd), ay0=dm(i1,m3);  // h++ next En(z). ay0:dummy
            sat mr;                             // round and sturate it.
innerLoop:              
            dm(i6+=m6)=mr1;                     // *(output++)=result

outerloop:      
        i1=my1;                          // reset h pointer 

フィルター処理が終わったら、呼び出し側にディレイラインの現在の処理位置を示すcurPtrを返します。m0は文法上必要なので添えているだけで、意味はありません。この処理でi0の値は破壊されますが、もう使いませんのでかまいません。

    //******************** finalization of pointer          
                                    // put back curPtr*
    get_arg( i0, curPtr );
    dm(i0+=m0)=i7;                  // curPtr* = current delayline pointer. 

最後にスタックフレームを廃棄します。MSTATを復帰することによって整数/固定小数点乗算モードが呼び出し時の設定に戻ります。

    //******************** Stack frame                                      
    l7=0;                   // restore l7
    l1=0;                   // restore l1
    pops(mstat);            // restore mstat
    and_pops(m0);           // restore m0
    and_pops(i7);           // restore i7
    exit;                   // dispose frame and say good bye

次は⇒デシメータの実装

2191空挺団 | プログラム | EZ-KIT | こぼれ話 | アーキテクチャー | 命令 | レジスタ | DSP掲示板 | FAQ |