-
Notifications
You must be signed in to change notification settings - Fork 0
/
tendonModelFitting.m
113 lines (92 loc) · 4.39 KB
/
tendonModelFitting.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
datapath = 'C:\Users\Nick\Documents\Courses\Mechanics of Biological Tissues\Final Project\tendon-model-fitting\data';
resultspath = 'C:\Users\Nick\Documents\Courses\Mechanics of Biological Tissues\Final Project\tendon-model-fitting\results';
boxpath = 'C:\Users\Nick\Box Sync\tendon-model-fitting';
% Original sample numbers from instron tests
mapping = struct();
mapping.control.S1 = 'Sample 18';
mapping.control.S2 = 'Sample 22';
mapping.control.S3 = 'Sample 26';
mapping.static.S1 = 'Sample 14';
mapping.static.S2 = 'Sample 19';
mapping.static.S3 = 'Sample 24';
mapping.dynamic.S1 = 'Sample 17';
mapping.dynamic.S2 = 'Sample 21';
mapping.dynamic.S3 = 'Sample 25';
% Set conditions, samples, and strain rates
conditions = { {'control',{'S1','S2','S3'}} , ...
{'static' ,{'S1','S2','S3'}} , ...
{'dynamic',{'S1','S2','S3'}} };
test = {'low','high'};
% Create data and results structs
sdat = sampledata(datapath);
data = ramp(datapath,conditions);
results = struct();
save(fullfile(boxpath,'data','sampledata.mat'),'sdat')
save(fullfile(boxpath,'data','trialdata.mat'),'data')
try load(fullfile(resultspath,'results.mat'))
catch ME
control = struct();
save(fullfile(resultspath,'results.mat'),'control')
end
% Set max isometric force
Fmax = 60; % N
for c = 1:length(conditions)
cond = conditions{c}{1};
sample = conditions{c}{2};
for s = 1:length(sample)
for t = 1:length(test)
results.(cond).(sample{s}).mapping = mapping.(cond).(sample{s});
% Sample data
lo = sdat.(cond).(sample{s}).lo;
A = sdat.(cond).(sample{s}).Aell;
% Trial data
l = data.(cond).(sample{s}).(test{t}).ext;
F = data.(cond).(sample{s}).(test{t}).load;
% Recover tendon slack length estimation
if ~strcmp(cond,'control') && strcmp(test{t},'low')
results.(cond).(sample{s}).lTs_est = l(1);
end
% Tendon measurements to fit
l = l-l(1);
lmda = (l/lo)+1; % strain
F_tilde = F/Fmax; % normalized force
T = F/A; % Piola stress
results.(cond).(sample{s}).(test{t}).lmda = lmda;
results.(cond).(sample{s}).(test{t}).F_tilde = F_tilde;
results.(cond).(sample{s}).(test{t}).T = T;
figure
plot(lmda,F_tilde)
title([cond ' / ' sample{s} ' / ' test{t}])
text(1.002,max(F_tilde)-0.1,'Select two points (Millard/Thelen models): ')
text(1.002,max(F_tilde)-0.167,' 1) Curve intersection at one normalized force')
text(1.002,max(F_tilde)-0.23,' 2) Curve intersection at toe region end')
[X,Y] = ginput(2);
% Millard tendon model fitting
[eIso,kIso,fToe,curviness,F_tilde_fit] = bestFitMillard(lmda,F_tilde,X,Y);
results.(cond).(sample{s}).(test{t}).Millard.eIso = eIso;
results.(cond).(sample{s}).(test{t}).Millard.kIso = kIso;
results.(cond).(sample{s}).(test{t}).Millard.fToe = fToe;
results.(cond).(sample{s}).(test{t}).Millard.curviness = curviness;
% Thelen tendon model fitting TODO
% [eToe,fToe,kToe,kLin] = bestFitThelen(lmda,F_tilde,X,Y)
% results.(cond).(sample{s}).(test{t}).Thelen.eToe = eToe
% results.(cond).(sample{s}).(test{t}).Thelen.fToe = fToe
% results.(cond).(sample{s}).(test{t}).Thelen.kToe = kToe
% results.(cond).(sample{s}).(test{t}).Thelen.kLin = kLin
% HW2 tendon model fitting
[alpha,beta,Etan0,Etan20,T_fit] = bestFitHW2(lmda,T);
results.(cond).(sample{s}).(test{t}).HW2.alpha = alpha;
results.(cond).(sample{s}).(test{t}).HW2.beta = beta;
results.(cond).(sample{s}).(test{t}).HW2.Etan0 = Etan0;
results.(cond).(sample{s}).(test{t}).HW2.Etan20 = Etan20;
hold on
plot(lmda,F_tilde_fit)
figure
plot(lmda,T)
hold on
plot(lmda,T_fit)
save(fullfile(resultspath,'results.mat'),'-append','-struct','results') %,cond,sample{s},test{t})
save(fullfile(boxpath,'results','results.mat'),'-append','-struct','results')
end
end
end