% calc_lag_MVPdata
% =========================================================================
%
% USAGE:
%  Alter the pathnames where indicated in the comment line and then run.
%  All required data is loaded.


% DESCRIPTION:
%  Determines the conductivity lag that should be applied to downcast MVP
%  conductivity data to correct for the differing response times of
%  temperature and conductivity sensors, which is further exacerbated by
%  the rapid descent rate of the MVP. This is achieved by calculating the
%  largest cross correlationa of dT/dt and dC/dt. The method was developed
%  by Barth et al. (1996) and implemented by Ullman & Hebert (2013).
%
% OUTPUT
%  A structure array of MVP downcasts containing metadata, pressure,
%  temperature, conductivity, salinity and sigmat.
%
% REFERENCES
%  Barth, J. A., R. O'Malley, J. Fleischbein, R. L. Smith and A. Huyer,
%   1996. SeaSoar and CTD observations during the Coastal Jet Separation
%   cruise W9408A  August to September 1994. College of Oceanic and
%   Atmospheric Sciences, Oregon State University, Corvallis. Reference
%   96-1, Data Report 162, November 1996.
%
%  Ullman, D. S. and D. Hebert (2013), Processing of underway CTD data, J.
%   Atmos. Oceanic Tech., (submitted).
%
% AUTHOR:
%  Dr. Susanne Craig, Bedford Institute of Oceanography (BIO)
%  Adapted from code written by Dr. Susanne Craig & Dr. Dave Hebert (BIO)
%
% VERSION NUMBER:
% 1.0 (25th June 2013)
% 1.1 (17th July 2013)
% 1.2 (16th Oct 2013)
%
% =========================================================================

clear
close all
clc

% Change pathnames as appropriate
addpath(genpath('/Users/scraig/Documents/Ross Sea/MVP/final MVP code')) % <==
load '/Users/scraig/Documents/Ross Sea/MVP/MVPrawData' % <==
path2save = '/Users/scraig/Documents/Ross Sea/MVP'; % <==

%% Identify the best lag

% Calculate the lag of the largest cross-correlation of DT/dt and DC/dt
% Use this section to step through each of the MVP profiles, plot the lag
% matrix and get a feel for which lag minimises the salinity spikes.
lag = zeros(length(MVPd),1);

for n = 1:length(MVPd)
    
    
    dtemp = diff(MVPd(n,1).temp);
    dcond = diff(MVPd(n,1).cond);
    [c, lags] = xcorr(dtemp,dcond,20);
    idx = find(c == max(c));
    lag(n,:) = lags(idx);
end

% Average the lag to obtain a mean
lagMean = mean(lag);



