4

逆離散コサイン変換の実装があり、このコードにどのように到達したかを理解しようとしています。これまでのところ、これはおそらく、DFT (離散フーリエ変換) ではなく DCTのCooley-Tukey radix-2 Decimation-in-timeの最適化された実装であることがわかりました。

ただし、各段階で正確に何が起こるかについては、まだ途方に暮れています。Wx 定数はおそらくひねり要因であると考えました。

誰でも説明への参照を提供したり、このコードに説明を提供したりできますか?

//Twiddle factors
#define W1 2841 /* 2048*sqrt(2)*cos(1*pi/16) */
#define W2 2676 /* 2048*sqrt(2)*cos(2*pi/16) */
#define W3 2408 /* 2048*sqrt(2)*cos(3*pi/16) */
#define W5 1609 /* 2048*sqrt(2)*cos(5*pi/16) */
#define W6 1108 /* 2048*sqrt(2)*cos(6*pi/16) */
#define W7 565 /* 2048*sqrt(2)*cos(7*pi/16) */

//Discrete Cosine Transform on a row of 8 DCT coefficients.
NJ_INLINE void njRowIDCT(int* blk) {
    int x0, x1, x2, x3, x4, x5, x6, x7, x8;
    int t;
    if (!((x1 = blk[4] << 11)
        | (x2 = blk[6])
        | (x3 = blk[2])
        | (x4 = blk[1])
        | (x5 = blk[7])
        | (x6 = blk[5])
        | (x7 = blk[3])))
    {
        blk[0] = blk[1] = blk[2] = blk[3] = blk[4] = blk[5] = blk[6] = blk[7] = blk[0] << 3;
        return;
    }
    x0 = (blk[0] << 11) + 128;  //For rounding at fourth stage

    //First stage
    /*What exactly are we doing here? Do the x values have a meaning?*/
    x8 = W7 * (x4 + x5);
    x4 = x8 + (W1 - W7) * x4;
    x5 = x8 - (W1 + W7) * x5;
    x8 = W3 * (x6 + x7);
    x6 = x8 - (W3 - W5) * x6;
    x7 = x8 - (W3 + W5) * x7;

    //Second stage
    x8 = x0 + x1;
    x0 -= x1;
    x1 = W6 * (x3 + x2);
    x2 = x1 - (W2 + W6) * x2;
    x3 = x1 + (W2 - W6) * x3;
    x1 = x4 + x6;
    x4 -= x6;
    x6 = x5 + x7;
    x5 -= x7;

    //Third stage
    x7 = x8 + x3;
    x8 -= x3;
    x3 = x0 + x2;
    x0 -= x2;
    x2 = (181 * (x4 + x5) + 128) >> 8;
    x4 = (181 * (x4 - x5) + 128) >> 8;

    //Fourth stage
    blk[0] = (x7 + x1) >> 8;  //bit shift is to emulate 8 bit fixed point precision
    blk[1] = (x3 + x2) >> 8;
    blk[2] = (x0 + x4) >> 8;
    blk[3] = (x8 + x6) >> 8;
    blk[4] = (x8 - x6) >> 8;
    blk[5] = (x0 - x4) >> 8;
    blk[6] = (x3 - x2) >> 8;
    blk[7] = (x7 - x1) >> 8;

}
4

2 に答える 2

4

私は DCT の専門家ではありませんが、これまでにFFT の実装をいくつか書いたことがあります。以下をひとつまみの塩でお召し上がりください。

void njRowIDCT(int* blk)

このアルゴリズムは、24:8 の精度で固定小数点演算を使用する長さ 8 の基数 2 の DCT のように見えると正しく言います。最後のステージが8だけ右にシフトして目的の値を取得するため、精度を推測しています(それとテルテールコメント;)

長さが 8 であるため、累乗は 3 (2^3 = 8) であり、DCT に 3 つのステージがあることを意味します。これまでのところ、これはすべて FFT と非常によく似ています。「第4段階」は、固定小数点演算後に元の精度を回復するためのスケーリングにすぎないようです。

私が見る限り、入力ステージは、入力配列 blk からローカル変数 x0-x7 へのビット反転です。x8 は一時変数のようです。申し訳ありませんが、それ以上説明することはできません。

ビット反転ステージ

入力のビット反転

アップデート

DSP For Scientists and Engineersをご覧ください。信号処理のトピックを明確かつ正確に説明します。この章は DCT に関するものです (p497 まで飛ばしてください)。

Wn (回転因子) はこの章の基底関数に対応しますが、これは 8x8 (2D) DCT を記述していることに注意してください。

私が言及した 3 つのステージに関して、8 ポイント FFT の説明と比較してください。

8ptFFTグラフ

FFTは、ビット反転した入力アレイ(本質的に複雑なマルチアドである)で蝶を実行しており、1つのパスにWNまたはTwiddle因子を途中で掛けています。FFTは段階的に実行されます。私はあなたのDCTコードが何をしているのかをまだ解決していませんが、このような図に分解するかもしれません。

それか、彼らが話していることを知っている人がステップアップします;-)

このページに再度アクセスして、さらにコードを解読しながら編集します

于 2012-01-04T21:15:03.540 に答える
1

DFT と DCT はどちらも単なる線形変換であり、単一の複素行列乗算として表すことができます (厳密に実数の入力に対して剪定される場合もあります)。したがって、上記の方程式を組み合わせて、各最終項の式を取得できます。これは、最終的に線形変換行列の 1 行に相当するはずです (丸めの問題は無視します)。次に、上記のコード シーケンスが、一般的な部分式の最適化またはリファクタリングを行計算間および/または行計算内で手動で行う方法を確認します。

于 2012-01-04T21:03:28.450 に答える