Tell me more ×
Mathematics Stack Exchange is a question and answer site for people studying math at any level and professionals in related fields. It's 100% free, no registration required.
up vote 1 down vote favorite

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 :)

share|improve this question
1  
Apply the FFT two both sequences, multiply the results pointwise, then transform back. – dls Mar 16 '12 at 16:11
Well that's the general idea. Except that to do it for each $m$ is no time gain. I was guessing the clever thing to do is to approximate the function $\sum_{j=1}^N K(s - x_j) y_j$ using FFT, and this bit I'm not sure how to do. Besides, the FFT theorem convolution applies directly to convolutions like $\sum a_{i-j} b_j$, so to do the straightforward thing I would need to recompute the sequence $K(x_m-k_j)$ for each $m$, which would eat away any time gain..? – Maciej Mar 16 '12 at 16:20
I've not heard of approximating a function by its Fourier transform. You have $$ (K \ast y) (x_m) = \sum_{j = 1}^N K(x_m - x_j) y_j$$ and $$\widehat{K \ast y} (x_m) = \hat{K}(x_m) \hat{y_m}$$ You can use FFT to compute these: do the multiplication and then do the inverse transform for each $m$. I don't understand how to use the Fourier transform to approximate a function. The other thing you could think about is to numerically approximate it: compute only a few of the function values and then use interpolation. – Matt N. Mar 17 '12 at 8:32
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.) – dbaupp May 23 '12 at 8:34

1 Answer

up vote 1 down vote

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!

share|improve this answer

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.