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\ProcessThese\'
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(i).name);
    dataTable{i} = readtable(strcat(path, sortedDataFileList(i).name), opts);

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

    titleString = "Vrs vs Time " + sortedDataFileList(i).name;
    figure('name', titleString, 'NumberTitle', 'off');
    plot(dataTable{i}.dtLV, dataTable{i}.vrsVolts)
    title(titleString)

    % Separate out cycles
    [cycle, numCycles] = getCyclesFromPressure(dataTable, i);
    cycleData{i} = cycle;

    % Waterfall plot
    plotPressureWaterfall(cycle, numCycles)
    title(sortedDataFileList(i).name)

    if(0)
        %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 j=2:2:50
            cycleIndex = j/2;
            windowText = sprintf("Cycle Num: %i", j/2);
            figure('name', windowText, 'NumberTitle', 'off', 'Unit', 'normalized','Position', [left, bottom, width, height]);
            hold on

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

            %Get rid of data from pressure holds
            dPressure = diff(x1);x1
            holdIndexes = find(dPressure < 1);
            x1(holdIndexes) = [];
            cycle{cycleIndex}.xIncr = x1;
            y1(holdIndexes) = [];
            cycle{cycleIndex}.yIncr = y1;

            lastwarn('')
            [cycle{cycleIndex}.polyFitValsIncr, cycle{cycleIndex}.SIncr] = polyfit(cycle{cycleIndex}.xIncr, cycle{cycleIndex}.yIncr, polyOrder);
            [p,S] = polyfit(cycle{cycleIndex}.xIncr, cycle{cycleIndex}.yIncr, polyOrder);
            [warnMsg, warnId] = lastwarn;
            if isempty(warnMsg)

                [cycle{cycleIndex}.yFitIncr, cycle{cycleIndex}.deltaIncr] = polyval(cycle{cycleIndex}.polyFitValsIncr, cycle{cycleIndex}.xIncr, cycle{cycleIndex}.SIncr);
                [yFit, delta] = polyval(p, cycle{cycleIndex}.xIncr, S)
                cycle{cycleIndex}.maxResidualIncr = max(abs(cycle{cycleIndex}.yFitIncr - cycle{cycleIndex}.yIncr));
                cycle{cycleIndex}.meanResidualIncr = mean(abs(cycle{cycleIndex}.yFitIncr - cycle{cycleIndex}.yIncr));
                plot(cycle{cycleIndex}.xIncr, cycle{cycleIndex}.yIncr, '+', cycle{cycleIndex}.xIncr, cycle{cycleIndex}.yFitIncr);
                titleString1 = sprintf('Incr P fit mean(abs(residual)): %.3e max(abs(residual)): %.3e\n', cycle{cycleIndex}.meanResidualIncr, cycle{cycleIndex}.maxResidualIncr);

                indexes = find(dataTable{1}.StepNum == i+1);
                x2 = dataTable{i}.pressureDBar(indexes);
                y2 = dataTable{i}.vrsVolts(indexes);

                %Get rid of data from pressure holds
                dPressure = abs(diff(x2));
                holdIndexes = find(dPressure < 1);
                x2(holdIndexes) = [];
                cycle{cycleIndex}.xDecr = x2;
                y2(holdIndexes) = [];
                cycle{cycleIndex}.yDecr = y2;

                [cycle{cycleIndex}.polyFitValsDecr, cycle{cycleIndex}.SDecr] = polyfit(cycle{cycleIndex}.xDecr, cycle{cycleIndex}.yDecr, polyOrder);
                [cycle{cycleIndex}.yFitDecr, cycle{cycleIndex}.deltaDecr] = polyval(cycle{cycleIndex}.polyFitValsDecr, cycle{cycleIndex}.xDecr, cycle{cycleIndex}.SDecr);
                cycle{cycleIndex}.maxResidualDecr = max(abs(cycle{cycleIndex}.yFitDecr - cycle{cycleIndex}.yDecr));
                cycle{cycleIndex}.meanResidualDecr = mean(abs(cycle{cycleIndex}.yFitDecr - cycle{cycleIndex}.yDecr));
                plot(cycle{cycleIndex}.xDecr, cycle{cycleIndex}.yDecr, '+', cycle{cycleIndex}.xDecr, cycle{cycleIndex}.yFitDecr);
                titleString2 = sprintf('Decr P fit mean(abs(residual)): %.3e max(abs(residual)): %.3e\n', cycle{cycleIndex}.meanResidualDecr, cycle{cycleIndex}.maxResidualDecr);

                cycle{cycleIndex}.xBoth = [cycle{cycleIndex}.xIncr; cycle{cycleIndex}.xDecr];
                cycle{cycleIndex}.yBoth = [cycle{cycleIndex}.yIncr; cycle{cycleIndex}.yDecr];
                cycle{cycleIndex}.polyFitValsBoth = polyfit(cycle{cycleIndex}.xBoth, cycle{cycleIndex}.yBoth, polyOrder);
                cycle{cycleIndex}.yFitBoth = polyval(cycle{cycleIndex}.polyFitValsBoth, cycle{cycleIndex}.xBoth);
                cycle{cycleIndex}.maxResidualBoth = max(abs(cycle{cycleIndex}.yFitBoth - cycle{cycleIndex}.yBoth));
                cycle{cycleIndex}.meanResidualBoth = mean(abs(cycle{cycleIndex}.yFitBoth - cycle{cycleIndex}.yBoth));
                plot(cycle{cycleIndex}.xBoth, cycle{cycleIndex}.yFitBoth);
                titleString3 = sprintf('All data fit mean(abs(residual)): %.3e max(abs(residual)): %.3e\n', cycle{cycleIndex}.meanResidualBoth, cycle{cycleIndex}.maxResidualBoth);

                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'
                cycle{cycleIndex}.yFitIncr = NaN;
                cycle{cycleIndex}.maxResidualIncr = NaN;
                cycle{cycleIndex}.meanResidualIncr = NaN;

                cycle{cycleIndex}.xDecr = NaN;
                cycle{cycleIndex}.yDecr = NaN;
                cycle{cycleIndex}.polyFitValsDecr = NaN;
                cycle{cycleIndex}.yFitDecr = NaN;
                cycle{cycleIndex}.maxResidualDecr = NaN;
                cycle{cycleIndex}.meanResidualDecr = NaN;

                cycle{cycleIndex}.xBoth = NaN;
                cycle{cycleIndex}.yBoth = NaN;
                cycle{cycleIndex}.polyFitValsBoth = NaN;
                cycle{cycleIndex}.yFitBoth = NaN;
                cycle{cycleIndex}.maxResidualBoth = NaN;
                cycle{cycleIndex}.meanResidualBoth = NaN;
            end
        end
    end
end
