提高FFT的執行速度

我是一名編程初學者，目前正在嘗試處理需要快速傅里葉變換實現的項目。提高FFT的執行速度

我迄今成功地實現以下：

沒有人有任何的替代品和建議，以提高程序的速度，而不對精度失去了。

short FFTMethod::FFTcalc(short int dir,long m,double *x,double *y) 
{ 
long n,i,i1,j,k,i2,l,l1,l2; 
double c1,c2,tx,ty,t1,t2,u1,u2,z; 

/* Calculate the number of points */ 
n = 1; 
for (i=0;i<m;i++) 
    n *= 2; 

/* Do the bit reversal */ 
i2 = n >> 1; 
j = 0; 
for (i=0;i<n-1;i++) { 
    if (i < j) { 
    tx = x[i]; 
    ty = y[i]; 
    x[i] = x[j]; 
    y[i] = y[j]; 
    x[j] = tx; 
    y[j] = ty; 
    } 
    k = i2; 
    while (k <= j) { 
    j -= k; 
    k >>= 1; 
    } 
    j += k; 
} 

/* Compute the FFT */ 
c1 = -1.0; 
c2 = 0.0; 
l2 = 1; 
for (l=0;l<m;l++) { 
    l1 = l2; 
    l2 <<= 1; 
    u1 = 1.0; 
    u2 = 0.0; 
    for (j=0;j<l1;j++) { 
    for (i=j;i<n;i+=l2) { 
     i1 = i + l1; 
     t1 = u1 * x[i1] - u2 * y[i1]; 
     t2 = u1 * y[i1] + u2 * x[i1]; 
     x[i1] = x[i] - t1; 
     y[i1] = y[i] - t2; 
     x[i] += t1; 
     y[i] += t2; 
    } 
    z = u1 * c1 - u2 * c2; 
    u2 = u1 * c2 + u2 * c1; 
    u1 = z; 
    } 
    c2 = sqrt((1.0 - c1)/2.0); 
    if (dir == 1) 
    c2 = -c2; 
    c1 = sqrt((1.0 + c1)/2.0); 
    } 

/* Scaling for forward transform */ 
if (dir == 1) { 
    for (i=0;i<n;i++) { 
     x[i] /= n; 
     y[i] /= n; 
    } 
} 


    return(1); 
}

來源

2011-12-21 sagarn

除非您爲了理解目的需要自己寫，否則FFTW（http://www.fftw.org/）是一個很棒的圖書館。這是一個自我調整，超快速和可靠的實現，你可以從C++調用它（請參閱http://www.fftw.org/faq/section2.html#cplusplus）。 – 2011-12-21 09:29:21

我非常喜歡FFTReal。 http://ldesoras.free.fr/prod.html – 2011-12-21 09:36:39

爲什麼你寫自己的實現，而不是使用其中一個無數的庫，它們可能都更快，更好的測試，更準確和更多的功能？ – PlasmaHH 2011-12-21 10:50:35

我最近在Eric Postpischil的Construction of a high performance FFTs上發現了這個優秀的PDF。我自己開發了幾個FFT，我知道與商業圖書館競爭是多麼困難。如果你的FFT只比Intel或FFTW慢4倍，相信我你做得很好！然而，你可以競爭，這是如何。

總結這篇文章，作者指出Radix2 FFT簡單但效率低下，最有效的結構是基4 FFT。更有效的方法是Radix8，但這通常不適合CPU的寄存器，因此Radix4是首選。您可以執行Radix2 FFT的10個階段（如2^10 - 1024）或5個階段的Radix4 FFT（4^5 = 1024），可以分階段構建FFT，因此要計算1024點FFT，。如果您願意，您甚至可以按8 * 4 * 4 * 4 * 2的級數計算1024點FFT。更少的階段意味着更少的讀寫內存（FFT性能的瓶頸是內存帶寬），因此動態選擇4,8或更高的基數是必須的。由於所有權重爲1 + 0i，0 + 1i，-1 + 0i，0-1i和Radix4蝶形碼可以被寫入以完全適合緩存，所以Radix4階段是特別有效的。其次，FFT中的每個階段是不一樣的。第一階段的權重都等於1 + 0i。計算這個權重並且甚至乘以它是沒有意義的，因爲它是一個複數乘以1，所以第一階段可以在沒有權重的情況下執行。最後一個階段也可以被區別對待，並且可以用來執行時間抽取（位反轉）。 Eric Postpischil的文件涵蓋了所有這些。

權重可以預先計算並存儲在表中。正弦/餘弦計算在x86硬件上每次需要大約100-150個週期，因此預計算它們可以節省總計算時間的10-20％，因爲在這種情況下，存儲器訪問速度比CPU計算速度快。使用快速算法一次性計算sincos特別有用（請注意，cos等於sqrt（1.0 - 正弦*正弦），或者使用查表，cos僅爲正弦相移）。

最後，一旦您擁有了超級簡化的FFT實現，您可以利用SIMD向量化來計算蝶形程序內每個週期的4倍浮點運算或2倍雙浮點運算，以進一步提高100-300％的速度。綜合以上所述，你可以擁有一個漂亮而又快速的FFT！

若要進一步您可以通過提供針對特定處理器體系結構的FFT階段的不同實現來進行優化。緩存大小，寄存器數量，SSE/SSE2/3/4指令集等在每臺機器上都不相同，因此選擇一種適合所有方法的方法往往會被有針對性的例程毆打。例如，在FFTW中，許多較小尺寸的FFT是針對特定體系結構高度優化的展開（無循環）實現。通過結合這些較小的結構（例如RadixN例程），您可以爲當前任務選擇最快，最好的例程。

來源

2011-12-21 14:12:35

非常感謝。你一直很有幫助。我會嘗試做出改變。 – sagarn 2011-12-23 09:16:19

性能調整是一種黑色藝術。我會建議創建一個測試應用程序，它運行多個不同FFT方法的迭代並對它們進行計時，並將結果的準確性和轉換速度與已知的FFT實現（例如FFTW）進行比較。不要完全改變一個實現，而應該保留它，但創建新的實現並進行比較。你會驚訝什麼會不會提高性能。例如。減少乘法的次數可能不會像確保按順序執行RAM讀取一樣有效，並且儘可能少的次數！ – 2011-12-23 09:27:21

如果評論對你有幫助，請投票。謝謝！ :-) – 2011-12-31 10:44:46

雖然我不能給你一個表現暗示，現在，我想給你一些優化建議，是一個評論太長：

如果您還沒有爲此，請立即爲您的代碼編寫一些正確性測試。簡單的測試，比如「對這個數組進行FFT並查看結果是否與我提供的結果相匹配」就足夠了，但是在優化代碼之前，您需要一個堅定而自動的單元測試來確認您的優化代碼是否正確。
然後分析您的代碼，看看實際的瓶頸在哪裏。雖然我懷疑最內層的循環for (i=j;i<n;i+=l2) {，但看到比相信更好。

來源

2011-12-21 09:44:56 thiton

這看起來是一個基本的基數radix-2 FFT的實現，從舊的教科書出來。有許多數十年的關於以各種方式優化FFT的論文，這取決於許多因素。例如，你的數據是否比CPU緩存小？例如，如果數據向量加上係數表將適合CPU dcache，並且/或者如果乘積比CPU上的存儲器訪問慢得多，則預先計算旋轉因子表可能會減少總週期計數重複使用FFT。但是，如果不是，預計算實際上可能會變慢。基準。因人而異。

來源

2011-12-21 10:07:37 hotpaw2

是的，你是對的@ hotpaw2，我提到了一本名爲C的數值食譜，因爲我發現它是最好的開始。然而，這只是第一次嘗試，在完成項目之前我有很多優化工作要做。是的，數據比CPU緩存小。 – sagarn 2011-12-21 10:31:36

有幾件事情，我可以建議嘗試：

不要交換輸入要素，而不是計算出位反轉指標。這將爲您節省大量內存讀取和寫入操作。
如果您正在執行許多相同大小的FFT，請預先計算係數。這將節省一些計算。
使用基數-4 FFT而不是基數-2。這將導致內部循環中更少的迭代。

最終答案當然可以通過分析代碼來找到。

來源

2011-12-21 10:07:50

謝謝@亞歷克斯。我會嘗試這樣做。 – sagarn 2011-12-21 10:29:56

如果我理解你是正確的，（1）是一個壞主意。你正在節省一些內存操作，但你也會隨機化更多的內存操作，因爲它破壞了主循環中CPU緩存的優點。 – 2012-04-12 13:24:38

@JonHarrop：不會交換「隨機化」嗎？無論如何，無論如何，在交換時間或以後如果沒有交換，您將不可避免地訪問相同的數據*和*。 – 2012-04-13 05:38:08

提高FFT的執行速度

回答

相關問題