close all
clear all

addpath(genpath('/Volumes/Chem/ARGO_PROCESSING/MFILES/FLOATS/'));
addpath(genpath('/Users/Mau/Documents/MBARIprojects/CPF/CPFPROCESSING'));
addpath(genpath('/Users/Mau/Documents/MATLABhelp_new/ARGO_PROCESSING-master/')); %
addpath(genpath('/Users/Mau/Documents/MATLABhelp_new/'));

%% PURPOSE: Loads all .mat data processed for a given mission (done in processCPFprofiles.m within the calc_BGC.m function)
% Makes plots to compare Optode temperature and CTD (SBE41CP) temperature.
% NOTE: Ken notes a lag of ~ 5m in the DOXY data, Josh points out there
% seem to be significant differences between the Optode and CTD temperature
% data

% Reprocessed all mission directories under /ParagonTests
CPF_ID='CPF003';
dirs.root=['/Volumes/CPF/',CPF_ID,'/Data/'];
dirs.mat=[dirs.root,'matfiles/'];
% folders=dir([dirs.mat,'ParagonTests/MICRORIDER/2021-11-16*']); % each folder corresponds to a mission identifiable by deployment date

dirs.root=['/Users/mau/Documents/MBARIprojects/CPF/',CPF_ID,'/Data/'];
dirs.mat=[dirs.root,'matfiles/'];
% folders=dir([dirs.mat,'ParagonTests/MICRORIDER/2021-11-16*']); % each folder corresponds to a mission identifiable by deployment date
% folders=dir([dirs.mat,'ParagonTests/2021-11-16*']); % each folder corresponds to a mission identifiable by deployment date
folders=dir([dirs.mat,'ParagonTests/2022-3-29*']); % each folder corresponds to a mission identifiable by deployment date

maxp=35; % maxp for plotting

nf=length(folders);

for f=1 %:nf;
    path_data=[folders(f).folder,'/',folders(f).name,'/'];
    dirs.temp=path_data;
    
    mat_list   = dir([dirs.temp,'*.mat']); % get list of file names to process
    nP=length(mat_list);
     
    HR=cell(nP,1);
    LR=cell(nP,1);
    INFO=cell(nP,1);
    cast=cell(1,nP);

    for p=1:nP
         mat_file=mat_list(p).name;
         cast{p}=regexp(mat_file,'\d+(?=\.mat)','once','match');

        D=load([path_data,mat_file]);
        HR{p}=D.HR;
        LR{p}=D.LR;
        INFO{p}=D.INFO;
    end
    
    % Extract date stamp from folder name
    ddate=regexp(path_data,['(?<year>\d{4})-(?<month>\d+)-(?<day>\d+)'],'names');
    dirs.date=[ddate.year,'-',ddate.month,'-',ddate.day];

end
% NOTE: should add a line to copy or save .mat files from atlas into my local directory, as
% plot_cpf_HRctd.m will attempt to read matfiles in my local directory

%% Plot coms for TEMP and optode TEMP
path_plots='/Users/mau/Documents/MBARIprojects/CPF/CPFPROCESSING/PLOTS/CTDoptode_coms/';
path_plots='/Volumes/CPF/CPF003/Data/data_processing/plots/CTDoptode_temp_coms/'

close all

cmap=cbrewer('div','Spectral',nP); % colorbar used to loop through profile data

for p=1:nP
plot(LR{p}.TEMP-LR{p}.TEMP_DOXY,LR{p}.PRES,'color',cmap(p,:)); hold on
axis ij
end

legend(cellstr(cast),'location','NE')

ylim([0 maxp])

ylabel('Pressure [dbar]')
xlabel('CTD - Toptode Temperature')
title(dirs.date)
print_setup(20,30)

print([path_plots,'CTD-Optode_Temperature_',dirs.date,'.jpg'],'-djpeg');

%% Plot gradient for CTD Temperature
close all

for p=1:nP
    deltaT_ctd=diff(LR{p}.TEMP);
    plot(deltaT_ctd,LR{p}.PRES(2:end),'color',cmap(p,:)); hold on
    axis ij
end
legend(cellstr(cast),'location','NE')

xlim([-0.2 1])
ylim([0 maxp])

ylabel('Pressure [dbar]')
xlabel('deltaT CTD')
title(dirs.date)
print_setup(20,30)

print([path_plots,'deltaT_CTD_',dirs.date,'.jpg'],'-djpeg');

%% Plot gradient for T Optode
close all

for p=1:nP
    deltaT_optode=diff(LR{p}.TEMP_DOXY);
    plot(deltaT_optode,LR{p}.PRES(2:end),'color',cmap(p,:)); hold on
    axis ij
end
legend(cellstr(cast),'location','NE')

xlim([-0.2 1])
ylim([0 maxp])

ylabel('Pressure [dbar]')
xlabel('deltaT optode')
title(dirs.date)
print_setup(20,30)

print([path_plots,'deltaT_OPTODE_',dirs.date,'.jpg'],'-djpeg');

%% Plot CTD-Optode Temp Difference against DeltaT/dz
% Plot analogous to Gene's
% solid line is the difference between the two
close all

cmap=cbrewer('div','Spectral',nP); % colorbar used to loop through profile data

for p=1:nP
plot(LR{p}.TEMP-LR{p}.TEMP_DOXY,LR{p}.PRES,'color',cmap(p,:)); hold on
end

deltaT_ctd=cell(nP,1);
for p=1:nP
deltaT_ctd{p}=diff(LR{p}.TEMP);
plot(deltaT_ctd{p},LR{p}.PRES(2:end),':','color',cmap(p,:)); hold on
axis ij
end

legend([cellstr(cast), cellstr(cast)],'location','SE')

ylim([0 maxp])

ylabel('Pressure [dbar]')
xlabel('CTD - Toptode Temperature')
title(dirs.date)
print_setup(20,30)

print([path_plots,'CTD-Optode_Temperature_Diff_and_DeltaTctd_',dirs.date,'.jpg'],'-djpeg');

%% Idem as above but the correlation
close all

cmap=cbrewer('div','Spectral',nP); % colorbar used to loop through profile data

for p=1:nP
    diffTs=LR{p}.TEMP-LR{p}.TEMP_DOXY;
    scatter(deltaT_ctd{p},diffTs(2:end),'color',cmap(p,:)); hold on
    [r{p},lags{p}]=xcorr(deltaT_ctd{p},diffTs(2:end));
end
legend([cellstr(cast)],'location','SE')

xlabel('deltaT/dz CTD')
ylabel('CTD - Toptode Temperature')
title(dirs.date)
print_setup(20,30)

% print([path_plots,'CTD-Optode_Temperature_Diff_and_DeltaTctd_',dirs.date,'.jpg'],'-djpeg');

%%
close all

subplot(121)
for p=1:nP
plot(LR{p}.TEMP,LR{p}.PRES,'-o','color',cmap(p,:)); hold on
axis ij
end
legend(cellstr(cast),'location','NW')

cmin=floor(nanmin(LR{p}.TEMP)); cmax=ceil(nanmax(LR{p}.TEMP));
xlim([cmin cmax])
ylim([0 maxp])

ylabel('Pressure [dbar]')
xlabel('CTD Temperature')

title(dirs.date)

subplot(122)
for p=1:nP
plot(LR{p}.TEMP_DOXY,LR{p}.PRES,'-o','color',cmap(p,:)); hold on
axis ij
end
ylabel('Pressure [dbar]')
xlabel('Optode Temperature')
legend(cellstr(cast),'location','NW')

xlim([cmin cmax])
ylim([0 maxp])

title(dirs.date)
print_setup(40,30)

print([path_plots,'CTDTemp_OptodeTemp_',dirs.date,'.jpg'],'-djpeg');

%% "Correcting" the missalignment
close all

for p=1:nP
plot(LR{p}.TEMP(1:end-1)-LR{p}.TEMP_DOXY(2:end),LR{p}.PRES(1:end-1),'color',cmap(p,:)); hold on
axis ij
end
legend(cellstr(cast),'location','NE')

ylim([0 maxp])

ylabel('Pressure [dbar]')
xlabel('CTD - Toptode Temperature')
title(dirs.date)

print_setup(20,30)

print([path_plots,'CTD-Optode_Temperature_',dirs.date,'_corrected.jpg'],'-djpeg');

%% Plot Oxygen data
% and grey out those  that do have anchor data
% NOTE: the "kink", increase in O2 in the first bottom sample only happens
% when the float was in "false anchor mode" during 2022-03-14 deployment

close all

for p=1:nP
plot(LR{p}.DOXY,LR{p}.PRES,'-o','color',cmap(p,:)); hold on
axis ij
cmin(p)=floor(nanmin(LR{p}.DOXY)); cmax(p)=ceil(nanmax(LR{p}.DOXY));
end

% cmin=floor(nanmin(cmin)); cmax=ceil(nanmax(cmax));
% xlim([cmin cmax])

legend([cellstr(cast(1:nP))],'location','SE')

ylim([0 maxp])

% % grey out profiles with anchor data
% panchor=[1 2 3 5 13];
% for p=1:5
% plot(LR{panchor(p)}.DOXY,LR{panchor(p)}.PRES,'-o','color',[0.4 0.4 0.4]); hold on
% axis ij
% end
% 
% legend([cellstr(cast(1:nP)),'001','002','003','005','013'],'location','SE')

ylabel('Pressure [dbar]')
xlabel('Optode DOXY')

title(dirs.date)
print_setup(30,40)

print([path_plots,'DOXY_Optode_',dirs.date,'.jpg'],'-djpeg');
