interpPhi

 INTERPPHI Interpolate (state variable, phase) pairs to obtain a 
 single-valued function of phase
   kin = interpPhi( phi0, kin0, phi )
 INPUT:
   phi0 - n x 1 - phase samples
   kin0 - n x p - corresponding state samples, p different variables
   phi  - m x 1 - desired phase samples

 OUTPUT:
   kin  - m x p - interpolated state value at desired phases

 Phase is not strictly increasing.  To compensate for this,
 we oversample, interpolate linearly, then smooth the result.

 The result should be extremely close to the original data in regions
 where phase is strictly increasing.

 Only considers non-NaN region of data for convenience.
 
 $Revision: $
 By Sam Burden, Berkeley 2009

0001 function kin = interpPhi( phi0, kin0, phi )
0002 % INTERPPHI Interpolate (state variable, phase) pairs to obtain a
0003 % single-valued function of phase
0004 %   kin = interpPhi( phi0, kin0, phi )
0005 % INPUT:
0006 %   phi0 - n x 1 - phase samples
0007 %   kin0 - n x p - corresponding state samples, p different variables
0008 %   phi  - m x 1 - desired phase samples
0009 %
0010 % OUTPUT:
0011 %   kin  - m x p - interpolated state value at desired phases
0012 %
0013 % Phase is not strictly increasing.  To compensate for this,
0014 % we oversample, interpolate linearly, then smooth the result.
0015 %
0016 % The result should be extremely close to the original data in regions
0017 % where phase is strictly increasing.
0018 %
0019 % Only considers non-NaN region of data for convenience.
0020 %
0021 % $Revision: $
0022 % By Sam Burden, Berkeley 2009
0023 
0024 phi0 = phi0(:);
0025 % kin0 = kin0(:);
0026 phi = phi(:);
0027 
0028 kin = nan*zeros(size(phi,1),size(kin0,2));
0029 
0030 pnn = ~isnan(phi0);
0031 
0032 for p = 1:size(kin0,2)
0033   
0034   [Bo,Ao] = butter(3, 0.01);
0035   
0036   knn = ~isnan(kin0(:,p));
0037   
0038   k = find(pnn & knn);  
0039 
0040   ovrsmp = sortrows(interp1(1:length(phi0(k)),[phi0(k),kin0(k,p)],1:0.025:length(phi0(k)),'linear'),1);
0041   ovrsmp(diff(ovrsmp(:,1)) == 0,:) = nan;
0042   notnans = find(~isnan(ovrsmp(:,2)));
0043 
0044   ovrsmp2 = nan*ovrsmp(:,2);
0045   ovrsmp2(notnans) = filtfilt(Bo, Ao, ovrsmp(notnans,2));
0046 
0047   kin(:,p) = interp1(ovrsmp(notnans,1),ovrsmp2(notnans),phi,'linear');
0048 end
0049