clear all
close all

opts = delimitedTextImportOptions("NumVariables", 30);

% Specify range and delimiter
opts.DataLines = [25, Inf];
opts.Delimiter = ",";

% Specify column names and types
%opts.VariableNames = ["LabviewDateTime", "ThermistorTemp", "BathTemp", "BathSetPoint", "PressuredBar", "MaxFlow", "FlowRate", "PumpVolume", "CRC", "ID", "DATETIME", "UTCSECS", "CTDDEPTH", "CTDTEMP", "CTDSALINITY", "SAMPLECOUNTER", "RESETCOUNTER", "ERRORCOUNTER", "HOUSINGTEMPERATURE", "HOUSINGHUMIDITY", "INPUTVOLTAGE", "INPUTCURRENT", "FOOBARPHVALUE", "BACKUPBATTERY", "VRSVOLTS", "VRSSTDDEV", "VKVOLTS", "VKSTDDEV", "IKAMPS", "IBAMPS"];
%opts.VariableNames = ["dtLV",   "thermTemp", "bathTemp", "bathSP", "pressureDBar", "maxFlow", "flowRate", "cycleNum", "pumpVolume", "CRC", "ID", "dateTime", "utcSecs", "ctdDepth", "ctdTemp", "ctdSalinity", "sampleCounter", "resetCounter", "errorCounter", "housingTemp", "housingHumidity", "inputVoltage", "inputCurrent", "foobarPHValue", "backupBattery", "vrsVolts", "vrsStdDev", "vkVolts", "vkStdDev", "ikAmps", "ibAmps"];
%opts.VariableTypes = ["datetime", "double", "double",   "double", "double",         "double", "double",     "double",   "double", "string", "categorical", "datetime", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double"];
%For the new VI with SetTemp, SetPressure and StepNum
opts.VariableNames = ["dtLV",   "thermTemp", "bathTemp", "bathSP", "pressureDBar", "maxFlow", "flowRate", "pumpVolume", "CRC", "ID", "dateTime", "utcSecs", "ctdDepth", "ctdTemp", "ctdSalinity", "sampleCounter", "resetCounter", "errorCounter", "housingTemp", "housingHumidity", "inputVoltage", "inputCurrent", "foobarPHValue", "backupBattery", "vrsVolts", "vrsStdDev", "vkVolts", "vkStdDev", "ikAmps", "ibAmps", "SetTemp", "SetPressure", "StepNum"];
opts.VariableTypes = ["datetime", "double", "double",   "double", "double",         "double", "double",    "double", "string", "categorical", "datetime", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double", "double"];

% Specify file level properties
opts.ExtraColumnsRule = "ignore";
opts.EmptyLineRule = "read";

% Specify variable properties
opts = setvaropts(opts, "CRC", "WhitespaceRule", "preserve");
opts = setvaropts(opts, ["CRC", "ID"], "EmptyFieldRule", "auto");
%opts = setvaropts(opts, "dtLV", "InputFormat", "MM/dd/uuuu hh:mm:ss aa");
%opts = setvaropts(opts, "dateTime", "InputFormat", "MM/dd/uuuu hh:mm:ss aa");

%% Process files
path = '\\atlas\ProjectLibrary\ChemSensors\GeneM\GDF pH Sensor\Data\Cal Data\Calcs\DataFiles\'
fullPath = strcat(path, '*.csv');

dataFiles = dir(strcat(path, '*.csv'));
fileDatenums = dataFiles.date;
[sortedTime, i] = sort([dataFiles.datenum]);
sortedDataFileList = dataFiles(i);

for i=1:length(dataFiles)
    % Import the data
    disp("Processing: " + dataFiles(1).name);
    dataTable{i} = readtable(strcat(path, sortedDataFileList(i).name), opts);
end

% Create overview plot
createOverviewPlot2(dataTable, sortedDataFileList, 1);

% Separate out cycles
[cycle, numCycles] = getCycles(dataTable, "pressure");

% Waterfall plot
plotPressureWaterfall(cycle, numCycles)

figure('name', "Vrs vs Time", 'NumberTitle', 'off');
plot(dataTable{1}.dtLV, dataTable{1}.vrsVolts)
title('time vs Vrs')

%Plot individual cycles
numBins = 100;
xIncr = 0.2;
yIncr = 0.3;
left = 0.05;
bottom = 0.65;
width = 0.25;
height = 0.25;

polyOrder = 3;

for i=2:2:50
    windowText = sprintf("Cycle Num: %i", i/2);
    figure('name', windowText, 'NumberTitle', 'off', 'Unit', 'normalized','Position', [left, bottom, width, height]);
    hold on

    indexes = find(dataTable{1}.StepNum == i);
    x1 = dataTable{1}.pressureDBar(indexes);
    y1 = dataTable{1}.vrsVolts(indexes);

    lastwarn('')
    polyFitVals = polyfit(x1, y1, polyOrder);
    [warnMsg, warnId] = lastwarn;
    if isempty(warnMsg)
    
    yFit = polyval(polyFitVals, x1);
    maxResidual = max(abs(yFit - y1));
    meanResidual = mean(abs(yFit - y1));
    plot(x1, y1, x1, yFit);
    titleString1 = sprintf('Incr P fit residual mean abs diff: %.3e max: %.3e\n', meanResidual, maxResidual);

    indexes = find(dataTable{1}.StepNum == i+1);
    x2 = dataTable{1}.pressureDBar(indexes);
    y2 = dataTable{1}.vrsVolts(indexes);
    polyFitVals = polyfit(x2, y2, polyOrder);
    yFit = polyval(polyFitVals, x2);
    maxResidual = max(abs(yFit - y2));
    meanResidual = mean(abs(yFit - y2));
    plot(x2, y2, x2, yFit);
    titleString2 = sprintf('Decr P fit residual mean abs diff: %.3e max: %.3e\n', meanResidual, maxResidual);

    x3 = [x1; x2];
    y3 = [y1; y2];
    polyFitVals = polyfit(x3, y3, polyOrder);
    yFit = polyval(polyFitVals, x3);
    maxResidual = max(abs(yFit - y3));
    meanResidual = mean(abs(yFit - y3));
    plot(x3, yFit);
    titleString3 = sprintf('All data fit residual mean abs diff: %.3e max: %.3e\n', meanResidual, maxResidual);

    legend(titleString1, 'incr fit', titleString2, 'decr fit', titleString3, 'FontSize', 16, 'location', 'northeast');
    xlabel('Pressure (dBar)');
    ylabel('Vrs (volts)');

    left = left + xIncr;
    if(mod(i, 10) == 0)
        bottom = bottom - yIncr;
        left = 0.05;
    end
    else
        disp 'skipping this cycle'
    end
end


% 
% 
% for i=1:2:length(cycle)-1
%     cycleNum = i;
% 
%     [indexes, edges] = discretize(cycle{1,cycleNum}.pressureDBar, numBins);
%     for j=1:numBins
%         pressureMeanIncr(j) = mean(cycle{1,cycleNum}.pressureDBar(indexes(j)));
%         vrsMeanIncr(j) = mean(cycle{1,cycleNum}.vrsVolts(indexes(j)));
%     end
% 
%     cycleNum = cycleNum + 1;
%     [indexes, edges] = discretize(cycle{1,cycleNum}.pressureDBar, numBins);
%     for j=1:numBins
%         pressureMeanDecr(j) = mean(cycle{1,cycleNum}.pressureDBar(indexes(j)));
%         vrsMeanDecr(j) = mean(cycle{1,cycleNum}.vrsVolts(indexes(j)));
%     end
%     hyst = vrsMeanIncr - vrsMeanDecr;
% 
%     windowText = sprintf("Cycle Num: %i", i/2 + 0.5);
%     figure('name', windowText, 'NumberTitle', 'off', 'Unit', 'normalized','Position', [left, bottom, width, height]);
%     left = left + xIncr;
%     bottom = bottom - yIncr;
%     plot(pressureMeanIncr, hyst, 'g');
%     ax = gca;
%     ax.YAxis.Exponent = -3;
%     ylabel('Vrs Incr - Vrs Decr')
%     hold on
%     
%     yyaxis right
%     plot(pressureMeanIncr, vrsMeanIncr, '-+b', pressureMeanDecr, vrsMeanDecr, '-or')
% %     plot(cycle{cycleNum - 1}.pressureDBar, cycle{cycleNum - 1}.vrsVolts, '-y')
% %     plot(cycle{cycleNum}.pressureDBar, cycle{cycleNum}.vrsVolts, '-c')
% 
%     ylabel('Vrs')
%     xlabel('Binned Pressure (dBar)')
%     legend('Hysteresis', 'Binned Increasing P', 'Binned Decreasing P')
% end
