FFT phase

Calculating the FFT of a waveform we have img and real values and is
possible to find magnitude and phase of sines.

magnitude = (sqrt(real^2 +img^2))/float(number_of_samples);

phase = (atan2(img/real) * 180) / PI;

I see that the phase is wrong. In a web page I see:

"Each of the harmonic frequencies is defined by a magnitude (amplitude)
and a phase. The phase indicates how to shift the harmonic before
adding it to the sum. The phase information can be difficult to
interpret and its use is restricted to a few specialized applications.
The rest of this document is concerned only with harmonic magnitudes"

WOW

It's so difficult to calculate the correct phase of a sine? What's the
right way?

john john wrote:
> Calculating the FFT of a waveform we have img and real values and is
> possible to find magnitude and phase of sines.
>
> magnitude = (sqrt(real^2 +img^2))/float(number_of_samples);
>
> phase = (atan2(img/real) * 180) / PI;
>
> I see that the phase is wrong. In a web page I see:
>
> "Each of the harmonic frequencies is defined by a magnitude (amplitude)
> and a phase. The phase indicates how to shift the harmonic before
> adding it to the sum. The phase information can be difficult to
> interpret and its use is restricted to a few specialized applications.
> The rest of this document is concerned only with harmonic magnitudes"
>
> WOW
>
> It's so difficult to calculate the correct phase of a sine? What's the
> right way?


Hello John,

I noticed you are doing an atan2() but with just one argument. Normaly
it is atan2(x,y). To help you track down your problems, try just
feeding a cosine's data to your FFT. Pick your frequency so the result
is in just one bin.  I.e., if you are using N samples then let your
data be  cos(2*pi*w*x/N) where x goes from 0,1,2,3,...,N-1 and w is an
integer. Then your FFT should just give you purely real data (the im
part will be essentially zero) and calculate your phase. Then use sin()
instead of cosine and get purely imag. data.

Clay

> I noticed you are doing an atan2() but with just one argument. Normaly
> it is atan2(x,y). To help you track down your problems, try just
> feeding a cosine's data to your FFT. Pick your frequency so the result
> is in just one bin.  I.e., if you are using N samples then let your
> data be  cos(2*pi*w*x/N) where x goes from 0,1,2,3,...,N-1 and w is an
> integer. Then your FFT should just give you purely real data (the im
> part will be essentially zero) and calculate your phase. Then use sin()
> instead of cosine and get purely imag. data.
>
> Clay

Yes I use atan2(x,y). Your suggestion make me in a wrong way becouse if
I do what you say the results from FFT function are erroneus and the
img data are always different than zero. I implement the FFT function
explained in this page:

http://gemic.e-technik.uni-ulm.de/lehre/basic_mathematics/fourier/node7.php3

have you any idea about how to calculate the right phase and amplitude
using real and img data from this code? Others codes are ok too.

For a real FFT :

magnitude = (sqrt(real^2 +img^2))/(FFT_Size/2);

Value will be phasor length i.e. Peak. Remember, this will be modified if
window function used.

phase = atan2(img/real);

Value will be +/- PI

Jeff

john john wrote:
> Yes I use atan2(x,y). Your suggestion make me in a wrong way becouse if
> I do what you say the results from FFT function are erroneus and the
> img data are always different than zero. I implement the FFT function
> explained in this page:
>
> http://gemic.e-technik.uni-ulm.de/lehre/basic_mathematics/fourier/node7.php3
>
> have you any idea about how to calculate the right phase and amplitude
> using real and img data from this code? Others codes are ok too.

one thing to think about is both the concepts of "negative frequency"
in the FFT output as well as "negative time" in the FFT input (this is
what the MATLAB function "fftshift()" is for).

if N is even (it should be a power of 2 to make your FFT truly "F"),
the N/2 bin is ambiguously +N/2 or -N/2 and usually gets assigned to
-N/2.  but the latter half, from N/2+1 up to N-1 are almost
consistently thought of as negative frequency in the frequency domain
or negative time in the time domain.  remember that the FFT or DFT is
an inherently circular operation and you want to treat x[N-1] as x[-1],
x[N-2] as x[-2], and so on until you get to x[N-N/2] or x[N/2] which is
the ambiguous bin but usually interpreted as x[-N/2].

also window your sinusoid nicely, center it at x[0] (putting the
negative indexed stuff in the upper half of the FFT buffer) and your
magnitude and phase will come out looking okay.  (if you don't do it
right, there will be a factor of +/- 180 degrees added to alternating
bins,which looks ugly.)

r b-j

Guys,

Much to my annoyance over the years, the convention for atan2()
actually seems to be atan2(y,x) rather than atan2(x,y). Makes a
difference.

Dirk

Dirk Bell
DSP Consultant


john john wrote:
> > I noticed you are doing an atan2() but with just one argument. Normaly
> > it is atan2(x,y). To help you track down your problems, try just
> > feeding a cosine's data to your FFT. Pick your frequency so the result
> > is in just one bin.  I.e., if you are using N samples then let your
> > data be  cos(2*pi*w*x/N) where x goes from 0,1,2,3,...,N-1 and w is an
> > integer. Then your FFT should just give you purely real data (the im
> > part will be essentially zero) and calculate your phase. Then use sin()
> > instead of cosine and get purely imag. data.
> >
> > Clay
>
> Yes I use atan2(x,y). Your suggestion make me in a wrong way becouse if
> I do what you say the results from FFT function are erroneus and the
> img data are always different than zero. I implement the FFT function
> explained in this page:
>
> http://gemic.e-technik.uni-ulm.de/lehre/basic_mathematics/fourier/node7.php3
>
> have you any idea about how to calculate the right phase and amplitude
> using real and img data from this code? Others codes are ok too.

dbell wrote:
> Guys,
> 
> Much to my annoyance over the years, the convention for atan2()
> actually seems to be atan2(y,x) rather than atan2(x,y). Makes a
> difference.

Why should that annoy you? Viven V = x + jy, |V| = sqrt(x^2 + y^2) and 
arg(V) = arctan(y/z). Opposite over adjacent.

Jerry
-- 
Engineering is the art of making what you want from things you can get.
�����������������������������������������������������������������������

robert bristow-johnson wrote:

> one thing to think about is both the concepts of "negative frequency"
> in the FFT output as well as "negative time" in the FFT input (this is
> what the MATLAB function "fftshift()" is for).

I don't know Matlab. Does fftshift() give the same result as inverting 
the sign of all the odd-numbered samples?

Jerry
-- 
Engineering is the art of making what you want from things you can get.
�����������������������������������������������������������������������

dbell wrote:
 > Guys,
 >
 > Much to my annoyance over the years, the convention for atan2()
 > actually seems to be atan2(y,x) rather than atan2(x,y). Makes a
 > difference.

Why should that annoy you? Given V = x + jy, |V| = sqrt(x2 + y2) and 
arg(V) = arctan(y/z). Opposite over adjacent.

Jerry
-- 
Engineering is the art of making what you want from things you can get.
�����������������������������������������������������������������������

Jerry Avins skrev:
> robert bristow-johnson wrote:
>
> > one thing to think about is both the concepts of "negative frequency"
> > in the FFT output as well as "negative time" in the FFT input (this is
> > what the MATLAB function "fftshift()" is for).
>
> I don't know Matlab. Does fftshift() give the same result as inverting
> the sign of all the odd-numbered samples?

I don't know what effects that sign inversion might have, but I do
know that FFTSHIFT re-arranges the output of the FFT function.

X = fft(x);
XX = fftshift(X);

After the first line, X(1) contains the DC component of the spectrum
of x (matlab start indexing at 1, not 0). After the call to FFTSHIFT
the DC coefficient is located near N/2, N being the length of x, and
the spectrum XX is conjugate symmetric around the DC coefficient.

Rune