% function [ t, pairs ] = location_to_timedelay(geometry, pairs, locations, fs, c)
%
% For each row of "pairs", compute the "delay of signal 2 relative to signal 1".
% -> signal 1: microphone in column 1 of "pairs".
% -> signal 2: microphone in column 2 of "pairs".
%
% geometry:  nDimensions-by-nChannels
% pairs:     nPairs-by-2    ( if [], we build it: all possible pairs )
% locations: nDimensions-by-nLocations
% fs:        frequency in Hz
% c:         speed of sound in the air in m/s
%
% t:         nPairs-by-nLocations    ( time-delays in samples )

function [ t, pairs ] = location_to_timedelay( geometry, pairs, locations, fs, c )

a = 5; 
if nargin < a
  error( sprintf( 'location_to_timedelay() needs at least %d parameters', a ) );
end

if isempty( pairs )
  nChannels = size( geometry, 2 );
  for a = 1:nChannels
    for b = (a+1):nChannels
      pairs = [ pairs; a b ];
    end
  end
end

nPairs      = size( pairs, 1 );
nLocations  = size( locations, 2 );
nDimensions = size( geometry, 1 );

maxdelay = pair_distance( geometry, pairs ) * fs / c;

% Where to store the result
t = zeros( nPairs, nLocations );

for pair = 1:nPairs
  
  % Location of microphone one
  x1 = geometry(:, pairs(pair, 1));
  
  % Location of microphone two
  x2 = geometry(:, pairs(pair, 2));
  
  d1 = zeros( 1, nLocations );
  d2 = zeros( 1, nLocations );
  for dim = 1:nDimensions
    d1 = d1 + ( locations( dim, : ) - x1( dim ) ) .^ 2;
    d2 = d2 + ( locations( dim, : ) - x2( dim ) ) .^ 2;
  end
  
  % Euclidian distance from each location to each of the two microphones
  d1 = sqrt( d1 );
  d2 = sqrt( d2 );
  
  t( pair, : ) = ( d2 - d1 ) * fs / c;
  
  % Make sure we do not exceed physical bounds
  % ( sometimes a small calculation error happens, about 1e-13 for example )
  t( pair, : ) = max( -maxdelay( pair ), t( pair, : ) );
  t( pair, : ) = min( +maxdelay( pair ), t( pair, : ) ); 
  
end