3
$\begingroup$

I have two sequences of the same length, $(x_i), i=1, 2, \ldots, N$ and $(y_i), i=1, 2, \ldots, N$ and a function $K(t) = -t \times \exp(-t^2 / 2)/ \sqrt{2 \pi}$.

I need to compute the following quantity for each $m=1, 2, \ldots, N$:

$\sum_{j=1}^N K(x_m - x_j) \times y_j$

which is a tad slow when done directly (I need it when $N = O(10^4)$ ). I know this can be much improved on with the use of FFT, however it is something I never really worked with. Could anyone suggest any links or a way to rewrite it in FFT form?

I'd be very grateful for any help :)

  • 1
    I'm not sure if this is precisely what you are talking about but naive convolution is $O(n^2)$, while the FFT method is only $O(n \log n)$ (it is 2 FFTs, which are $O(n \log n)$, and one pointwise multiply, which is $O(n)$, so the process is $O(n \log n)$ overall.)2012-05-23

2 Answers 2

2

Here you have a good link.

Basically you need to compute the FFT of each signal individually, multiply the spectrums and the do the inverse FFT of the resulting sectrum. That will be the result of the convolution.

Good luck!

0

Need to account for wrap around issues

1) Shift the K sequence such that it is centered in the middle of the array call K1 of lenght n1 2) Store the K1 array in a new array K2 with double the size of K1 n2 = 2Xn1 such that K2(1 to n1 ) = complex(K1,0) and K2( n1+1 to n2 ) = 0 where K2 is complex 3) Next shift the K2 array over to the left by n1/2-1 this will center the maximum of the K function in the first bin of the K2 array 4) Store the y data in an array y2 of of lenght n2 such that y2(1 to n1 ) = complex(y,0) and y2(n1+1 to n2) = 0 where y2 is complex

5) calculate the FFT of K2, i.e., FFT_K2 = FFT( K2, 1 ) 1 => exp(iwt) 6) calculate the FFT of y2, i.e., FFT_y2 = FFT( y2, 1 ) 1 => exp(iwt) 7) Let FFT_d2 = CONJUGATE(FFT_K2) * FFT_y2 8) Inverse FFT, i.e., d2 = FFT( FFT_d2, -1 ) -1 => exp(-iwt) 9) convolution = REAL( ds(1 to n1) ) 10) this assumes the exp(-iwt) has the 1/n factor