Dismiss Notice
Join Physics Forums Today!
The friendliest, high quality science and math community on the planet! Everyone who loves science is here!

Matlab - calculate phase shift using fft

  1. Apr 15, 2015 #1
    hello,

    I have 3 signals in the form of sampled values. When I plot them using

    plot (t,vPa,t,vPb,t,vPc)

    where vPa, vPb, vPc contains the values and t contains the sampling istants I get this:



    tttt.jpg


    when I calculate phase shift using fft I get phase angle = 0.
    I have tried using my own code ( which works normally on other vectors) and also a function downloaded from Mathworks and the result is always zero.

    The mathworks function is at the end.

    what am I doing wrong?

    Thanks
    Mattia


    Code (Matlab M):

    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    %              Phase Difference Measurement            %
    %               with MATLAB Implementation             %
    %                                                      %
    % Author: M.Sc. Eng. Hristo Zhivomirov        12/01/14 %
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


    function PhDiff = phdiffmeasure(x, y)

    % function: PhDiff = phdiffmeasure(x, y)
    % x - first signal in the time domain
    % y - second signal in the time domain
    % PhDiff - phase difference Y -> X, degrees

    % represent x as column-vector if it is not
    if size(x, 2) > 1
        x = x';
    end

    % represent y as column-vector if it is not
    if size(y, 2) > 1
        y = y';
    end

    % signals length
    xlen = length(x);
    ylen = length(y);

    % window preparation
    xwin = hanning(xlen, 'periodic');
    ywin = hanning(ylen, 'periodic');

    % fft of the first signal
    X = fft(x.*xwin);

    % fft of the second signal
    Y = fft(y.*ywin);

    % phase difference calculation
    [~, indx] = max(abs(X));
    [~, indy] = max(abs(Y));
    PhDiff = angle(Y(indy)) - angle(X(indx));
    PhDiff = PhDiff*180/pi;

    end
     
     
  2. jcsd
  3. Apr 15, 2015 #2

    jedishrfu

    Staff: Mentor

    I don't see what you did wrong.

    However, have you looked at the intermediate values for X and Y and the angle(X(indx)) values that are used to compute the phase difference?
     
  4. Apr 15, 2015 #3
    I have now, angle(Y(indy)) and angle(X(indx)) are both zero.

    I have attached the original 3 vectors.
     

    Attached Files:

    Last edited: Apr 15, 2015
  5. Apr 15, 2015 #4
    I have found a different way to calculate phase shift that can be used in this case:

    h1 = hilbert(vPa);

    h2 = hilbert(vPb);

    phase_rad = max(angle(h1 ./ h2))

    which gives back the correct value of 2/3 * pi.

    However I would be couriois to know why for this specific case fft didnt work. I thought it was a pretty generic method
     
  6. Apr 15, 2015 #5

    jedishrfu

    Staff: Mentor

    I meant the X and Y in the function you showed.

    Aren't they supposed to be mostly zeros with a few spikes?

    I think you need to look at each line of the function to decide if its doing the right thing.

    The code looks like it applies an FFT to the two input signals to produce a frequency spectrum. From that it determines the highest spike and from that subtracts them to determine the phase difference.
     
  7. Apr 15, 2015 #6
    yes that is what it does. Also my own code did that. ffts are mostly zero with one big spike on 50 Hz.
    I am positive both my own code and the one I posted are correct. becaouse I have used them many times and they returned correct results
     
  8. Apr 15, 2015 #7

    jedishrfu

    Staff: Mentor

    Then I guess there's something unique about your input data and the function you listed. Perhaps if you step through the function interactively you'll be able to spot the difference. It could be a number of elements limit or datatype mismatch...
     
  9. Apr 15, 2015 #8

    jedishrfu

    Staff: Mentor

    Perhaps your data doesn't have an FFT with a significant frequency that can be picked out i.e. no one max value

    Try plotting the FFT's of both signals to see what you get.
     
  10. Apr 16, 2015 #9
    thanks for the spectrum plot tips jedishrfu
    I haven't solved it yes but I think you pointed me in the right direction.

    this is the spectrum of the 3 signals:

    spectrum.jpg


    So there is a huge spike at zero Hz.

    I have composed each signal multiplying 2 previous signals that where at 50Hz. So I expect the new signals to be at 100 Hz. And I see a spike there.
    But I dont understand why the spike in zero.
     
    Last edited: Apr 16, 2015
  11. Apr 16, 2015 #10
    You have a DC component in you signal. max function should detected the peak of the periodic signal ?
     
  12. Apr 16, 2015 #11

    when I use max function it gives back the index of frequency = 0 HZ
     
  13. Apr 16, 2015 #12
    Then the answer is correct, the DC components have zero phase shift :)
     
  14. Apr 16, 2015 #13
    spectrum1.jpg
    well, ok thats correct.

    Since I need the 100Hz I subtracted the mean value of the signal to itself.



    thanks a lot!
     
  15. Apr 16, 2015 #14
    The zero value of the fourier transform will always just be the mean give a or take a factor of pi and volume.
     
  16. Apr 16, 2015 #15

    boneh3ad

    User Avatar
    Science Advisor
    Gold Member

    In general, it's usually a good idea to make sure the signal you wish to pass through your FFT has zero mean, that way you are only dealing with nonzero frequency components.
     
  17. Apr 17, 2015 #16
    the original 2 signals had no DC component. But by multiplying them value by value the resulting signal was no longer centered on the zero for the vertical axis.
    Pretty obvious after you get it..
     
  18. May 7, 2016 #17
    what did you change ad your coding to make that output? thanks in advance :)
     
Know someone interested in this topic? Share this thread via Reddit, Google+, Twitter, or Facebook

Have something to add?
Draft saved Draft deleted



Similar Discussions: Matlab - calculate phase shift using fft
  1. Phase shift using FFT (Replies: 11)

Loading...