%   This script computes salinity in a new time base.  It is the first function
%   to call in the most accurate CTD computation I have been able to create.
%   It has been used with the MA Bay 1998 data set.
%
%   This function 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.108627273 6.03127	0.085341];		% Uses lag from T-S linearization
%
%   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.
%
%   Generates the following variables:
%
%      ctd_time:  	uniformly sampled time, 0.2s period
%      ctd_temp:   	water temperature sampled on time base
%      ctd_pres:  	pressure (bar) sampled on time base
%      ctd_depth:		depth (m) sampled on time base
%      ctd_salt:   	salinity sampled on time base
%      ctd_dens: 		density sampled on time base
%
%    James G. Bellingham
%    March 28, 1999
%

% k = [-0.05788 6.03127	0.085341];  		% Uses lag from despiking
k = [-0.108627273 6.03127	0.085341];		% Uses lag from T-S linearization

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

% If the run is obviously a yoyo run, use only the yoyo data...
if (length(going_up)+length(going_down)) / length(time) > 0.8
   t1 = time(going_down(1));
   if t1 < 45 % don't use first 45 seconds under any circumstances...
      t1 = 45;
   end;
   
   t2 = time(going_up(length(going_up)));
end;

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

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

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

%  Convert bar to decibars and make depth correction for AUV attitude
ok = find(~isnan(data(:,m_pressure)));
ctd_pres = interp1q(data(ok,m_present_time), ...
   cCTD_pres(data(ok, m_pressure), data(ok, m_pitch)), ctd_time);

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

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

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

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

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

% -----------------------------------------------------------------------
% Compute the salinity, density, depth, and sos...

ctd_salt = sw_salt(10*cond/sw_c3515,cond_temp,10*ctd_pres);
ctd_dens = sw_dens(ctd_salt, ctd_temp, 10*ctd_pres);
ctd_depth = sw_dpth(10*ctd_pres, 42);   % Use latitude of MA Bay....
ctd_sos = sw_svel(ctd_salt, ctd_temp, 10*ctd_pres);

clear cond; clear cond_temp; clear ok; clear fcond; clear dt; clear n;
clear cell_temp; clear cell_tau; clear g; clear t1; clear t2;