Phase

Phase - Compute instantaneous phase in signal.

 Compute instantaneous phase (wrapped and unwrapped) and amplitude
 in signal, using the Hilbert transform of the signal.

  USAGE

    [phase,amplitude,unwrapped] = Phase(samples,times)

    samples        signal (e.g. filtered local field potential) <a href="matlab:help samples">samples</a>
    times          optional timestamps where phase should be interpolated

  NOTE

    Angles are returned in radians.

  SEE

    See also FilterLFP, PhasePrecession, PhaseMap.

0001 function [phase,amplitude,unwrapped] = Phase(samples,times)
0002 
0003 %Phase - Compute instantaneous phase in signal.
0004 %
0005 % Compute instantaneous phase (wrapped and unwrapped) and amplitude
0006 % in signal, using the Hilbert transform of the signal.
0007 %
0008 %  USAGE
0009 %
0010 %    [phase,amplitude,unwrapped] = Phase(samples,times)
0011 %
0012 %    samples        signal (e.g. filtered local field potential) <a href="matlab:help samples">samples</a>
0013 %    times          optional timestamps where phase should be interpolated
0014 %
0015 %  NOTE
0016 %
0017 %    Angles are returned in radians.
0018 %
0019 %  SEE
0020 %
0021 %    See also FilterLFP, PhasePrecession, PhaseMap.
0022 
0023 % Copyright (C) 2004-2011 by Michaël Zugaro
0024 %
0025 % This program is free software; you can redistribute it and/or modify
0026 % it under the terms of the GNU General Public License as published by
0027 % the Free Software Foundation; either version 3 of the License, or
0028 % (at your option) any later version.
0029 
0030 % Check number of parameters
0031 if nargin < 1,
0032     error('Incorrect number of parameters (type ''help <a href="matlab:help Phase">Phase</a>'' for details).');
0033 end
0034 
0035 % Check parameter sizes
0036 if size(samples,2) < 2,
0037     error('Parameter ''samples'' is not a matrix (type ''help <a href="matlab:help Phase">Phase</a>'' for details).');
0038 end
0039 if nargin == 2,
0040     if ~isdvector(times),
0041         error('Parameter ''times'' is not a vector (type ''help <a href="matlab:help Phase">Phase</a>'' for details).');
0042     end
0043     if size(times,2) ~= 1, times = times'; end
0044 end
0045 
0046 phase = zeros(size(samples));
0047 phase(:,1) = samples(:,1);
0048 unwrapped = zeros(size(samples));
0049 unwrapped(:,1) = samples(:,1);
0050 amplitude = zeros(size(samples));
0051 amplitude(:,1) = samples(:,1);
0052 
0053 for j = 2:size(samples,2),
0054     % Compute phase and amplitude using Hilbert transform
0055     h = hilbert(samples(:,j));
0056     phase(:,j) = mod(angle(h),2*pi);
0057     amplitude(:,j) = abs(h);
0058     unwrapped(:,j) = unwrap(phase(:,j));
0059 end
0060 
0061 if nargin == 2,
0062 %      phase = Interpolate([phase(:,1) exp(i*phase(:,2))],times(:,1));
0063 %      phase(:,2) = angle(phase(:,2));
0064     phase = Interpolate(phase,times(:,1),'type','circular');
0065     amplitude = Interpolate(amplitude,times(:,1));
0066     unwrapped = Interpolate(unwrapped,times(:,1));
0067 end