Thread Subject: Understanding phase in the FFT

Hi,

Here the simple code of three series, which differ only by phase. I don't understand why I am not getting phases,which I expect to. Shouldn't I get 0 (or 360), 90 and 45 ? Thanks a lot!

fs = 10;
t = 0:1/fs:50-1/fs;

x1 = sin(2*pi*t);
x2 = sin(2*pi*t+pi/2);
x3 = sin(2*pi*t+pi/4);

m = length(x1); % Window length
n = pow2(nextpow2(m)); % Transform length

y1 = fft(x1,n);
y2 = fft(x2,n);
y3 = fft(x3,n);

FT_power1 = abs(y1(1:floor(n/2))).^2;
FT_phase1= (unwrap(angle(y1(1:floor(n/2))))) * 180/pi;
[c1,i1] = max(FT_power1);
phase(1) = FT_phase1(i1);

FT_power2 = abs(y2(1:floor(n/2))).^2;
FT_phase2= (unwrap(angle(y2(1:floor(n/2))))) * 180/pi;
[c2,i2] = max(FT_power2);
phase(2) = FT_phase2(i2);

FT_power3 = abs(y3(1:floor(n/2))).^2;
FT_phase3= (unwrap(angle(y3(1:floor(n/2))))) * 180/pi;
[c3,i3] = max(FT_power3);
phase(3) = FT_phase3(i3);

phase =

  305.1562 35.0154 -10.0109

On Jul 30, 9:35 am, "aaaa Vadimov" <a...@hotmail.com> wrote:
> Hi,
>
> Here the simple code of three series, which differ only by phase. I don't understand why I am not getting phases,which I expect to. Shouldn't I get 0 (or 360), 90 and 45 ? Thanks a lot!
>
> fs = 10;                        
> t = 0:1/fs:50-1/fs;              
>
> x1 = sin(2*pi*t);
> x2 = sin(2*pi*t+pi/2);
> x3 = sin(2*pi*t+pi/4);
>
> m = length(x1);          % Window length
> n = pow2(nextpow2(m));  % Transform length
>
> y1 = fft(x1,n);
> y2 = fft(x2,n);
> y3 = fft(x3,n);
>
> FT_power1 = abs(y1(1:floor(n/2))).^2;
> FT_phase1= (unwrap(angle(y1(1:floor(n/2))))) * 180/pi;
> [c1,i1] = max(FT_power1);
> phase(1) = FT_phase1(i1);
>
> FT_power2 = abs(y2(1:floor(n/2))).^2;
> FT_phase2= (unwrap(angle(y2(1:floor(n/2))))) * 180/pi;
> [c2,i2] = max(FT_power2);
> phase(2) = FT_phase2(i2);
>
> FT_power3 = abs(y3(1:floor(n/2))).^2;
> FT_phase3= (unwrap(angle(y3(1:floor(n/2))))) * 180/pi;
> [c3,i3] = max(FT_power3);
> phase(3) = FT_phase3(i3);
>
> phase =
>
>   305.1562   35.0154  -10.0109

It's called leakage.
Because n is not an integral multiple of 1/fs, the energy of each sine
wave is not in a single component, but spread over adjacent
components, so when you find the maximum, it is not at exactly fs.
You can fix this simply by using m instead of n in the FFT.
These days, there's no need to use dyadic numbers for FFT.

On Jul 30, 9:45 am, TideMan <mul...@gmail.com> wrote:
> On Jul 30, 9:35 am, "aaaa Vadimov" <a...@hotmail.com> wrote:
>
>
>
> > Hi,
>
> > Here the simple code of three series, which differ only by phase. I don't understand why I am not getting phases,which I expect to. Shouldn't I get 0 (or 360), 90 and 45 ? Thanks a lot!
>
> > fs = 10;                        
> > t = 0:1/fs:50-1/fs;              
>
> > x1 = sin(2*pi*t);
> > x2 = sin(2*pi*t+pi/2);
> > x3 = sin(2*pi*t+pi/4);
>
> > m = length(x1);          % Window length
> > n = pow2(nextpow2(m));  % Transform length
>
> > y1 = fft(x1,n);
> > y2 = fft(x2,n);
> > y3 = fft(x3,n);
>
> > FT_power1 = abs(y1(1:floor(n/2))).^2;
> > FT_phase1= (unwrap(angle(y1(1:floor(n/2))))) * 180/pi;
> > [c1,i1] = max(FT_power1);
> > phase(1) = FT_phase1(i1);
>
> > FT_power2 = abs(y2(1:floor(n/2))).^2;
> > FT_phase2= (unwrap(angle(y2(1:floor(n/2))))) * 180/pi;
> > [c2,i2] = max(FT_power2);
> > phase(2) = FT_phase2(i2);
>
> > FT_power3 = abs(y3(1:floor(n/2))).^2;
> > FT_phase3= (unwrap(angle(y3(1:floor(n/2))))) * 180/pi;
> > [c3,i3] = max(FT_power3);
> > phase(3) = FT_phase3(i3);
>
> > phase =
>
> >   305.1562   35.0154  -10.0109
>
> It's called leakage.
> Because n is not an integral multiple of 1/fs, the energy of each sine
> wave is not in a single component, but spread over adjacent
> components, so when you find the maximum, it is not at exactly fs.
> You can fix this simply by using m instead of n in the FFT.
> These days, there's no need to use dyadic numbers for FFT.

I should add that if you're interested in calculating the phase of a
sinusoidal signal of a particular frequency, using FFT is completely
the wrong way to go about it.
We do this all the time in tidal analysis - calculating the amplitude
and phase of various tidal constituents - and without almost exception
it is done in the time domain using Matlab's \ facility, or similar.

On Jul 29, 2:35 pm, "aaaa Vadimov" <a...@hotmail.com> wrote:
> Hi,
>
> Here the simple code of three series, which differ only by phase. I don't understand why I am not getting phases,which I expect to. Shouldn't I get 0 (or 360), 90 and 45 ? Thanks a lot!
>
> fs = 10;
> t = 0:1/fs:50-1/fs;
>
> x1 = sin(2*pi*t);
> x2 = sin(2*pi*t+pi/2);
> x3 = sin(2*pi*t+pi/4);
> ...

Why do you expect a sine wave to give a phase of 0 with an fft()? The
basis vectors are cos()+i*sin().

Dale B. Dalrymple

dbd <dbd@ieee.org> wrote in message <8ea06202-cfd1-44cf-8656-c6c03daa210d@q40g2000prh.googlegroups.com>...
> On Jul 29, 2:35 pm, "aaaa Vadimov" <a...@hotmail.com> wrote:
> > Hi,
> >
> > Here the simple code of three series, which differ only by phase. I don't understand why I am not getting phases,which I expect to. Shouldn't I get 0 (or 360), 90 and 45 ? Thanks a lot!
> >
> > fs = 10;
> > t = 0:1/fs:50-1/fs;
> >
> > x1 = sin(2*pi*t);
> > x2 = sin(2*pi*t+pi/2);
> > x3 = sin(2*pi*t+pi/4);
> > ...
>
> Why do you expect a sine wave to give a phase of 0 with an fft()? The
> basis vectors are cos()+i*sin().

Phase is relative. The phase between two signals is obtained from the cross spectrum. The phase of a single signal isn't that meaningful or useful.

On Jul 30, 4:50 pm, "Steve Amphlett" <Firstname.Lastn...@Where-I-
Work.com> wrote:
> dbd <d...@ieee.org> wrote in message <8ea06202-cfd1-44cf-8656-c6c03daa2...@q40g2000prh.googlegroups.com>...
> > On Jul 29, 2:35 pm, "aaaa Vadimov" <a...@hotmail.com> wrote:
> > > Hi,
>
> > > Here the simple code of three series, which differ only by phase. I don't understand why I am not getting phases,which I expect to. Shouldn't I get 0 (or 360), 90 and 45 ? Thanks a lot!
>
> > > fs = 10;
> > > t = 0:1/fs:50-1/fs;
>
> > > x1 = sin(2*pi*t);
> > > x2 = sin(2*pi*t+pi/2);
> > > x3 = sin(2*pi*t+pi/4);
> > > ...
>
> > Why do you expect a sine wave to give a phase of 0 with an fft()? The
> > basis vectors are cos()+i*sin().
>
> Phase is relative.  The phase between two signals is obtained from the cross spectrum.  The phase of a single signal isn't that meaningful or useful.

No it's not. It's absolute.
In the time domain, it defines the time when the cosine is a maximum.
In the frequency domain, it defines what proportion of the Fourier
transform is imaginary.

Thanks for all your answers.
Actually, I don't need to make FFT of sine function, it was just a test.
What I actually have are two signals, which have the maximal power at some frequency. I preselected those signals and this is the only significant power in them. However, when I calculate the correlation between these two signals the correlation is around "0". Would it be correct to expect to get difference between phases around 180 degrees? Would it be just abs(phase(2)-phase(1)), if I am using my previous example? If somebody can post some simple example of correct phase usage, I would really appreciate it.

TideMan <mulgor@gmail.com> wrote in message <a4004cfa-186e-43c7-ac6f-559fd9ef78c8@12g2000pri.googlegroups.com>...
> On Jul 30, 4:50?pm, "Steve Amphlett" <Firstname.Lastn...@Where-I-
> Work.com> wrote:
> > dbd <d...@ieee.org> wrote in message <8ea06202-cfd1-44cf-8656-c6c03daa2...@q40g2000prh.googlegroups.com>...
> > > On Jul 29, 2:35 pm, "aaaa Vadimov" <a...@hotmail.com> wrote:
> > > > Hi,
> >
> > > > Here the simple code of three series, which differ only by phase. I don't understand why I am not getting phases,which I expect to. Shouldn't I get 0 (or 360), 90 and 45 ? Thanks a lot!
> >
> > > > fs = 10;
> > > > t = 0:1/fs:50-1/fs;
> >
> > > > x1 = sin(2*pi*t);
> > > > x2 = sin(2*pi*t+pi/2);
> > > > x3 = sin(2*pi*t+pi/4);
> > > > ...
> >
> > > Why do you expect a sine wave to give a phase of 0 with an fft()? The
> > > basis vectors are cos()+i*sin().
> >
> > Phase is relative. ?The phase between two signals is obtained from the cross spectrum. ?The phase of a single signal isn't that meaningful or useful.
>
> No it's not. It's absolute.
> In the time domain, it defines the time when the cosine is a maximum.
> In the frequency domain, it defines what proportion of the Fourier
> transform is imaginary.

TideMan <mulgor@gmail.com> wrote in message <a4004cfa-186e-43c7-ac6f-559fd9ef78c8@12g2000pri.googlegroups.com>...
> On Jul 30, 4:50?pm, "Steve Amphlett" <Firstname.Lastn...@Where-I-
> Work.com> wrote:
> > dbd <d...@ieee.org> wrote in message <8ea06202-cfd1-44cf-8656-c6c03daa2...@q40g2000prh.googlegroups.com>...
> > > On Jul 29, 2:35 pm, "aaaa Vadimov" <a...@hotmail.com> wrote:
> > > > Hi,
> >
> > > > Here the simple code of three series, which differ only by phase. I don't understand why I am not getting phases,which I expect to. Shouldn't I get 0 (or 360), 90 and 45 ? Thanks a lot!
> >
> > > > fs = 10;
> > > > t = 0:1/fs:50-1/fs;
> >
> > > > x1 = sin(2*pi*t);
> > > > x2 = sin(2*pi*t+pi/2);
> > > > x3 = sin(2*pi*t+pi/4);
> > > > ...
> >
> > > Why do you expect a sine wave to give a phase of 0 with an fft()? The
> > > basis vectors are cos()+i*sin().
> >
> > Phase is relative. ?The phase between two signals is obtained from the cross spectrum. ?The phase of a single signal isn't that meaningful or useful.
>
> No it's not. It's absolute.
> In the time domain, it defines the time when the cosine is a maximum.
> In the frequency domain, it defines what proportion of the Fourier
> transform is imaginary.

Ok, phase is useful when it's relative, for example, when there is a common reference involved.

TideMan <mulgor@gmail.com> wrote in message <917fec82-e97e-4228-9349-a1053f19bc3b@z4g2000prh.googlegroups.com>...
> On Jul 30, 9:45 am, TideMan <mul...@gmail.com> wrote:
> > On Jul 30, 9:35 am, "aaaa Vadimov" <a...@hotmail.com> wrote:
> >
> >
> >
> > > Hi,
> >
> > > Here the simple code of three series, which differ only by phase. I don't understand why I am not getting phases,which I expect to. Shouldn't I get 0 (or 360), 90 and 45 ? Thanks a lot!
> >
> > > fs = 10;                        
> > > t = 0:1/fs:50-1/fs;              
> >
> > > x1 = sin(2*pi*t);
> > > x2 = sin(2*pi*t+pi/2);
> > > x3 = sin(2*pi*t+pi/4);
> >
> > > m = length(x1);          % Window length
> > > n = pow2(nextpow2(m));  % Transform length
> >
> > > y1 = fft(x1,n);
> > > y2 = fft(x2,n);
> > > y3 = fft(x3,n);
> >
> > > FT_power1 = abs(y1(1:floor(n/2))).^2;
> > > FT_phase1= (unwrap(angle(y1(1:floor(n/2))))) * 180/pi;
> > > [c1,i1] = max(FT_power1);
> > > phase(1) = FT_phase1(i1);
> >
> > > FT_power2 = abs(y2(1:floor(n/2))).^2;
> > > FT_phase2= (unwrap(angle(y2(1:floor(n/2))))) * 180/pi;
> > > [c2,i2] = max(FT_power2);
> > > phase(2) = FT_phase2(i2);
> >
> > > FT_power3 = abs(y3(1:floor(n/2))).^2;
> > > FT_phase3= (unwrap(angle(y3(1:floor(n/2))))) * 180/pi;
> > > [c3,i3] = max(FT_power3);
> > > phase(3) = FT_phase3(i3);
> >
> > > phase =
> >
> > >   305.1562   35.0154  -10.0109
> >
> > It's called leakage.
> > Because n is not an integral multiple of 1/fs, the energy of each sine
> > wave is not in a single component, but spread over adjacent
> > components, so when you find the maximum, it is not at exactly fs.
> > You can fix this simply by using m instead of n in the FFT.
> > These days, there's no need to use dyadic numbers for FFT.
>
> I should add that if you're interested in calculating the phase of a
> sinusoidal signal of a particular frequency, using FFT is completely
> the wrong way to go about it.
> We do this all the time in tidal analysis - calculating the amplitude
> and phase of various tidal constituents - and without almost exception
> it is done in the time domain using Matlab's \ facility, or similar.

Dear TideMan,
I'm in wave energy and desperately want to calculate the phase and amplitude of each frequency in the spectrum of measured wave elevation. This sounds like exactly what you've done with your tidal analysis. As you imply, the fft approach gives rubbish phase information. Would you be kind enough to provide some more detailed information on how you actually do this? Maybe a link to some documentation or even the code itself :)
Thanks in advance,
One grateful WaveMan...

On Aug 17, 10:32 pm, "Chris " <chris.signore...@wavebob.com> wrote:
> TideMan <mul...@gmail.com> wrote in message <917fec82-e97e-4228-9349-a1053f19b...@z4g2000prh.googlegroups.com>...
> > On Jul 30, 9:45 am, TideMan <mul...@gmail.com> wrote:
> > > On Jul 30, 9:35 am, "aaaa Vadimov" <a...@hotmail.com> wrote:
>
> > > > Hi,
>
> > > > Here the simple code of three series, which differ only by phase. I don't understand why I am not getting phases,which I expect to. Shouldn't I get 0 (or 360), 90 and 45 ? Thanks a lot!
>
> > > > fs = 10;                        
> > > > t = 0:1/fs:50-1/fs;              
>
> > > > x1 = sin(2*pi*t);
> > > > x2 = sin(2*pi*t+pi/2);
> > > > x3 = sin(2*pi*t+pi/4);
>
> > > > m = length(x1);          % Window length
> > > > n = pow2(nextpow2(m));  % Transform length
>
> > > > y1 = fft(x1,n);
> > > > y2 = fft(x2,n);
> > > > y3 = fft(x3,n);
>
> > > > FT_power1 = abs(y1(1:floor(n/2))).^2;
> > > > FT_phase1= (unwrap(angle(y1(1:floor(n/2))))) * 180/pi;
> > > > [c1,i1] = max(FT_power1);
> > > > phase(1) = FT_phase1(i1);
>
> > > > FT_power2 = abs(y2(1:floor(n/2))).^2;
> > > > FT_phase2= (unwrap(angle(y2(1:floor(n/2))))) * 180/pi;
> > > > [c2,i2] = max(FT_power2);
> > > > phase(2) = FT_phase2(i2);
>
> > > > FT_power3 = abs(y3(1:floor(n/2))).^2;
> > > > FT_phase3= (unwrap(angle(y3(1:floor(n/2))))) * 180/pi;
> > > > [c3,i3] = max(FT_power3);
> > > > phase(3) = FT_phase3(i3);
>
> > > > phase =
>
> > > >   305.1562   35.0154  -10.0109
>
> > > It's called leakage.
> > > Because n is not an integral multiple of 1/fs, the energy of each sine
> > > wave is not in a single component, but spread over adjacent
> > > components, so when you find the maximum, it is not at exactly fs.
> > > You can fix this simply by using m instead of n in the FFT.
> > > These days, there's no need to use dyadic numbers for FFT.
>
> > I should add that if you're interested in calculating the phase of a
> > sinusoidal signal of a particular frequency, using FFT is completely
> > the wrong way to go about it.
> > We do this all the time in tidal analysis - calculating the amplitude
> > and phase of various tidal constituents - and without almost exception
> > it is done in the time domain using Matlab's \ facility, or similar.
>
> Dear TideMan,
> I'm in wave energy and desperately want to calculate the phase and amplitude of each frequency in the spectrum of measured wave elevation. This sounds like exactly what you've done with your tidal analysis. As you imply, the fft approach gives rubbish phase information. Would you be kind enough to provide some more detailed information on how you actually do this? Maybe a link to some documentation or even the code itself :)
> Thanks in advance,
> One grateful WaveMan...

Google t_tide