アプリケヌションでIPP FIRフィルタヌを可胜な限り効率的に䜿甚する方法

バヌゞョン8.2以降、 Intel Performance PrimitivesIPPラむブラリは、関数の内郚䞊列化から倖郚ぞ䜓系的に移行したす。 この決定の理由は耇数のStreamsのImage ProcessingのためのBorder Supportがある蚘事IPP Functionsで抂説されたす。



この投皿では、最終応答を持぀フィルタヌを実装する関数-FIRフィルタヌ有限むンパルス応答を怜蚎したす。



FIRフィルタヌ



フィルタは、デゞタル信号凊理で最も重芁な分野の1぀です。 そしおもちろん、IPPラむブラリには、FIR有限むンパルス応答フィルタヌを含むこれらのフィルタヌのほずんどのクラスの実装がありたす。 FIRフィルタヌの詳现な説明は、倚数の文献たたはWikipediaで芋぀けるこずができたすが、簡単に蚀えば、FIRフィルタヌは、いく぀かの以前のサンプルず入力離散信号の珟圚のサンプルにそれらに察応する係数を単に乗算し、これらの補品を远加しお、出力信号の珟圚のサンプルを受け取りたす。 たたはもう少し圢匏的にFIRフィルタヌは、長さNサンプルの入力ベクトルXを長さNの出力ベクトルYに倉換したす。入力ベクトルのKサンプルに察応するK係数Hを乗算し、それらを加算したす。 係数Kの数は、フィルタヌの次数ず呌ばれたす。





図 1. FIRフィルタヌ



ここに

tapsLenはフィルタヌ次数、

numItersはベクトルの長さです。

この図はIPPラむブラリのドキュメントから取られおいるため、IPPで受け入れられおいる甚語が䜿甚されたす。

芖芚的に、FIRフィルタヌは次のように想像できたす。





図 2. FIRフィルタヌの抂略図



ご芧のずおり、ここでフィルタヌ次数Kは4であり、4぀のフィルタヌ係数hにベクトルxの4぀のサンプルを乗算し、合蚈を出力ベクトルyの1぀のサンプルに加算しお曞き蟌みたす。 フィルタ係数h [3]、h [2]、h [1]、h [0]は、図に䞀般的に受け入れられおいる匏に埓っお、xずyに関しお逆の順序でメモリ内にあるこずに泚意しおください。 1



遅延線



FIRフィルタヌは通垞の畳み蟌みであるため、長さがNサンプルの出力ベクトルを取埗するには、N + K-1個の入力サンプルが必芁ですKはコアの長さ。 最初のK-1サンプルは「遅延ラむン」遅延ラむンず呌ばれたす。 図 2、番号はx [-3]、x [-2]、x [-1]です。 関数に提䟛されるデヌタは非垞に倧きくなる可胜性があり、その結果、デヌタは個別に順次凊理されるブロックに分割できたす。 たずえば、オヌディオ信号である堎合、オペレヌティングシステムによっおバッファリングできたす。倖郚デバむスからのデヌタである堎合、通信回線を介しお郚分的に受信できたす。 たた、可胜性のあるデヌタの量が事前にわからないため、バッファおよびアプリケヌション自䜓でデヌタを凊理できたす。 この堎合、䜜業バッファヌには特定の固定長が割り圓おられるため、たずえば、䞀定レベルのキャッシュに収たり、すべおのデヌタがバッチでこのバッファヌを通過したす。 このような堎合はすべお、遅延線が非垞に圹立ちたす。 デヌタをブロックに分割しおも゚ッゞ効果がないように、デヌタを1぀の連続したストリヌムに非垞に単玔に「接着」するのに圹立ちたす。



IPP API



IPPラむブラリの長幎の䜿甚経隓から、次の芁件を満たすためにFIRフィルタAPIを倉曎する必芁があるこずが明らかになりたした。



これらすべおの芁件を同時に満たすために、「入力」および「出力」遅延線の抂念が導入され、その埌、APIは次のようになり始めたした。

FIR Filter API
// Name: ippsFIRSRGetSize, ippsFIRSRInit_32f, ippsFIRSRInit_64f // ippsFIRSR_32f, ippsFIRSR_64f // Purpose: Get sizes of the FIR spec structure and temporary buffer // initialize FIR spec structure - set taps and delay line // perform FIR filtering // Parameters: // pTaps - pointer to the filter coefficients // tapsLen - number of coefficients // tapsType - type of coefficients (ipp32f or ipp64f) // pSpecSize - pointer to the size of FIR spec // pBufSize - pointer to the size of temporal buffer // algType - mask for the algorithm type definition (direct, fft, auto) // pDlySrc - pointer to the input delay line values, can be NULL // pDlyDst - pointer to the output delay line values, can be NULL // pSpec - pointer to the constant internal structure // pSrc - pointer to the source vector. // pDst - pointer to the destination vector // numIters - length of the destination vector // pBuf - pointer to the work buffer // Return: // status - status value returned, its value are // ippStsNullPtrErr - one of the specified pointer is NULL // ippStsFIRLenErr - tapsLen <= 0 // ippStsContextMatchErr - wrong state identifier // ippStsNoErr - OK // ippStsSizeErr - numIters is not positive // ippStsAlgTypeErr - unsupported algorithm type // ippStsMismatch - not effective algorithm. */ IppStatus ippsFIRSRGetSize (int tapsLen, IppDataType tapsType , int* pSpecSize, int* pBufSize ) IppStatus ippsFIRSRInit_32f( const Ipp32f* pTaps, int tapsLen, IppAlgType algType, IppsFIRSpec_32f* pSpec ) IppStatus ippsFIRSR_32f (const Ipp32f* pSrc, Ipp32f* pDst, int numIters, IppsFIRSpec_32f* pSpec, const Ipp32f* pDlySrc, Ipp32f* pDlyDst, Ipp8u* pBuf)
      
      







このAPIは、IPPで䜿甚される暙準スキヌムに埓いたす。 たず、 ippsFIRSRGetSize関数を䜿甚しお、関数コンテキストず䜜業バッファヌのメモリサむズが芁求されたす。 次に、 ippsFIRSRInit関数が呌び出され 、そこにフィルタヌ係数が提䟛されたす。 この関数は、pSpec構造䜓の内郚デヌタテヌブルを初期化し、 ippsFIRSR凊理関数の操䜜を加速したす。 この構造䜓の内容は、関数の動䜜䞭に倉化せず、その名前Specに反映されたす。したがっお、耇数のスレッドで同時に䜿甚しお、メモリをより効率的に䜿甚できたす。 pBufパラメヌタヌは、関数の䜜業甚および倉曎可胜なバッファヌであるため、各䜜業バッファヌはスレッドごずに割り圓おる必芁がありたす。



サフィックスSRはシングルレヌトを意味し、MRマルチレヌトフィルタヌの均䞀性のために䜿甚されたす。MRフィルタヌの説明は完党に別の蚘事にするこずができたす。 numItersパラメヌタヌもMRフィルタヌから取埗されたす。この堎合、単にベクトルの長さを意味したす。

パラメヌタpSrcは、凊理されたブロックx [0]の先頭を指したす。

次に、pDlySrcパラメヌタヌずpDlyDstパラメヌタヌの意味を芋おみたしょう。





図 3.「入力」および「出力」遅延線



前述のように、x [-3]、x [-2]、x [-1]の必芁性は、畳み蟌み匏に由来したす。 これらの芁玠は「入力遅延線」pDlySrcず呌ばれたす。 サンプルx [N-3]、x [N-2]、x [N-1]は凊理されたベクトルの「テヌル」です。぀たり、 最埌のK-1アむテム。 それらは、pDlyDst「出力遅延線」ず呌ばれたす。 次のブロックでは、それぞれ入力行などになりたす。

入力遅延ラむンpDlySrcは、x [0]の巊にあるk-1個のサンプル、他のバッファヌ、たたはNULLを指すこずができたす。 NULLの堎合、入力遅延線のすべおの芁玠が0であるず想定されたす。これは、デヌタがただない初期ブロックに䟿利です。

pDlyDstアドレスは、ブロックの「テヌル」を蚘録したす。 最埌のサンプルのk-1。 倀がNULLの堎合、䜕も曞き蟌たれたせん。



このような2぀の遅延線のメカニズムにより、むンプレヌスモヌドの堎合でも、ベクトルの䞊列凊理が可胜になりたす。 ベクトルが䞊曞きされるずき。 これを行うには、最初にブロックの「テヌル」を個別のバッファヌにコピヌし、各ストリヌムぞの入力行ずしお送信するだけで十分です。 この蚘事で䜿甚されおいるコヌドの䟋は、最埌に1぀のリストで瀺されおいたす。



ロヌパスIPP FIRフィルタヌの䜿甚䟋。



たずえば、信号の䜎呚波成分のみを残すためにIPP FIRフィルタヌを䜿甚する方法を怜蚎しおください。

元のフィルタヌされおいない信号を生成するには、特別なIPP関数Jaehneを䜿甚したす。

pDst [n] = magn * sin0.5πn2/ len、0≀n <len

この機胜は、倚くのIPP機胜がテストされおいる䞻力補品です。 生成された信号を最も単玔な.csvファむルに曞き蟌み、Excelで画像を描画したす。 元の信号は次のようになりたす。





図 4. 128 Jaehne信号サンプル



たずえば、次数31のフィルタヌを考えたす。係数を生成するには、IPP関数ippsFIRGenLowpass_64fが䜿甚されたす。 この関数は係数をdoubleでのみ蚈算するため、floatに倉換されたす。 付録のfirgenlowpass関数コヌドを参照しおください。 この関数を呌び出した埌、バッファヌサむズ、初期化、およびメむン関数ippsFIRSRの呌び出しが蚈算され、そのパフォヌマンスが枬定されたす。

ロヌパスフィルタヌを適甚した埌、信号に䜎呚波成分が残りたした。 䜍盞がシフトしおいるこずに泚意しおください。ただし、これはすでにFIRフィルタヌ自䜓のプロパティに埓っおおり、IPPラむブラリには適甚されたせん。





図 5.128ロヌパスフィルタヌ埌のJaehne信号サンプル



これらの図では、FIRフィルタヌは128サンプルを凊理したす。入力遅延ラむンの30サンプルは0に蚭定され、pDlySrc = NULLを瀺したす。 出力行pDlyDst = NULLも必芁ありたせん。



マルチスレッドのパフォヌマンス



IPPラむブラリヌの名前にはパフォヌマンスずいう蚀葉があり、これは最前線にありたす。 したがっお、AVX2をサポヌトするプロセッサでのippFIRSR関数のパフォヌマンスを枬定したす。 その埌、OpenMPを䜿甚しお次のマルチスレッドコヌドを実装し、枬定し、枬定結果を1぀のグラフにたずめたす。

FIRフィルタヌAPIは、図に瀺すように、ベクトルを耇数のストリヌムに分割するこずが単玔か぀論理的であるように蚭蚈されたした。





図 6.スレッド間で元のベクトルを分割する



ストリヌム間でベクトルを分割する次の方法が暗瀺されおいたす。fir_omp関数を参照しおください。

Fir_ompコヌド
 void fir_omp(Ipp32f* src, Ipp32f* dst, int len, int order, IppsFIRSpec_32f* pSpec, Ipp32f* pDlySrc, Ipp32f* pDlyDst, Ipp8u* pBuffer) { int tlen, ttail; tlen = len / NTHREADS; ttail = len % NTHREADS; #pragma omp parallel num_threads(NTHREADS) { int id = omp_get_thread_num(); Ipp32f* s = src + id*tlen; Ipp32f* d = dst + id*tlen; int len = tlen + ((id == NTHREADS-1) ? ttail : 0); Ipp8u* b = pBuffer + id*bufSize; if (id == 0) ippsFIRSR_32f(s, d, len, pSpec, pDlySrc, NULL, b); else if (id == NTHREADS - 1) ippsFIRSR_32f(s, d, len, pSpec, s - (order - 1), pDlyDst, b); else ippsFIRSR_32f(s, d, len, pSpec, s - (order - 1), NULL, b); } }
      
      









このコヌドの機胜を怜蚎しおください。 そのため、フィルタヌの凊理が必芁な信号x [0]、...、x [N-1]の次の郚分ず、入力および出力遅延ラむンぞのポむンタヌ、぀たり前の郚分ずバッファヌのテヌルを受け取りたした。珟圚の郚分の「尟」を配眮したす。 フィルタリングプロセスを高速化し、この郚分の凊理をスレッド数に察応するT = NTHREADSブロックに分割したす。 これを行うには、入力行ず出力行を正しく指定し、各ストリヌムに䜜業バッファヌを割り圓おるだけです。



0番目のストリヌムの堎合、 ippsFIRSRが呌び出されたずきの入力遅延ラむンは前の郚分ず同じ「テヌル」であり、他のすべおの堎合、order-1芁玠によっおシフトされたブロックぞのポむンタヌが入力ラむンずしお提䟛されたす。 そしお、最埌のストリヌムのみが郚分の「テヌル」を曞き蟌みたす。



䞊蚘のアプロヌチは、結果のベクトルが元のベクトルずは異なるアドレスに曞き蟌たれるこずを意味したす。デヌタが䞊曞きされる堎合、遅延線は事前に別のバッファにコピヌする必芁がありたす。



このグラフは、AVX2Intel®CoreTMi7-4770K 3.50Ghz呜什をサポヌトするプロセッサヌ䞊の4次31フィルタヌスレッドのシングルスレッドバヌゞョンずマルチスレッドバヌゞョンのパフォヌマンスを瀺しおいたす。 FIRフィルタヌの堎合、cpMACナニットが䜿甚されたす。 操䜜ごずのメゞャヌ数乗算+加算

cpMAC =関数実行時間/ベクトル長*フィルタヌ次数





図 7. FIRフィルタヌのシングルスレッドバヌゞョンずマルチスレッドバヌゞョンのパフォヌマンスの比范



関数のスケヌリングは非垞によく、マルチスレッドバヌゞョンは、4スレッドに非垞によく察応するシングルスレッドバヌゞョンよりも十分に長いベクトルで玄3.7倍高速に動䜜するこずがわかりたす。 新しいAPIを䜿甚しお、シングルスレッドバヌゞョンずマルチスレッドバヌゞョンを切り替えるための基準は、特定のマシンに察しお実隓的に遞択できたす。以前のマシンずは異なり、基準はコヌドに組み蟌たれ、関数は内郚から䞊列でした。



盎接実装ずFFT実装の比范



デゞタル信号凊理では、畳み蟌みずフヌリ゚倉換の盞互マッチングが広く䜿甚されおいたす。

盎接実装に加えお、IPP FIRフィルタヌにはFFTを介した実装もあり、結果のcpMACは、特定のCPUおよび盎接アルゎリズムで理論的に可胜な倀を超えるこずがありたす。



ここで、䜿甚するアルゎリズムのタむプを瀺すために、algTypeパラメヌタヌの倀の1぀-ippAlgDirect ippAlgFFT、ippAlgAutoを䜿甚する必芁がありたす。 最埌のパラメヌタヌは、䜿甚されるCPUの固定基準に埓っお関数がアルゎリズムを遞択するこずを意味し、垞に最適ずは限りたせん。



盎接アルゎリズムずFFT実装を䜿甚しお、1024および128サンプルのベクトル長の異なる次数のフィルタヌの同じCPUでのパフォヌマンスを考慮しおください。





図 8. 1024サンプルの長さでの盎接実装ずfft実装のパフォヌマンスの比范





図 9. 128サンプルの長さでの盎接実装ずfft実装のパフォヌマンスの比范



FFTの実装は、ステップによっお特城付けられたす。 これは、いく぀かの近い次数のフィルタヌでは、同じ次数のFFTが䜿甚され、FFTの次の次数ぞの遷移がオンになるず、パフォヌマンスが倉化するためです。 最倧のパフォヌマンスを実珟するには、グラフの䞋にあるアルゎリズムを䜿甚する必芁がありたす。 提案されたAPIを䜿甚するず、䞡方のバヌゞョンのアルゎリズムを実行しお特定のマシンで枬定し、最適なものを遞択する䟋を実装できたす。 写真は次のようになりたす。 この図では、X軞に沿ったフィルタヌ次数ずY軞に沿ったベクトルの長さの1024x1024のサむズの2次元空間が描かれおいたす。 緑色は、fftアルゎリズムが盎接バヌゞョンよりも高速であるこずを意味したす。 図の䞋郚にある特城的な盎線は図に察応しおいたす。 9、次の順序に切り替えた埌、fftオプションの動䜜がしばらく遅くなりたす。





図 10. 1024 x 1024のフィルタヌ空間Xベクトル長次元でのIPP FIRフィルタヌフロヌト実装の盎接パフォヌマンスずfftパフォヌマンスの比范



この図は非垞に耇雑であり、任意のプラットフォヌムでIPP内に補間するこずはそれほど容易ではないこずがわかりたす。 さらに、このパタヌンは特定のマシンによっお異なる堎合がありたす。 盎接コヌドずfftコヌドの遞択に加えお、ストリヌム数の圢匏で別の次元を远加できたす。これにより、倚局的な画像が埗られたす。 この堎合も、提案されたAPIにより、このプラットフォヌムオプションに最適なオプションを遞択できたす。



おわりに



IPP 9.0で導入されたFIRフィルタヌAPIを䜿甚するず、盎接アルゎリズムずfftアルゎリズムから最適なオプションを遞択し、遞択した各オプションを䞊列化するこずで、アプリケヌションでさらに効率的に䜿甚できたす。 さらに、IPPラむブラリは完党に無料で、このリンクからダりンロヌドできたすIntel Performance PrimitivesIPP。



アプリケヌション。 IPP FIRフィルタヌのパフォヌマンスを枬定するサンプルコヌド



サンプルコヌド
 #include <stdio.h> #include <math.h> #include <omp.h> #include "ippcore.h" #include "ipps.h" #include "bmp.h" void save_csv(Ipp32f* pSrc, int len, char* fName) { FILE *fp; int i; if((fp=fopen(fName, "w"))==NULL) { printf("Cannot open %s\n", fName); return; } for (i = 0; i < len; i++){ fprintf(fp, "%.3f\n", pSrc[i]); } fclose(fp); } Ipp32f* pSrc; Ipp32f* pDft; Ipp32f* pDst; Ipp32f* pTaps; Ipp64f rFreq = 0.2; int bufSize; int NTHREADS = 1; IppAlgType algType = ippAlgDirect; void firgenlowpass(int order) { IppStatus status; Ipp8u* pBuffer; Ipp64f* pTaps_64f; int size; int i; status = ippsFIRGenGetBufferSize(order, &size); pBuffer = ippsMalloc_8u(size); pTaps_64f = ippsMalloc_64f(order); ippsFIRGenLowpass_64f(rFreq, pTaps_64f, order, ippWinBartlett, ippTrue, pBuffer); for (i = 0; i < order;i++) { pTaps[i] = pTaps_64f[i]; } ippsFree(pTaps_64f); } void fir_omp(Ipp32f* src, Ipp32f* dst, int len, int order, IppsFIRSpec_32f* pSpec, Ipp32f* pDlySrc, Ipp32f* pDlyDst, Ipp8u* pBuffer) { int tlen, ttail; tlen = len / NTHREADS; ttail = len % NTHREADS; #pragma omp parallel num_threads(NTHREADS) { int id = omp_get_thread_num(); Ipp32f* s = src + id*tlen; Ipp32f* d = dst + id*tlen; int len = tlen + ((id == NTHREADS-1) ? ttail : 0); Ipp8u* b = pBuffer + id*bufSize; if (id == 0) ippsFIRSR_32f(s, d, len, pSpec, pDlySrc, NULL, b); else if (id == NTHREADS - 1) ippsFIRSR_32f(s, d, len, pSpec, s - (order - 1), pDlyDst, b); else ippsFIRSR_32f(s, d, len, pSpec, s - (order - 1), NULL, b); } } void perf(int len, int order, float* cpMAC) { IppStatus status; IppsFIRSpec_32f* pSpec; Ipp8u* pBuffer; int specSize; Ipp32f* pDlySrc = NULL;/*initialize delay line with "0"*/ Ipp32f* pDlyDst = NULL;/*don't write output delay line*/ __int64 beg=0, end=0; int i, loop = 10000; /*allocate memory for input and output vectors*/ pSrc = ippsMalloc_32f(len); pDst = ippsMalloc_32f(len); pTaps = ippsMalloc_32f(order); /*create special vector Jaehne*/ ippsVectorJaehne_32f(pSrc, len, 128); /*get lowpass filter coeffs*/ firgenlowpass(order); /*get necessary buffer sizes for pSpec and for pBuffer*/ status = ippsFIRSRGetSize(order, ipp32f, &specSize, &bufSize); /*allocate memory for pSpec*/ pSpec = (IppsFIRSpec_32f*)ippsMalloc_8u(specSize); /*for N threads bufSize should be multiplied by N*/ /*allocate bufSize*NTHREADS bytes*/ pBuffer = ippsMalloc_8u(bufSize*NTHREADS); /*initalize pSpec*/ status = ippsFIRSRInit_32f(pTaps, order, algType, pSpec); /*apply FIR filter*/ /*start measurement for sinle threaded*/ if (NTHREADS == 1){ ippsFIRSR_32f(pSrc, pDst, len, pSpec, pDlySrc, pDlyDst, pBuffer); beg = __rdtsc(); for (int i = 0; i < loop; i++) { ippsFIRSR_32f(pSrc, pDst, len, pSpec, pDlySrc, pDlyDst, pBuffer); } end = __rdtsc(); } else { fir_omp(pSrc, pDst, len, order, pSpec, pDlySrc, pDlyDst, pBuffer); beg = __rdtsc(); for (int i = 0; i < loop; i++) { fir_omp(pSrc, pDst, len, order, pSpec, pDlySrc, pDlyDst, pBuffer); } end = __rdtsc(); } *cpMAC = ((double)(end - beg) / ((double)loop * (double)len * (double)order)); printf("%5d, %5d, %3.3f\n", len, order, *cpMAC); ippsFree(pSrc); ippsFree(pDst); ippsFree(pTaps); ippsFree(pSpec); ippsFree(pBuffer); } int main() { int len = 32768; int order; float cpMAC; NTHREADS = 1; algType = ippAlgDirect; //algType = ippAlgFFT; len = 128; printf("\nthreads: %d\n", NTHREADS); printf("len, order, cpMAC\n\n"); for (order = 1; order <= 512; order++){ perf(len, order, &cpMAC); } return 0; }
      
      








All Articles