2011-08-30 111 views
5

我需要對2個大數組進行一維卷積。我在C#中使用這段代碼,但它需要很長時間才能運行。1D無FFT快速卷積

我知道,我知道! FFT卷積非常快。但在這個項目中,我無法使用它。 這是該項目不使用FFT的一個約束(請不要問爲什麼:/)。

這是我在C#代碼(從MATLAB移植,順便說一句):

var result = new double[input.Length + filter.Length - 1]; 
for (var i = 0; i < input.Length; i++) 
{ 
    for (var j = 0; j < filter.Length; j++) 
    { 
     result[i + j] += input[i] * filter[j]; 
    } 
} 

所以任何人都知道任何快速卷積算法WIDTHOUT FFT?

+5

雖然你說過不要問,爲什麼你不能使用FFT?如果這是針對明確禁止的類項目,則應該將其標記爲家庭作業。 – templatetypedef

+0

C#可以調用CUDA嗎?如果可以的話,你可以使用並行指令,這相當大地加速了幼稚卷積。或者你可以使用Winograd變換或其他東西(不是Cooley-Tukey經典FFT,如果這足夠滿足你的「不FFT」規則)。或者,如果您知道輸入或濾波器的某些信息(例如只有某些頻率存在或某種情況),則可以使用該知識。你將不得不更加具體地瞭解你的限制以及你可能擁有的任何外部知識。 – mtrw

回答

5

卷積是數字相同的與額外的纏繞步驟的多項式乘法。因此,可以使用所有的多項式和大整數乘法算法來執行卷積。

FFT是獲得快速O(n log(n))運行時的唯一方法。但是你仍然可以使用分治法(如Karatsuba's algorithm)獲得亞二次運行時間。

Karatsuba的算法很容易實現,一旦你明白它是如何工作的。它運行在O(n^1.585),並且可能會比嘗試超級優化經典的O(n^2)方法更快。

0

這裏有兩種可能會導致較小的加速,但你需要測試以確保。

  1. 展開內部循環以刪除一些測試。如果知道過濾器長度將總是一個倍數(如果有的話),這將更容易。
  2. 顛倒循環的順序。 filter.length是否通過整個陣列。這在內部循環中的解引用較少,但可能具有較差的緩存行爲。
4

你可以減少索引訪問次數result,還有Length屬性:

int inputLength = filter.Length; 
int filterLength = filter.Length; 
var result = new double[inputLength + filterLength - 1]; 
for (int i = resultLength; i >= 0; i--) 
{ 
    double sum = 0; 
    // max(i - input.Length + 1,0) 
    int n1 = i < inputLength ? 0 : i - inputLength + 1; 
    // min(i, filter.Length - 1) 
    int n2 = i < filterLength ? i : filterLength - 1; 
    for (int j = n1; j <= n2; j++) 
    { 
     sum += input[i - j] * filter[j]; 
    } 
    result[i] = sum; 
} 

如果進一步分裂外循環,你可以擺脫一些重複的條件句。 (這裏假定0 < filterLengthinputLengthresultLength ≤)

int inputLength = filter.Length; 
int filterLength = filter.Length; 
int resultLength = inputLength + filterLength - 1; 

var result = new double[resultLength]; 

for (int i = 0; i < filterLength; i++) 
{ 
    double sum = 0; 
    for (int j = i; j >= 0; j--) 
    { 
     sum += input[i - j] * filter[j]; 
    } 
    result[i] = sum; 
} 
for (int i = filterLength; i < inputLength; i++) 
{ 
    double sum = 0; 
    for (int j = filterLength - 1; j >= 0; j--) 
    { 
     sum += input[i - j] * filter[j]; 
    } 
    result[i] = sum; 
} 
for (int i = inputLength; i < resultLength; i++) 
{ 
    double sum = 0; 
    for (int j = i - inputLength + 1; j < filterLength; j++) 
    { 
     sum += input[i - j] * filter[j]; 
    } 
    result[i] = sum; 
} 
0

您可以使用特殊的IIR濾波器。 然後處理,如:

y(n)= a1*y(n-1)+b1*y(n-2)...+a2*x(n-1)+b2*x(n-2)...... 

我認爲它更快。