0

I'm trying to recover the phase of a simple (audio) signal in matlab: In matlab I do the following:

% This wave is perfectly periodic in the sample.  That is, 
% there are exactly 1000 periods.
swave = sin( 2 * pi * (0:10000) * 441/44100);

% Find the fft
sFFT = fft(swave);
% Remove the duplicate data in the FFT
sFFT = sFFT(1:length(sFFT)/2);

% Take a look a the amplitudes from the FFT and it checks out
freqs = 44100/ 2*linspace(0,1,length(sFFT);
plot(freqs, abs(sFFT));

% Now to get the phase
plot(freqs, angle(sFFT));

This result makes almost no sense to me. Because this is a sin wave (not a cos wave). I expect to see 1/2*pi = 1.57079 for the value of the 441hz bin. Instead I see a nearly discontinuous jump from (441, -1.539) to (445, 1.603). Why is 441 so far from the correct value? Why is 445 so close?

The value for all of the bins besides 441 hz are a mystery to me. I've also tried several other methods of recovering the phase including unwrap(angle(sFFT)) and atan2(imag(sFFT), real(sFFT)); These change the output but also do not make any sense to me. Why are bins besides 441 any value but 0 (like the abs(FFT) shows?). Why is the 441 bin close but not the correct value?

Thanks for the help!

user2109
  • 1
  • 1
  • 1
  • 3
    did you try `unwrap` ? it corrects phase angles to produce smoother phase plots, see http://www.mathworks.com/help/matlab/ref/unwrap.html – bla Mar 20 '13 at 19:56

1 Answers1

0

To better estimate the phase, change the 0 phase reference to the center of the window by doing an fftshift, which should eliminate most of the alternating phase discontinuities, and then interpolate, but noting that you have changed your 0 point.

Note the the phase of a vector of length zero is meaningless in most cases, and thus anything close is likely just numerical noise (near zero vectors in various random directions) due to finite precision math.

Thus some people just clamp the plotted phase to zero for any magnitude below some noise floor.

hotpaw2
  • 70,107
  • 14
  • 90
  • 153
  • I have shifted the fft so that 0hz is centered on the graph; the left 1/2 of the graph is -frequencies and the right half is +frequencies. How does this help anything? What do you mean when you say that this "should eliminate most of the alternating phase discontinuities." Also, what interpolation are you suggesting I do? You suggestion to zero out phase measurements that are not above a certain threshold is helpful in removing most of the "noise" phase data. However, it does not help me understand the phase reading at 441hz which is what I'm most interested in. – user2109 Mar 25 '13 at 19:02
  • An fftshift shifts the data before the FFT or rotates every other complex bin after. That moves the reference for a phase of zero to the middle data sample where there is no discontinuity. Interpolation with a Sinc kernel (windowed) would be best. – hotpaw2 Mar 25 '13 at 20:03