function [ YI ] = interp1Complex( X, Y, XI, varargin )
% YI = interp1Complex( X, Y, XI, dim, ADDL_ARGS )
% Interpolates to find YI, the values of the underlying function Y at the 
% points in the array XI under the assumption that Y is complex and the
% X values are sufficiently close so that the angle between adjacent Y 
% values is less than pi.
%
% In contrast to conventional interp1, interp1Complex interpolates the
% magnitude and phase information independently.  This is appropriate for 
% some applications (e.g. spectral interpolation), but not all, so use with
% caution.
%
% Complex interpolation is built upon the MATLAB function interp1 and
% if X and Y are real, the interp1Complex will return exactly the same
% results as interp1.  
%
% X must be a vector of length N.   
% 
% If Y is a vector, then it must also have length N, and YI is the same 
% size as XI.  If Y is an array of size [N,D1,D2,...,Dk], then the 
% interpolation is performed for each D1-by-D2-by-...-Dk value in 
% Y(i,:,:,...,:). 
%
% If XI is a vector of length M, then YI has size [M,D1,D2,...,Dk]. If XI 
% is an array of size [M1,M2,...,Mj], then YI is of size 
% [M1,M2,...,Mj,D1,D2,...,Dk].
% 
%
% Inputs:   X         - [Nx1]        Points at which function is evaluated
%           Y         - [NxD1xD2...] Function values
%           XI        - [M1xM2x...]  New points to evaluate function
%           ADDL_ARGS -              Additional arguments to pass to interp1
% 
% Outputs:  YI        - [M1xM2x...xD1xD2x...] Interpolated data
%
% Copyright (c) 2011, Hidden Solutions, LLC
% All rights reserved.
%
% Redistribution and use in source and binary forms, with or without
% modification, are permitted provided that the following conditions are met:
%
%   * Redistributions of source code must retain the above copyright
%      notice, this list of conditions and the following disclaimer.
%    * Redistributions in binary form must reproduce the above copyright
%      notice, this list of conditions and the following disclaimer in the
%      documentation and/or other materials provided with the distribution.
%    * Neither the name Hidden Solutions, LLC nor the names of any
%      contributors may be used to endorse or promote products derived from 
%      this software without specific prior written permission.
%
% THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
% ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
% WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
% DISCLAIMED. IN NO EVENT SHALL HIDDEN SOLUTIONS, LLC BE LIABLE FOR ANY
% DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
% (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
% LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
% ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
% (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
% SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
%
% Author: James S. Hall (Hidden Solutions, LLC)
% Date:   10 August 2011
%

% Perform interpolation on magnitude and phase separately
YI =        interp1(X,      abs(Y),     XI, varargin{:}) .* ...
     exp(1i*interp1(X,unwrap(angle(Y)), XI, varargin{:}));

end

