%%% This is my test file for making sure all functions work together nicely
clear all; close all;
%% Read in data from sensors and specto pH

% sensorfile = '20240221.txt';
% specpHfile = '20240222_Nano_mFET_array.txt';
% CTDfile = '2024-02.csv';
%% 
sensorfile = '20240306.txt';
specpHfile = '20240222_Nano_mFET_array.txt';
CTDfile = '2024-02-03_combined.csv';
%% 
DataStart(sensorfile); % edits the file to be readable by matlab
MetaStart(sensorfile);
SensorData = readdata('datafile.txt'); % reads xfets and creates struct
MetaData = readmeta('metafile.txt');
delete('datafile.txt') % removes file created by DataStart
delete('metafile.txt') % removes file created by DataStart
SpecpHData = readspecpH(specpHfile); 
CTD = readCTD(CTDfile);
%% 
[~, iCTD, idata] = intersect(CTD.utc_time,SensorData.SampTime1);
SensorData = SensorData(idata,:);
CTD = CTD(iCTD,:);
[~ ,~ , ispec] = intersect(CTD.utc_time,SpecpHData.DateTime);
SpecpHData = SpecpHData(ispec,:);
%% Explore sensor data
EDA(SensorData)
WarmUpTime(SensorData)
Compare(SensorData)

%savefigs % Save nice figures as pdfs


%% Calculate pH Insitu

pH_insitu = spec_pH_insitu(SpecpHData);

%% Calculate sensor pH measurements
% TC_Cont1 is controller temp not insitu. Pull data from 
% http://pubdata.mlml.calstate.edu/mlml_last/seawater/
spH = pH_insitu(SpecpHData.QC);

k0_1 = nan(length(SensorData.Vrse1),2);
for i = 1:length(SpecpHData.pH_tot(SpecpHData.QC))
    [k0_1(:,i), ~] = calc_k0_ext(SensorData.Vrse1, CTD.temp,...
        CTD.sal, spH(i), -0.0010626); %mFET
end
k0_1_std = std(k0_1);
k0_1 = mean(k0_1,'all');

[pHext_tot1, pHext_free1] = calc_dfet_pHext_k0(SensorData.Vrse1,...
    CTD.temp, CTD.sal, k0_1, -0.0010626); % 
% I need a better way to handle this...

%Maybe do a line where you compare the insitu sal from MLML vs our measured
%densiy meter salinity we took in lab. Would help with QF assesment once
%you have more samples.
%%
k0_2 = nan(length(SensorData.Vrse2),2);
for i = 1:length(SpecpHData.pH_tot(SpecpHData.QC))
    [k0_2(:,i), ~] = calc_k0_ext(SensorData.Vrse2, CTD.temp,...
        CTD.sal, spH(i), -0.000963044); %Nano
end
k0_2_std = std(k0_2);
k0_2 = mean(k0_2,'all');

[pHext_tot2, pHext_free2] = calc_dfet_pHext_k0(SensorData.Vrse2,...
    CTD.temp, CTD.sal, k0_2, -0.000963044); % 

[~ , isame, ~] = intersect(SensorData.SampTime1,SpecpHData.DateTime);
meandiffmfet = mean(SpecpHData.pH_tot - pHext_tot1(isame));
meandiffnano = mean(SpecpHData.pH_tot - pHext_tot2(isame));
%% Plot sensor pH and specto pH

% pH = pHPlots(SensorData, SpecpHData, pHext_tot1, pHext_tot2);

% print('pHPlots','-dpdf','-fillpage')

pH = figure;
clf

tiledlayout(3,1)

nexttile
plot(SensorData.SampTime1,pHext_tot1)
hold on
scatter(SpecpHData.DateTime(SpecpHData.QC),SpecpHData.pH_tot(SpecpHData.QC),'filled')
ylabel('pH')
legend('mfet pH','Lab Spec pH',Location='southeast')
title('mfet')
subtitle(['k0: ',num2str(round(k0_1,3)),' +/- ', num2str(round(k0_1_std(1),3)),...
    ' Mean Diff: ', num2str(meandiffmfet)])
% add mean diff, k0, k0 std, error bars
nexttile
plot(SensorData.SampTime2,pHext_tot2)
hold on
scatter(SpecpHData.DateTime(SpecpHData.QC),SpecpHData.pH_tot(SpecpHData.QC),'filled')
ylabel('pH')
legend('nano pH', 'Lab Spec pH',Location='southeast')
title('nano')
subtitle(['k0: ',num2str(round(k0_2,3)),' +/- ', num2str(round(k0_2_std(1),3)),...
    ' Mean Diff: ', num2str(meandiffnano)])
nexttile
scatter(SensorData.SampNum1,pHext_tot1-pHext_tot2)
title('pH Residuals (mfet-nano)')

%%
pHScatter = figure;

tiledlayout(1,1)
axis('square')
nexttile
yyaxis left
scatter(SpecpHData.pH_tot,pHext_tot1(isame))
ylabel('mfet pH')
hold on
yyaxis right
scatter(SpecpHData.pH_tot,pHext_tot2(isame))
ylabel('nano pH')
grid on
xlabel('Spec pH')

title('mfet & nano pH vs spec pH')

legend(['mfet, k0: ',num2str(round(k0_1,3)),' +/- ', num2str(round(k0_1_std(1),3)),...
    ', Mean Diff: ', num2str(meandiffmfet), ' Sensor: T51'], ['nano, k0: ',num2str(round(k0_2,3)),' +/- ', num2str(round(k0_2_std(1),3)),...
    ', Mean Diff: ', num2str(meandiffnano),', Sensor: DF402'],'Location','northwest');