function [toa_date,toa]=doppcorr(toa_date,toa,sampl_int);

% doppcorr.m
% 10-NOV-2000, HFurey.  This version based on the centered differencing 
% method recommended by Phil Richardson and Amy Bower.  
% Constants in both the old code, and Olaf's code is based
% on a 250 Hz transmission frequency and 80 second signal length, which should
% not be assumed constant.  There are currently no input variable that
% can be used for signal strength and transmission frequency, perhaps
% add to rfb-file header?
% Also, if sound velocity is to be used, please make sure chosen sound velocity
% is used, not the default (usually del grosso).

% ALGORITHM:
% deltaT = (delta * velocity_towards_soso * signal_start_frequency)/
%          (signal sweep * sound speed);
% dTOA = ( diff(TOA) * sv )./(diff(TOA date) * secperday);
% doppler corrected TOAs = TOAs - (deltaT * dTOAs);
% Note: sound velocity does not enter into code below because it is cancelled 
% out in algorithm above. 
%------------------------------------------------------------------------------

% Doppler Constant:
% signal_length = 80 seconds
% signal_start_frequency = 259.375 Hz
% signal_sweep = 1.523 Hz;
% sec_per_day = 86400;
%
% doppcon = (signal_length*signal_start_frequency)/(signal_sweep*sec_per_day);
doppcon = 0.151990;

% innocuous.
[toa_date,ix]=sort(toa_date);
toa=toa(ix);

% I really need to know sampling interval to correctly apply the time.
% I am not going to use the diff routine because will not do centered difference.
% this makes for a long loop - but ebbene.
dctoa = toa;
diffdate = diff(toa_date); % n-1 long
% doppcorrTOA = toa - doppcon*(d(toa)/d(t)), provided d(t) is sampling interval
% I want to avoid calculating the correction over a period of time > samplerate.
sampl_int = sampl_int./24; % convert to days

for m = 1:length(toa_date)
   % if first toa
	if m==1
		if diffdate(m)==sampl_int % else a single point, and leave TOA unchanged
		dctoa(m) = toa(m) - doppcon*((toa(m+1)-toa(m))/sampl_int);
		end
	
	% if last toa
	elseif m==length(toa_date)
		if diffdate(m-1)==sampl_int % else a single point, and leave TOA unchanged
		dctoa(m) = toa(m) - doppcon*((toa(m)-toa(m-1))/sampl_int);
		end

	% if middle toas
   else
		if diffdate(m)==sampl_int & diffdate(m-1)==sampl_int % centered difference
		dctoa(m) = toa(m) - doppcon*((toa(m+1)-toa(m-1))/(2*sampl_int));
		elseif diffdate(m)==sampl_int & diffdate(m-1)~=sampl_int % forward difference
		dctoa(m) = toa(m) - doppcon*((toa(m+1)-toa(m))/sampl_int);
		elseif diffdate(m)~=sampl_int & diffdate(m-1)==sampl_int % backwards difference
		dctoa(m) = toa(m) - doppcon*((toa(m)-toa(m-1))/sampl_int);
		end
		% if point does not qualify for date restrictions then a stand-alone point,
		% and should remain unchanged.
	end
end

toa=dctoa;