MATLAB Central Newsreader - Understanding phase in the FFT

Understanding phase in the FFT

aaaa Vadimov — Wed, 29 Jul 2009 21:35:03 -0400

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

Re: Understanding phase in the FFT

TideMan — Wed, 29 Jul 2009 21:45:45 -0400

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.

Re: Understanding phase in the FFT

TideMan — Wed, 29 Jul 2009 22:16:12 -0400

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.

Re: Understanding phase in the FFT

dbd — Thu, 30 Jul 2009 03:37:31 -0400

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

Re: Understanding phase in the FFT

Steve Amphlett — Thu, 30 Jul 2009 04:50:04 -0400

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.

Re: Understanding phase in the FFT

TideMan — Thu, 30 Jul 2009 05:31:30 -0400

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.

Re: Understanding phase in the FFT

aaaa Vadimov — Thu, 30 Jul 2009 06:57:01 -0400

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.

Re: Understanding phase in the FFT

Steve Amphlett — Thu, 30 Jul 2009 07:02:01 -0400

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.

Re: Understanding phase in the FFT

Chris — Wed, 17 Aug 2011 10:32:09 -0400

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...

Re: Understanding phase in the FFT

TideMan — Wed, 17 Aug 2011 10:52:37 -0400

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