%% Cleaning Messy Tension Data
% This example demonstrates the typical steps for cleaning and averaging
% messy uniaxial tension data that might occur when measuring at the low
% end of a load cell and/or when measuring challenging materials like a 
% non-woven sheet.
% 
% Please note: The comments in this example are truly kept to a minimum, 
% unlike the other Kornucopia Workshops which have large blocks of 
% descriptive comments. The intent with this workshop example is to more
% accurately mimic a typical script of a user. 
% (Even the comments here are probably excessive!)

%% Set-up
k_cleanup();
adv = k_adv();
k_unitsPreferenceActivate('mm_N_s');
% k_unitsVariables('N');
k_set('inheritPropsAutoWarning', 'off');

figH = figure('Name', 'Cleaning Messy Tension Data', ...
    'pos', k_figZoom(2));

tabStuff = k_figTabsCreate(10, figH);

%% Overview Tab

currentTabH = tabStuff.TabsH(1);
k_figTabsClear(currentTabH);
currentTabH.Title = 'Overview';

dataDir = '../data';
picFileName = 'test_setup.png';

k_picDisplay(fullfile(dataDir, picFileName), 'north', ...
    'parent', currentTabH, ...
    'nudge', [0, -0.03], ...
    'TitleText', 'Test Setup');

txt = {
    'This is a tension test of Material polyP.'
    'The test was run with initial slack in each sample.'
    ' '
    'Note: The test used a laser for displacement measurement,'
    'but the laser is not depicted in this picture. '
    ' '
    'This script reads the data and cleans it up.'
    };

k_displayOnFigure(0.5, 0.3, txt, ...
    'parent', currentTabH, ...
    adv.k_displayOnFigure.box.on, ...
    adv.k_displayOnFigure.alignVertical.top);
    
k_figTabsDisplay(currentTabH);


%% Find and load the raw data files

fNames = k_dirList(fullfile(dataDir, '*.txt'), 'files');
nFiles = numel(fNames);

raw = cell(nFiles,1);

% For file importing 
delim = ',';
startStr = 'Specimen Dimensions';
rs1 = 1;
colHeadStr = 'Time';
rs2 = 0;
endStr = [];
rs3 = [];

for z = 1 : nFiles
    fName = fullfile(dataDir, fNames{z});
    tmp = k_readText(fName,delim,startStr,rs1,colHeadStr,rs2,endStr,rs3);
    raw{z} = tmp;
end

%% Plot all the data on 1 plot

currentTabH = tabStuff.TabsH(2);
k_figTabsClear(currentTabH);
currentTabH.Title = 'Raw Data';

k_plot(raw, [], ...
    'parent', currentTabH, ...
    'independ', 'Laser', ...
    'depend', 'Load', ...
    'grid', 'both', ...
    'lineWidth', 1);

k_figTabsDisplay(currentTabH);

%% Initial Cleanup of the Raw data

cleaned = cell(nFiles,1);
smthFac = 0.95;
fLow = '1*N';
for z = 1 : nFiles
    tmp = raw{z};
    
    % Trim at fLow and the max for each curve
    idx = k_trimIndices(tmp('Load'), fLow, max(tmp('Load')));
    trimmed = tmp(idx,:);
    
    % Smooth the curve
    smth = k_smooth(trimmed, [], smthFac); 
    
    % Heal the curve back to zero Load
    %{
    healed = k_tweakXY(smth, [], 'linear', 0.1, 'start', smth('Load'), '0*N', ...
        adv.k_tweakXY.return.plotAndTweaked); 
    %}
    healed = k_tweakXY(smth, [], 'linear', 0.1, 'start', smth('Load'), '0*N');    
    
    % Record the displacement offset, then reset Time and Laser to zero values.
    disp_offset = healed(1, 'Laser');
    healed = healed - healed(1, :);
    healed.Props.UserData.disp_offset = disp_offset;
    
    cleaned{z} = healed;
end

% Plot the cleaned data
currentTabH = tabStuff.TabsH(3);
k_figTabsClear(currentTabH);
currentTabH.Title = 'Cleaned Data';

aH = subplot(2,1,1, 'parent', currentTabH);
k_plot(raw, [], ...
    'parent', aH, ...
    'independ', 'Laser', ...
    'depend', 'Load', ...
    'grid', 'both', ...
    'lineWidth', 1, ...
    'title', {'user', 'Raw Data'});

aH = subplot(2,1,2, 'parent', currentTabH);
k_plot(cleaned, [], ...
    'parent', aH, ...
    'independ', 'Laser', ...
    'depend', 'Load', ...
    'grid', 'both', ...
    'lineWidth', 1, ...
    adv.k_plot.legendLocation.SouthEast, ...
    'title', {'user', 'Cleaned Data'});

k_figTabsDisplay(currentTabH);

%% Convert Force vs. Disp to Stress vs Strain

nomSS = cell(nFiles,1);
for z = 1 : nFiles
    tmp = cleaned{z};
    
    % unpack the header to get sample geom
    geom = k_unpackHeader(tmp.Props.UserData.DatasetHeader, ':');
    
    % Compute an adjusted gage length
    % geom('adj gage length') = geom('gage length') + tmp.Props.UserData.disp_offset;
    geom{'adj gage length'} = geom('gage length') + tmp.Props.UserData.disp_offset;
    area = geom('thickness') * geom('width');
    
    nomSS_tmp = [tmp('Laser') / geom('adj gage length'), tmp('Load') / area];
    nomSS_tmp.Props.ColNames = 'nomStrain, nomStress';
    nomSS_tmp.Props.Description = tmp.Props.Description;
    nomSS_tmp('nomStrain') = nomSS_tmp('nomStrain').convert('%');
    
    nomSS{z} = nomSS_tmp;    
end

currentTabH = tabStuff.TabsH(4);
k_figTabsClear(currentTabH);
currentTabH.Title = 'Nom SS Curves';

k_plot(nomSS, [], ...
    'parent', currentTabH, ...
    'grid', 'both', ...
    'lineWidth', 1);

txt = {
    'In a more sophisticated analysis, this data would be'
    'converted to True Stress vs Strain via a constant volumed'
    'deformation assumption or some other appropriate transform.'
    };

k_displayOnFigure(0.6, 0.25, txt, ...
    'box', 'on');

% Note for above, without a 'parent' ADV , the function puts
% the text in the last "hot" plot axes. 

k_figTabsDisplay(currentTabH);

%% Compute Average Curve

% Compute ave curve with default ADV options
ave_a = k_averageXY(nomSS, [], 'x');
ave_a.Props.Description = 'Default ADV options';

% Use advanced extension methods on the ave calc
ave_b = k_averageXY(nomSS, [], 'x', ...
    adv.k_averageXY.extend_to.aveLastVal, ...
    adv.k_averageXY.extend_using.p2, ...
    adv.k_averageXY.extend_refLengthRegion.extendRegion, ...
    adv.k_averageXY.extend_refLengthFactor, 0.25);
ave_b.Props.Description = 'Extended Ave via ADV options';

nomSS_all = k_stackForPlot(nomSS);

currentTabH = tabStuff.TabsH(5);
k_figTabsClear(currentTabH);
currentTabH.Title = 'Ave of Nom SS Curves';


toPlt = {
    nomSS_all
    ave_a
    ave_b
    };

legendTxt = {
    'All the nominal curves'
    'Ave with default ADV options'
    'More elaborate Extended Ave'
    };

k_plot(toPlt, [], ...
    'parent', currentTabH, ...
    'grid', 'both', ...
    'lineWidth', [0.5,2,1.5], ...
    'markerSymb', 'b,o,x', ...
    'markerOccur', {'total', 8}, ...
    adv.k_plot.markersAtEndPoints.on, ...    
    'legendText', {'user', legendTxt});

k_figTabsDisplay(currentTabH);

%% Isolate the Initial Slope to Estimate a Modulus

currentTabH = tabStuff.TabsH(6);
k_figTabsClear(currentTabH);
currentTabH.Title = 'Modulus Determination';

idx = k_trimIndices(ave_b('nomStress'), [], '4*MPa');
tmp = ave_b(idx, :);

fitFunc = k_fitFunc(tmp, [], 'CP1', ...
    'parent', currentTabH);

Eo = fitFunc.FitConstants('C1');

txt = ['The modulus is Eo = ', num2str(Eo)];
    
k_displayOnFigure(0.75, 0.15, txt, ...
    'parent', currentTabH, ...
    'box', 'on', ...
    'fontSize', 12);

k_figTabsDisplay(currentTabH);

