Any FFT expert around?

classic Classic list List threaded Threaded
12 messages Options
Reply | Threaded
Open this post in threaded view
|

Any FFT expert around?

oxyopes
hey guys,

the (simple) code bellow is how i ve learned to do FFT according to
several docs online. However i observe a dependency of the signal
frequency in the spectrum on the constant FsFactor. In other words,
the signal frequency depends on the sampling rate. I should'nt be,
right? So what is wrong here?

  #---------- start code ------------
  clear all
  nu=10;                    % signal frequency
  at=5;                       % acquisition time
  T1=2                       % signal decay constant in time domain
  FsFactor=8;             % the ratio (sampling rate)/(signal
frequency), or Fs/nu

  clf
  Fs = FsFactor*nu;                             % Sampling rate
  t = ((0:(at*Fs-1))/Fs)';                         % Time vector
  freq=linspace(-1,1,Fs*at)' * Fs/2;         % frequency vector
  tdsig=exp(-i*2*pi*nu*t).*exp(-t./(T1));    % time domain signal
  freqsig=fft(tdsig);                                % freq. domain signal
  subplot(1,2,1)
  plot(t,tdsig)
  axis([ 0 0.4])       % zooming time domain to see that period=1/nu
  subplot(1,2,2)
  plot(freq, freqsig)
  #---------- end code ------------

Do it yourself. Just rerun the code trying different values of
FsFactor (eg: 2, 4, 8).
Thx a lot for any hint!!!

oxy

ps: cross post http://www.mathworks.com/matlabcentral/answers/115700-fft-why-signal-frequency-depends-on-sampling-rate
_______________________________________________
Help-octave mailing list
[hidden email]
https://mailman.cae.wisc.edu/listinfo/help-octave
Reply | Threaded
Open this post in threaded view
|

Re: Any FFT expert around?

briankaz
Hi Oxy,

The frequency vector for FFT should go from DC as bin 1 to slightly less than Fs in the last bin (Fs - Fs/N where N is the number of FFT points).

If you want to plot this as a range from -Fs/2 to +Fs/2 (minus the same fraction) you can do this, but in this case you should use the fftshift function as follows:

freqsig = fftshift(fft(tdsig));

Also, make sure you are intending to plot the REAL part of the FFT only (I believe that is the default behavior for how you are executing the plot).

-Brian


2014-02-12 14:22 GMT+01:00 oxy <[hidden email]>:
hey guys,

the (simple) code bellow is how i ve learned to do FFT according to
several docs online. However i observe a dependency of the signal
frequency in the spectrum on the constant FsFactor. In other words,
the signal frequency depends on the sampling rate. I should'nt be,
right? So what is wrong here?

  #---------- start code ------------
  clear all
  nu=10;                    % signal frequency
  at=5;                       % acquisition time
  T1=2                       % signal decay constant in time domain
  FsFactor=8;             % the ratio (sampling rate)/(signal
frequency), or Fs/nu

  clf
  Fs = FsFactor*nu;                             % Sampling rate
  t = ((0:(at*Fs-1))/Fs)';                         % Time vector
  freq=linspace(-1,1,Fs*at)' * Fs/2;         % frequency vector
  tdsig=exp(-i*2*pi*nu*t).*exp(-t./(T1));    % time domain signal
  freqsig=fft(tdsig);                                % freq. domain signal
  subplot(1,2,1)
  plot(t,tdsig)
  axis([ 0 0.4])       % zooming time domain to see that period=1/nu
  subplot(1,2,2)
  plot(freq, freqsig)
  #---------- end code ------------

Do it yourself. Just rerun the code trying different values of
FsFactor (eg: 2, 4, 8).
Thx a lot for any hint!!!

oxy

ps: cross post http://www.mathworks.com/matlabcentral/answers/115700-fft-why-signal-frequency-depends-on-sampling-rate
_______________________________________________
Help-octave mailing list
[hidden email]
https://mailman.cae.wisc.edu/listinfo/help-octave


_______________________________________________
Help-octave mailing list
[hidden email]
https://mailman.cae.wisc.edu/listinfo/help-octave
Reply | Threaded
Open this post in threaded view
|

Re: Any FFT expert around?

Peter L. Soendergaard
In reply to this post by oxyopes
> Hi Oxy,
>
> The frequency vector for FFT should go from DC as bin 1 to slightly less
> than Fs in the last bin (Fs - Fs/N where N is the number of FFT points).
>
> If you want to plot this as a range from -Fs/2 to +Fs/2 (minus the same
> fraction) you can do this, but in this case you should use the fftshift
> function as follows:
>
> freqsig = fftshift(fft(tdsig));
>
> Also, make sure you are intending to plot the REAL part of the FFT only (I
> believe that is the default behavior for how you are executing the plot).
>

More blatant promotion: LTFAT has the "plotfft" and "plotfftreal"
functions for visualing coefficients coming from fft operations:

http://ltfat.sourceforge.net/doc/fourier/plotfft.php
http://ltfat.sourceforge.net/doc/fourier/plotfftreal.php

Cheers,
Peter.


_______________________________________________
Help-octave mailing list
[hidden email]
https://mailman.cae.wisc.edu/listinfo/help-octave
Reply | Threaded
Open this post in threaded view
|

Re: Any FFT expert around?

oxyopes
In reply to this post by briankaz
Hi Brian,

> The frequency vector for FFT should go from DC as bin 1 to slightly less
> than Fs in the last bin (Fs - Fs/N where N is the number of FFT points).

If i understand u correctly, the frequency vector (freq) must be written
as in this modified version of code below. However, according to Nyquist we
cannot measure a frequency higher than Fs/2. Thus i do not see
the meaning of plotting up to ~Fs.

Also, if I rerun the code below changing FsFactor, I again see this
dependency of the signal frequency on FsFactor.
I cannot understand it. Looks like basics, yet not that obvious.

  #---------- start code ------------
  clear all
  nu=10;                    % signal frequency
  at=5;                       % acquisition time
  T1=2                       % signal decay constant in time domain
  FsFactor=16;             % the ratio (sampling rate)/(signal
frequency), or Fs/nu

  clf
  Fs = FsFactor*nu;                             % Sampling rate
  t = ((0:(at*Fs-1))/Fs)';                         % Time vector
  % freq=linspace(-1,1,Fs*at)' * Fs/2;         % frequency vector, first version
  freq=linspace(0,1,Fs*at)' * Fs- Fs/(Fs*at);  % frequency vector
suggested by Brian
  tdsig=exp(-i*2*pi*nu*t).*exp(-t./(T1));    % time domain signal
  freqsig=fft(tdsig);                                % freq. domain signal
  subplot(1,2,1)
  plot(t,tdsig)
  axis([ 0 0.4])       % zooming time domain to see that period=1/nu
  subplot(1,2,2)
  plot(freq, freqsig)
  #---------- end code ------------

thx guys ...


> 2014-02-12 14:22 GMT+01:00 oxy <[hidden email]>:
>
>> hey guys,
>>
>> the (simple) code bellow is how i ve learned to do FFT according to
>> several docs online. However i observe a dependency of the signal
>> frequency in the spectrum on the constant FsFactor. In other words,
>> the signal frequency depends on the sampling rate. I should'nt be,
>> right? So what is wrong here?
>>
>>   #---------- start code ------------
>>   clear all
>>   nu=10;                    % signal frequency
>>   at=5;                       % acquisition time
>>   T1=2                       % signal decay constant in time domain
>>   FsFactor=8;             % the ratio (sampling rate)/(signal
>> frequency), or Fs/nu
>>
>>   clf
>>   Fs = FsFactor*nu;                             % Sampling rate
>>   t = ((0:(at*Fs-1))/Fs)';                         % Time vector
>>   freq=linspace(-1,1,Fs*at)' * Fs/2;         % frequency vector
>>   tdsig=exp(-i*2*pi*nu*t).*exp(-t./(T1));    % time domain signal
>>   freqsig=fft(tdsig);                                % freq. domain
>> signal
>>   subplot(1,2,1)
>>   plot(t,tdsig)
>>   axis([ 0 0.4])       % zooming time domain to see that period=1/nu
>>   subplot(1,2,2)
>>   plot(freq, freqsig)
>>   #---------- end code ------------
>>
>> Do it yourself. Just rerun the code trying different values of
>> FsFactor (eg: 2, 4, 8).
>> Thx a lot for any hint!!!
>>
>> oxy
>>
>> ps: cross post
>> http://www.mathworks.com/matlabcentral/answers/115700-fft-why-signal-frequency-depends-on-sampling-rate
>> _______________________________________________
>> Help-octave mailing list
>> [hidden email]
>> https://mailman.cae.wisc.edu/listinfo/help-octave
>>
>
_______________________________________________
Help-octave mailing list
[hidden email]
https://mailman.cae.wisc.edu/listinfo/help-octave
Reply | Threaded
Open this post in threaded view
|

Re: Any FFT expert around?

briankaz
Hi Oxy,

Since you're taking the FFT of a complex vector it will have unique terms for negative vs. positive frequency.  You are correct that it makes more sense to think of the FFT from -Fs/2 to +Fs/2.  It's just that the Octave fft function doesn't compute the bins in that order.  Due to aliasing, any frequency Fs/2 + df is equivalent to (-Fs/2 + df) so it's more a matter of choice how you want to display the data.

If you prefer plotting from -Fs/2 to +Fs/2 I would change these two lines of your code as follows:

freq=linspace(-Fs/2,Fs/2-Fs/(Fs*at),Fs*at)'
freqsig=fftshift(fft(tdsig);

Let us know if that gives you what you want!

-Brian


2014-02-12 17:15 GMT+01:00 oxy <[hidden email]>:
Hi Brian,

> The frequency vector for FFT should go from DC as bin 1 to slightly less
> than Fs in the last bin (Fs - Fs/N where N is the number of FFT points).

If i understand u correctly, the frequency vector (freq) must be written
as in this modified version of code below. However, according to Nyquist we
cannot measure a frequency higher than Fs/2. Thus i do not see
the meaning of plotting up to ~Fs.

Also, if I rerun the code below changing FsFactor, I again see this
dependency of the signal frequency on FsFactor.
I cannot understand it. Looks like basics, yet not that obvious.

  #---------- start code ------------
  clear all
  nu=10;                    % signal frequency
  at=5;                       % acquisition time
  T1=2                       % signal decay constant in time domain
  FsFactor=16;             % the ratio (sampling rate)/(signal
frequency), or Fs/nu

  clf
  Fs = FsFactor*nu;                             % Sampling rate
  t = ((0:(at*Fs-1))/Fs)';                         % Time vector
  % freq=linspace(-1,1,Fs*at)' * Fs/2;         % frequency vector, first version
  freq=linspace(0,1,Fs*at)' * Fs- Fs/(Fs*at);  % frequency vector
suggested by Brian
  tdsig=exp(-i*2*pi*nu*t).*exp(-t./(T1));    % time domain signal
  freqsig=fft(tdsig);                                % freq. domain signal
  subplot(1,2,1)
  plot(t,tdsig)
  axis([ 0 0.4])       % zooming time domain to see that period=1/nu
  subplot(1,2,2)
  plot(freq, freqsig)
  #---------- end code ------------

thx guys ...


> 2014-02-12 14:22 GMT+01:00 oxy <[hidden email]>:
>
>> hey guys,
>>
>> the (simple) code bellow is how i ve learned to do FFT according to
>> several docs online. However i observe a dependency of the signal
>> frequency in the spectrum on the constant FsFactor. In other words,
>> the signal frequency depends on the sampling rate. I should'nt be,
>> right? So what is wrong here?
>>
>>   #---------- start code ------------
>>   clear all
>>   nu=10;                    % signal frequency
>>   at=5;                       % acquisition time
>>   T1=2                       % signal decay constant in time domain
>>   FsFactor=8;             % the ratio (sampling rate)/(signal
>> frequency), or Fs/nu
>>
>>   clf
>>   Fs = FsFactor*nu;                             % Sampling rate
>>   t = ((0:(at*Fs-1))/Fs)';                         % Time vector
>>   freq=linspace(-1,1,Fs*at)' * Fs/2;         % frequency vector
>>   tdsig=exp(-i*2*pi*nu*t).*exp(-t./(T1));    % time domain signal
>>   freqsig=fft(tdsig);                                % freq. domain
>> signal
>>   subplot(1,2,1)
>>   plot(t,tdsig)
>>   axis([ 0 0.4])       % zooming time domain to see that period=1/nu
>>   subplot(1,2,2)
>>   plot(freq, freqsig)
>>   #---------- end code ------------
>>
>> Do it yourself. Just rerun the code trying different values of
>> FsFactor (eg: 2, 4, 8).
>> Thx a lot for any hint!!!
>>
>> oxy
>>
>> ps: cross post
>> http://www.mathworks.com/matlabcentral/answers/115700-fft-why-signal-frequency-depends-on-sampling-rate
>> _______________________________________________
>> Help-octave mailing list
>> [hidden email]
>> https://mailman.cae.wisc.edu/listinfo/help-octave
>>
>
_______________________________________________
Help-octave mailing list
[hidden email]
https://mailman.cae.wisc.edu/listinfo/help-octave


_______________________________________________
Help-octave mailing list
[hidden email]
https://mailman.cae.wisc.edu/listinfo/help-octave
Reply | Threaded
Open this post in threaded view
|

Re: Any FFT expert around?

oxyopes
On 2/12/14, Brian Kaczynski <[hidden email]> wrote:

> Hi Oxy,
>
> Since you're taking the FFT of a complex vector it will have unique terms
> for negative vs. positive frequency.  You are correct that it makes more
> sense to think of the FFT from -Fs/2 to +Fs/2.  It's just that the Octave
> fft function doesn't compute the bins in that order.  Due to aliasing, any
> frequency Fs/2 + df is equivalent to (-Fs/2 + df) so it's more a matter of
> choice how you want to display the data.
>
> If you prefer plotting from -Fs/2 to +Fs/2 I would change these two lines
> of your code as follows:
>
> freq=linspace(-Fs/2,Fs/2-Fs/(Fs*at),Fs*at)'
> freqsig=fftshift(fft(tdsig);
>
> Let us know if that gives you what you want!
>
> -Brian
>
>
> 2014-02-12 17:15 GMT+01:00 oxy <[hidden email]>:
>
>> Hi Brian,
>>
>> > The frequency vector for FFT should go from DC as bin 1 to slightly
>> > less
>> > than Fs in the last bin (Fs - Fs/N where N is the number of FFT
>> > points).
>>
>> If i understand u correctly, the frequency vector (freq) must be written
>> as in this modified version of code below. However, according to Nyquist
>> we
>> cannot measure a frequency higher than Fs/2. Thus i do not see
>> the meaning of plotting up to ~Fs.
>>
>> Also, if I rerun the code below changing FsFactor, I again see this
>> dependency of the signal frequency on FsFactor.
>> I cannot understand it. Looks like basics, yet not that obvious.
>>
>>   #---------- start code ------------
>>   clear all
>>   nu=10;                    % signal frequency
>>   at=5;                       % acquisition time
>>   T1=2                       % signal decay constant in time domain
>>   FsFactor=16;             % the ratio (sampling rate)/(signal
>> frequency), or Fs/nu
>>
>>   clf
>>   Fs = FsFactor*nu;                             % Sampling rate
>>   t = ((0:(at*Fs-1))/Fs)';                         % Time vector
>>   % freq=linspace(-1,1,Fs*at)' * Fs/2;         % frequency vector, first
>> version
>>   freq=linspace(0,1,Fs*at)' * Fs- Fs/(Fs*at);  % frequency vector
>> suggested by Brian
>>   tdsig=exp(-i*2*pi*nu*t).*exp(-t./(T1));    % time domain signal
>>   freqsig=fft(tdsig);                                % freq. domain
>> signal
>>   subplot(1,2,1)
>>   plot(t,tdsig)
>>   axis([ 0 0.4])       % zooming time domain to see that period=1/nu
>>   subplot(1,2,2)
>>   plot(freq, freqsig)
>>   #---------- end code ------------
>>
>> thx guys ...
>>
>>
>> > 2014-02-12 14:22 GMT+01:00 oxy <[hidden email]>:
>> >
>> >> hey guys,
>> >>
>> >> the (simple) code bellow is how i ve learned to do FFT according to
>> >> several docs online. However i observe a dependency of the signal
>> >> frequency in the spectrum on the constant FsFactor. In other words,
>> >> the signal frequency depends on the sampling rate. I should'nt be,
>> >> right? So what is wrong here?
>> >>
>> >>   #---------- start code ------------
>> >>   clear all
>> >>   nu=10;                    % signal frequency
>> >>   at=5;                       % acquisition time
>> >>   T1=2                       % signal decay constant in time domain
>> >>   FsFactor=8;             % the ratio (sampling rate)/(signal
>> >> frequency), or Fs/nu
>> >>
>> >>   clf
>> >>   Fs = FsFactor*nu;                             % Sampling rate
>> >>   t = ((0:(at*Fs-1))/Fs)';                         % Time vector
>> >>   freq=linspace(-1,1,Fs*at)' * Fs/2;         % frequency vector
>> >>   tdsig=exp(-i*2*pi*nu*t).*exp(-t./(T1));    % time domain signal
>> >>   freqsig=fft(tdsig);                                % freq. domain
>> >> signal
>> >>   subplot(1,2,1)
>> >>   plot(t,tdsig)
>> >>   axis([ 0 0.4])       % zooming time domain to see that period=1/nu
>> >>   subplot(1,2,2)
>> >>   plot(freq, freqsig)
>> >>   #---------- end code ------------
>> >>
>> >> Do it yourself. Just rerun the code trying different values of
>> >> FsFactor (eg: 2, 4, 8).
>> >> Thx a lot for any hint!!!
>> >>
>> >> oxy
>> >>
>> >> ps: cross post
>> >>
>> http://www.mathworks.com/matlabcentral/answers/115700-fft-why-signal-frequency-depends-on-sampling-rate
>> >> _______________________________________________
>> >> Help-octave mailing list
>> >> [hidden email]
>> >> https://mailman.cae.wisc.edu/listinfo/help-octave
>> >>
>> >
>> _______________________________________________
>> Help-octave mailing list
>> [hidden email]
>> https://mailman.cae.wisc.edu/listinfo/help-octave
>>
>
_______________________________________________
Help-octave mailing list
[hidden email]
https://mailman.cae.wisc.edu/listinfo/help-octave
Reply | Threaded
Open this post in threaded view
|

Re: Any FFT expert around?

oxyopes
hi brian,

now it works ... great!!!!!!!

Look at the code below ajusted according to your hints.

Just one question more: you are calculating the frequency
vector as:

         freq1=-linspace(-Fs1/2,   Fs1/2-1/at,   Fs1*at)';

that makes it non-symetric. Isn't it better ...

         freq1=-linspace(-Fs1/2+1/at,   Fs1/2-1/at,   Fs1*at)';

...? Then it's symetric. Hey, thanks a lot! Great help!!!!!

oxy

%======== DEMONSTRATION CODE =========
% TTF Demo: this code shows how to calculate the
% frequency domain when doing FFT of td-signals.

%--- parameter setup: same signal, 2 sampling rates ---
clear all
nu=10;            % signal frequency
at=5;             % acquisition time
T1=2              % signal decay constant in time domain
FsFactor1=4;      % sampling rate factor 1, Fs1/nu
FsFactor2=16;     % sampling rate factor 2, Fs2/nu

%----------- calculating -------------------
Fs1 = FsFactor1*nu;            % Sampling rate 1
Fs2 = FsFactor2*nu;            % Sampling rate 2

t1 = ((0:(at*Fs1-1))/Fs1)';    % Time vector 1
t2 = ((0:(at*Fs2-1))/Fs2)';    % Time vector 2

tdsig1=exp(-i*2*pi*nu*t1).*exp(-t1./(T1)); % time domain signal 1
tdsig2=exp(-i*2*pi*nu*t2).*exp(-t2./(T1)); % time domain signal 2

% I realized, the signal was appearing negative, thus minus...
freq1=-linspace(-Fs1/2,Fs1/2-Fs1/(Fs1*at),Fs1*at)'; % freq. vector 1
freq2=-linspace(-Fs2/2,Fs2/2-Fs2/(Fs2*at),Fs2*at)'; % freq. vector 2

freqsig1=fftshift(fft(tdsig1));  % signal vector 1
freqsig2=fftshift(fft(tdsig2));  % signal vector 2

%----------- end ------------------

=================================================
On 2/19/14, oxy <[hidden email]> wrote:

> On 2/12/14, Brian Kaczynski <[hidden email]> wrote:
>> Hi Oxy,
>>
>> Since you're taking the FFT of a complex vector it will have unique terms
>> for negative vs. positive frequency.  You are correct that it makes more
>> sense to think of the FFT from -Fs/2 to +Fs/2.  It's just that the Octave
>> fft function doesn't compute the bins in that order.  Due to aliasing,
>> any
>> frequency Fs/2 + df is equivalent to (-Fs/2 + df) so it's more a matter
>> of
>> choice how you want to display the data.
>>
>> If you prefer plotting from -Fs/2 to +Fs/2 I would change these two lines
>> of your code as follows:
>>
>> freq=linspace(-Fs/2,Fs/2-Fs/(Fs*at),Fs*at)'
>> freqsig=fftshift(fft(tdsig);
>>
>> Let us know if that gives you what you want!
>>
>> -Brian
>>
>>
>> 2014-02-12 17:15 GMT+01:00 oxy <[hidden email]>:
>>
>>> Hi Brian,
>>>
>>> > The frequency vector for FFT should go from DC as bin 1 to slightly
>>> > less
>>> > than Fs in the last bin (Fs - Fs/N where N is the number of FFT
>>> > points).
>>>
>>> If i understand u correctly, the frequency vector (freq) must be written
>>> as in this modified version of code below. However, according to Nyquist
>>> we
>>> cannot measure a frequency higher than Fs/2. Thus i do not see
>>> the meaning of plotting up to ~Fs.
>>>
>>> Also, if I rerun the code below changing FsFactor, I again see this
>>> dependency of the signal frequency on FsFactor.
>>> I cannot understand it. Looks like basics, yet not that obvious.
>>>
>>>   #---------- start code ------------
>>>   clear all
>>>   nu=10;                    % signal frequency
>>>   at=5;                       % acquisition time
>>>   T1=2                       % signal decay constant in time domain
>>>   FsFactor=16;             % the ratio (sampling rate)/(signal
>>> frequency), or Fs/nu
>>>
>>>   clf
>>>   Fs = FsFactor*nu;                             % Sampling rate
>>>   t = ((0:(at*Fs-1))/Fs)';                         % Time vector
>>>   % freq=linspace(-1,1,Fs*at)' * Fs/2;         % frequency vector, first
>>> version
>>>   freq=linspace(0,1,Fs*at)' * Fs- Fs/(Fs*at);  % frequency vector
>>> suggested by Brian
>>>   tdsig=exp(-i*2*pi*nu*t).*exp(-t./(T1));    % time domain signal
>>>   freqsig=fft(tdsig);                                % freq. domain
>>> signal
>>>   subplot(1,2,1)
>>>   plot(t,tdsig)
>>>   axis([ 0 0.4])       % zooming time domain to see that period=1/nu
>>>   subplot(1,2,2)
>>>   plot(freq, freqsig)
>>>   #---------- end code ------------
>>>
>>> thx guys ...
>>>
>>>
>>> > 2014-02-12 14:22 GMT+01:00 oxy <[hidden email]>:
>>> >
>>> >> hey guys,
>>> >>
>>> >> the (simple) code bellow is how i ve learned to do FFT according to
>>> >> several docs online. However i observe a dependency of the signal
>>> >> frequency in the spectrum on the constant FsFactor. In other words,
>>> >> the signal frequency depends on the sampling rate. I should'nt be,
>>> >> right? So what is wrong here?
>>> >>
>>> >>   #---------- start code ------------
>>> >>   clear all
>>> >>   nu=10;                    % signal frequency
>>> >>   at=5;                       % acquisition time
>>> >>   T1=2                       % signal decay constant in time domain
>>> >>   FsFactor=8;             % the ratio (sampling rate)/(signal
>>> >> frequency), or Fs/nu
>>> >>
>>> >>   clf
>>> >>   Fs = FsFactor*nu;                             % Sampling rate
>>> >>   t = ((0:(at*Fs-1))/Fs)';                         % Time vector
>>> >>   freq=linspace(-1,1,Fs*at)' * Fs/2;         % frequency vector
>>> >>   tdsig=exp(-i*2*pi*nu*t).*exp(-t./(T1));    % time domain signal
>>> >>   freqsig=fft(tdsig);                                % freq. domain
>>> >> signal
>>> >>   subplot(1,2,1)
>>> >>   plot(t,tdsig)
>>> >>   axis([ 0 0.4])       % zooming time domain to see that period=1/nu
>>> >>   subplot(1,2,2)
>>> >>   plot(freq, freqsig)
>>> >>   #---------- end code ------------
>>> >>
>>> >> Do it yourself. Just rerun the code trying different values of
>>> >> FsFactor (eg: 2, 4, 8).
>>> >> Thx a lot for any hint!!!
>>> >>
>>> >> oxy
>>> >>
>>> >> ps: cross post
>>> >>
>>> http://www.mathworks.com/matlabcentral/answers/115700-fft-why-signal-frequency-depends-on-sampling-rate
>>> >> _______________________________________________
>>> >> Help-octave mailing list
>>> >> [hidden email]
>>> >> https://mailman.cae.wisc.edu/listinfo/help-octave
>>> >>
>>> >
>>> _______________________________________________
>>> Help-octave mailing list
>>> [hidden email]
>>> https://mailman.cae.wisc.edu/listinfo/help-octave
>>>
>>
>
_______________________________________________
Help-octave mailing list
[hidden email]
https://mailman.cae.wisc.edu/listinfo/help-octave
Reply | Threaded
Open this post in threaded view
|

Re: Any FFT expert around?

briankaz
2014-02-19 10:49 GMT+01:00, oxy <[hidden email]>:

> hi brian,
>
> now it works ... great!!!!!!!
>
> Look at the code below ajusted according to your hints.
>
> Just one question more: you are calculating the frequency
> vector as:
>
>          freq1=-linspace(-Fs1/2,   Fs1/2-1/at,   Fs1*at)';
>
> that makes it non-symetric. Isn't it better ...
>
>          freq1=-linspace(-Fs1/2+1/at,   Fs1/2-1/at,   Fs1*at)';
>
> ...? Then it's symetric. Hey, thanks a lot! Great help!!!!!

That would be "better" if you really want to insist on having a
symmetrical frequency vector at the expense of accuracy.
Unfortunately, it's inaccurate.  For one thing, your symmetrical
frequency vector does not contain a DC term if you have an even number
of FFT points (which you should based on your examples).

If you want the FFT value at the nth point to really correspond to
freq1(n) you should use the slightly offset frequency vector that I
provided.

If it is more important to whatever you're doing that the frequency
vector be exactly symmetrical and you can live with sacrificing the
true 1-to-1 correspondence between frequency points and FFT points,
you can use your symmetrical frequency vector.

> oxy
>
> %======== DEMONSTRATION CODE =========
> % TTF Demo: this code shows how to calculate the
> % frequency domain when doing FFT of td-signals.
>
> %--- parameter setup: same signal, 2 sampling rates ---
> clear all
> nu=10;            % signal frequency
> at=5;             % acquisition time
> T1=2              % signal decay constant in time domain
> FsFactor1=4;      % sampling rate factor 1, Fs1/nu
> FsFactor2=16;     % sampling rate factor 2, Fs2/nu
>
> %----------- calculating -------------------
> Fs1 = FsFactor1*nu;            % Sampling rate 1
> Fs2 = FsFactor2*nu;            % Sampling rate 2
>
> t1 = ((0:(at*Fs1-1))/Fs1)';    % Time vector 1
> t2 = ((0:(at*Fs2-1))/Fs2)';    % Time vector 2
>
> tdsig1=exp(-i*2*pi*nu*t1).*exp(-t1./(T1)); % time domain signal 1
> tdsig2=exp(-i*2*pi*nu*t2).*exp(-t2./(T1)); % time domain signal 2
>
> % I realized, the signal was appearing negative, thus minus...
> freq1=-linspace(-Fs1/2,Fs1/2-Fs1/(Fs1*at),Fs1*at)'; % freq. vector 1
> freq2=-linspace(-Fs2/2,Fs2/2-Fs2/(Fs2*at),Fs2*at)'; % freq. vector 2
>
> freqsig1=fftshift(fft(tdsig1));  % signal vector 1
> freqsig2=fftshift(fft(tdsig2));  % signal vector 2
>
> %----------- end ------------------
>
> =================================================
> On 2/19/14, oxy <[hidden email]> wrote:
>> On 2/12/14, Brian Kaczynski <[hidden email]> wrote:
>>> Hi Oxy,
>>>
>>> Since you're taking the FFT of a complex vector it will have unique
>>> terms
>>> for negative vs. positive frequency.  You are correct that it makes more
>>> sense to think of the FFT from -Fs/2 to +Fs/2.  It's just that the
>>> Octave
>>> fft function doesn't compute the bins in that order.  Due to aliasing,
>>> any
>>> frequency Fs/2 + df is equivalent to (-Fs/2 + df) so it's more a matter
>>> of
>>> choice how you want to display the data.
>>>
>>> If you prefer plotting from -Fs/2 to +Fs/2 I would change these two
>>> lines
>>> of your code as follows:
>>>
>>> freq=linspace(-Fs/2,Fs/2-Fs/(Fs*at),Fs*at)'
>>> freqsig=fftshift(fft(tdsig);
>>>
>>> Let us know if that gives you what you want!
>>>
>>> -Brian
>>>
>>>
>>> 2014-02-12 17:15 GMT+01:00 oxy <[hidden email]>:
>>>
>>>> Hi Brian,
>>>>
>>>> > The frequency vector for FFT should go from DC as bin 1 to slightly
>>>> > less
>>>> > than Fs in the last bin (Fs - Fs/N where N is the number of FFT
>>>> > points).
>>>>
>>>> If i understand u correctly, the frequency vector (freq) must be
>>>> written
>>>> as in this modified version of code below. However, according to
>>>> Nyquist
>>>> we
>>>> cannot measure a frequency higher than Fs/2. Thus i do not see
>>>> the meaning of plotting up to ~Fs.
>>>>
>>>> Also, if I rerun the code below changing FsFactor, I again see this
>>>> dependency of the signal frequency on FsFactor.
>>>> I cannot understand it. Looks like basics, yet not that obvious.
>>>>
>>>>   #---------- start code ------------
>>>>   clear all
>>>>   nu=10;                    % signal frequency
>>>>   at=5;                       % acquisition time
>>>>   T1=2                       % signal decay constant in time domain
>>>>   FsFactor=16;             % the ratio (sampling rate)/(signal
>>>> frequency), or Fs/nu
>>>>
>>>>   clf
>>>>   Fs = FsFactor*nu;                             % Sampling rate
>>>>   t = ((0:(at*Fs-1))/Fs)';                         % Time vector
>>>>   % freq=linspace(-1,1,Fs*at)' * Fs/2;         % frequency vector,
>>>> first
>>>> version
>>>>   freq=linspace(0,1,Fs*at)' * Fs- Fs/(Fs*at);  % frequency vector
>>>> suggested by Brian
>>>>   tdsig=exp(-i*2*pi*nu*t).*exp(-t./(T1));    % time domain signal
>>>>   freqsig=fft(tdsig);                                % freq. domain
>>>> signal
>>>>   subplot(1,2,1)
>>>>   plot(t,tdsig)
>>>>   axis([ 0 0.4])       % zooming time domain to see that period=1/nu
>>>>   subplot(1,2,2)
>>>>   plot(freq, freqsig)
>>>>   #---------- end code ------------
>>>>
>>>> thx guys ...
>>>>
>>>>
>>>> > 2014-02-12 14:22 GMT+01:00 oxy <[hidden email]>:
>>>> >
>>>> >> hey guys,
>>>> >>
>>>> >> the (simple) code bellow is how i ve learned to do FFT according to
>>>> >> several docs online. However i observe a dependency of the signal
>>>> >> frequency in the spectrum on the constant FsFactor. In other words,
>>>> >> the signal frequency depends on the sampling rate. I should'nt be,
>>>> >> right? So what is wrong here?
>>>> >>
>>>> >>   #---------- start code ------------
>>>> >>   clear all
>>>> >>   nu=10;                    % signal frequency
>>>> >>   at=5;                       % acquisition time
>>>> >>   T1=2                       % signal decay constant in time domain
>>>> >>   FsFactor=8;             % the ratio (sampling rate)/(signal
>>>> >> frequency), or Fs/nu
>>>> >>
>>>> >>   clf
>>>> >>   Fs = FsFactor*nu;                             % Sampling rate
>>>> >>   t = ((0:(at*Fs-1))/Fs)';                         % Time vector
>>>> >>   freq=linspace(-1,1,Fs*at)' * Fs/2;         % frequency vector
>>>> >>   tdsig=exp(-i*2*pi*nu*t).*exp(-t./(T1));    % time domain signal
>>>> >>   freqsig=fft(tdsig);                                % freq. domain
>>>> >> signal
>>>> >>   subplot(1,2,1)
>>>> >>   plot(t,tdsig)
>>>> >>   axis([ 0 0.4])       % zooming time domain to see that period=1/nu
>>>> >>   subplot(1,2,2)
>>>> >>   plot(freq, freqsig)
>>>> >>   #---------- end code ------------
>>>> >>
>>>> >> Do it yourself. Just rerun the code trying different values of
>>>> >> FsFactor (eg: 2, 4, 8).
>>>> >> Thx a lot for any hint!!!
>>>> >>
>>>> >> oxy
>>>> >>
>>>> >> ps: cross post
>>>> >>
>>>> http://www.mathworks.com/matlabcentral/answers/115700-fft-why-signal-frequency-depends-on-sampling-rate
>>>> >> _______________________________________________
>>>> >> Help-octave mailing list
>>>> >> [hidden email]
>>>> >> https://mailman.cae.wisc.edu/listinfo/help-octave
>>>> >>
>>>> >
>>>> _______________________________________________
>>>> Help-octave mailing list
>>>> [hidden email]
>>>> https://mailman.cae.wisc.edu/listinfo/help-octave
>>>>
>>>
>>
>
_______________________________________________
Help-octave mailing list
[hidden email]
https://mailman.cae.wisc.edu/listinfo/help-octave
Reply | Threaded
Open this post in threaded view
|

Re: Any FFT expert around?

oxyopes
Hi Brian,

you won again  :-)
Yeah, i prefer an even number of points on the freq vector.
Thus following your suggestion. Below is my final demo code.

Today is maybe time for a new small donation to octave proj.

Thank you very much!
oxy


% TTF Demo: this code shows how to calculate the
% frequency domain when doing FFT of td-signals.

%--- parameter setup: same signal, 2 sampling rates ---
clear all
nu=10;            % signal frequency
at=5;             % acquisition time
T1=2              % signal decay constant in time domain
FsFactor1=4;      % sampling rate factor 1, Fs1/nu
FsFactor2=16;     % sampling rate factor 2, Fs2/nu

%----------- calculating -------------------
Fs1 = FsFactor1*nu;            % Sampling rate 1
Fs2 = FsFactor2*nu;            % Sampling rate 2

t1 = ((0:(at*Fs1-1))/Fs1)';    % Time vector 1
t2 = ((0:(at*Fs2-1))/Fs2)';    % Time vector 2

tdsig1=exp(-i*2*pi*nu*t1).*exp(-t1./(T1)); % time domain signal 1
tdsig2=exp(-i*2*pi*nu*t2).*exp(-t2./(T1)); % time domain signal 2

% I realized, the signal was appearing negative, thus minus...
freq1=linspace(-Fs1/2,Fs1/2-Fs1/(Fs1*at),Fs1*at)'; % freq. vector 1
freq2=linspace(-Fs2/2,Fs2/2-Fs2/(Fs2*at),Fs2*at)'; % freq. vector 2

freqsig1=fftshift(fft(tdsig1));  % signal vector 1
freqsig2=fftshift(fft(tdsig2));  % signal vector 2


%----------- ploting -------------------
clf
  subplot(1,2,1)
  plot(freq1, freqsig1)
  axis([-20 0]), title('FsFactor=4')

  subplot(1,2,2)
  plot(freq2, freqsig2)
  axis([-20  0]), title('FsFactor=16')


On 2/19/14, Brian Kaczynski <[hidden email]> wrote:

> 2014-02-19 10:49 GMT+01:00, oxy <[hidden email]>:
>> hi brian,
>>
>> now it works ... great!!!!!!!
>>
>> Look at the code below ajusted according to your hints.
>>
>> Just one question more: you are calculating the frequency
>> vector as:
>>
>>          freq1=-linspace(-Fs1/2,   Fs1/2-1/at,   Fs1*at)';
>>
>> that makes it non-symetric. Isn't it better ...
>>
>>          freq1=-linspace(-Fs1/2+1/at,   Fs1/2-1/at,   Fs1*at)';
>>
>> ...? Then it's symetric. Hey, thanks a lot! Great help!!!!!
>
> That would be "better" if you really want to insist on having a
> symmetrical frequency vector at the expense of accuracy.
> Unfortunately, it's inaccurate.  For one thing, your symmetrical
> frequency vector does not contain a DC term if you have an even number
> of FFT points (which you should based on your examples).
>
> If you want the FFT value at the nth point to really correspond to
> freq1(n) you should use the slightly offset frequency vector that I
> provided.
>
> If it is more important to whatever you're doing that the frequency
> vector be exactly symmetrical and you can live with sacrificing the
> true 1-to-1 correspondence between frequency points and FFT points,
> you can use your symmetrical frequency vector.
>
>> oxy
>>
>> %======== DEMONSTRATION CODE =========
>> % TTF Demo: this code shows how to calculate the
>> % frequency domain when doing FFT of td-signals.
>>
>> %--- parameter setup: same signal, 2 sampling rates ---
>> clear all
>> nu=10;            % signal frequency
>> at=5;             % acquisition time
>> T1=2              % signal decay constant in time domain
>> FsFactor1=4;      % sampling rate factor 1, Fs1/nu
>> FsFactor2=16;     % sampling rate factor 2, Fs2/nu
>>
>> %----------- calculating -------------------
>> Fs1 = FsFactor1*nu;            % Sampling rate 1
>> Fs2 = FsFactor2*nu;            % Sampling rate 2
>>
>> t1 = ((0:(at*Fs1-1))/Fs1)';    % Time vector 1
>> t2 = ((0:(at*Fs2-1))/Fs2)';    % Time vector 2
>>
>> tdsig1=exp(-i*2*pi*nu*t1).*exp(-t1./(T1)); % time domain signal 1
>> tdsig2=exp(-i*2*pi*nu*t2).*exp(-t2./(T1)); % time domain signal 2
>>
>> % I realized, the signal was appearing negative, thus minus...
>> freq1=-linspace(-Fs1/2,Fs1/2-Fs1/(Fs1*at),Fs1*at)'; % freq. vector 1
>> freq2=-linspace(-Fs2/2,Fs2/2-Fs2/(Fs2*at),Fs2*at)'; % freq. vector 2
>>
>> freqsig1=fftshift(fft(tdsig1));  % signal vector 1
>> freqsig2=fftshift(fft(tdsig2));  % signal vector 2
>>
>> %----------- end ------------------
>>
>> =================================================
>> On 2/19/14, oxy <[hidden email]> wrote:
>>> On 2/12/14, Brian Kaczynski <[hidden email]> wrote:
>>>> Hi Oxy,
>>>>
>>>> Since you're taking the FFT of a complex vector it will have unique
>>>> terms
>>>> for negative vs. positive frequency.  You are correct that it makes
>>>> more
>>>> sense to think of the FFT from -Fs/2 to +Fs/2.  It's just that the
>>>> Octave
>>>> fft function doesn't compute the bins in that order.  Due to aliasing,
>>>> any
>>>> frequency Fs/2 + df is equivalent to (-Fs/2 + df) so it's more a matter
>>>> of
>>>> choice how you want to display the data.
>>>>
>>>> If you prefer plotting from -Fs/2 to +Fs/2 I would change these two
>>>> lines
>>>> of your code as follows:
>>>>
>>>> freq=linspace(-Fs/2,Fs/2-Fs/(Fs*at),Fs*at)'
>>>> freqsig=fftshift(fft(tdsig);
>>>>
>>>> Let us know if that gives you what you want!
>>>>
>>>> -Brian
>>>>
>>>>
>>>> 2014-02-12 17:15 GMT+01:00 oxy <[hidden email]>:
>>>>
>>>>> Hi Brian,
>>>>>
>>>>> > The frequency vector for FFT should go from DC as bin 1 to slightly
>>>>> > less
>>>>> > than Fs in the last bin (Fs - Fs/N where N is the number of FFT
>>>>> > points).
>>>>>
>>>>> If i understand u correctly, the frequency vector (freq) must be
>>>>> written
>>>>> as in this modified version of code below. However, according to
>>>>> Nyquist
>>>>> we
>>>>> cannot measure a frequency higher than Fs/2. Thus i do not see
>>>>> the meaning of plotting up to ~Fs.
>>>>>
>>>>> Also, if I rerun the code below changing FsFactor, I again see this
>>>>> dependency of the signal frequency on FsFactor.
>>>>> I cannot understand it. Looks like basics, yet not that obvious.
>>>>>
>>>>>   #---------- start code ------------
>>>>>   clear all
>>>>>   nu=10;                    % signal frequency
>>>>>   at=5;                       % acquisition time
>>>>>   T1=2                       % signal decay constant in time domain
>>>>>   FsFactor=16;             % the ratio (sampling rate)/(signal
>>>>> frequency), or Fs/nu
>>>>>
>>>>>   clf
>>>>>   Fs = FsFactor*nu;                             % Sampling rate
>>>>>   t = ((0:(at*Fs-1))/Fs)';                         % Time vector
>>>>>   % freq=linspace(-1,1,Fs*at)' * Fs/2;         % frequency vector,
>>>>> first
>>>>> version
>>>>>   freq=linspace(0,1,Fs*at)' * Fs- Fs/(Fs*at);  % frequency vector
>>>>> suggested by Brian
>>>>>   tdsig=exp(-i*2*pi*nu*t).*exp(-t./(T1));    % time domain signal
>>>>>   freqsig=fft(tdsig);                                % freq. domain
>>>>> signal
>>>>>   subplot(1,2,1)
>>>>>   plot(t,tdsig)
>>>>>   axis([ 0 0.4])       % zooming time domain to see that period=1/nu
>>>>>   subplot(1,2,2)
>>>>>   plot(freq, freqsig)
>>>>>   #---------- end code ------------
>>>>>
>>>>> thx guys ...
>>>>>
>>>>>
>>>>> > 2014-02-12 14:22 GMT+01:00 oxy <[hidden email]>:
>>>>> >
>>>>> >> hey guys,
>>>>> >>
>>>>> >> the (simple) code bellow is how i ve learned to do FFT according to
>>>>> >> several docs online. However i observe a dependency of the signal
>>>>> >> frequency in the spectrum on the constant FsFactor. In other words,
>>>>> >> the signal frequency depends on the sampling rate. I should'nt be,
>>>>> >> right? So what is wrong here?
>>>>> >>
>>>>> >>   #---------- start code ------------
>>>>> >>   clear all
>>>>> >>   nu=10;                    % signal frequency
>>>>> >>   at=5;                       % acquisition time
>>>>> >>   T1=2                       % signal decay constant in time domain
>>>>> >>   FsFactor=8;             % the ratio (sampling rate)/(signal
>>>>> >> frequency), or Fs/nu
>>>>> >>
>>>>> >>   clf
>>>>> >>   Fs = FsFactor*nu;                             % Sampling rate
>>>>> >>   t = ((0:(at*Fs-1))/Fs)';                         % Time vector
>>>>> >>   freq=linspace(-1,1,Fs*at)' * Fs/2;         % frequency vector
>>>>> >>   tdsig=exp(-i*2*pi*nu*t).*exp(-t./(T1));    % time domain signal
>>>>> >>   freqsig=fft(tdsig);                                % freq. domain
>>>>> >> signal
>>>>> >>   subplot(1,2,1)
>>>>> >>   plot(t,tdsig)
>>>>> >>   axis([ 0 0.4])       % zooming time domain to see that
>>>>> >> period=1/nu
>>>>> >>   subplot(1,2,2)
>>>>> >>   plot(freq, freqsig)
>>>>> >>   #---------- end code ------------
>>>>> >>
>>>>> >> Do it yourself. Just rerun the code trying different values of
>>>>> >> FsFactor (eg: 2, 4, 8).
>>>>> >> Thx a lot for any hint!!!
>>>>> >>
>>>>> >> oxy
>>>>> >>
>>>>> >> ps: cross post
>>>>> >>
>>>>> http://www.mathworks.com/matlabcentral/answers/115700-fft-why-signal-frequency-depends-on-sampling-rate
>>>>> >> _______________________________________________
>>>>> >> Help-octave mailing list
>>>>> >> [hidden email]
>>>>> >> https://mailman.cae.wisc.edu/listinfo/help-octave
>>>>> >>
>>>>> >
>>>>> _______________________________________________
>>>>> Help-octave mailing list
>>>>> [hidden email]
>>>>> https://mailman.cae.wisc.edu/listinfo/help-octave
>>>>>
>>>>
>>>
>>
>
_______________________________________________
Help-octave mailing list
[hidden email]
https://mailman.cae.wisc.edu/listinfo/help-octave
Reply | Threaded
Open this post in threaded view
|

Re: Any FFT expert around?

oxyopes
Hi Brian,

> freq=linspace(-Fs1*(N-1)/(2*N),  Fs1*(N-1)/(2*N),  Fs1*at)';

ok, i will rewrite it with an example:

If i sample 11 points (N=11) at the rate Fs=4 and fft it,
my first point on the frequency vector will be zero.
The next will be df=Fs/N=0.36364.
The 3rd will be 2*df=Fs/N=0.72727
...
...

In other words:

point k, frequency
point 1:   0  * 0.36364 =   0.00000
point 2:   1   *  0.36364 =   0.36364
point 2:   2   *  0.36364 =   0.72727
...
point 6:   5   *  0.36364 =   1.81818
point 7:   -1   *  0.36364 =   -0.36364
point 8:   -2   *  0.36364 =   -0.72727
...
point 11:   -5   *  0.36364 =   -1.81818


Again in other words, here is the code:


Fs=4                          % sampling rate
N=11                          % nr. of points
df=Fs/N                       % freq. increment per point
MaxFreq=df*(N-1)/2            % max measurable freq
freqIdx=[(0:N/2)'; -(1:N/2)'] % index vector to get the freqs
freq=freqIdx*df;              % frequency vector
[freqIdx freq]                % index beside frequency
fftshift(freq)                % shifted freq vector


Thanks very much again for your very helpful wisdom!!!
All the best  ...



On 3/13/14, oxy <[hidden email]> wrote:

> Hi Brian,
>
>> freq=linspace(-Fs1*(N-1)/(2*N),  Fs1*(N-1)/(2*N),  Fs1*at)';
>
> ok, i will rewrite it with an example:
>
> If i sample 11 points (N=11) at the rate Fs=4 and fft it,
> my first point on the frequency vector will be zero.
> The next will be df=Fs/N=0.36364.
> The 3rd will be 2*df=Fs/N=0.72727
> ...
> ...
>
> In other words:
>
> point k, frequency
> point 1:   0  * 0.36364 =   0.00000
> point 2:   1   *  0.36364 =   0.36364
> point 2:   2   *  0.36364 =   0.72727
> ...
> point 6:   5   *  0.36364 =   1.81818
> point 7:   -1   *  0.36364 =   -0.36364
> point 8:   -2   *  0.36364 =   -0.72727
> ...
> point 11:   -5   *  0.36364 =   -1.81818
>
>
> Again in other words, here is the code:
>
>
> Fs=4                          % sampling rate
> N=11                          % nr. of points
> df=Fs/N                       % freq. increment per point
> MaxFreq=df*(N-1)/2            % max measurable freq
> freqIdx=[(0:N/2)'; -(1:N/2)'] % index vector to get the freqs
> freq=freqIdx*df;              % frequency vector
> [freqIdx freq]                % index beside frequency
> fftshift(freq)                % shifted freq vector
>
>
> Thanks very much again for your very helpful wisdom!!!
> All the best  ...
>
>
>> 2014-03-12 11:45 GMT+01:00 oxy <[hidden email]>:
>>
>>> hi all,
>>>
>>> again in the question about FFT frequency axis.
>>>
>>> Pls confirm if I'm doing it right here:
>>>
>>> If I do fftshift(fft(time_domain_signal))
>>> of a certain time_domain_signal with even number of points,
>>> then the frequency axis becomes
>>>
>>> freq=linspace(-Fs1/2,  Fs1/2-Fs1/(Fs1*at),  Fs1*at)';
>>>
>>> where Fs= sampling rate, at=acquisition time.
>>>
>>> If however my time domain signal has an odd number of points, then
>>>
>>> freq=linspace(-Fs1/2,  Fs1/2,  Fs1*at)';
>>>
>>> Am id doing it right?
>>>
>>> We had a recent thread on that. See some msg below.
>>>
>>> Thx for your wisdom!
>>>
>>> oxy
>>>
>>>
>>> > % TTF Demo: this code shows how to calculate the
>>> > % frequency domain when doing FFT of td-signals.
>>> >
>>> > %--- parameter setup: same signal, 2 sampling rates ---
>>> > clear all
>>> > nu=10;            % signal frequency
>>> > at=5;             % acquisition time
>>> > T1=2              % signal decay constant in time domain
>>> > FsFactor1=4;      % sampling rate factor 1, Fs1/nu
>>> > FsFactor2=16;     % sampling rate factor 2, Fs2/nu
>>> >
>>> > %----------- calculating -------------------
>>> > Fs1 = FsFactor1*nu;            % Sampling rate 1
>>> > Fs2 = FsFactor2*nu;            % Sampling rate 2
>>> >
>>> > t1 = ((0:(at*Fs1-1))/Fs1)';    % Time vector 1
>>> > t2 = ((0:(at*Fs2-1))/Fs2)';    % Time vector 2
>>> >
>>> > tdsig1=exp(-i*2*pi*nu*t1).*exp(-t1./(T1)); % time domain signal 1
>>> > tdsig2=exp(-i*2*pi*nu*t2).*exp(-t2./(T1)); % time domain signal 2
>>> >
>>> > % I realized, the signal was appearing negative, thus minus...
>>> > freq1=linspace(-Fs1/2,Fs1/2-Fs1/(Fs1*at),Fs1*at)'; % freq. vector 1
>>> > freq2=linspace(-Fs2/2,Fs2/2-Fs2/(Fs2*at),Fs2*at)'; % freq. vector 2
>>> >
>>> > freqsig1=fftshift(fft(tdsig1));  % signal vector 1
>>> > freqsig2=fftshift(fft(tdsig2));  % signal vector 2
>>> >
>>> >
>>> > %----------- ploting -------------------
>>> > clf
>>> >   subplot(1,2,1)
>>> >   plot(freq1, freqsig1)
>>> >   axis([-20 0]), title('FsFactor=4')
>>> >
>>> >   subplot(1,2,2)
>>> >   plot(freq2, freqsig2)
>>> >   axis([-20  0]), title('FsFactor=16')
>>> > % --------- end ---------------------------
>>> >
>>> >>>================================================
>>> >>> Just one question more: you are calculating the frequency
>>> >>> vector as:
>>> >>>
>>> >>>          freq1=-linspace(-Fs1/2,   Fs1/2-1/at,   Fs1*at)';
>>> >>>
>>> >>> that makes it non-symetric. Isn't it better ...
>>> >>>
>>> >>>          freq1=-linspace(-Fs1/2+1/at,   Fs1/2-1/at,   Fs1*at)';
>>> >>>
>>> >>> ...? Then it's symetric. Hey, thanks a lot! Great help!!!!!
>>> >>
>>> >> That would be "better" if you really want to insist on having a
>>> >> symmetrical frequency vector at the expense of accuracy.
>>> >> Unfortunately, it's inaccurate.  For one thing, your symmetrical
>>> >> frequency vector does not contain a DC term if you have an even
>>> >> number
>>> >> of FFT points (which you should based on your examples).
>>> >>
>>> >> If you want the FFT value at the nth point to really correspond to
>>> >> freq1(n) you should use the slightly offset frequency vector that I
>>> >> provided.
>>> >>
>>> >> If it is more important to whatever you're doing that the frequency
>>> >> vector be exactly symmetrical and you can live with sacrificing the
>>> >> true 1-to-1 correspondence between frequency points and FFT points,
>>> >> you can use your symmetrical frequency vector.
>>> >>
>>> >>> oxy
>>> >>>
>>> >>> =================================================
>>> >>> On 2/19/14, oxy <[hidden email]> wrote:
>>> >>>> On 2/12/14, Brian Kaczynski <[hidden email]> wrote:
>>> >>>>> Hi Oxy,
>>> >>>>>
>>> >>>>> Since you're taking the FFT of a complex vector it will have
>>> >>>>> unique
>>> >>>>> terms
>>> >>>>> for negative vs. positive frequency.  You are correct that it
>>> >>>>> makes
>>> >>>>> more
>>> >>>>> sense to think of the FFT from -Fs/2 to +Fs/2.  It's just that the
>>> >>>>> Octave
>>> >>>>> fft function doesn't compute the bins in that order.  Due to
>>> aliasing,
>>> >>>>> any
>>> >>>>> frequency Fs/2 + df is equivalent to (-Fs/2 + df) so it's more a
>>> >>>>> matter
>>> >>>>> of
>>> >>>>> choice how you want to display the data.
>>> >>>>>
>>> >>>>> If you prefer plotting from -Fs/2 to +Fs/2 I would change these
>>> >>>>> two
>>> >>>>> lines
>>> >>>>> of your code as follows:
>>> >>>>>
>>> >>>>> freq=linspace(-Fs/2,Fs/2-Fs/(Fs*at),Fs*at)'
>>> >>>>> freqsig=fftshift(fft(tdsig);
>>> >>>>>
>>> >>>>> Let us know if that gives you what you want!
>>> >>>>>
>>> >>>>> -Brian
>>> >>>>>
>>> >>>>>
>>> >>>>> 2014-02-12 17:15 GMT+01:00 oxy <[hidden email]>:
>>> >>>>>
>>> >>>>>> Hi Brian,
>>> >>>>>>
>>> >>>>>> > The frequency vector for FFT should go from DC as bin 1 to
>>> slightly
>>> >>>>>> > less
>>> >>>>>> > than Fs in the last bin (Fs - Fs/N where N is the number of FFT
>>> >>>>>> > points).
>>> >>>>>>
>>> >>>>>> If i understand u correctly, the frequency vector (freq) must be
>>> >>>>>> written
>>> >>>>>> as in this modified version of code below. However, according to
>>> >>>>>> Nyquist
>>> >>>>>> we
>>> >>>>>> cannot measure a frequency higher than Fs/2. Thus i do not see
>>> >>>>>> the meaning of plotting up to ~Fs.
>>> >>>>>>
>>> >>>>>> Also, if I rerun the code below changing FsFactor, I again see
>>> >>>>>> this
>>> >>>>>> dependency of the signal frequency on FsFactor.
>>> >>>>>> I cannot understand it. Looks like basics, yet not that obvious.
>>> >>>>>>
>>> >>>>>>   #---------- start code ------------
>>> >>>>>>   clear all
>>> >>>>>>   nu=10;                    % signal frequency
>>> >>>>>>   at=5;                       % acquisition time
>>> >>>>>>   T1=2                       % signal decay constant in time
>>> >>>>>> domain
>>> >>>>>>   FsFactor=16;             % the ratio (sampling rate)/(signal
>>> >>>>>> frequency), or Fs/nu
>>> >>>>>>
>>> >>>>>>   clf
>>> >>>>>>   Fs = FsFactor*nu;                             % Sampling rate
>>> >>>>>>   t = ((0:(at*Fs-1))/Fs)';                         % Time vector
>>> >>>>>>   % freq=linspace(-1,1,Fs*at)' * Fs/2;         % frequency
>>> >>>>>> vector,
>>> >>>>>> first
>>> >>>>>> version
>>> >>>>>>   freq=linspace(0,1,Fs*at)' * Fs- Fs/(Fs*at);  % frequency vector
>>> >>>>>> suggested by Brian
>>> >>>>>>   tdsig=exp(-i*2*pi*nu*t).*exp(-t./(T1));    % time domain signal
>>> >>>>>>   freqsig=fft(tdsig);                                % freq.
>>> >>>>>> domain
>>> >>>>>> signal
>>> >>>>>>   subplot(1,2,1)
>>> >>>>>>   plot(t,tdsig)
>>> >>>>>>   axis([ 0 0.4])       % zooming time domain to see that
>>> >>>>>> period=1/nu
>>> >>>>>>   subplot(1,2,2)
>>> >>>>>>   plot(freq, freqsig)
>>> >>>>>>   #---------- end code ------------
>>> >>>>>>
>>> >>>>>> thx guys ...
>>> >>>>>>
>>> >>>>>>
>>> >>>>>> > 2014-02-12 14:22 GMT+01:00 oxy <[hidden email]>:
>>> >>>>>> >
>>> >>>>>> >> hey guys,
>>> >>>>>> >>
>>> >>>>>> >> the (simple) code bellow is how i ve learned to do FFT
>>> >>>>>> >> according
>>> >>>>>> >> to
>>> >>>>>> >> several docs online. However i observe a dependency of the
>>> >>>>>> >> signal
>>> >>>>>> >> frequency in the spectrum on the constant FsFactor. In other
>>> >>>>>> >> words,
>>> >>>>>> >> the signal frequency depends on the sampling rate. I should'nt
>>> be,
>>> >>>>>> >> right? So what is wrong here?
>>> >>>>>> >>
>>> >>>>>> >>   #---------- start code ------------
>>> >>>>>> >>   clear all
>>> >>>>>> >>   nu=10;                    % signal frequency
>>> >>>>>> >>   at=5;                       % acquisition time
>>> >>>>>> >>   T1=2                       % signal decay constant in time
>>> >>>>>> >> domain
>>> >>>>>> >>   FsFactor=8;             % the ratio (sampling rate)/(signal
>>> >>>>>> >> frequency), or Fs/nu
>>> >>>>>> >>
>>> >>>>>> >>   clf
>>> >>>>>> >>   Fs = FsFactor*nu;                             % Sampling
>>> >>>>>> >> rate
>>> >>>>>> >>   t = ((0:(at*Fs-1))/Fs)';                         % Time
>>> >>>>>> >> vector
>>> >>>>>> >>   freq=linspace(-1,1,Fs*at)' * Fs/2;         % frequency
>>> >>>>>> >> vector
>>> >>>>>> >>   tdsig=exp(-i*2*pi*nu*t).*exp(-t./(T1));    % time domain
>>> >>>>>> >> signal
>>> >>>>>> >>   freqsig=fft(tdsig);                                % freq.
>>> >>>>>> >> domain
>>> >>>>>> >> signal
>>> >>>>>> >>   subplot(1,2,1)
>>> >>>>>> >>   plot(t,tdsig)
>>> >>>>>> >>   axis([ 0 0.4])       % zooming time domain to see that
>>> >>>>>> >> period=1/nu
>>> >>>>>> >>   subplot(1,2,2)
>>> >>>>>> >>   plot(freq, freqsig)
>>> >>>>>> >>   #---------- end code ------------
>>> >>>>>> >>
>>> >>>>>> >> Do it yourself. Just rerun the code trying different values of
>>> >>>>>> >> FsFactor (eg: 2, 4, 8).
>>> >>>>>> >> Thx a lot for any hint!!!
>>> >>>>>> >>
>>> >>>>>> >> oxy
>>> >>>>>> >>
>>> >>>>>> >> ps: cross post
>>> >>>>>> >>
>>> >>>>>>
>>> http://www.mathworks.com/matlabcentral/answers/115700-fft-why-signal-frequency-depends-on-sampling-rate
>>> >>>>>> >> _______________________________________________
>>> >>>>>> >> Help-octave mailing list
>>> >>>>>> >> [hidden email]
>>> >>>>>> >> https://mailman.cae.wisc.edu/listinfo/help-octave
>>> >>>>>> >>
>>> >>>>>> >
>>> >>>>>> _______________________________________________
>>> >>>>>> Help-octave mailing list
>>> >>>>>> [hidden email]
>>> >>>>>> https://mailman.cae.wisc.edu/listinfo/help-octave
>>> >>>>>>
>>> >>>>>
>>> >>>>
>>> >>>
>>> >>
>>> >
>>>
>>
>
_______________________________________________
Help-octave mailing list
[hidden email]
https://mailman.cae.wisc.edu/listinfo/help-octave
Reply | Threaded
Open this post in threaded view
|

Re: Any FFT expert around?

briankaz
Hi Oxy,

Not quite.  See below.

2014-03-13 14:57 GMT+01:00 oxy <[hidden email]>:
Hi Brian,

> freq=linspace(-Fs1*(N-1)/(2*N),  Fs1*(N-1)/(2*N),  Fs1*at)';

ok, i will rewrite it with an example:

If i sample 11 points (N=11) at the rate Fs=4 and fft it,
my first point on the frequency vector will be zero.
The next will be df=Fs/N=0.36364.
The 3rd will be 2*df=Fs/N=0.72727
...
...

In other words:

point k, frequency
point 1:   0  * 0.36364 =   0.00000
point 2:   1   *  0.36364 =   0.36364
point 2:   2   *  0.36364 =   0.72727
...
point 6:   5   *  0.36364 =   1.81818
point 7:   -1   *  0.36364 =   -0.36364
point 8:   -2   *  0.36364 =   -0.72727
...
point 11:   -5   *  0.36364 =   -1.81818
 
 
Points 1-6 are correct but 7-11 are backwards.  Point 7 should be -5 * 0.36364 and point 11 should be -1 * 0.36364.

Let's say you place the original points in order as follows:
[1 2 3 4 5 6 7 8 9 10 11]
When you do fftshift() on the fft output, the points get rearranged as follows:
[7 8 9 10 11 1 2 3 4 5 6]

So after this rearrangement, the frequencies must be in this order:
(-5:5)*0.36364

I hope this makes sense!
 

Again in other words, here is the code:


Fs=4                          % sampling rate
N=11                          % nr. of points
df=Fs/N                       % freq. increment per point
MaxFreq=df*(N-1)/2            % max measurable freq
freqIdx=[(0:N/2)'; -(1:N/2)'] % index vector to get the freqs
freq=freqIdx*df;              % frequency vector
[freqIdx freq]                % index beside frequency
fftshift(freq)                % shifted freq vector


Thanks very much again for your very helpful wisdom!!!
All the best  ...



On 3/13/14, oxy <[hidden email]> wrote:
> Hi Brian,
>
>> freq=linspace(-Fs1*(N-1)/(2*N),  Fs1*(N-1)/(2*N),  Fs1*at)';
>
> ok, i will rewrite it with an example:
>
> If i sample 11 points (N=11) at the rate Fs=4 and fft it,
> my first point on the frequency vector will be zero.
> The next will be df=Fs/N=0.36364.
> The 3rd will be 2*df=Fs/N=0.72727
> ...
> ...
>
> In other words:
>
> point k, frequency
> point 1:   0  * 0.36364 =   0.00000
> point 2:   1   *  0.36364 =   0.36364
> point 2:   2   *  0.36364 =   0.72727
> ...
> point 6:   5   *  0.36364 =   1.81818
> point 7:   -1   *  0.36364 =   -0.36364
> point 8:   -2   *  0.36364 =   -0.72727
> ...
> point 11:   -5   *  0.36364 =   -1.81818
>
>
> Again in other words, here is the code:
>
>
> Fs=4                          % sampling rate
> N=11                          % nr. of points
> df=Fs/N                       % freq. increment per point
> MaxFreq=df*(N-1)/2            % max measurable freq
> freqIdx=[(0:N/2)'; -(1:N/2)'] % index vector to get the freqs
> freq=freqIdx*df;              % frequency vector
> [freqIdx freq]                % index beside frequency
> fftshift(freq)                % shifted freq vector
>
>
> Thanks very much again for your very helpful wisdom!!!
> All the best  ...
>
>
>> 2014-03-12 11:45 GMT+01:00 oxy <[hidden email]>:
>>
>>> hi all,
>>>
>>> again in the question about FFT frequency axis.
>>>
>>> Pls confirm if I'm doing it right here:
>>>
>>> If I do fftshift(fft(time_domain_signal))
>>> of a certain time_domain_signal with even number of points,
>>> then the frequency axis becomes
>>>
>>> freq=linspace(-Fs1/2,  Fs1/2-Fs1/(Fs1*at),  Fs1*at)';
>>>
>>> where Fs= sampling rate, at=acquisition time.
>>>
>>> If however my time domain signal has an odd number of points, then
>>>
>>> freq=linspace(-Fs1/2,  Fs1/2,  Fs1*at)';
>>>
>>> Am id doing it right?
>>>
>>> We had a recent thread on that. See some msg below.
>>>
>>> Thx for your wisdom!
>>>
>>> oxy
>>>
>>>
>>> > % TTF Demo: this code shows how to calculate the
>>> > % frequency domain when doing FFT of td-signals.
>>> >
>>> > %--- parameter setup: same signal, 2 sampling rates ---
>>> > clear all
>>> > nu=10;            % signal frequency
>>> > at=5;             % acquisition time
>>> > T1=2              % signal decay constant in time domain
>>> > FsFactor1=4;      % sampling rate factor 1, Fs1/nu
>>> > FsFactor2=16;     % sampling rate factor 2, Fs2/nu
>>> >
>>> > %----------- calculating -------------------
>>> > Fs1 = FsFactor1*nu;            % Sampling rate 1
>>> > Fs2 = FsFactor2*nu;            % Sampling rate 2
>>> >
>>> > t1 = ((0:(at*Fs1-1))/Fs1)';    % Time vector 1
>>> > t2 = ((0:(at*Fs2-1))/Fs2)';    % Time vector 2
>>> >
>>> > tdsig1=exp(-i*2*pi*nu*t1).*exp(-t1./(T1)); % time domain signal 1
>>> > tdsig2=exp(-i*2*pi*nu*t2).*exp(-t2./(T1)); % time domain signal 2
>>> >
>>> > % I realized, the signal was appearing negative, thus minus...
>>> > freq1=linspace(-Fs1/2,Fs1/2-Fs1/(Fs1*at),Fs1*at)'; % freq. vector 1
>>> > freq2=linspace(-Fs2/2,Fs2/2-Fs2/(Fs2*at),Fs2*at)'; % freq. vector 2
>>> >
>>> > freqsig1=fftshift(fft(tdsig1));  % signal vector 1
>>> > freqsig2=fftshift(fft(tdsig2));  % signal vector 2
>>> >
>>> >
>>> > %----------- ploting -------------------
>>> > clf
>>> >   subplot(1,2,1)
>>> >   plot(freq1, freqsig1)
>>> >   axis([-20 0]), title('FsFactor=4')
>>> >
>>> >   subplot(1,2,2)
>>> >   plot(freq2, freqsig2)
>>> >   axis([-20  0]), title('FsFactor=16')
>>> > % --------- end ---------------------------
>>> >
>>> >>>================================================
>>> >>> Just one question more: you are calculating the frequency
>>> >>> vector as:
>>> >>>
>>> >>>          freq1=-linspace(-Fs1/2,   Fs1/2-1/at,   Fs1*at)';
>>> >>>
>>> >>> that makes it non-symetric. Isn't it better ...
>>> >>>
>>> >>>          freq1=-linspace(-Fs1/2+1/at,   Fs1/2-1/at,   Fs1*at)';
>>> >>>
>>> >>> ...? Then it's symetric. Hey, thanks a lot! Great help!!!!!
>>> >>
>>> >> That would be "better" if you really want to insist on having a
>>> >> symmetrical frequency vector at the expense of accuracy.
>>> >> Unfortunately, it's inaccurate.  For one thing, your symmetrical
>>> >> frequency vector does not contain a DC term if you have an even
>>> >> number
>>> >> of FFT points (which you should based on your examples).
>>> >>
>>> >> If you want the FFT value at the nth point to really correspond to
>>> >> freq1(n) you should use the slightly offset frequency vector that I
>>> >> provided.
>>> >>
>>> >> If it is more important to whatever you're doing that the frequency
>>> >> vector be exactly symmetrical and you can live with sacrificing the
>>> >> true 1-to-1 correspondence between frequency points and FFT points,
>>> >> you can use your symmetrical frequency vector.
>>> >>
>>> >>> oxy
>>> >>>
>>> >>> =================================================
>>> >>> On 2/19/14, oxy <[hidden email]> wrote:
>>> >>>> On 2/12/14, Brian Kaczynski <[hidden email]> wrote:
>>> >>>>> Hi Oxy,
>>> >>>>>
>>> >>>>> Since you're taking the FFT of a complex vector it will have
>>> >>>>> unique
>>> >>>>> terms
>>> >>>>> for negative vs. positive frequency.  You are correct that it
>>> >>>>> makes
>>> >>>>> more
>>> >>>>> sense to think of the FFT from -Fs/2 to +Fs/2.  It's just that the
>>> >>>>> Octave
>>> >>>>> fft function doesn't compute the bins in that order.  Due to
>>> aliasing,
>>> >>>>> any
>>> >>>>> frequency Fs/2 + df is equivalent to (-Fs/2 + df) so it's more a
>>> >>>>> matter
>>> >>>>> of
>>> >>>>> choice how you want to display the data.
>>> >>>>>
>>> >>>>> If you prefer plotting from -Fs/2 to +Fs/2 I would change these
>>> >>>>> two
>>> >>>>> lines
>>> >>>>> of your code as follows:
>>> >>>>>
>>> >>>>> freq=linspace(-Fs/2,Fs/2-Fs/(Fs*at),Fs*at)'
>>> >>>>> freqsig=fftshift(fft(tdsig);
>>> >>>>>
>>> >>>>> Let us know if that gives you what you want!
>>> >>>>>
>>> >>>>> -Brian
>>> >>>>>
>>> >>>>>
>>> >>>>> 2014-02-12 17:15 GMT+01:00 oxy <[hidden email]>:
>>> >>>>>
>>> >>>>>> Hi Brian,
>>> >>>>>>
>>> >>>>>> > The frequency vector for FFT should go from DC as bin 1 to
>>> slightly
>>> >>>>>> > less
>>> >>>>>> > than Fs in the last bin (Fs - Fs/N where N is the number of FFT
>>> >>>>>> > points).
>>> >>>>>>
>>> >>>>>> If i understand u correctly, the frequency vector (freq) must be
>>> >>>>>> written
>>> >>>>>> as in this modified version of code below. However, according to
>>> >>>>>> Nyquist
>>> >>>>>> we
>>> >>>>>> cannot measure a frequency higher than Fs/2. Thus i do not see
>>> >>>>>> the meaning of plotting up to ~Fs.
>>> >>>>>>
>>> >>>>>> Also, if I rerun the code below changing FsFactor, I again see
>>> >>>>>> this
>>> >>>>>> dependency of the signal frequency on FsFactor.
>>> >>>>>> I cannot understand it. Looks like basics, yet not that obvious.
>>> >>>>>>
>>> >>>>>>   #---------- start code ------------
>>> >>>>>>   clear all
>>> >>>>>>   nu=10;                    % signal frequency
>>> >>>>>>   at=5;                       % acquisition time
>>> >>>>>>   T1=2                       % signal decay constant in time
>>> >>>>>> domain
>>> >>>>>>   FsFactor=16;             % the ratio (sampling rate)/(signal
>>> >>>>>> frequency), or Fs/nu
>>> >>>>>>
>>> >>>>>>   clf
>>> >>>>>>   Fs = FsFactor*nu;                             % Sampling rate
>>> >>>>>>   t = ((0:(at*Fs-1))/Fs)';                         % Time vector
>>> >>>>>>   % freq=linspace(-1,1,Fs*at)' * Fs/2;         % frequency
>>> >>>>>> vector,
>>> >>>>>> first
>>> >>>>>> version
>>> >>>>>>   freq=linspace(0,1,Fs*at)' * Fs- Fs/(Fs*at);  % frequency vector
>>> >>>>>> suggested by Brian
>>> >>>>>>   tdsig=exp(-i*2*pi*nu*t).*exp(-t./(T1));    % time domain signal
>>> >>>>>>   freqsig=fft(tdsig);                                % freq.
>>> >>>>>> domain
>>> >>>>>> signal
>>> >>>>>>   subplot(1,2,1)
>>> >>>>>>   plot(t,tdsig)
>>> >>>>>>   axis([ 0 0.4])       % zooming time domain to see that
>>> >>>>>> period=1/nu
>>> >>>>>>   subplot(1,2,2)
>>> >>>>>>   plot(freq, freqsig)
>>> >>>>>>   #---------- end code ------------
>>> >>>>>>
>>> >>>>>> thx guys ...
>>> >>>>>>
>>> >>>>>>
>>> >>>>>> > 2014-02-12 14:22 GMT+01:00 oxy <[hidden email]>:
>>> >>>>>> >
>>> >>>>>> >> hey guys,
>>> >>>>>> >>
>>> >>>>>> >> the (simple) code bellow is how i ve learned to do FFT
>>> >>>>>> >> according
>>> >>>>>> >> to
>>> >>>>>> >> several docs online. However i observe a dependency of the
>>> >>>>>> >> signal
>>> >>>>>> >> frequency in the spectrum on the constant FsFactor. In other
>>> >>>>>> >> words,
>>> >>>>>> >> the signal frequency depends on the sampling rate. I should'nt
>>> be,
>>> >>>>>> >> right? So what is wrong here?
>>> >>>>>> >>
>>> >>>>>> >>   #---------- start code ------------
>>> >>>>>> >>   clear all
>>> >>>>>> >>   nu=10;                    % signal frequency
>>> >>>>>> >>   at=5;                       % acquisition time
>>> >>>>>> >>   T1=2                       % signal decay constant in time
>>> >>>>>> >> domain
>>> >>>>>> >>   FsFactor=8;             % the ratio (sampling rate)/(signal
>>> >>>>>> >> frequency), or Fs/nu
>>> >>>>>> >>
>>> >>>>>> >>   clf
>>> >>>>>> >>   Fs = FsFactor*nu;                             % Sampling
>>> >>>>>> >> rate
>>> >>>>>> >>   t = ((0:(at*Fs-1))/Fs)';                         % Time
>>> >>>>>> >> vector
>>> >>>>>> >>   freq=linspace(-1,1,Fs*at)' * Fs/2;         % frequency
>>> >>>>>> >> vector
>>> >>>>>> >>   tdsig=exp(-i*2*pi*nu*t).*exp(-t./(T1));    % time domain
>>> >>>>>> >> signal
>>> >>>>>> >>   freqsig=fft(tdsig);                                % freq.
>>> >>>>>> >> domain
>>> >>>>>> >> signal
>>> >>>>>> >>   subplot(1,2,1)
>>> >>>>>> >>   plot(t,tdsig)
>>> >>>>>> >>   axis([ 0 0.4])       % zooming time domain to see that
>>> >>>>>> >> period=1/nu
>>> >>>>>> >>   subplot(1,2,2)
>>> >>>>>> >>   plot(freq, freqsig)
>>> >>>>>> >>   #---------- end code ------------
>>> >>>>>> >>
>>> >>>>>> >> Do it yourself. Just rerun the code trying different values of
>>> >>>>>> >> FsFactor (eg: 2, 4, 8).
>>> >>>>>> >> Thx a lot for any hint!!!
>>> >>>>>> >>
>>> >>>>>> >> oxy
>>> >>>>>> >>
>>> >>>>>> >> ps: cross post
>>> >>>>>> >>
>>> >>>>>>
>>> http://www.mathworks.com/matlabcentral/answers/115700-fft-why-signal-frequency-depends-on-sampling-rate
>>> >>>>>> >> _______________________________________________
>>> >>>>>> >> Help-octave mailing list
>>> >>>>>> >> [hidden email]
>>> >>>>>> >> https://mailman.cae.wisc.edu/listinfo/help-octave
>>> >>>>>> >>
>>> >>>>>> >
>>> >>>>>> _______________________________________________
>>> >>>>>> Help-octave mailing list
>>> >>>>>> [hidden email]
>>> >>>>>> https://mailman.cae.wisc.edu/listinfo/help-octave
>>> >>>>>>
>>> >>>>>
>>> >>>>
>>> >>>
>>> >>
>>> >
>>>
>>
>
_______________________________________________
Help-octave mailing list
[hidden email]
https://mailman.cae.wisc.edu/listinfo/help-octave


_______________________________________________
Help-octave mailing list
[hidden email]
https://mailman.cae.wisc.edu/listinfo/help-octave
Reply | Threaded
Open this post in threaded view
|

Re: Any FFT expert around?

oxyopes
Hi Brian,

>> point k, frequency
>> point 1:   0  * 0.36364 =   0.00000
>> point 2:   1   *  0.36364 =   0.36364
>> point 2:   2   *  0.36364 =   0.72727
>> ...
>> point 6:   5   *  0.36364 =   1.81818
>> point 7:   -1   *  0.36364 =   -0.36364
>> point 8:   -2   *  0.36364 =   -0.72727
>> ...
>> point 11:   -5   *  0.36364 =   -1.81818
>
> Points 1-6 are correct but 7-11 are backwards.  Point 7 should be -5 *
> 0.36364 and point 11 should be -1 * 0.36364.
>
> Let's say you place the original points in order as follows:
> [1 2 3 4 5 6 7 8 9 10 11]
> When you do fftshift() on the fft output, the points get rearranged as
> follows:
> [7 8 9 10 11 1 2 3 4 5 6]
> I hope this makes sense!

Yes that makes a lot of sense!
It was actually straight forward to deduce from you previous email.
And I did not pay attention on that.

Thank you a lot for your very precise answers!
A precious email correspondent!

All the best ...




>> On 3/13/14, oxy <[hidden email]> wrote:
>> > Hi Brian,
>> >
>> >> freq=linspace(-Fs1*(N-1)/(2*N),  Fs1*(N-1)/(2*N),  Fs1*at)';
>> >
>> > ok, i will rewrite it with an example:
>> >
>> > If i sample 11 points (N=11) at the rate Fs=4 and fft it,
>> > my first point on the frequency vector will be zero.
>> > The next will be df=Fs/N=0.36364.
>> > The 3rd will be 2*df=Fs/N=0.72727
>> > ...
>> > ...
>> >
>> > In other words:
>> >
>> > point k, frequency
>> > point 1:   0  * 0.36364 =   0.00000
>> > point 2:   1   *  0.36364 =   0.36364
>> > point 2:   2   *  0.36364 =   0.72727
>> > ...
>> > point 6:   5   *  0.36364 =   1.81818
>> > point 7:   -1   *  0.36364 =   -0.36364
>> > point 8:   -2   *  0.36364 =   -0.72727
>> > ...
>> > point 11:   -5   *  0.36364 =   -1.81818
>> >
>> >
>> > Again in other words, here is the code:
>> >
>> >
>> > Fs=4                          % sampling rate
>> > N=11                          % nr. of points
>> > df=Fs/N                       % freq. increment per point
>> > MaxFreq=df*(N-1)/2            % max measurable freq
>> > freqIdx=[(0:N/2)'; -(1:N/2)'] % index vector to get the freqs
>> > freq=freqIdx*df;              % frequency vector
>> > [freqIdx freq]                % index beside frequency
>> > fftshift(freq)                % shifted freq vector
>> >
>> >
>> > Thanks very much again for your very helpful wisdom!!!
>> > All the best  ...
>> >
>> >
>> >> 2014-03-12 11:45 GMT+01:00 oxy <[hidden email]>:
>> >>
>> >>> hi all,
>> >>>
>> >>> again in the question about FFT frequency axis.
>> >>>
>> >>> Pls confirm if I'm doing it right here:
>> >>>
>> >>> If I do fftshift(fft(time_domain_signal))
>> >>> of a certain time_domain_signal with even number of points,
>> >>> then the frequency axis becomes
>> >>>
>> >>> freq=linspace(-Fs1/2,  Fs1/2-Fs1/(Fs1*at),  Fs1*at)';
>> >>>
>> >>> where Fs= sampling rate, at=acquisition time.
>> >>>
>> >>> If however my time domain signal has an odd number of points, then
>> >>>
>> >>> freq=linspace(-Fs1/2,  Fs1/2,  Fs1*at)';
>> >>>
>> >>> Am id doing it right?
>> >>>
>> >>> We had a recent thread on that. See some msg below.
>> >>>
>> >>> Thx for your wisdom!
>> >>>
>> >>> oxy
>> >>>
>> >>>
>> >>> > % TTF Demo: this code shows how to calculate the
>> >>> > % frequency domain when doing FFT of td-signals.
>> >>> >
>> >>> > %--- parameter setup: same signal, 2 sampling rates ---
>> >>> > clear all
>> >>> > nu=10;            % signal frequency
>> >>> > at=5;             % acquisition time
>> >>> > T1=2              % signal decay constant in time domain
>> >>> > FsFactor1=4;      % sampling rate factor 1, Fs1/nu
>> >>> > FsFactor2=16;     % sampling rate factor 2, Fs2/nu
>> >>> >
>> >>> > %----------- calculating -------------------
>> >>> > Fs1 = FsFactor1*nu;            % Sampling rate 1
>> >>> > Fs2 = FsFactor2*nu;            % Sampling rate 2
>> >>> >
>> >>> > t1 = ((0:(at*Fs1-1))/Fs1)';    % Time vector 1
>> >>> > t2 = ((0:(at*Fs2-1))/Fs2)';    % Time vector 2
>> >>> >
>> >>> > tdsig1=exp(-i*2*pi*nu*t1).*exp(-t1./(T1)); % time domain signal 1
>> >>> > tdsig2=exp(-i*2*pi*nu*t2).*exp(-t2./(T1)); % time domain signal 2
>> >>> >
>> >>> > % I realized, the signal was appearing negative, thus minus...
>> >>> > freq1=linspace(-Fs1/2,Fs1/2-Fs1/(Fs1*at),Fs1*at)'; % freq. vector 1
>> >>> > freq2=linspace(-Fs2/2,Fs2/2-Fs2/(Fs2*at),Fs2*at)'; % freq. vector 2
>> >>> >
>> >>> > freqsig1=fftshift(fft(tdsig1));  % signal vector 1
>> >>> > freqsig2=fftshift(fft(tdsig2));  % signal vector 2
>> >>> >
>> >>> >
>> >>> > %----------- ploting -------------------
>> >>> > clf
>> >>> >   subplot(1,2,1)
>> >>> >   plot(freq1, freqsig1)
>> >>> >   axis([-20 0]), title('FsFactor=4')
>> >>> >
>> >>> >   subplot(1,2,2)
>> >>> >   plot(freq2, freqsig2)
>> >>> >   axis([-20  0]), title('FsFactor=16')
>> >>> > % --------- end ---------------------------
>> >>> >
>> >>> >>>================================================
>> >>> >>> Just one question more: you are calculating the frequency
>> >>> >>> vector as:
>> >>> >>>
>> >>> >>>          freq1=-linspace(-Fs1/2,   Fs1/2-1/at,   Fs1*at)';
>> >>> >>>
>> >>> >>> that makes it non-symetric. Isn't it better ...
>> >>> >>>
>> >>> >>>          freq1=-linspace(-Fs1/2+1/at,   Fs1/2-1/at,   Fs1*at)';
>> >>> >>>
>> >>> >>> ...? Then it's symetric. Hey, thanks a lot! Great help!!!!!
>> >>> >>
>> >>> >> That would be "better" if you really want to insist on having a
>> >>> >> symmetrical frequency vector at the expense of accuracy.
>> >>> >> Unfortunately, it's inaccurate.  For one thing, your symmetrical
>> >>> >> frequency vector does not contain a DC term if you have an even
>> >>> >> number
>> >>> >> of FFT points (which you should based on your examples).
>> >>> >>
>> >>> >> If you want the FFT value at the nth point to really correspond to
>> >>> >> freq1(n) you should use the slightly offset frequency vector that
>> >>> >> I
>> >>> >> provided.
>> >>> >>
>> >>> >> If it is more important to whatever you're doing that the
>> >>> >> frequency
>> >>> >> vector be exactly symmetrical and you can live with sacrificing
>> >>> >> the
>> >>> >> true 1-to-1 correspondence between frequency points and FFT
>> >>> >> points,
>> >>> >> you can use your symmetrical frequency vector.
>> >>> >>
>> >>> >>> oxy
>> >>> >>>
>> >>> >>> =================================================
>> >>> >>> On 2/19/14, oxy <[hidden email]> wrote:
>> >>> >>>> On 2/12/14, Brian Kaczynski <[hidden email]> wrote:
>> >>> >>>>> Hi Oxy,
>> >>> >>>>>
>> >>> >>>>> Since you're taking the FFT of a complex vector it will have
>> >>> >>>>> unique
>> >>> >>>>> terms
>> >>> >>>>> for negative vs. positive frequency.  You are correct that it
>> >>> >>>>> makes
>> >>> >>>>> more
>> >>> >>>>> sense to think of the FFT from -Fs/2 to +Fs/2.  It's just that
>> the
>> >>> >>>>> Octave
>> >>> >>>>> fft function doesn't compute the bins in that order.  Due to
>> >>> aliasing,
>> >>> >>>>> any
>> >>> >>>>> frequency Fs/2 + df is equivalent to (-Fs/2 + df) so it's more
>> >>> >>>>> a
>> >>> >>>>> matter
>> >>> >>>>> of
>> >>> >>>>> choice how you want to display the data.
>> >>> >>>>>
>> >>> >>>>> If you prefer plotting from -Fs/2 to +Fs/2 I would change these
>> >>> >>>>> two
>> >>> >>>>> lines
>> >>> >>>>> of your code as follows:
>> >>> >>>>>
>> >>> >>>>> freq=linspace(-Fs/2,Fs/2-Fs/(Fs*at),Fs*at)'
>> >>> >>>>> freqsig=fftshift(fft(tdsig);
>> >>> >>>>>
>> >>> >>>>> Let us know if that gives you what you want!
>> >>> >>>>>
>> >>> >>>>> -Brian
>> >>> >>>>>
>> >>> >>>>>
>> >>> >>>>> 2014-02-12 17:15 GMT+01:00 oxy <[hidden email]>:
>> >>> >>>>>
>> >>> >>>>>> Hi Brian,
>> >>> >>>>>>
>> >>> >>>>>> > The frequency vector for FFT should go from DC as bin 1 to
>> >>> slightly
>> >>> >>>>>> > less
>> >>> >>>>>> > than Fs in the last bin (Fs - Fs/N where N is the number of
>> FFT
>> >>> >>>>>> > points).
>> >>> >>>>>>
>> >>> >>>>>> If i understand u correctly, the frequency vector (freq) must
>> >>> >>>>>> be
>> >>> >>>>>> written
>> >>> >>>>>> as in this modified version of code below. However, according
>> >>> >>>>>> to
>> >>> >>>>>> Nyquist
>> >>> >>>>>> we
>> >>> >>>>>> cannot measure a frequency higher than Fs/2. Thus i do not see
>> >>> >>>>>> the meaning of plotting up to ~Fs.
>> >>> >>>>>>
>> >>> >>>>>> Also, if I rerun the code below changing FsFactor, I again see
>> >>> >>>>>> this
>> >>> >>>>>> dependency of the signal frequency on FsFactor.
>> >>> >>>>>> I cannot understand it. Looks like basics, yet not that
>> >>> >>>>>> obvious.
>> >>> >>>>>>
>> >>> >>>>>>   #---------- start code ------------
>> >>> >>>>>>   clear all
>> >>> >>>>>>   nu=10;                    % signal frequency
>> >>> >>>>>>   at=5;                       % acquisition time
>> >>> >>>>>>   T1=2                       % signal decay constant in time
>> >>> >>>>>> domain
>> >>> >>>>>>   FsFactor=16;             % the ratio (sampling rate)/(signal
>> >>> >>>>>> frequency), or Fs/nu
>> >>> >>>>>>
>> >>> >>>>>>   clf
>> >>> >>>>>>   Fs = FsFactor*nu;                             % Sampling
>> >>> >>>>>> rate
>> >>> >>>>>>   t = ((0:(at*Fs-1))/Fs)';                         % Time
>> >>> >>>>>> vector
>> >>> >>>>>>   % freq=linspace(-1,1,Fs*at)' * Fs/2;         % frequency
>> >>> >>>>>> vector,
>> >>> >>>>>> first
>> >>> >>>>>> version
>> >>> >>>>>>   freq=linspace(0,1,Fs*at)' * Fs- Fs/(Fs*at);  % frequency
>> vector
>> >>> >>>>>> suggested by Brian
>> >>> >>>>>>   tdsig=exp(-i*2*pi*nu*t).*exp(-t./(T1));    % time domain
>> signal
>> >>> >>>>>>   freqsig=fft(tdsig);                                % freq.
>> >>> >>>>>> domain
>> >>> >>>>>> signal
>> >>> >>>>>>   subplot(1,2,1)
>> >>> >>>>>>   plot(t,tdsig)
>> >>> >>>>>>   axis([ 0 0.4])       % zooming time domain to see that
>> >>> >>>>>> period=1/nu
>> >>> >>>>>>   subplot(1,2,2)
>> >>> >>>>>>   plot(freq, freqsig)
>> >>> >>>>>>   #---------- end code ------------
>> >>> >>>>>>
>> >>> >>>>>> thx guys ...
>> >>> >>>>>>
>> >>> >>>>>>
>> >>> >>>>>> > 2014-02-12 14:22 GMT+01:00 oxy <[hidden email]>:
>> >>> >>>>>> >
>> >>> >>>>>> >> hey guys,
>> >>> >>>>>> >>
>> >>> >>>>>> >> the (simple) code bellow is how i ve learned to do FFT
>> >>> >>>>>> >> according
>> >>> >>>>>> >> to
>> >>> >>>>>> >> several docs online. However i observe a dependency of the
>> >>> >>>>>> >> signal
>> >>> >>>>>> >> frequency in the spectrum on the constant FsFactor. In
>> >>> >>>>>> >> other
>> >>> >>>>>> >> words,
>> >>> >>>>>> >> the signal frequency depends on the sampling rate. I
>> should'nt
>> >>> be,
>> >>> >>>>>> >> right? So what is wrong here?
>> >>> >>>>>> >>
>> >>> >>>>>> >>   #---------- start code ------------
>> >>> >>>>>> >>   clear all
>> >>> >>>>>> >>   nu=10;                    % signal frequency
>> >>> >>>>>> >>   at=5;                       % acquisition time
>> >>> >>>>>> >>   T1=2                       % signal decay constant in
>> >>> >>>>>> >> time
>> >>> >>>>>> >> domain
>> >>> >>>>>> >>   FsFactor=8;             % the ratio (sampling
>> >>> >>>>>> >> rate)/(signal
>> >>> >>>>>> >> frequency), or Fs/nu
>> >>> >>>>>> >>
>> >>> >>>>>> >>   clf
>> >>> >>>>>> >>   Fs = FsFactor*nu;                             % Sampling
>> >>> >>>>>> >> rate
>> >>> >>>>>> >>   t = ((0:(at*Fs-1))/Fs)';                         % Time
>> >>> >>>>>> >> vector
>> >>> >>>>>> >>   freq=linspace(-1,1,Fs*at)' * Fs/2;         % frequency
>> >>> >>>>>> >> vector
>> >>> >>>>>> >>   tdsig=exp(-i*2*pi*nu*t).*exp(-t./(T1));    % time domain
>> >>> >>>>>> >> signal
>> >>> >>>>>> >>   freqsig=fft(tdsig);                                %
>> >>> >>>>>> >> freq.
>> >>> >>>>>> >> domain
>> >>> >>>>>> >> signal
>> >>> >>>>>> >>   subplot(1,2,1)
>> >>> >>>>>> >>   plot(t,tdsig)
>> >>> >>>>>> >>   axis([ 0 0.4])       % zooming time domain to see that
>> >>> >>>>>> >> period=1/nu
>> >>> >>>>>> >>   subplot(1,2,2)
>> >>> >>>>>> >>   plot(freq, freqsig)
>> >>> >>>>>> >>   #---------- end code ------------
>> >>> >>>>>> >>
>> >>> >>>>>> >> Do it yourself. Just rerun the code trying different values
>> of
>> >>> >>>>>> >> FsFactor (eg: 2, 4, 8).
>> >>> >>>>>> >> Thx a lot for any hint!!!
>> >>> >>>>>> >>
>> >>> >>>>>> >> oxy
>> >>> >>>>>> >>
>> >>> >>>>>> >> ps: cross post
>> >>> >>>>>> >>
>> >>> >>>>>>
>> >>>
>> http://www.mathworks.com/matlabcentral/answers/115700-fft-why-signal-frequency-depends-on-sampling-rate
>> >>> >>>>>> >> _______________________________________________
>> >>> >>>>>> >> Help-octave mailing list
>> >>> >>>>>> >> [hidden email]
>> >>> >>>>>> >> https://mailman.cae.wisc.edu/listinfo/help-octave
>> >>> >>>>>> >>
>> >>> >>>>>> >
>> >>> >>>>>> _______________________________________________
>> >>> >>>>>> Help-octave mailing list
>> >>> >>>>>> [hidden email]
>> >>> >>>>>> https://mailman.cae.wisc.edu/listinfo/help-octave
>> >>> >>>>>>
>> >>> >>>>>
>> >>> >>>>
>> >>> >>>
>> >>> >>
>> >>> >
>> >>>
>> >>
>> >
>> _______________________________________________
>> Help-octave mailing list
>> [hidden email]
>> https://mailman.cae.wisc.edu/listinfo/help-octave
>>
>
_______________________________________________
Help-octave mailing list
[hidden email]
https://mailman.cae.wisc.edu/listinfo/help-octave