>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
>
A few questions:
1) Are you using the same arrays for fft_in and fft_out?
2) If yes to 1), did you fftw_malloc with sizeof(fftw_complex) *
(n+zeropadsize)?
3) Are you using the same plan for each iteration?
Mark
Reply by Rune Allnor●June 3, 20112011-06-03
On Jun 2, 7:13�pm, steve <bungalow_st...@yahoo.com> 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.
>
> > --
>
> > 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)
The *logic* is the same with DC as with any other Fourier coefficent.
The *psychology* is not, though:
- Novices and newbies don't always consider DC a 'frequency' since
it representa constant.
- Seasoned DSP practitioners often forget about DC since we usually
work under the assumption that the signal is zero mean.
- Often enough the DC component in plots coincide with the vertical
plot axis, so one does not see it quite as clearly as the other
Fourier coefficients.
Rune
Reply by Fred Marshall●June 2, 20112011-06-02
On 6/2/2011 1:15 PM, steve wrote:
>
> fair enough, I suppose I look at DC like any other bin, and design
> such that it can' t interfere with bins of interest to begin with
My inclination is to agree with Steve - that any bin is any bin.
Also, recall that zero padding doesn't add information and expecting it
to isn't such a good thing. Reasonably done, the original samples can
still be there, unchanged. So what we're talking about are the
interpolated samples. And the argument seems to be that "they shouldn't
be there or should be suppressed". When, in fact one is looking at a
finer resolution version of the results of windowing which might be
viewed as part of the signal under analysis. Better to use a narrow
bandpass filter if it's going to be an issue. Better windowing can be
nice but it still applies to every bin, DC or not.
I think the rationale is, as somewhat stated:
- the DC term *can* be big.
- it's easy to pick on.
- it's often of no interest.
- it's easy to get rid of.
So, why not?
That's not the same as saying it's necessary. Maybe advisable? I can
see where it would be advisable in some cases. And, as above, why
bother figuring it out a priori? Just do it - to Tim's point. That's
kinda like "necessary" .. almost.
What's "perturbed" by it would be the interpolation values. So, one has
to examine how those are going to be used.
Fred
Reply by steve●June 2, 20112011-06-02
On Jun 2, 3:52�pm, Tim Wescott <t...@seemywebsite.com> wrote:
> 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.
I suppose if that was the case I make it easy on myself and get it out
before it hit the A/D (with an simple 2 pole analog HPF), otherwise a
overwhelmingly large DC would waste all my dynamic range with my
signals of interest riding on top of it.
Reply by steve●June 2, 20112011-06-02
On Jun 2, 3:16�pm, glen herrmannsfeldt <g...@ugcs.caltech.edu> wrote:
> steve <bungalow_st...@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
fair enough, I suppose I look at DC like any other bin, and design
such that it can' t interfere with bins of interest to begin with
Reply by Tim Wescott●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
Reply by glen herrmannsfeldt●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 steve●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 Andrew Holme●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?
>
> Fred
The 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 Tim Wescott●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