FFT phase

Started by john john December 31, 2006
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
Ron N. wrote:
> Jerry Avins wrote:
...
>> Why are you busting my chops? > > Not you in particular, just the idea in general that sampling and > bandlimiting are related by *necessity*. They aren't. It's just > that samples are usually worthless (or worth less) unless one > does some sort of suitable low pass (or bandpass) filtering, > which a beginner, perhaps such as the OP, doesn't know about > or has forgotten.
Actually, my point was different. Without prior bandpass filtering, aliasing corrupts the signal; agreed. My point was that after sampling, there is only this signal (and its images) within the band allowed by the sample rate. For a baseband signal, only components up to Fs/2 exist in the samples. These are replicated with or without inversion on up the spectrum. Whichever image of an aliased signal one chooses, it has all of the aliases that the baseband image has. The act of sampling converts out-of-band components to in-band aliases. After sampling, there are no out-of-band components; they've been converted to in-band components. I call that bandlimiting. If you have a better term, I'd like to adopt it. Jerry -- Engineering is the art of making what you want from things you can get. ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯
Jerry Avins wrote:
> Ron N. wrote: > > Jerry Avins wrote: > >> Ron N. wrote: > >>> Jerry Avins wrote: > >>>> The spectra are easy to compute *for continuous waveforms*. They are > >>>> tabulated in many books. You cannot deal with continuous waveforms using > >>>> a computer, and if you sample the waveform you bandlimit it *whether you > >>>> intend to or not*. > >>> Sampling does not automatically bandlimit. You can undersample > >>> for some demodulation schemes, for instance. The sampler does > >>> not know which band you want. You have to pick a band filter > >>> (low pass, or band pass for your selected frequency range). > >> Sampling aliases all the frequencies above Fs/2 to other frequencies > >> below Fs/2. Any proper reconstruction procedure will generate a > >> continuous waveform from those samples all of whose components are below > >> Fs/2. I call that bandlimiting. What do you call it? > > > > Improper reconstruction. If your signal is above Fs/2 (but narrow- > > band enough not to be aliased with itself) then you should use a > > reconstruction formula for that higher frequency band. > > > > The reconstruction below Fs/2 is only valid if you have out-of-band > > information telling you that that is where your signal belongs. > > Otherwise any reconstruction is ambiguous. > > Do you really think the OP was into subband sampling? It was clearly a > question about baseband and I think you know that too.
Neither. I wouldn't be surprised if the OP was dealing with a non-bandlimited signal (which means a baseband reconstruction will be incorrect, or inexact). Sampling doesn't fix that.
> Why are you busting my chops?
Not you in particular, just the idea in general that sampling and bandlimiting are related by *necessity*. They aren't. It's just that samples are usually worthless (or worth less) unless one does some sort of suitable low pass (or bandpass) filtering, which a beginner, perhaps such as the OP, doesn't know about or has forgotten. IMHO. YMMV. -- rhn A.T nicholson d.0.t C-o-M
Ron N. wrote:
> Jerry Avins wrote: >> Ron N. wrote: >>> Jerry Avins wrote: >>>> The spectra are easy to compute *for continuous waveforms*. They are >>>> tabulated in many books. You cannot deal with continuous waveforms using >>>> a computer, and if you sample the waveform you bandlimit it *whether you >>>> intend to or not*. >>> Sampling does not automatically bandlimit. You can undersample >>> for some demodulation schemes, for instance. The sampler does >>> not know which band you want. You have to pick a band filter >>> (low pass, or band pass for your selected frequency range). >> Sampling aliases all the frequencies above Fs/2 to other frequencies >> below Fs/2. Any proper reconstruction procedure will generate a >> continuous waveform from those samples all of whose components are below >> Fs/2. I call that bandlimiting. What do you call it? > > Improper reconstruction. If your signal is above Fs/2 (but narrow- > band enough not to be aliased with itself) then you should use a > reconstruction formula for that higher frequency band. > > The reconstruction below Fs/2 is only valid if you have out-of-band > information telling you that that is where your signal belongs. > Otherwise any reconstruction is ambiguous.
Do you really think the OP was into subband sampling? It was clearly a question about baseband and I think you know that too. Why are you busting my chops? Jerry -- Engineering is the art of making what you want from things you can get. ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯
Jerry Avins wrote:
> Ron N. wrote: > > Jerry Avins wrote: > >> The spectra are easy to compute *for continuous waveforms*. They are > >> tabulated in many books. You cannot deal with continuous waveforms using > >> a computer, and if you sample the waveform you bandlimit it *whether you > >> intend to or not*. > > > > Sampling does not automatically bandlimit. You can undersample > > for some demodulation schemes, for instance. The sampler does > > not know which band you want. You have to pick a band filter > > (low pass, or band pass for your selected frequency range). > > Sampling aliases all the frequencies above Fs/2 to other frequencies > below Fs/2. Any proper reconstruction procedure will generate a > continuous waveform from those samples all of whose components are below > Fs/2. I call that bandlimiting. What do you call it?
Improper reconstruction. If your signal is above Fs/2 (but narrow- band enough not to be aliased with itself) then you should use a reconstruction formula for that higher frequency band. The reconstruction below Fs/2 is only valid if you have out-of-band information telling you that that is where your signal belongs. Otherwise any reconstruction is ambiguous. IMHO. YMMV. -- rhn A.T nicholson d.0.t C-o-M
Ron N. wrote:
> Jerry Avins wrote: > > > No frequencies above half the sample rate exist in a > > set of samples. > > Of course they exist. It's just that if they are there *and* the > baseband spectrum is also there, then they are aliased together > and you can't tell them apart (without some other external > information).
dunno if i agree with this. we have to be specific about what we mean. by set of samples, i mean { x[n]: n integer, -inf < n < +inf } in this set there is no concept of inbetween integer n. x[r] is not defined unless r is an integer. this can be represented with a fourier transform (the DTFT) which is the Z transform of x[n] X(z) = SUM x[n] z^(-n) n with z = e^(jw) for real w . this X(e^(jw)) is periodic with period 2*pi. normally we deal with -pi < w < +pi. the repeated *images* are not *aliases*. they are images and simply from X(z) or from x[n] we don't know whether or not there was aliasing. strictly speaking, that information is gone and lacking any outside info, i think you normally have to assume that a frequency component that exists was originally between -pi and +pi. that's why they're called "aliases", they were originally outside that range and masquarade as if they were always between -pi and +pi. in this sense, Jerry is right, i think. and if you make the hard-core assumption that x[n] is strictly associated with the bandlimited reconstruction: x(t) = SUM x[n] sinc(t-n) n then there are explicitly no (radian) frequencies above pi in magnitude. if you instead say x(t) = SUM x[n] delta(t-n) n there are images, but no aliases (unless you have outside information that tells you something about x(t) before it was sampled). r b-j
Ron N. wrote:
> Jerry Avins wrote: >> The spectra are easy to compute *for continuous waveforms*. They are >> tabulated in many books. You cannot deal with continuous waveforms using >> a computer, > > You forget about symbolic math programs, which can deal with some > continuous waveforms if they are functions of a form that the symbolic > programs can handle (symbolic integration, convolution, etc.) > A computer can also use the book tabulations.
I didn't forget. If I have to touch all bases and enumerate all possible exceptions in every discussion, I won't be able to discuss much. Jerry -- Engineering is the art of making what you want from things you can get. ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯
Ron N. wrote:
> Jerry Avins wrote: >> The spectra are easy to compute *for continuous waveforms*. They are >> tabulated in many books. You cannot deal with continuous waveforms using >> a computer, and if you sample the waveform you bandlimit it *whether you >> intend to or not*. > > Sampling does not automatically bandlimit. You can undersample > for some demodulation schemes, for instance. The sampler does > not know which band you want. You have to pick a band filter > (low pass, or band pass for your selected frequency range).
Sampling aliases all the frequencies above Fs/2 to other frequencies below Fs/2. Any proper reconstruction procedure will generate a continuous waveform from those samples all of whose components are below Fs/2. I call that bandlimiting. What do you call it?
>> No frequencies above half the sample rate exist in a >> set of samples. > > Of course they exist. It's just that if they are there *and* the > baseband spectrum is also there, then they are aliased together > and you can't tell them apart (without some other external > information).
In general, they can't be separated even with special information. It's hard to unscramble eggs or sounds.
> Actually, in practice, since low pass filters only approximate > true bandlimiting filters of infinite width, they are always there, > just hopefully below the desired noise floor.
So is thermal noise. We live with what we have to. Jerry -- Engineering is the art of making what you want from things you can get. ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯
Jerry Avins wrote:
> The spectra are easy to compute *for continuous waveforms*. They are > tabulated in many books. You cannot deal with continuous waveforms using > a computer,
You forget about symbolic math programs, which can deal with some continuous waveforms if they are functions of a form that the symbolic programs can handle (symbolic integration, convolution, etc.) A computer can also use the book tabulations. IMHO. YMMV. -- rhn A.T nicholson d.0.t C-o-M
Jerry Avins wrote:
> The spectra are easy to compute *for continuous waveforms*. They are > tabulated in many books. You cannot deal with continuous waveforms using > a computer, and if you sample the waveform you bandlimit it *whether you > intend to or not*.
Sampling does not automatically bandlimit. You can undersample for some demodulation schemes, for instance. The sampler does not know which band you want. You have to pick a band filter (low pass, or band pass for your selected frequency range).
> No frequencies above half the sample rate exist in a > set of samples.
Of course they exist. It's just that if they are there *and* the baseband spectrum is also there, then they are aliased together and you can't tell them apart (without some other external information). Actually, in practice, since low pass filters only approximate true bandlimiting filters of infinite width, they are always there, just hopefully below the desired noise floor. IMHO. YMMV. -- rhn A.T nicholson d.0.t C-o-M
john john wrote:
> Jerry Avins ha scritto: > >> No. The spectra of square and triangular waves are infinite series. >> Computers can't deal with infinities. If you truncate the series, the >> waveforms you get might surprise you. >> >> What are you trying to see?
> I know about infinite series and I don't cut the series. I'm just > trying to see the frequency spectrum of a waveform but every frequency > with its correct phase. I won't the infinite but just the firsts ten > frequency.
The spectra are easy to compute *for continuous waveforms*. They are tabulated in many books. You cannot deal with continuous waveforms using a computer, and if you sample the waveform you bandlimit it *whether you intend to or not*. No frequencies above half the sample rate exist in a set of samples. Jerry -- Engineering is the art of making what you want from things you can get. ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯