Hi.
I use FFTW for FFT analysis. Padded my data with zeros.
Before padding my data, my fft looked like what is expected with some
frequencies showing stronger power than others. After I padded my data, the
resulting plot showed a decreasing exponential trend interleaved with the
actual expected data. I spent hours but couldnt figure out a reason except
to conclude that there is something I'm missing about fftw.
My codes are as follow:
// initialized the usual fftw variables
// fftw_complex *fft_in, *fft_out;
// fftw_plan p;
// copy the data
for (i = 0; i < n; i++)
{fft_in[i][0] = mydata[i];
fft_in[i][1] = 0;}
// for the part where there are padded zeros, the codes are
for (i = 0; i < n+zeropadsize; i++)
{
if (i < n) fft_in[i][0] = mydata[i];
else fft_in[i][0] = 0;
fft_in[i][1] = 0;}
// follow by the usual compute the plan and execute the fft.
// of course when there is no zeropadding, the zeropadsize is = 0
p = fftw_plan_dft_1d((n+zeropadsize), fft_in, fft_out, FFTW_FORWARD,
FFTW_ESTIMATE);
fftw_execute(p);
Is my way of padding the data correct for use with fftw?
Any suggestions are much appreciated.
Many thanks!
Marc
Strange results from zero padded FFT
Started by ●June 2, 2011
Reply by ●June 2, 20112011-06-02
On Jun 2, 8:46�am, "Marc2050" <maarcc@n_o_s_p_a_m.gmail.com> wrote:> Hi. > > I use FFTW for FFT analysis. Padded my data with zeros. > Before padding my data, my fft looked like what is expected with some > frequencies showing stronger power than others. After I padded my data, the > resulting plot showed a decreasing exponential trend interleaved with the > actual expected data. I spent hours but couldnt figure out a reason except > to conclude that there is something I'm missing about fftw. > > My codes are as follow: > // initialized the usual fftw variables > // fftw_complex *fft_in, *fft_out; > // fftw_plan p; > > // copy the data > for (i = 0; i < n; i++) > {fft_in[i][0] = mydata[i]; > fft_in[i][1] = 0;} > > // for the part where there are padded zeros, the codes are > for (i = 0; i < n+zeropadsize; i++) > { > if (i < n) fft_in[i][0] = mydata[i]; > else fft_in[i][0] = 0; > fft_in[i][1] = 0;} > > // follow by the usual compute the plan and execute the fft. > // of course when there is no zeropadding, the zeropadsize is = 0 > p = fftw_plan_dft_1d((n+zeropadsize), fft_in, fft_out, FFTW_FORWARD, > FFTW_ESTIMATE); > fftw_execute(p); > > Is my way of padding the data correct for use with fftw? > Any suggestions are much appreciated. > > Many thanks! > MarcWould have to look at your plot, but zero padding does smear the energy around the bins quite a bit, you can no longer have a bin centered results which you might of had before you zero padded, the more % padding the worse the FFT will look.
Reply by ●June 2, 20112011-06-02
On Jun 2, 8:46 am, "Marc2050" <maarcc@n_o_s_p_a_m.gmail.com> wrote:> Hi. > > I use FFTW for FFT analysis. Padded my data with zeros. > Before padding my data, my fft looked like what is expected with some > frequencies showing stronger power than others. After I padded my data, > the > resulting plot showed a decreasing exponential trend interleaved with the > actual expected data. I spent hours but couldnt figure out a reason except > to conclude that there is something I'm missing about fftw. > > My codes are as follow: > // initialized the usual fftw variables > // fftw_complex *fft_in, *fft_out; > // fftw_plan p; > > // copy the data > for (i = 0; i < n; i++) > {fft_in[i][0] = mydata[i]; > fft_in[i][1] = 0;} > > // for the part where there are padded zeros, the codes are > for (i = 0; i < n+zeropadsize; i++) > { > if (i < n) fft_in[i][0] = mydata[i]; > else fft_in[i][0] = 0; > fft_in[i][1] = 0;} > > // follow by the usual compute the plan and execute the fft. > // of course when there is no zeropadding, the zeropadsize is = 0 > p = fftw_plan_dft_1d((n+zeropadsize), fft_in, fft_out, FFTW_FORWARD, > FFTW_ESTIMATE); > fftw_execute(p); > > Is my way of padding the data correct for use with fftw? > Any suggestions are much appreciated. > > Many thanks! > MarcWhat is the mean (DC) value of mydata[0...n-1] ? Is it non-zero? Try padding with that.
Reply by ●June 2, 20112011-06-02
On 06/02/2011 05:46 AM, Marc2050 wrote:> Hi. > > I use FFTW for FFT analysis. Padded my data with zeros. > Before padding my data, my fft looked like what is expected with some > frequencies showing stronger power than others. After I padded my data, the > resulting plot showed a decreasing exponential trend interleaved with the > actual expected data. I spent hours but couldnt figure out a reason except > to conclude that there is something I'm missing about fftw. > > My codes are as follow: > // initialized the usual fftw variables > // fftw_complex *fft_in, *fft_out; > // fftw_plan p; > > // copy the data > for (i = 0; i< n; i++) > {fft_in[i][0] = mydata[i]; > fft_in[i][1] = 0;} > > // for the part where there are padded zeros, the codes are > for (i = 0; i< n+zeropadsize; i++) > { > if (i< n) fft_in[i][0] = mydata[i]; > else fft_in[i][0] = 0; > fft_in[i][1] = 0;} > > // follow by the usual compute the plan and execute the fft. > // of course when there is no zeropadding, the zeropadsize is = 0 > p = fftw_plan_dft_1d((n+zeropadsize), fft_in, fft_out, FFTW_FORWARD, > FFTW_ESTIMATE); > fftw_execute(p); > > Is my way of padding the data correct for use with fftw? > Any suggestions are much appreciated.It's almost certainly not fftw, but the properties of the discrete Fourier transform itself that are giving you grief. If you're going to zero pad then you certainly need to remove the data average, and may even want to detrend it as well. Then zero pad the results. While you're at it, using something other than a rectangular window may not be a bad thing. -- Tim Wescott Wescott Design Services http://www.wescottdesign.com Do you need to implement control loops in software? "Applied Control Theory for Embedded Systems" was written for you. See details at http://www.wescottdesign.com/actfes/actfes.html
Reply by ●June 2, 20112011-06-02
On 6/2/2011 7:49 AM, Tim Wescott wrote:> On 06/02/2011 05:46 AM, Marc2050 wrote: >> Hi. >> >> I use FFTW for FFT analysis. Padded my data with zeros. >> Before padding my data, my fft looked like what is expected with some >> frequencies showing stronger power than others. After I padded my >> data, the >> resulting plot showed a decreasing exponential trend interleaved with the >> actual expected data. I spent hours but couldnt figure out a reason >> except >> to conclude that there is something I'm missing about fftw. >> >> My codes are as follow: >> // initialized the usual fftw variables >> // fftw_complex *fft_in, *fft_out; >> // fftw_plan p; >> >> // copy the data >> for (i = 0; i< n; i++) >> {fft_in[i][0] = mydata[i]; >> fft_in[i][1] = 0;} >> >> // for the part where there are padded zeros, the codes are >> for (i = 0; i< n+zeropadsize; i++) >> { >> if (i< n) fft_in[i][0] = mydata[i]; >> else fft_in[i][0] = 0; >> fft_in[i][1] = 0;} >> >> // follow by the usual compute the plan and execute the fft. >> // of course when there is no zeropadding, the zeropadsize is = 0 >> p = fftw_plan_dft_1d((n+zeropadsize), fft_in, fft_out, FFTW_FORWARD, >> FFTW_ESTIMATE); >> fftw_execute(p); >> >> Is my way of padding the data correct for use with fftw? >> Any suggestions are much appreciated. > > It's almost certainly not fftw, but the properties of the discrete > Fourier transform itself that are giving you grief. > > If you're going to zero pad then you certainly need to remove the data > average, and may even want to detrend it as well. Then zero pad the > results. While you're at it, using something other than a rectangular > window may not be a bad thing. >Tim, ??? Removing the average only affects one output sample and that one is obvious. Am I missing something? Fred
Reply by ●June 2, 20112011-06-02
On 06/02/2011 09:42 AM, Fred Marshall wrote:> On 6/2/2011 7:49 AM, Tim Wescott wrote: >> On 06/02/2011 05:46 AM, Marc2050 wrote: >>> Hi. >>> >>> I use FFTW for FFT analysis. Padded my data with zeros. >>> Before padding my data, my fft looked like what is expected with some >>> frequencies showing stronger power than others. After I padded my >>> data, the >>> resulting plot showed a decreasing exponential trend interleaved with >>> the >>> actual expected data. I spent hours but couldnt figure out a reason >>> except >>> to conclude that there is something I'm missing about fftw. >>> >>> My codes are as follow: >>> // initialized the usual fftw variables >>> // fftw_complex *fft_in, *fft_out; >>> // fftw_plan p; >>> >>> // copy the data >>> for (i = 0; i< n; i++) >>> {fft_in[i][0] = mydata[i]; >>> fft_in[i][1] = 0;} >>> >>> // for the part where there are padded zeros, the codes are >>> for (i = 0; i< n+zeropadsize; i++) >>> { >>> if (i< n) fft_in[i][0] = mydata[i]; >>> else fft_in[i][0] = 0; >>> fft_in[i][1] = 0;} >>> >>> // follow by the usual compute the plan and execute the fft. >>> // of course when there is no zeropadding, the zeropadsize is = 0 >>> p = fftw_plan_dft_1d((n+zeropadsize), fft_in, fft_out, FFTW_FORWARD, >>> FFTW_ESTIMATE); >>> fftw_execute(p); >>> >>> Is my way of padding the data correct for use with fftw? >>> Any suggestions are much appreciated. >> >> It's almost certainly not fftw, but the properties of the discrete >> Fourier transform itself that are giving you grief. >> >> If you're going to zero pad then you certainly need to remove the data >> average, and may even want to detrend it as well. Then zero pad the >> results. While you're at it, using something other than a rectangular >> window may not be a bad thing. >> > > Tim, > > ??? Removing the average only affects one output sample and that one is > obvious. Am I missing something?Yes, you are. This is a graph of a function with just DC offset, no other (or little other) content: -------------------- .................... What is it's Fourier transform? (Hint -- you're correct so far) This is a graph of a function with just DC offset that's been zero-padded: -------------------- -----------------....................--------------------- What's it's Fourier transform? What is the Fourier transform of this zero-padded DC offset going to do when it's added to the Fourier transform of the function that you're trying to extract from the un-padded signal, particularly if that signal's spectrum is well distributed in frequency? If you window (even zero padding is 'windowing' in this context), you should _always_ take out the DC offset -- and very possibly any linear trend -- before you window. Otherwise your result will be a great big spectrum of the windowing function itself, swamping out the information that you want to see. -- Tim Wescott Wescott Design Services http://www.wescottdesign.com Do you need to implement control loops in software? "Applied Control Theory for Embedded Systems" was written for you. See details at http://www.wescottdesign.com/actfes/actfes.html
Reply by ●June 2, 20112011-06-02
"Fred Marshall" <fmarshallxremove_the_x@acm.org> wrote in message news:j4PFp.16088$Ky6.9168@en-nntp-16.dc1.easynews.com...> On 6/2/2011 7:49 AM, Tim Wescott wrote: >> On 06/02/2011 05:46 AM, Marc2050 wrote: >>> Hi. >>> >>> I use FFTW for FFT analysis. Padded my data with zeros. >>> Before padding my data, my fft looked like what is expected with some >>> frequencies showing stronger power than others. After I padded my >>> data, the >>> resulting plot showed a decreasing exponential trend interleaved with >>> the >>> actual expected data. I spent hours but couldnt figure out a reason >>> except >>> to conclude that there is something I'm missing about fftw. >>> >>> My codes are as follow: >>> // initialized the usual fftw variables >>> // fftw_complex *fft_in, *fft_out; >>> // fftw_plan p; >>> >>> // copy the data >>> for (i = 0; i< n; i++) >>> {fft_in[i][0] = mydata[i]; >>> fft_in[i][1] = 0;} >>> >>> // for the part where there are padded zeros, the codes are >>> for (i = 0; i< n+zeropadsize; i++) >>> { >>> if (i< n) fft_in[i][0] = mydata[i]; >>> else fft_in[i][0] = 0; >>> fft_in[i][1] = 0;} >>> >>> // follow by the usual compute the plan and execute the fft. >>> // of course when there is no zeropadding, the zeropadsize is = 0 >>> p = fftw_plan_dft_1d((n+zeropadsize), fft_in, fft_out, FFTW_FORWARD, >>> FFTW_ESTIMATE); >>> fftw_execute(p); >>> >>> Is my way of padding the data correct for use with fftw? >>> Any suggestions are much appreciated. >> >> It's almost certainly not fftw, but the properties of the discrete >> Fourier transform itself that are giving you grief. >> >> If you're going to zero pad then you certainly need to remove the data >> average, and may even want to detrend it as well. Then zero pad the >> results. While you're at it, using something other than a rectangular >> window may not be a bad thing. >> > > Tim, > > ??? Removing the average only affects one output sample and that one is > obvious. Am I missing something? > > FredThe FFT length is n+zeropadsize Samples x[0] to x[n-1] have a different mean value from samples x[n] to x[n+zeropadsize-1] There is a rectangular waveform (with duty cycle n : n+zeropadsize) superimposed on the data.
Reply by ●June 2, 20112011-06-02
On Jun 2, 12:55�pm, Tim Wescott <t...@seemywebsite.com> wrote:> On 06/02/2011 09:42 AM, Fred Marshall wrote: > > > > > > > On 6/2/2011 7:49 AM, Tim Wescott wrote: > >> On 06/02/2011 05:46 AM, Marc2050 wrote: > >>> Hi. > > >>> I use FFTW for FFT analysis. Padded my data with zeros. > >>> Before padding my data, my fft looked like what is expected with some > >>> frequencies showing stronger power than others. After I padded my > >>> data, the > >>> resulting plot showed a decreasing exponential trend interleaved with > >>> the > >>> actual expected data. I spent hours but couldnt figure out a reason > >>> except > >>> to conclude that there is something I'm missing about fftw. > > >>> My codes are as follow: > >>> // initialized the usual fftw variables > >>> // fftw_complex *fft_in, *fft_out; > >>> // fftw_plan p; > > >>> // copy the data > >>> for (i = 0; i< n; i++) > >>> {fft_in[i][0] = mydata[i]; > >>> fft_in[i][1] = 0;} > > >>> // for the part where there are padded zeros, the codes are > >>> for (i = 0; i< n+zeropadsize; i++) > >>> { > >>> if (i< n) fft_in[i][0] = mydata[i]; > >>> else fft_in[i][0] = 0; > >>> fft_in[i][1] = 0;} > > >>> // follow by the usual compute the plan and execute the fft. > >>> // of course when there is no zeropadding, the zeropadsize is = 0 > >>> p = fftw_plan_dft_1d((n+zeropadsize), fft_in, fft_out, FFTW_FORWARD, > >>> FFTW_ESTIMATE); > >>> fftw_execute(p); > > >>> Is my way of padding the data correct for use with fftw? > >>> Any suggestions are much appreciated. > > >> It's almost certainly not fftw, but the properties of the discrete > >> Fourier transform itself that are giving you grief. > > >> If you're going to zero pad then you certainly need to remove the data > >> average, and may even want to detrend it as well. Then zero pad the > >> results. While you're at it, using something other than a rectangular > >> window may not be a bad thing. > > > Tim, > > > ??? Removing the average only affects one output sample and that one is > > obvious. Am I missing something? > > Yes, you are. > > This is a graph of a function with just DC offset, no other (or little > other) content: > > -------------------- > > .................... > > What is it's Fourier transform? �(Hint -- you're correct so far) > > This is a graph of a function with just DC offset that's been zero-padded: > > � � � � � � � � � -------------------- > > -----------------....................--------------------- > > What's it's Fourier transform? �What is the Fourier transform of this > zero-padded DC offset going to do when it's added to the Fourier > transform of the function that you're trying to extract from the > un-padded signal, particularly if that signal's spectrum is well > distributed in frequency? > > If you window (even zero padding is 'windowing' in this context), you > should _always_ take out the DC offset -- and very possibly any linear > trend -- before you window. �Otherwise your result will be a great big > spectrum of the windowing function itself, swamping out the information > that you want to see. > > -- > > Tim Wescott > Wescott Design Serviceshttp://www.wescottdesign.com > > Do you need to implement control loops in software? > "Applied Control Theory for Embedded Systems" was written for you. > See details athttp://www.wescottdesign.com/actfes/actfes.html- Hide quoted text - > > - Show quoted text -that doesn't make any sense since that logic applies to any input of any frequency, not just dc (any sine wave of frequency x is chopped off when you pad, not just the DC, and the transform of a chopped off sinewave isn't pretty)
Reply by ●June 2, 20112011-06-02
steve <bungalow_steve@yahoo.com> wrote:> On Jun 2, 12:55�pm, Tim Wescott <t...@seemywebsite.com> wrote:>> Yes, you are.(snip)>> What's it's Fourier transform? �What is the Fourier transform of this >> zero-padded DC offset going to do when it's added to the Fourier >> transform of the function that you're trying to extract from the >> un-padded signal, particularly if that signal's spectrum is well >> distributed in frequency?>> If you window (even zero padding is 'windowing' in this context), you >> should _always_ take out the DC offset -- and very possibly any linear >> trend -- before you window. �Otherwise your result will be a great big >> spectrum of the windowing function itself, swamping out the information >> that you want to see.(snip)> that doesn't make any sense since that logic applies to any input of > any frequency, not just dc (any sine wave of frequency x is chopped > off when you pad, not just the DC, and the transform of a chopped off > sinewave isn't pretty)He also said: "While you're at it, using something other than a rectangular window may not be a bad thing." But in many cases the DC offset is large relative to the others. If you use a rounded rectangular window, that will reduce the high-frequency terms due to the chopped-off sine. -- glen
Reply by ●June 2, 20112011-06-02
On 06/02/2011 10:13 AM, steve wrote:> On Jun 2, 12:55 pm, Tim Wescott<t...@seemywebsite.com> wrote: >> On 06/02/2011 09:42 AM, Fred Marshall wrote: >> >> >> >> >> >>> On 6/2/2011 7:49 AM, Tim Wescott wrote: >>>> On 06/02/2011 05:46 AM, Marc2050 wrote: >>>>> Hi. >> >>>>> I use FFTW for FFT analysis. Padded my data with zeros. >>>>> Before padding my data, my fft looked like what is expected with some >>>>> frequencies showing stronger power than others. After I padded my >>>>> data, the >>>>> resulting plot showed a decreasing exponential trend interleaved with >>>>> the >>>>> actual expected data. I spent hours but couldnt figure out a reason >>>>> except >>>>> to conclude that there is something I'm missing about fftw. >> >>>>> My codes are as follow: >>>>> // initialized the usual fftw variables >>>>> // fftw_complex *fft_in, *fft_out; >>>>> // fftw_plan p; >> >>>>> // copy the data >>>>> for (i = 0; i< n; i++) >>>>> {fft_in[i][0] = mydata[i]; >>>>> fft_in[i][1] = 0;} >> >>>>> // for the part where there are padded zeros, the codes are >>>>> for (i = 0; i< n+zeropadsize; i++) >>>>> { >>>>> if (i< n) fft_in[i][0] = mydata[i]; >>>>> else fft_in[i][0] = 0; >>>>> fft_in[i][1] = 0;} >> >>>>> // follow by the usual compute the plan and execute the fft. >>>>> // of course when there is no zeropadding, the zeropadsize is = 0 >>>>> p = fftw_plan_dft_1d((n+zeropadsize), fft_in, fft_out, FFTW_FORWARD, >>>>> FFTW_ESTIMATE); >>>>> fftw_execute(p); >> >>>>> Is my way of padding the data correct for use with fftw? >>>>> Any suggestions are much appreciated. >> >>>> It's almost certainly not fftw, but the properties of the discrete >>>> Fourier transform itself that are giving you grief. >> >>>> If you're going to zero pad then you certainly need to remove the data >>>> average, and may even want to detrend it as well. Then zero pad the >>>> results. While you're at it, using something other than a rectangular >>>> window may not be a bad thing. >> >>> Tim, >> >>> ??? Removing the average only affects one output sample and that one is >>> obvious. Am I missing something? >> >> Yes, you are. >> >> This is a graph of a function with just DC offset, no other (or little >> other) content: >> >> -------------------- >> >> .................... >> >> What is it's Fourier transform? (Hint -- you're correct so far) >> >> This is a graph of a function with just DC offset that's been zero-padded: >> >> -------------------- >> >> -----------------....................--------------------- >> >> What's it's Fourier transform? What is the Fourier transform of this >> zero-padded DC offset going to do when it's added to the Fourier >> transform of the function that you're trying to extract from the >> un-padded signal, particularly if that signal's spectrum is well >> distributed in frequency? >> >> If you window (even zero padding is 'windowing' in this context), you >> should _always_ take out the DC offset -- and very possibly any linear >> trend -- before you window. Otherwise your result will be a great big >> spectrum of the windowing function itself, swamping out the information >> that you want to see. >> > > that doesn't make any sense since that logic applies to any input of > any frequency, not just dc (any sine wave of frequency x is chopped > off when you pad, not just the DC, and the transform of a chopped off > sinewave isn't pretty)The DC component of any real-world signal is often overwhelmingly larger than any other. If you had a signal (e.g. samples of power line voltage) with a strong sine component you'd do well to remove that, too (or instead). And the transform of a chopped-off sine wave _isn't_ pretty, which is why a rectangular window is rarely the best choice. -- Tim Wescott Wescott Design Services http://www.wescottdesign.com Do you need to implement control loops in software? "Applied Control Theory for Embedded Systems" was written for you. See details at http://www.wescottdesign.com/actfes/actfes.html
