%% Glider Mission 21903401
gdata = readtable("Diel_GPP_log.csv");
i = 6;

filename = char(gdata.Filename(i));
load(filename, 's'); % .mat file

idiel = gdata.StartProf(i):gdata.EndProf(i);
iupper = 1:5;
pH = mean(s.pH25atm(iupper,idiel),1,"omitmissing");
TempC = mean(s.tc(iupper,idiel),1,"omitmissing");
Sal = mean(s.psal(iupper,idiel),1,"omitmissing");
s.o2satconc = reshape(calcO2sat(s.tc(:), s.psal(:)), size(s.pres)); % Is this the correct method for this calc?
s.o2anom = s.doxy - s.o2satconc;
meano2anom = mean(s.o2anom(iupper,idiel),1);
sdn = s.sdn(idiel)-7/24; % GMT to local

[diel_pH, SDNs] = parse_to_diel(sdn, pH, 24);
[diel_TempC, ~] = parse_to_diel(sdn, TempC, 24);
[diel_Sal, ~] = parse_to_diel(sdn, Sal, 24);
[diel_meano2anom, ~] = parse_to_diel(sdn, meano2anom, 24);
LegendNames = cellstr(datestr(SDNs, 'mm/dd/yy'));

fig1 = figure('Visible','on');
clf
tl = tiledlayout(4,1);
title(tl,['Glider Mission ID: ', num2str(gdata.MissionID(i))])
subtitle(tl,{['pH Diel Cycle ', char(gdata.StartDate(i)), ' to ', char(gdata.EndDate(i))],...
    ['Mean Layer: ', num2str(abs(s.depth(iupper(1)))), ' m - ', num2str(abs(s.depth(iupper(end)))), ' m']})

nexttile
plot(sdn,pH,'.')
datetick('x',6,'keeplimits')
ylabel('pH 25C')
xlabel('date')
title('pH full time series')

nexttile
plot(sdn,TempC,'.')
datetick('x',6,'keeplimits')
ylabel('Temp (C)')
xlabel('date')
title('Temperature full time series')

nexttile
plot(sdn,Sal,'.')
datetick('x',6,'keeplimits')
ylabel('Sal (PSU)')
xlabel('date')
title('Salinity full time series')

nexttile
plot(sdn,meano2anom,'.')
datetick('x',6,'keeplimits')
ylabel('Sal (PSU)')
xlabel('date')
title('mean o2 anom full time series')

%% Diels
clf
plot(hour(sdn),mean(pH,1,"omitmissing"),'.')

%% Do some parsing

idiel_QC = gdata.StartProf(i)+8:gdata.EndProf(i)-40; 
pH = mean(s.pH25atm(iupper,idiel_QC),1,"omitmissing");
mn = movmean(pH,5);
sd = std(pH);
% ii_QC = pH > mn - 2*sd & pH < mn + 2* sd;
pH_QC = pH; %(ii_QC); % Filter eroneous measurements
sdn_QC = s.sdn(idiel_QC)-7/24; % GMT to local


[diel, SDNs] = parse_to_diel(sdn_QC, pH_QC, 24);
LegendNames = cellstr(datestr(SDNs, 'mm/dd/yy'));



diel_detrend = detrend(diel,1,'omitmissing');
d = reshape(diel_detrend,1,[]);




fig1 = figure('Visible','on');
clf
tl = tiledlayout(3,1);
title(tl,['Glider Mission ID: ', num2str(gdata.MissionID(i))])
% subtitle(tl,{['pH Diel Cycle ', char(gdata.StartDate(i)), ' to ', char(gdata.EndDate(i))],...
    % ['Mean Layer: ', num2str(abs(s.depth(iupper(1)))), ' m - ', num2str(abs(s.depth(iupper(end)))), ' m']})

nexttile
plot(sdn_QC,pH_QC,'.')
datetick('x',6,'keeplimits')
ylabel('pH 25C')
xlabel('date')
title('pH full time series')

nexttile
plot(0:23,diel_detrend,'.')
legend(LegendNames,Location='eastoutside')
ylabel('pH 25C')
xlabel('Local Time hrs')
xticks([0 6 12 18 23])
xticklabels({'00:00', '6:00', '12:00', '18:00', '00:00'});
xlim([0 23])
title('diel cycle by day')

nexttile
plot(0:23,mean(diel_detrend,1,"omitmissing"),'.')
ylabel('pH 25C')
xlabel('Local Time hrs')
title('mean diel cycle')
xlim([0 23])
xticks([0 6 12 18 23])
xticklabels({'00:00', '6:00', '12:00', '18:00', '00:00'});