% This function calculates a 2-channel stereo matrix S for a given mono wave
% file and the desired parameters for the difference of time and amplitude
% the 2 resulting channels shall have.
% The intention is to move a (mono) sound from one virtual position between two
% loudspeakers to another, so that it appears clearly as a "phanton sound source"
%
% it is recommended to use loudspeakers instead of headphones
%
% usage: [S, FS] = stereomove(amp_diff_start,time_diff_start,amp_diff_dest,time_diff_dest);
%        e.g. [S, FS] = stereomove(-6,-500,+6,+500);  
%
%        use sound(S, FS) to play the matrix S
%        or wavwrite(S, FS, 'stereo.wav'); to store it on your disc
%
% PARAMETERS:
%
% amp_diff_...   : difference between signal amplitude in dB ref. to left channel
%                  pos. values -> the right speaker is louder
%                  neg. values -> the left speaker is louder
%                  amp_diff = 0 -> no difference between left and right speaker
%         _start : the start value
%         _dest  : the destination value
%                  (just use identical values for static amplitude differences)

%
% time_diff_     : desired time difference in microseconds between L and R   
%                  use positive values to play the sound earlier(!) on the
%                  right speaker
%         _start : the start value
%         _dest  : the destination value
%                  (just use identical values for static time differences)
%
% OUTPUT:
%
% S : [Nx2]matrix where N is the number of sound samples
%
% authors: Frithjof Hummes, Hendrik Buschmeier 
% date   : 2005-10-03


function [S,FS] = stereomove(amp_diff_start, time_diff_start, amp_diff_dest, time_diff_dest)

wavname = 'move';
              
% load the mono wave file 
[Y, FS, NBITS] = wavread(wavname);

% initialization of the stereo matrix
S = zeros(length(Y),2);

% calculate the "factor-vectors" for the channels amplitudes 
% for a given difference between R and L in dB 

   volume_steps = 20; % number of volume steps

   ratio_start = 10^(amp_diff_start/10);
   rf_start = ratio_start / (ratio_start + 1);   % start factor for right channel
   lf_start = 1 / (ratio_start + 1);             % start factor for left channel
   ratio_dest = 10^(amp_diff_dest/10);
   rf_dest = ratio_dest / (ratio_dest + 1);      % last factor for right channel
  
   
   delta_f = rf_dest - rf_start;           % range between start and end amplitude
   stepsize = delta_f / volume_steps;      % size of one volume step
   
   RF = zeros(length(Y)+volume_steps,1);   % vectors of the amplification factors
   LF = zeros(length(Y)+volume_steps,1);   % including space for errors due to rounding
   
   for i = 1:volume_steps           
       A = (i-1)*(ceil(length(RF)/volume_steps)) + 1;
       B = i*(ceil(length(RF)/volume_steps));
       RF(A:B) = (rf_start + ((i-1)*stepsize)) * ones(ceil(length(RF)/volume_steps),1);
       LF(A:B) = (lf_start - ((i-1)*stepsize)) * ones(ceil(length(RF)/volume_steps),1);
   end
           
         S(:,2) = RF(1:length(Y)) .* Y;
         S(:,1) = LF(1:length(Y)) .* Y;
         
% calculate the number of samples to shift the channel
% FS/1,000,000 samples per microsecond

   shft_start = round(time_diff_start * (FS / 1000000));
   shft_dest  = round(time_diff_dest * (FS / 1000000));
   
   delta_shft = shft_dest - shft_start;
   
% shift the right channel to the starting position

   shft = shft_start;
   T = S(:,2);                         % the temporary right channel
   R = zeros(length(S(:,2)),1);        % the new right channel
   if (shft < 0)  
        T = T(1:length(T)+shft);       % take the first samples 
        R(1-shft:length(R)) = T;       % place it at the desired position               

   elseif (shft > 0) 
        T = T(shft+1:length(T));       % take the last samples 
        R(1:length(T)) = T;            % place it at the desired position                    

   else R = T;                         % no shifting
   end  
   S(:,2) = R;      % copy to right channel in stereo martix

% now begins the continous shifting

if (delta_shft > 0)         % kick samples of the right channel
      kick = delta_shft
      T = zeros(length(S(:,2)),1);
      for i = 1:(kick+1)
          A = (i-1)*floor(length(T)/(kick+1)) + 1;
          B = i*floor(length(T)/(kick+1));
          if (i < (kick+1)) 
              T(A:B) = S(A+(i-1):B+(i-1) ,2);
          else 
              T(A:length(S(:,2))-kick) = S(A+(i-1):length(S(:,2)),2); % just the rest 
          end
      end
      S(:,2) = T;    % copy into the stereo matrix          
      
elseif (delta_shft < 0)     % kick samples of the left channel
      kick = (-1)*delta_shft;
      T = zeros(length(S(:,1)),1);
      for i = 1:(kick+1)
          A = (i-1)*floor(length(T)/(kick+1)) + 1;
          B = i*floor(length(T)/(kick+1));
          if (i < (kick+1)) 
              T(A:B) = S(A+(i-1):B+(i-1) ,1);
          else 
              T(A:length(S(:,1))-kick) = S(A+(i-1):length(S(:,1)),1); % just the rest 
          end
      end
      S(:,1) = T;      % copy into the stereo matrix          
   
else             % No further changings
   
end

return % stereomove