Hi, I'm using FFT convolution for smoothing. Basically FFT kernel is exp(-x^2) and the signal is something positive. All of it is done in a single pass, so no need for overlap-add and stuff like that. Seems working fine, except now I tried it on quite long data - the FFT length was 262144, kernel being 32768 samples (zero-padded to 262144), calculated in 32-bit floating points, FFT implementation from Intel IPP. Also the kernel obviously goes very low on both ends, approximately 1E-6, and the signal being processed reaches 1E-10 on the right end. The results are fine until some point, where it goes below say 1E-8 and then it goes even negative! Or more precisely the values are jumping around 0. I mean it's convolution of 2 nonnegative signals and the result is still negative. Could it be an accuracy problem or should I look for some trouble in the algorithm? (it's not a big deal, I'm just worried if there isn't some actual algorithmic trouble) jungledmnc
Bad accuracy of very long FFT convolution
Started by ●January 8, 2011
Reply by ●January 8, 20112011-01-08
On Jan 8, 7:11�am, "jungledmnc" <jungledmnc@n_o_s_p_a_m.gmail.com> wrote:> Hi, > > I'm using FFT convolution for smoothing. Basically FFT kernel is exp(-x^2) > and the signal is something positive. All of it is done in a single pass, > so no need for overlap-add and stuff like that. > > Seems working fine, except now I tried it on quite long data - the FFT > length was 262144, kernel being 32768 samples (zero-padded to 262144), > calculated in 32-bit floating points, FFT implementation from Intel IPP. > Also the kernel obviously goes very low on both ends, approximately 1E-6, > and the signal being processed reaches 1E-10 on the right end. The results > are fine until some point, where it goes below say 1E-8 and then it goes > even negative! Or more precisely the values are jumping around 0. > > I mean it's convolution of 2 nonnegative signals and the result is still > negative. Could it be an accuracy problem or should I look for some trouble > in the algorithm? (it's not a big deal, I'm just worried if there isn't > some actual algorithmic trouble)you have to define *specifically* what the *finite* impulse response is (and it cannot be as long as the FFT) and then do overlap-add or overlap-save with that. r b-j
Reply by ●January 8, 20112011-01-08
On Jan 8, 4:11�am, "jungledmnc" <jungledmnc@n_o_s_p_a_m.gmail.com> wrote:> Hi, > > I'm using FFT convolution for smoothing. Basically FFT kernel is exp(-x^2) > and the signal is something positive. All of it is done in a single pass, > so no need for overlap-add and stuff like that. > ...> I mean it's convolution of 2 nonnegative signals and the result is still > negative. Could it be an accuracy problem or should I look for some trouble > in the algorithm? (it's not a big deal, I'm just worried if there isn't > some actual algorithmic trouble) > > jungledmncHow do the magnitudes of the negative values compare to the quantization step size of the most positive convolution value in 32 bit float? That is, how do your 'errors' compare to the size of the expected error sources in your calculations? Dale B. Dalrymple
Reply by ●January 8, 20112011-01-08
>How do the magnitudes of the negative values compare to the >quantization step size of the most positive convolution value in 32 >bit float? That is, how do your 'errors' compare to the size of the >expected error sources in your calculations?I'm not really sure how to answer that, I think it kind of is my question :). So the thing is I'm don't really know how big error I can assume. I usually use this to smoothen out anything (since the kernel is basically a bell curve) and since often it is like 1024 samples kernel and 4096 samples signal, then I just put it together into 8192 samples and convolve it without overlap-add, because I measured 8192 point FFT is faster. In these cases I never ran into problems with accuracy. But I never tried such a huge blocks. It's an exception really, but the kernel is big as well, so I don't think it makes sense to bother with overlap-add (moreover this time it is not such a performance thing). I think with the FFT size of 2^18 and O(log(n)) per sample, it becomes 18 multiplications per FFT per sample + IFFT + complex multiplication. So I'd say about 50 multiplications. And the error is about +- 1E-10. Seems a little bit too much, but who knows. Due to the nature of FFT the sample error is probably propagated between samples (I mean every sample affects every other sample), so maybe that's the problem.
Reply by ●January 8, 20112011-01-08
On Jan 8, 1:33�pm, "jungledmnc" <jungledmnc@n_o_s_p_a_m.gmail.com> wrote:> >How do the magnitudes of the negative values compare to the > >quantization step size of the most positive convolution value in 32 > >bit float? That is, how do your 'errors' compare to the size of the > >expected error sources in your calculations? > > I'm not really sure how to answer that, I think it kind of is my question > :). > ...Since you haven't told us about your signal, no one else can answer it. What is the maximum value of the convolution? In 32 bit floating point there are about 23 bits used in a linear representation of magnitude. How does 23 bits below the convolution maximum compare with the magnitude of the negative values? That is, how do the expected noise sources (quantization in 32 bit float) compare to the 'errors' that concern you? Dale B. Dalrymple
Reply by ●January 9, 20112011-01-09
On 08/01/2011 21:33, jungledmnc wrote: ...> I think with the FFT size of 2^18 and O(log(n)) per sample, it becomes 18 > multiplications per FFT per sample + IFFT + complex multiplication. So I'd > say about 50 multiplications. And the error is about +- 1E-10. Seems a > little bit too much, but who knows. Due to the nature of FFT the sample > error is probably propagated between samples (I mean every sample affects > every other sample), so maybe that's the problem. >A simple way to look at 32bit floats is to understand that it supports approximately seven decimal digits of precision (doubles gives you more like 17 decimal digits). The issue is typically in adding a small number to a large one, which is where things start to break down accuracy-wise - the exponent runs out of space, so the mantissa loses precision. An FFT of that size is not so unusual in audio when performing vanilla convolution reverb (when of course you do need overlap add/save), but I would expect to need to use double precision; i.e. I never bothered to try a 32bit floats version. 262144 works out (at the CD rate) as about 6 seconds. The doubles version runs fine on an "impulse response" of half a minute or so. Believe it or not, that was the first thing my tester tried, being the usual mad musician type. Richard Dobson
Reply by ●January 9, 20112011-01-09
On Jan 8, 4:33�pm, "jungledmnc" <jungledmnc@n_o_s_p_a_m.gmail.com> wrote: ...> I think with the FFT size of 2^18 and O(log(n)) per sample, it becomes 18 > multiplications per FFT per sample + IFFT + complex multiplication. So I'd > say about 50 multiplications. And the error is about +- 1E-10. Seems a > little bit too much, but who knows. Due to the nature of FFT the sample > error is probably propagated between samples (I mean every sample affects > every other sample), so maybe that's the problem.jungle, whether you are doing it with overlap-whatever or it's one big FFT, the question is the same. what is the definition of the FIR? how long is it? it better be shorter than 2^18 or whatever because the FFT will be doing circular convolution in any case. if there is wrapping in the convolution, you will have error. r b-j
