function error = LagTTError(k_in, data, m_present_time, water_temp, water_cond, m_pressure, type);

% error = LagTTError(k, data, m_present_time, water_temp, water_cond, m_pressure);
%
%   This function computes the error for (data(t,a) - k(1)*data(t-k(2),b))^2
%   It computes the temperature of the  conductivity cell, and then assumes
%   that some fraction of the water's heat content goes to the thermal mass
%   of the conductivity cell.
%
%   For k:
%     k(1) = time lag due to plumbing (s)
%     k(2) = thermal time constant of conductivity cell (s)
%     k(3) = thermal coupling between water and cell
%     k = [-0.35 10 0.1]  (best first guess)
%
%   This function does NOT assume the parameters have been sampled on a
%   uniform time interval.  It resamples the parameters onto a uniform
%   time-base using MATLABs INTERP1Q function.  This makes the function
%   a little slow, but hopefully more accurate.
%
%   This function high-pass filters the data before looking for the error
%   since we are looking for a short lag time.  It uses the fft function
%   for the filtering, by simply truncating the low frequency component.
%
%   The variable 'type' determines what sort of error is computed, as below:
%
%		type = 1   tries to get the most linear T-S curve 
%      type = 2   minimizes the time derivative of S to eliminate spiking...
%
%   J. G. Bellingham
%   December 27, 1998
%

% Check and see if all the parameters were passed across...
k = k_in;

if length(k) == 1
   k(2) = 6.2940;
   k(3) = 0.0644;
elseif length(k) == 2
   k(2:3) = k(1:2);
   k(1) = -0.05788;
end;

% Determine bounds of data.  Don't use first or last 45 seconds...
lastn = length(data(:, m_present_time));
t1 = data(5*45, m_present_time);
t2 = data(lastn-5*45, m_present_time);

% Establish new time-base with uniform interval
dt = 0.2;
time = (t1:dt:t2)';
n = length(time);

% -----------------------------------------------------------------------
% Resample parameters on new (uniform) time-base
ok = find(~isnan(data(:,water_temp)));
temp = interp1q(data(ok,m_present_time), data(ok,water_temp), time);

% Add lag to conductivity
ok = find(~isnan(data(:,water_cond)));
cond = interp1q(data(ok,m_present_time), data(ok,water_cond), (time-k(1)));

%  Convert bar to decibars
ok = find(~isnan(data(:,m_pressure)));
pres = interp1q(data(ok,m_present_time), data(ok, m_pressure), time);

% -----------------------------------------------------------------------
% Compute the temperature of the conductivity cell...

cell_tau = k(2);   % time constant of cell...
cell_temp = temp;
g = exp(-dt/cell_tau);  % integration kernel

for indx = 2:length(time),
   cell_temp(indx) = g*cell_temp(indx-1) + (1-g)*temp(indx);
end

% -----------------------------------------------------------------------
% Compute the corrected water temperature at the conductivity cell

fcond = k(3);
cond_temp = temp - fcond * (temp - cell_temp);

% -----------------------------------------------------------------------
% Compute the salinity

salt = sw_salt(10*cond/sw_c3515,cond_temp,pres);
dens = sw_dens(salt, temp, pres);

% close all;

% figure;
% plot(temp, salt, '+')
% zoom on;

% figure;
% plot(time, temp, time, cond_temp);
% zoom on;

% up = find(diff(pres) < -0.01);
% dwn = find(diff(pres) > 0.01);
% figure;
% plot(pres(up), salt(up), '+', pres(dwn), salt(dwn), '.');
% zoom on;

% -----------------------------------------------------------------------
% Compute deviation from linear temp vs salinity relation...

if type == 1
   param = polyfit(temp, salt, 1);
   error = sum((salt - polyval(param, temp)).^2);
elseif type == 2
   error = sum((diff(salt)).^2)
end;

return;   