4

ほぼ完成したプロジェクトのプロファイリングを行ったところ、CPU 時間の約 4 分の 3 がこの IIR フィルター関数に費やされていることがわかりました (現在、ターゲット ハードウェアでは約 1 秒間に数十万回呼び出されています)。他のすべてがうまく機能しているので、特定のハードウェアとソフトウェアのターゲットに合わせて最適化できるかどうか疑問に思っています。私のターゲットは、iPhone 4 以降、iOS 4.3 以降、LLVM 4.x のみです。利益が得られる場合は、多少の不正確さはおそらく問題ありません。

static float filter(const float a, const float *b, const float c, float *d, const int e, const float x)
{
    float return_value = 0;

    d[0] = x;
    d[1] = c * d[0] + a * d[1];

    int j;

    for (j = 2; j <= e; j++) {
        return_value += (d[j] += a * (d[j + 1] - d[j - 1])) * b[j];
    }

    for (j = e + 1; j > 1; j--) {
        d[j] = d[j - 1];
    }
    return (return_value);
}

デフォルトのコンパイラ最適化を超えて最適化することが可能であれば、あなたの意見にも興味があります。NEON SIMD が役立つものなのか (これは私にとって新しい分野です)、VFP を活用できるのか、LLVM の自動ベクトル化が役立つのか疑問に思っています。

次のLLVMフラグを試しました：

-ffast-math (顕著な違いはありませんでした)

-O4 (iPhone 4S では 25% の時間短縮で大きな違いがありましたが、最小限のターゲット デバイスである iPhone 4 では顕著な違いはありませんでした。これを改善することが私の主な目標です)

-O3 -mllvm -unroll-allow-partial -mllvm -unroll-runtime -funsafe-math-optimizations -ffast-math -mllvm -vectorize -mllvm -bb-vectorize-aligned-only (Hal Finkel のスライドからの LLVM 自動ベクトル化フラグ: http://llvm.org/devmtg/2012-04-12/Slides/Hal_Finkel.pdf、Xcode リリース ターゲットのデフォルトの LLVM 最適化よりも遅くなります)

他のフラグ、異なるアプローチ、関数への変更を受け入れます。入力と戻り値の型と値はそのままにしておくことをお勧めします。ここで実際にFIRにNEON組み込み関数を使用することについての議論があります: https://pixhawk.ethz.ch/_media/software/optimization/neon_support_in_the_arm_compiler.pdfしかし、私は情報をうまく適用するのに十分な経験がありません。自分のケースに。明確にしていただきありがとうございます。

編集これに早く気付かなかったことをお詫びします。aka.nice の提案を調査した後、e、a、および c に渡される値は常に同じ値であり、実行前にそれらを知っているため、この情報を組み込むアプローチがオプションであることに気付きました。

4

3 に答える 3

7

以下は、vDSP ルーチンを使用するためにコードに加えることができるいくつかの変換です。これらの変換では、T0、T1、および T2 という名前のさまざまな一時バッファーが使用されます。これらはそれぞれ、e-1 要素に十分なスペースを持つ float の配列です。

まず、一時バッファを使用して を計算しa * b[j]ます。これにより、元のコードが変更されます。

for (j = 2; j <= e; j++) {
    return_value += (d[j] += a * (d[j + 1] - d[j - 1])) * b[j];
}

に：

vDSP_vsmul(b+2, 1, &a, T0, 1, e-1);
for (j = 2; j <= e; j++)
    return_value += (d[j] += (d[j+1] - d[j-1])) * T0[j-2];

次に、vDSP_vmul を使用して計算しd[j+1] * T0[j-2]ます。

vDSP_vsmul(b+2, 1, &a, T0, 1, e-1);
vDSP_vmul(d+3, 1, T0, 1, T1, 1, e-1);
for (j = 2; j <= e; j++)
    return_value += (d[j] += T1[j-2] - d[j-1] * T0[j-2];

次に、vDSP_vmul を vDSP_vma (ベクトル乗算加算) に昇格させて計算しd[j] + d[j+1] * T0[j-2]ます。

vDSP_vsmul(b+2, 1, &a, T0, 1, e-1);
vDSP_vma(d+3, 1, T0, 1, d+2, 1, T1, 1, e-1);
for (j = 2; j <= e; j++)
    return_value += (d[j] = T1[j-2] - d[j-1] * T0[j-2];

時間を計って、改善があるかどうかを確認すると思います。いくつかの問題があります:

  • SIMD コードは、データが 16 バイトでアラインされている場合に最適に機能します。やなどの配列インデックスを使用すると、これj-1をj+1防ぐことができます。携帯電話の ARM プロセッサは、他の一部のプロセッサほど整列されていないデータで悪くはありませんが、パフォーマンスはモデルによって異なります。
  • e が大きい (数千を超える) 場合、vDSP_vma 操作中に T0 と d がキャッシュから削除される可能性があり、次のループでそれらを再読み込みする必要があります。この影響を軽減するために、ストリップマイニングという手法があります。ここでは詳しく説明しませんが、基本的に、操作は配列の小さなストリップに分割されます。
  • 最終ループの IIR は、依然としてプロセッサのボトルネックになる可能性があります。一部の IIR (vDSP_deq22 など) を実行するための vDSP のルーチンがありますが、変換によって失われる可能性があるよりも多くのパフォーマンスを得るために、vDSP ルーチンに十分に一致する方法でこのフィルターを表現できるかどうかは明らかではありません。 .
  • return_value を計算するための最終ループの合計も、ループから削除して vDSP ルーチン (おそらく vDSP_sve) に置き換えることができますが、IIR によって引き起こされるスラックにより、実行時間を大幅に追加することなく追加を行うことができるのではないかと思います。ループ。

上記は私の頭の上から外れています。私はコードをテストしていません。変換を 1 つずつ行うことをお勧めします。これにより、各変更後にコードをテストし、続行する前にエラーを特定できます。

IIR ではない満足のいくフィルターを見つけることができれば、さらにパフォーマンスを最適化できる可能性があります。

于 2012-10-15T19:06:22.830 に答える
0

遅延演算子 z を使用してフィルターを書き直すために、この小さな演習を試みました

たとえば、e=4 の場合、入力 u と出力 y の名前を変更しました

d1*z= u
d2*z= c*u + a*d1
d3*z= d2 + a*(d3-d1*z)
d4*z= d3 + a*(d4-d2*z)
d5*z= d4 + a*(d5-d3*z)
y = (b2*d3*z + b3*d4*z + b4*d5*z)

di はフィルター状態であることに注意してください。
d3*z は d3 の次の値です (コードでは変数 d2 のように見えます)。
その後、di を削除して、伝達関数 y/u を z に書き込むことができます。
次に、上記の伝達関数を因数分解/単純化することにより、最小の表現には e 状態のみが必要であることがわかります。
分母は ですz*(z-a)^3。つまり、0 の極と、多重度 (e-1) の a の極です。

次に、フィルターを標準の状態空間行列表現に入れることができます。

z*X = A*X + B*u
y = C*X + d*u

極の特定の形式を使用すると、部分分数展開で転送を分解し、この特別な形式で行列 A と B を取得できます (matlab のような表記法)。

A = [0 1 0 0;    B=[0;
     0 a 1 0;       0;
     0 0 a 1;       0;
     0 0 0 a]       1]

C と d は少し簡単ではありませんが...
これらは部分分数展開の分子と直接項から抽出されます それらは
bi、c (次数 1) および a (次数 e) の多項式です
e=4 の場合、

C=[   a^3*b2            - a^2*b3                      + a*b4 ,
     -a^2*b2              + a*b3                + (c-a^2)*b4 ,
        a*b2        + (c-a^2)*b3 + (-2*a^2*c-2*a^2-a+a^4)*b4 ,
  (c-a^2)*b2 + (-a^2-a^2*c-a)*b3 +         (-2*a*c+2*a^3)*b4 ]
d=     -a*b2            - a*c*b3                    + a^2*b4

C と d を支配する e で再帰を見つけ、それらを事前に計算できる
場合、フィルターを単純なベクトル演算に減らすことができます。

z*X = a*[ 0; x2 ; x3 ; x4 ... xe ] + [x2 ; x3 ; x4 ... xe ; u ];
Y = C*[ x1 ; x2 ; x3 ; x4 ... xe ] + d*u

または、関数の(Xnext,y)=filter(X,u,a,C,d,e)擬似コードとして表現されます。

y = dot_product( C , X) + d*u; // (like BLAS _DOT)
Xnext(1:e-1) = X(2:e); // this is a memcopy (like BLAS _COPY)
Xnext(e)=u;
X(1)=0;
Xnext=a*X+Xnext; // this is an inplace vector muladd (like BLAS _AXPY)

X=Xnext; // another memcopy outside the function (can be moved inside).

BLAS 関数を使用すると、コードは Applecentric だけでなく、多くのハードウェアに移植できることに注意してください。パフォーマンスはそれほど変わらないと思います。

EDIT：部分分数展開について

純粋な部分分数展開では、対角状態空間表現と、1 でいっぱいの行列 B が得られます。これも興味深い変種になる可能性があります。(並列フィルター)
上記で使用した私のバリアントは、カスケードまたははしご (一連のフィルター) に似ています。

于 2012-10-16T22:48:50.593 に答える
0

入力と戻り値の型と値はそのままにしておきたいと思います…</p>

それにもかかわらず、レンダリングを float から integer に移動すると、かなり役立ちます。

その変更を提示する実装にローカライズしても役に立ちません。しかし、FIR だけを整数として再実装するように拡張すると、すぐに成果が得られます (ただし、サイズが常に信じられないほど小さいことが保証されている場合を除きます。その場合、変換/移動時間はより多くのコストがかかります)。もちろん、レンダー グラフの大部分を整数に移動すると、ゲインが大きくなり、必要な変換がさらに少なくなります。

もう 1 つの考慮事項は、Accelerate.framework のユーティリティを調べることです (独自の asm を作成する手間を省ける可能性があります)。

于 2012-10-15T18:43:56.840 に答える