%
%  gliderCTD_CorrectThermalLag.m
%
%  Purpose: Apply thermal lag correction to glider CTD salinity and density
%           and generate mat file
%
%  NOTE: At top of file, set glider index, deployment number and desired correction parameters
%
%  Requires:
%  correctThermalLag.m
%  MATLAB folder - contains util, matutil, seawater, plots, strfun, opnml
%
%
%  Author:  William Stark
%           Marine Sciences Department
%           UNC-Chapel Hill
%
%  Created: 21 May 2012
%
%//////////////////////////////////////////////////////////////////////////

clear all;


%!!!!!!  SET PRIOR TO RUNNING CODE  !!!!!!!!!!!!!!!!!!!!!!!!
%!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
%
% SET THE GLIDER INDEX (Pelagia = 1, Ramses = 2) ...
gliderIndex = 1;

% SET THE DEPLOYMENT NUMBER (1, 2 or 3) ...
deploymentNumber = 1;

% SET CORRECTION PARAMETERS STRUCTURE...
correctionParams = [0.1587 0.0214 6.5316 1.5969];

%
%!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
%!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!



% add paths for required files...
addpath('MATLAB/util/');
addpath('MATLAB/matutil/');
addpath('MATLAB/seawater/');
addpath('MATLAB/plots/');
addpath('MATLAB/strfun/');
addpath('MATLAB/opnml/');
addpath('MATLAB/opnml/FEM/');



% glider name string...
if (gliderIndex==1)
    strGliderName = 'Pelagia';
else
    strGliderName = 'Ramses';
end

% populate arrays for the deployment start and end dates...
% ex. strStart(2, 3) is start date for Ramses, Deployment 3
strStart = {'26-Jan-2012', '16-Feb-2012', '16-Mar-2012'; '26-Jan-2012', '16-Feb-2012', '16-Mar-2012'};
strEnd   = {'14-Feb-2012', '08-Mar-2012', '04-Apr-2012'; '14-Feb-2012', '12-Mar-2012', '03-Apr-2012'};

% deployment number string...
strDeploymentNumber = num2str(deploymentNumber);

% deployment start date string...
strStartDate = strStart(gliderIndex, deploymentNumber);

% deployment end date string...
strEndDate = strEnd(gliderIndex, deploymentNumber);

% define the path to the glider ascii files...
datadir = strcat('/Users/haloboy/Documents/MASC/MATLAB/CTD_data_correction/GLIDER_DATA_LEVEL0/',...
                 strGliderName, '_Deployment', strDeploymentNumber, '/');

% define default bounds for use in plots...
switch gliderIndex
    case 1  % Pelagia
        switch deploymentNumber
            case 1  % Deployment 1
                tempBounds =  [17.0 24.0];
                salinBounds = [36.0 36.4];
                densBounds =  [1025.0 1026.6];
                chlorBounds = [0.0 4.0];
                
            case 2  % Deployment 2
                tempBounds =  [17.0 24.0];
                salinBounds = [36.0 36.5];
                densBounds =  [1024.5 1026.8];
                chlorBounds = [0.0 4.0];
                
            case 3  % Deployment 3
                tempBounds =  [17.0 24.0];
                salinBounds = [35.9 36.7];
                densBounds =  [1024.4 1026.4];
                chlorBounds = [0.0 4.0];
        end
    case 2  % Ramses
        switch deploymentNumber
            case 1  % Deployment 1
                tempBounds =  [8.0 23.0];
                salinBounds = [35.0 36.4];
                densBounds =  [1024.5 1027.5];
                chlorBounds = [0.0 4.0];
            
            case 2  % Deployment 2
                tempBounds =  [9.0 25.0];
                salinBounds = [35.2 36.6];
                densBounds =  [1024.0 1027.5];
                chlorBounds = [0.0 4.0];
            
            case 3  % Deployment 3
                tempBounds =  [10.0 24.5];
                salinBounds = [35.3 36.7];
                densBounds =  [1024.4 1027.4];
                chlorBounds = [0.0 4.0];
        end
end


%##########################################################################################
    
%*** READ IN EBD DATA *****************************************************
% declare variables for storing data...
temp=[];
cond=[];
pres=[];
chlor=[];
ptime=[];
% mtime=[]; scioxy=[]; scibb=[]; scicdom=[]; scichlor=[]; scibbam=[];

% try to load all *.ebdasc files at once...
[files, Dstruct] = wilddir(datadir, '.ebdasc');
nfile = size(files, 1);

for i=1:nfile-1
    % protect against empty ebd file
    if(Dstruct(i).bytes>0)
        data = read_gliderasc2([datadir, files(i,:)]);
        %data = read_gliderasc3([datadir, files(i,:)]);

        % if the number of values (in data.data) is less than the number of 
        % vars (in data.vars), this means that the data were not completely read
        % in.  To correct this, pad data.data with NaNs until its length
        % equals that of data.vars...
        if (length(data.data) < length(data.vars))
            data.data = padarray(data.data, [0 length(data.vars)-length(data.data)], NaN, 'post');
        end

        % populate variables with data...
        if(~isempty(data.data))
            temp = [temp;data.data(:,strmatch('sci_water_temp', data.vars, 'exact'))];             % temperature
            cond = [cond;data.data(:,strmatch('sci_water_cond', data.vars, 'exact'))];             % conductivity
            pres = [pres;data.data(:,strmatch('sci_water_pressure', data.vars, 'exact'))];         % pressure (measure of depth) science bay
            switch gliderIndex
                case 1  % this is Pelagia...
                    chlor = [chlor;data.data(:,strmatch('sci_bbfl2s_chlor_scaled', data.vars, 'exact'))];  % chlorophyll
                case 2  % this is Ramses...
                    chlor = [chlor;data.data(:,strmatch('sci_flbbcd_chlor_units', data.vars, 'exact'))];  % chlorophyll
            end
            ptime = [ptime;data.data(:,strmatch('sci_m_present_time', data.vars, 'exact'))];       % present time
        end

        data = [];
    end  
end
%**************************************************************************

%*** READ IN DBD DATA *****************************************************
% declare variables for storing data...
ptime_dbd=[];
horizontalVelocity=[];
depth = [];
pitch=[];
avgDepthRate = [];
angleOfAttack = [];
lat = [];
lon = [];
gpsLat = [];
gpsLon = [];
wptLat = [];
wptLon = [];

% try to load all *.dbdasc files at once...
[files, Dstruct] = wilddir(datadir, '.dbdasc');
nfile = size(files, 1);

clear data;

for i=1:nfile
    % protect against empty dbd file
    if(Dstruct(i).bytes>0)
        data = read_gliderasc2([datadir, files(i,:)]);
        %data = read_gliderasc3([datadir, files(i,:)]);

        % if the number of values (in data.data) is less than the number of 
        % vars (in data.vars), this means that the data were not completely read
        % in.  To correct this, pad data.data with NaNs until its length
        % equals that of data.vars...
        if (length(data.data) < length(data.vars))
            data.data = padarray(data.data, [0 length(data.vars)-length(data.data)], NaN, 'post');
        end

        % populate variables with data...
        if(~isempty(data.data))
            ptime_dbd = [ptime_dbd; data.data(:,strmatch('m_present_time', data.vars, 'exact'))];               % present time
            horizontalVelocity = [horizontalVelocity; data.data(:,strmatch('m_speed', data.vars, 'exact'))];    % horizontal glider velocity      
            depth = [depth; data.data(:,strmatch('m_depth', data.vars, 'exact'))];                              % depth      
            pitch = [pitch; data.data(:,strmatch('m_pitch', data.vars, 'exact'))];                              % pitch (radians)
            avgDepthRate = [avgDepthRate; data.data(:,strmatch('m_avg_depth_rate', data.vars, 'exact'))];       % avg depth rate
            angleOfAttack = [angleOfAttack; data.data(:,strmatch('u_angle_of_attack', data.vars, 'exact'))];    % angle of attack (radians)
            wptLat = [wptLat; data.data(:,strmatch('c_wpt_lat', data.vars, 'exact'))];                          % Waypoint latitude
            wptLon = [wptLon; data.data(:,strmatch('c_wpt_lon', data.vars, 'exact'))];                          % Waypoint longitude
            gpsLat = [gpsLat; data.data(:,strmatch('m_gps_lat', data.vars, 'exact'))];                          % GPS latitude
            gpsLon = [gpsLon; data.data(:,strmatch('m_gps_lon', data.vars, 'exact'))];                          % GPS longitude
            lat = [lat; data.data(:,strmatch('m_lat', data.vars, 'exact'))];                                    % latitude
            lon = [lon; data.data(:,strmatch('m_lon', data.vars, 'exact'))];                                    % longitude
        end

        data = [];
    end
end
%**************************************************************************


% first, apply the sort() function to make sure that values in the time vectors
% (ptime and ptime_dbd) increase monotonically...
[Y,I] = sort(ptime);
ptime = Y;
temp = temp(I);
cond = cond(I);
pres = pres(I);
chlor = chlor(I);

[Y,I] = sort(ptime_dbd);
ptime_dbd = Y;
horizontalVelocity = horizontalVelocity(I);
depth = depth(I);
pitch = pitch(I);
avgDepthRate = avgDepthRate(I);
angleOfAttack = angleOfAttack(I);
wptLat = wptLat(I);
wptLon = wptLon(I);
gpsLat = gpsLat(I);
gpsLon = gpsLon(I);
lat = lat(I);
lon = lon(I);


% remove nans from EBD data...
i = find(~isnan(temp));
ptime = ptime(i); temp = temp(i); cond = cond(i); pres = pres(i); chlor = chlor(i);
i = find(~isnan(pres));
ptime = ptime(i); temp = temp(i); cond = cond(i); pres = pres(i); chlor = chlor(i);
i = find(~isnan(cond));
ptime = ptime(i); temp = temp(i); cond = cond(i); pres = pres(i); chlor = chlor(i);
i = find(~isnan(chlor));
ptime = ptime(i); temp = temp(i); cond = cond(i); pres = pres(i); chlor = chlor(i);

% remove conductivity values less than 1...
i = find(cond>=1);
ptime = ptime(i);  temp = temp(i);  pres = pres(i);  cond = cond(i);

% remove pressure values less than 0...
i = find(pres>=0);
ptime = ptime(i);  temp = temp(i);  pres = pres(i);  cond = cond(i);

% convert pitch and angle of attack from radians to degrees...
pitch = pitch*180/pi;
angleOfAttack = angleOfAttack*180/pi;

% compute actual glide angle = pitch + angle of attack...
glideAngle = pitch + angleOfAttack;

% make copy of dbd time stamp vector for use in salinity/density correction...
ptime1_dbd = ptime_dbd;


% remove nans from DBD data...
i = find(~isnan(horizontalVelocity));
ptime1_dbd = ptime1_dbd(i); horizontalVelocity = horizontalVelocity(i);
depth = depth(i); pitch = pitch(i); avgDepthRate = avgDepthRate(i);
glideAngle = glideAngle(i); lat = lat(i); lon = lon(i);
i = find(~isnan(depth));
ptime1_dbd = ptime1_dbd(i); horizontalVelocity = horizontalVelocity(i);
depth = depth(i); pitch = pitch(i); avgDepthRate = avgDepthRate(i);
glideAngle = glideAngle(i); lat = lat(i); lon = lon(i);
i = find(~isnan(pitch));
ptime1_dbd = ptime1_dbd(i); horizontalVelocity = horizontalVelocity(i);
depth = depth(i); pitch = pitch(i); avgDepthRate = avgDepthRate(i);
glideAngle = glideAngle(i); lat = lat(i); lon = lon(i);
i = find(~isnan(avgDepthRate));
ptime1_dbd = ptime1_dbd(i); horizontalVelocity = horizontalVelocity(i);
depth = depth(i); pitch = pitch(i); avgDepthRate = avgDepthRate(i);
glideAngle = glideAngle(i); lat = lat(i); lon = lon(i);
i = find(~isnan(glideAngle));
ptime1_dbd = ptime1_dbd(i); horizontalVelocity = horizontalVelocity(i);
depth = depth(i); pitch = pitch(i); avgDepthRate = avgDepthRate(i);
glideAngle = glideAngle(i); lat = lat(i); lon = lon(i);
i = find(~isnan(lat));
ptime1_dbd = ptime1_dbd(i); horizontalVelocity = horizontalVelocity(i);
depth = depth(i); pitch = pitch(i); avgDepthRate = avgDepthRate(i);
glideAngle = glideAngle(i); lat = lat(i); lon = lon(i);
i = find(~isnan(lon));
ptime1_dbd = ptime1_dbd(i); horizontalVelocity = horizontalVelocity(i);
depth = depth(i); pitch = pitch(i); avgDepthRate = avgDepthRate(i);
glideAngle = glideAngle(i); lat = lat(i); lon = lon(i);

% clip horizontalVelocity data at 0.6...
%i = find(horizontalVelocity~=0 & horizontalVelocity<0.6);
i = find(horizontalVelocity<0.6);
ptime1_dbd = ptime1_dbd(i); horizontalVelocity = horizontalVelocity(i);
depth = depth(i); pitch = pitch(i); avgDepthRate = avgDepthRate(i);
glideAngle = glideAngle(i); lat = lat(i); lon = lon(i);

% interpolate DBD data to align with EBD data...
horizontalVelocity = interp1(ptime1_dbd, horizontalVelocity, ptime);
depth = interp1(ptime1_dbd, depth, ptime);
pitch = interp1(ptime1_dbd, pitch, ptime);
avgDepthRate = interp1(ptime1_dbd, avgDepthRate, ptime);
glideAngle = interp1(ptime1_dbd, glideAngle, ptime);
lat = interp1(ptime1_dbd, lat, ptime);
lon = interp1(ptime1_dbd, lon, ptime);

% make sure there are no NaNs in the final set of data...
i = find(~isnan(horizontalVelocity));
horizontalVelocity = horizontalVelocity(i); depth = depth(i); pitch = pitch(i);
avgDepthRate = avgDepthRate(i); glideAngle = glideAngle(i); lat = lat(i); lon = lon(i);
ptime = ptime(i); temp = temp(i); cond = cond(i); pres = pres(i); chlor = chlor(i);
i = find(~isnan(depth));
horizontalVelocity = horizontalVelocity(i); depth = depth(i); pitch = pitch(i);
avgDepthRate = avgDepthRate(i); glideAngle = glideAngle(i); lat = lat(i); lon = lon(i);
ptime = ptime(i); temp = temp(i); cond = cond(i); pres = pres(i); chlor = chlor(i);
i = find(~isnan(pitch));
horizontalVelocity = horizontalVelocity(i); depth = depth(i); pitch = pitch(i);
avgDepthRate = avgDepthRate(i); glideAngle = glideAngle(i); lat = lat(i); lon = lon(i);
ptime = ptime(i); temp = temp(i); cond = cond(i); pres = pres(i); chlor = chlor(i);
i = find(~isnan(avgDepthRate));
horizontalVelocity = horizontalVelocity(i); depth = depth(i); pitch = pitch(i);
avgDepthRate = avgDepthRate(i); glideAngle = glideAngle(i); lat = lat(i); lon = lon(i);
ptime = ptime(i); temp = temp(i); cond = cond(i); pres = pres(i); chlor = chlor(i);
i = find(~isnan(glideAngle));
horizontalVelocity = horizontalVelocity(i); depth = depth(i); pitch = pitch(i);
avgDepthRate = avgDepthRate(i); glideAngle = glideAngle(i); lat = lat(i); lon = lon(i);
ptime = ptime(i); temp = temp(i); cond = cond(i); pres = pres(i); chlor = chlor(i);
i = find(~isnan(lat));
horizontalVelocity = horizontalVelocity(i); depth = depth(i); pitch = pitch(i);
avgDepthRate = avgDepthRate(i); glideAngle = glideAngle(i); lat = lat(i); lon = lon(i);
ptime = ptime(i); temp = temp(i); cond = cond(i); pres = pres(i); chlor = chlor(i);
i = find(~isnan(lon));
horizontalVelocity = horizontalVelocity(i); depth = depth(i); pitch = pitch(i);
avgDepthRate = avgDepthRate(i); glideAngle = glideAngle(i); lat = lat(i); lon = lon(i);
ptime = ptime(i); temp = temp(i); cond = cond(i); pres = pres(i); chlor = chlor(i);
i = find(~isnan(temp));
horizontalVelocity = horizontalVelocity(i); depth = depth(i); pitch = pitch(i);
avgDepthRate = avgDepthRate(i); glideAngle = glideAngle(i); lat = lat(i); lon = lon(i);
ptime = ptime(i); temp = temp(i); cond = cond(i); pres = pres(i); chlor = chlor(i);
i = find(~isnan(cond));
horizontalVelocity = horizontalVelocity(i); depth = depth(i); pitch = pitch(i);
avgDepthRate = avgDepthRate(i); glideAngle = glideAngle(i); lat = lat(i); lon = lon(i);
ptime = ptime(i); temp = temp(i); cond = cond(i); pres = pres(i); chlor = chlor(i);
i = find(~isnan(pres));
horizontalVelocity = horizontalVelocity(i); depth = depth(i); pitch = pitch(i);
avgDepthRate = avgDepthRate(i); glideAngle = glideAngle(i); lat = lat(i); lon = lon(i);
ptime = ptime(i); temp = temp(i); cond = cond(i); pres = pres(i); chlor = chlor(i);
i = find(~isnan(chlor));
horizontalVelocity = horizontalVelocity(i); depth = depth(i); pitch = pitch(i);
avgDepthRate = avgDepthRate(i); glideAngle = glideAngle(i); lat = lat(i); lon = lon(i);
ptime = ptime(i); temp = temp(i); cond = cond(i); pres = pres(i); chlor = chlor(i);

% eliminate any zeros in temp, cond, pres, chlor, depth, lat and lon
i = find(temp~=0);
horizontalVelocity = horizontalVelocity(i); depth = depth(i); pitch = pitch(i);
avgDepthRate = avgDepthRate(i); glideAngle = glideAngle(i); lat = lat(i); lon = lon(i);
ptime = ptime(i); temp = temp(i); cond = cond(i); pres = pres(i); chlor = chlor(i);
i = find(cond~=0);
horizontalVelocity = horizontalVelocity(i); depth = depth(i); pitch = pitch(i);
avgDepthRate = avgDepthRate(i); glideAngle = glideAngle(i); lat = lat(i); lon = lon(i);
ptime = ptime(i); temp = temp(i); cond = cond(i); pres = pres(i); chlor = chlor(i);
i = find(pres~=0);
horizontalVelocity = horizontalVelocity(i); depth = depth(i); pitch = pitch(i);
avgDepthRate = avgDepthRate(i); glideAngle = glideAngle(i); lat = lat(i); lon = lon(i);
ptime = ptime(i); temp = temp(i); cond = cond(i); pres = pres(i); chlor = chlor(i);
i = find(chlor~=0);
horizontalVelocity = horizontalVelocity(i); depth = depth(i); pitch = pitch(i);
avgDepthRate = avgDepthRate(i); glideAngle = glideAngle(i); lat = lat(i); lon = lon(i);
ptime = ptime(i); temp = temp(i); cond = cond(i); pres = pres(i); chlor = chlor(i);
i = find(depth~=0);
horizontalVelocity = horizontalVelocity(i); depth = depth(i); pitch = pitch(i);
avgDepthRate = avgDepthRate(i); glideAngle = glideAngle(i); lat = lat(i); lon = lon(i);
ptime = ptime(i); temp = temp(i); cond = cond(i); pres = pres(i); chlor = chlor(i);
i = find(lat~=0);
horizontalVelocity = horizontalVelocity(i); depth = depth(i); pitch = pitch(i);
avgDepthRate = avgDepthRate(i); glideAngle = glideAngle(i); lat = lat(i); lon = lon(i);
ptime = ptime(i); temp = temp(i); cond = cond(i); pres = pres(i); chlor = chlor(i);
i = find(lon~=0);
horizontalVelocity = horizontalVelocity(i); depth = depth(i); pitch = pitch(i);
avgDepthRate = avgDepthRate(i); glideAngle = glideAngle(i); lat = lat(i); lon = lon(i);
ptime = ptime(i); temp = temp(i); cond = cond(i); pres = pres(i); chlor = chlor(i);

% eliminate any negative chlorophyll values
i = find(chlor>0);
horizontalVelocity = horizontalVelocity(i); depth = depth(i); pitch = pitch(i);
avgDepthRate = avgDepthRate(i); glideAngle = glideAngle(i); lat = lat(i); lon = lon(i);
ptime = ptime(i); temp = temp(i); cond = cond(i); pres = pres(i); chlor = chlor(i);



% scale up the pressure...
pres = pres*10;

% calculate salinity (without correction)...
salin = sw_salt(10*cond/sw_c3515, temp, pres);

% calculate density (without correction)...
dens = sw_pden(salin, temp, pres, 0);

% calculate glider velocity using horizontal velocity (m_speed) and average depth rate (m_avg_depth_rate)...
gliderVelocity = sqrt(horizontalVelocity.^2 + avgDepthRate.^2);



% pass the correction parameters into the correctThermalLags function, which
% returns the corrected profile structure (with corrected temp and cond added)...
%
% profileStructure:  ptime: Present time instant at which this row was collected
%                    depth: Depth (pressure in decibars) measured by the CTD
%                    temp: Temperature measured by the CTD
%                    cond: Conductivity measured by the CTD
%                    pitch: Pitch angle of the glider (optional)
%
%

% CORRECTION SCHEME:  Pass in a glider velocity vector, so that correctThermalLag() will return
% profile data that has been corrected using flow speed equal to this glider velocity.
% Correction parameters are calculated using passed-in glider velocity.
profileStructure = struct('ptime', ptime, 'depth', pres, 'temp', temp, 'cond', cond, 'pitch', glideAngle);
[correctedProfileData, varargout] = correctThermalLag(profileStructure, correctionParams, gliderVelocity);

% get the corrected temperature...
tempCorrected = correctedProfileData.tempInCell;

% calculate the corrected salinity using the corrected temperature...
salinCorrected = sw_salt(10*cond/sw_c3515, tempCorrected, pres);

% calculate density...
densCorrected = sw_pden(salinCorrected, temp, pres, 0);


% convert ptime into datenum style...
ptime_datenum = ptime/3600/24+datenum(1970, 1, 1, 0, 0, 0);
ptime_datenum_dbd = ptime_dbd/3600/24+datenum(1970, 1, 1, 0, 0, 0);
ptimehrly = fix(ptime_datenum*24)/24;
ptimedaily = fix(ptime_datenum);
ptimedaily2 = unique(ptimedaily);
ptimedaily2 = ptimedaily2(1:2:end);


% make copies of dbd time stamp vector for use in lat/lon interpolation...
ptime_dbd_gps = ptime_dbd;
ptime_dbd_wpt = ptime_dbd;

% convert lats and lons to digital degrees...
gpsLat = ddmm2decdeg(gpsLat);
gpsLon = ddmm2decdeg(gpsLon);
wptLat = ddmm2decdeg(wptLat);
wptLon = ddmm2decdeg(wptLon);

% eliminate outliers in gpsLat, gpsLon...
i = find(abs(gpsLat) <= 90.0);
gpsLat = gpsLat(i);  gpsLon = gpsLon(i);  ptime_dbd_gps = ptime_dbd_gps(i);
i = find(abs(gpsLon) <= 180.0);
gpsLat = gpsLat(i);  gpsLon = gpsLon(i);  ptime_dbd_gps = ptime_dbd_gps(i);

% eliminate outliers in wptLat, wptLon...
i = find(abs(wptLat) <= 90.0);
wptLat = wptLat(i);  wptLon = wptLon(i);  ptime_dbd_wpt = ptime_dbd_wpt(i);
i = find(abs(wptLon) <= 180.0);
wptLat = wptLat(i);  wptLon = wptLon(i);  ptime_dbd_wpt = ptime_dbd_wpt(i);

% eliminate entries where wptLat==0
i = find(wptLat~=0);
wptLat = wptLat(i);  wptLon = wptLon(i);  ptime_dbd_wpt = ptime_dbd_wpt(i);

% eliminate nans before interpolating...
i = find(~isnan(gpsLat));
gpsLat = gpsLat(i);  gpsLon = gpsLon(i);  ptime_dbd_gps = ptime_dbd_gps(i);
i = find(~isnan(wptLat));
wptLat = wptLat(i);  wptLon = wptLon(i);  ptime_dbd_wpt = ptime_dbd_wpt(i);

% interpolate DBD lat/lon data to align with EBD data...
gpsLat = interp1(ptime_dbd_gps, gpsLat, ptime);
gpsLon = interp1(ptime_dbd_gps, gpsLon, ptime);
wptLat = interp1(ptime_dbd_wpt, wptLat, ptime);
wptLon = interp1(ptime_dbd_wpt, wptLon, ptime);

% use sw_dpth() to calculate depth from pres...
depth = sw_dpth(pres, gpsLat);


% create configuration struct...
units = struct('chlor', 'micrograms liter-1',...
               'dens', 'kg m-3',...
               'densCorrected', 'kg m-3',...
               'depth', 'm',...
               'gpsLat', 'decimal degrees',...
               'gpsLon', 'decimal degrees',...
               'pres', 'decibars',...
               'ptime', 'seconds since 0000-01-01T00:00',...
               'ptime_datenum', 'seconds since 0000-01-01T00:00',...
               'salin', 'psu',...
               'salinCorrected', 'psu',...
               'temp', 'deg C');

variable_description = struct('chlor', 'chlorophyll measured by glider',...
                              'chlorBounds', 'default chlorophyll bounds for plots',...
                              'dens', 'density measured by glider',...
                              'densBounds', 'default density bounds for plots',...
                              'densCorrected', 'density corrected for thermal lag',...
                              'depth', 'depth calculated as function of pressure and position latitude',...
                              'gpsLat', 'position latitude measured by glider GPS',...
                              'gpsLon', 'position longitude measured by glider GPS',...
                              'pres', 'pressure measured by glider',...
                              'ptime', 'time vector reported by glider',...
                              'ptime_datenum', 'Serial Date Number string',...
                              'salin', 'salinity measured by glider',...
                              'salinBounds', 'default salinity bounds for plots',...
                              'salinCorrected', 'salinity corrected for thermal lag',...
                              'temp', 'temperature measured by glider',...
                              'tempBounds', 'default temperature bounds for plots');

correction_parameters = struct('alpha_offset', correctionParams(1),...
                               'alpha_slope', correctionParams(2),...
                               'tau_offset', correctionParams(3),...
                               'tau_slope', correctionParams(4));
           
config = struct('glider_name', strGliderName,...
                'deployment_number', strDeploymentNumber,...
                'start_date', strStartDate,...
                'end_date', strEndDate,...
                'thermal_lag_correction_parameters', correction_parameters,...
                'var_descriptions', variable_description,...
                'var_units', units);

% set mat file name...
strMatFileName = strcat(strGliderName, '_CTD_Deployment', strDeploymentNumber, '.mat');



% save glider/deployment data to mat file...
save(strMatFileName,...
     'config',...
     'chlor',...
     'chlorBounds',...
     'dens',...
     'densBounds',...
     'densCorrected',...
     'depth',...
     'gpsLat',...
     'gpsLon',...
     'pres',...
     'ptime',...
     'ptime_datenum',...
     'salin',...
     'salinBounds',...
     'salinCorrected',...
     'temp',...
     'tempBounds');



