function samp = ifodelay(detId, cosTheta, phi, sampFreq)

% IFODELAY - find time delay for IFO from source sky location
%
% samp = ifodelay(detID, theta, phi, sampFreq);
%
% detId         detector identifier string
% cosTheta      cosine of source co-declination angle
%               N.B. only cos(theta) for 0<=theta<=pi supported
% phi           source longitude in rad
% sampFreq      sampling frequency of IFO
%
% samp          number of samples to ADD for delay
%
% SEE ALSO: GRAVEN, GETIFO

% $Id: ifodelay.m,v 1.9 2004/05/12 00:47:22 stuver Exp $

% Theory: Find the projection of the vector pointing from the center of the
% Earth to the detector onto the GW's unit vector (pointing from the center
% of the Earth to the source's sky location) and convert this into time of
% flight in units of samples.

% CHECK for errors
% GET the WGS-84 X, Y & Z coordinates for detId
% DEFINE the gravitational wave's unit vector from sky location
% DETERMINE the number of samples to add to the start time

samp = [];

errorcheck(detId, cosTheta, phi, sampFreq);

% Get the WGS-84 coordinates for detId
Cordnt = getifo(detId);
% Make sure that Cordnt is not empty
if isempty(Cordnt)
    Cordnt
    error('ifodelay: getifo returned empty Cordnt')
end


% Define the X, Y & Z WGS-84 coordinates for detId
XIFO = Cordnt.WGS84(1); % m
YIFO = Cordnt.WGS84(2); % m
ZIFO = Cordnt.WGS84(3); % m

% From the given cosine of the source's co-dec to line-of-sight,
% define the sine
sinTheta = sqrt(1-cosTheta.^2);
% N.B. only cos(theta) for 0<=theta<=pi supported

% Define the X, Y & Z coordinates of the GW's sky location
xgw = sinTheta*cos(phi);
ygw = sinTheta*sin(phi);
zgw = cosTheta;

% Define the GW's unit vector
gwhat = sqrt((xgw^2)+(ygw^2)+(zgw^2)); 

% Define the speed of light
C0=299792458; % m/s

% Determine the number of samples to add to the start time from the
% definition of dot product, speed of light and sampling frequency
%
% t = -(XIFO*xgw + YIFO*ygw + ZIFO*zgw)/c*gwhat
% samp = t*sampFreq

samp = -(((XIFO*xgw)+(YIFO*ygw)+(ZIFO*zgw))/(C0*gwhat))*sampFreq;

return

% ~~~ END MAIN ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

function errorcheck(detId, cosTheta, phi, sampFreq);

% ERRORCHECK - checks ifodelay's input data for errors and stops the code if
%              an error is found
%
% errorcheck(detId, cosTheta, phi, sampFreq)
%
% detId         detector identifier string
% cosTheta      cosine of source co-declination angle
% phi           source longitude in rad
% sampFreq      sampling frequency of IFO

% CHECK variables for error

% detId must be a string
if ~strcmp(class(detId), 'char')
    detId
    error('ifodelay: detId must be a string')
end

% cosTheta must be real and bewteen -1 and 1
if abs(cosTheta) > 1 | ~isreal(cosTheta)
    cosTheta
    error('ifodelay: cosTheta must be real and bewteen -1 and 1')
end

% phi must be real and between -pi and pi
if phi < -pi | phi > pi | ~isreal(phi)
    phi
    error('ifodelay: phi must be real and between -pi and pi')
end
    
% sampFreq must be a real positive number
if sampFreq <= 0 | isinf(sampFreq) | isnan(sampFreq) | ~isreal(sampFreq)
    sampFreq
    error('ifodelay: sampFreq must be a real finite number > 0')
end

return