Commit 2a8a6308 authored by Turnhout, M.C. van's avatar Turnhout, M.C. van
Browse files

Merge branch 'inge'

parents a30791f1 51eb9367
......@@ -61,6 +61,12 @@ slprj/
# Abaqus gitignore
*.rpy*
*.cae
*.sta
*.dat
*.msg
*.odb*
*.com
# Octave gitignore
octave-workspace
......
clear all;
close all;
clc;
[filename,path]=uigetfile({'*.tif';'*.*'},'MultiSelect','off');
path_and_file=[path filename];
info = imfinfo([path filename]);
filenames=regexp(filename,'\.','split');
frno = length(info);
xpix = info(1).Height;
ypix = info(1).Width;
% threshold=5;
% repeat=0;
stack = uint8(zeros(xpix,ypix,1,frno));
for frame=1:frno
color=imread(path_and_file,frame); % open color image
stack(:,:,:,frame)=rgb2gray(color); % convert from color to grayscale
end
for frame=1:frno
realstack(:,:,:,frame)=double(stack(:,:,:,frame));
end
%% scale
disp('click+enter at two of the ruler lines with 5 mm distance inbetween ')
for frame=1:frno
scrsz = get(0,'ScreenSize');
figure('Position',[1 scrsz(4)/2 scrsz(3)/2 scrsz(4)/2])
imagesc(stack(:,:,:,frame));
colormap gray
[ymin,xmin] = ginput;
[ymax,xmax] = ginput;
nrPix = abs(ymin-ymax);
cm = nrPix*2;
cmList(frame)=cm;
end
%% determine dimensions of flat films
flat = realstack(:,:,:,frno);
disp(['Choose basepoints for flat film'])
figure('Position',[0 (scrsz(4)/2-10) scrsz(3)/2 scrsz(4)/2])
imagesc(flat)
colormap gray
hold on
for film = 1:8
disp(['analysing film: ',num2str(film)])
fitindex=0;
film_ext = 2;
film_disc=100;
%% ================================================
% obtain total length of each MTF - use flat image
%==================================================
[a,b]=ginput;
basepoint_flat1=[ceil(a(1)),ceil(b(1))];
basepoint_flat2=[ceil(a(2)),ceil(b(2))];
plot(basepoint_flat1(1),basepoint_flat1(2),'yo',basepoint_flat2(1),basepoint_flat2(2),'yo')
bp_f1(film,1)=basepoint_flat1(1) - basepoint_flat1(1);
bp_f1(film,2)=basepoint_flat1(2) - basepoint_flat1(2);
bp_f2(film,1)=basepoint_flat2(1) - basepoint_flat1(1);
bp_f2(film,2)=basepoint_flat2(2) - basepoint_flat1(2);
disp(['Choose tippoints for flat film ' int2str(film)])
[a,b]=ginput;
tippoint_flat1=[ceil(a(1)),ceil(b(1))];
tippoint_flat2=[ceil(a(2)),ceil(b(2))];
plot(tippoint_flat1(1),tippoint_flat1(2),'ro',tippoint_flat2(1),tippoint_flat2(2),'ro')
tp_f1(film,1)=tippoint_flat1(1) - basepoint_flat1(1);
tp_f1(film,2)=tippoint_flat1(2) - basepoint_flat1(2);
tp_f2(film,1)=tippoint_flat2(1) - basepoint_flat1(1);
tp_f2(film,2)=tippoint_flat2(2) - basepoint_flat1(2);
%Calculate initial length of each MTF
L_a1 = (bp_f1(film,1) - tp_f1(film,1));
L_b1 = (bp_f1(film,2) - tp_f1(film,2));
init_L1(film) = sqrt(L_a1^2 + L_b1^2); % left side length
L_a2 = (bp_f2(film,1) - tp_f2(film,1));
L_b2 = (bp_f2(film,2) - tp_f2(film,2));
init_L2(film) = sqrt(L_a2^2 + L_b2^2); % right side length
init_L(film) = (init_L1(film) + init_L2(film))/2; % average length
%Calculate initial width of each MTF
W_a1 = (bp_f1(film,1) - bp_f2(film,1));
W_b1 = (bp_f1(film,2) - bp_f2(film,2));
init_W1(film) = sqrt(W_a1^2 + W_b1^2); % width at base
W_a2 = (tp_f1(film,1) - tp_f2(film,1));
W_b2 = (tp_f1(film,2) - tp_f2(film,2));
init_W2(film) = sqrt(W_a2^2 + W_b2^2); % width at top
init_W(film) = (init_W1(film) + init_W2(film))/2; % average width
%conversion of pixels to m
init_L_cm(film) = init_L(film)./cmList(frno);
L0(film) = init_L_cm(film)./100; %initial length in m
init_W_cm(film) = init_W(film)./cmList(frno);
W0(film) = init_W_cm(film)./100; %initial length in m
end
%% loop to measure multiple/all films of 1 sample
for film = 1:8
disp(['analysing film: ',num2str(film)])
fitindex=0;
film_ext = 2;
film_disc=100;
%% dimensions mTF
for frame=1:frno-1
disp(['analysing frame: ',num2str(frame)])
% Identify the corners of the film where the PIPAAm boundary was located
err=1;
while err==1
disp(['Choose basepoints for film ' int2str(film)])
h0=figure('Position',[0 (scrsz(4)/2-10) scrsz(3)/2 scrsz(4)/2]);
imagesc(realstack(:,:,:,frame))
colormap gray
hold on
[a,b]=ginput;
basepoint1=[ceil(a(1)),ceil(b(1))];
basepoint2=[ceil(a(2)),ceil(b(2))];
h1=plot(basepoint1(1),basepoint1(2),'yo',basepoint2(1),basepoint2(2),'yo');
bp1(film,1)=basepoint1(1);
bp1(film,2)=basepoint1(2);
bp2(film,1)=basepoint2(1);
bp2(film,2)=basepoint2(2);
% Identify the edges of the film corresponding with the corners above
disp(['Choose tippoints for film ' int2str(film)])
[a,b]=ginput;
tippoint1=[ceil(a(1)),ceil(b(1))];
tippoint2=[ceil(a(2)),ceil(b(2))];
h2=plot(tippoint1(1),tippoint1(2),'ro',tippoint2(1),tippoint2(2),'ro');
tp1(film,1)=tippoint1(1);
tp1(film,2)=tippoint1(2);
tp2(film,1)=tippoint2(1);
tp2(film,2)=tippoint2(2);
%% determine center of base and make center line for determining projection length
% Calculate the center of the base, the angle of film orienation, and the initial film length
centerbase(film,1)=(bp1(film,1)+bp2(film,1))./2;
centerbase(film,2)=(bp1(film,2)+bp2(film,2))./2;
basewidth(film)=sqrt((bp1(film,2)-bp2(film,2)).^2+(bp1(film,1)-bp2(film,1)).^2);
slope1=(tp1(film,2)-bp1(film,2))./(tp1(film,1)-bp1(film,1));
slope2=(tp2(film,2)-bp2(film,2))./(tp2(film,1)-bp2(film,1));
slope(film)=(slope1+slope2)./2;
theta(film)=atan(slope(film));
length1=sqrt((tp1(film,2)-bp1(film,2)).^2+(tp1(film,1)-bp1(film,1)).^2);
length2=sqrt((tp2(film,2)-bp2(film,2)).^2+(tp1(film,1)-bp1(film,1)).^2);
lengtht(film)=(length1+length2)./2; % average length between left and right
disp(['Click in direction of film relaxation and then hit return ' int2str(film)])
[filmcheck(1),filmcheck(2)]=ginput;
% for the vertical films, determine if they are curling up or down
if filmcheck(2)<=tp1(film,2) && theta(film)>0
dirline(film)=-1;
elseif filmcheck(2)>=tp1(film,2) && theta(film)<0
dirline(film)=-1;
else
dirline(film)=1;
end
% Calculate line of film orienation on which to calculate image intensity
filmlength=0;
for linecount=1:film_ext.*film_disc
if linecount==1
filmline(film,linecount,1)=centerbase(film,1);
filmline(film,linecount,2)=centerbase(film,2);
else
filmline(film,linecount,1)=filmline(film,linecount-1,1)+dirline(film).*cos(theta(film)).*lengtht(film)./film_disc;
filmline(film,linecount,2)=filmline(film,linecount-1,2)+dirline(film).*sin(theta(film)).*lengtht(film)./film_disc;
end
end
h3=plot(centerbase(film,1),centerbase(film,2),'go',filmline(film,:,1),filmline(film,:,2),'g');
% Establish 8 lines parallel to the one above on which to do parallel calculations
for linecount=1:9
if linecount==1
calclines(film,linecount,:,1)=filmline(film,:,1);
calclines(film,linecount,:,2)=filmline(film,:,2);
elseif linecount<=5
calclines(film,linecount,:,1)=filmline(film,:,1)+(linecount-1).*basewidth(film).*.0625.*sin(pi-theta(film));
calclines(film,linecount,:,2)=filmline(film,:,2)+(linecount-1).*basewidth(film).*.0625.*cos(pi-theta(film));
else
calclines(film,linecount,:,1)=filmline(film,:,1)-(linecount-5).*basewidth(film).*.0625.*sin(pi-theta(film));
calclines(film,linecount,:,2)=filmline(film,:,2)-(linecount-5).*basewidth(film).*.0625.*cos(pi-theta(film));
end
end
relax_check=input('Are you satisfied with the line path? (yes=0, no=1): ');
if relax_check==0
err=err+1;
elseif relax_check==1
delete(h1);delete(h2);delete(h3);
end
end
%% Treshold for bw image
hold off
close all
index = 0;
if film == 1
thresh(frame) = input('Enter Motion Threshold (0-130): ');
end
while index < 1;
thimag2 = realstack(:,:,:,frame)>thresh(frame);
thimag=imcomplement(thimag2);
imagesc(thimag)
colormap('gray')
w = input('Accept Threshold (yes = 0, no = 1): ');
if w == 0
index = 1;
else
thresh(frame) = input('Enter Motion Threshold (0-70): ');
end
close all
end
%% Threshold for opening and closing the image
index2=0;
if film == 1
disksize(frame) = input('Enter disk size for image close(2-15)');
end
while index2<1
cimag1 = imopen(thimag,strel('disk',2));
cimag2 = imclose(cimag1,strel('disk',disksize(frame)));
cimag(:,:,frame) = imfill(cimag2,'holes');
imagesc(cimag(:,:,frame))
colormap('gray')
w = input('Accept new image? (yes = 0, no = 1): ');
if w == 0
index2 = 1;
else
disksize(frame) = input('Enter disk size for image close(2-15)');
end
close all
end
close all
%% determine end of film
hold off
imagesc(realstack(:,:,:,frame))
colormap gray
hold on
for lineno=1:9
for lineloc=1:film_ext.*film_disc
lineval(lineno,lineloc)=cimag(ceil(calclines(film,lineno,lineloc,2)),ceil(calclines(film,lineno,lineloc,1)),frame);
end
end
lineavg=mean(lineval);
linecount=1;
linestop=0;
while linestop<1
if lineavg(linecount)+lineavg(linecount+1)+lineavg(linecount+2)<1
linestop=1;
line_end_x=filmline(film,linecount,1);
line_end_y=filmline(film,linecount,2);
plot(line_end_x,line_end_y,'go')
end
linecount=linecount+1;
end
film_projection(film,frame)=sqrt((line_end_x-centerbase(film,1)).^2+(line_end_y-centerbase(film,2)).^2);
%convert Length from pixels to m
film_projection_m(film,frame) = (film_projection(film,frame)./cmList(frame))./100;% MTF length in m
pause(.5)
end
%% ==================================================
%Calculate Radius of Curvature from film_projection
%==================================================
%conversion from length to Radius of curvature
for f = 1:frno-1
if film_projection_m(film,f) < (L0(film)/(2*pi()))
radius_c(film,f) = 0;
elseif (film_projection_m(film,f) >= (L0(film)/(2*pi()))) && (film_projection_m(film,f) <= (2*L0(film)/pi()))
radius_c(film,f) = film_projection_m(film,f);
elseif (film_projection_m(film,f) > (2*L0(film)/pi())) && (film_projection_m(film,f) <= L0(film))
options2=optimset('TolX',1e-100, 'maxFunEvals', 1e30, 'TolFun', 1e-50);
bound1 = 0.0;
bound2 = pi/(2*L0(film));
Roar = fminbnd(@(kappa) ((1/kappa)*sin(L0(film)*kappa)-film_projection_m(film,f)).^2, bound1, bound2,'options2');
kappa_track(film,f)=Roar;
radius_c(film,f) = 1/Roar;
elseif (film_projection_m(film,f) >= L0(film))
radius_c(film,f) = 10000*L0(film)*180/pi; % relative infinity
end
end
save([path,filenames{1},'_',num2str(film)])
end
%%
disp('All images are analyzed')
% write data of all timepoints (0 h = first frame)
for frame=1:frame
xlswrite([path,filenames{1},'_radius.xlsx'],{['frame ',num2str(frame)],'radius (mm)','curvature (1/mm)','pr length (mm)'},1,['A',num2str(10*frame-9)]);
xlswrite([path,filenames{1},'_radius.xlsx'],[radius_c(:,frame)*1000,1./(radius_c(:,frame)*1000),film_projection_m(:,frame)*1000],1,['B',num2str(1+10*frame-9)]);
end
%% write initial dimensions and corner coordinates of the films
xlswrite([path,filenames{1},'_dimensions.xlsx'],{'length (mm)','width (mm)'},1,'A1');
xlswrite([path,filenames{1},'_dimensions.xlsx'],[L0'*1000,W0'*1000],1,'A2'); % L0 and W0 in m
xlswrite([path,filenames{1},'_dimensions.xlsx'],{'coordinates (m)','bp_x1','bp_y1','bp_x2','bp_y2','tp_x1','tp_y1','tp_x2','tp_y2'},1,'D1');
xlswrite([path,filenames{1},'_dimensions.xlsx'],[bp_f1/(cmList(frno)*100),bp_f2/(cmList(frno)*100),tp_f1/(cmList(frno)*100),tp_f2/(cmList(frno)*100)],1,'E2'); % corner coordinates in m
clear all; close all
%% import nuclear orientation on films
exp1_0=xlsread('compleet overzicht experiment 1 - 4.xlsx','overview','I4:O11');
exp2_0=xlsread('compleet overzicht experiment 1 - 4.xlsx','overview','Y4:AE11');
exp3_0=xlsread('compleet overzicht experiment 1 - 4.xlsx','overview','AO4:AU11');
exp4_0=xlsread('compleet overzicht experiment 1 - 4.xlsx','overview','BE4:BK11');
exp1_15=xlsread('compleet overzicht experiment 1 - 4.xlsx','overview','I15:O21');
exp2_15=xlsread('compleet overzicht experiment 1 - 4.xlsx','overview','Y14:AE21');
exp3_15=xlsread('compleet overzicht experiment 1 - 4.xlsx','overview','AO14:AU21');
exp4_15=xlsread('compleet overzicht experiment 1 - 4.xlsx','overview','BE14:BK21');
exp1_30=xlsread('compleet overzicht experiment 1 - 4.xlsx','overview','I24:O31');
exp2_30=xlsread('compleet overzicht experiment 1 - 4.xlsx','overview','Y24:AE31');
exp3_30=xlsread('compleet overzicht experiment 1 - 4.xlsx','overview','AO24:AU28');
exp4_30=xlsread('compleet overzicht experiment 1 - 4.xlsx','overview','BE24:BK31');
exp1_45=xlsread('compleet overzicht experiment 1 - 4.xlsx','overview','I34:O38');
exp2_45=xlsread('compleet overzicht experiment 1 - 4.xlsx','overview','Y34:AE39');
exp3_45=xlsread('compleet overzicht experiment 1 - 4.xlsx','overview','AO34:AU40');
exp4_45=xlsread('compleet overzicht experiment 1 - 4.xlsx','overview','BE34:BK41');
%% get data
density_0 =[ exp1_0(:,1); exp2_0(:,1); exp3_0(:,1); exp4_0(:,1)];
density_15=[exp1_15(:,1);exp2_15(:,1);exp3_15(:,1);exp4_15(:,1)];
density_30=[exp1_30(:,1);exp2_30(:,1);exp3_30(:,1);exp4_30(:,1)];
density_45=[exp1_45(:,1);exp2_45(:,1);exp3_45(:,1);exp4_45(:,1)];
curv0_0 =[ exp1_0(:,2); exp2_0(:,2); exp3_0(:,2); exp4_0(:,2)];
curv0_15=[exp1_15(:,2);exp2_15(:,2);exp3_15(:,2);exp4_15(:,2)];
curv0_30=[exp1_30(:,2);exp2_30(:,2);exp3_30(:,2);exp4_30(:,2)];
curv0_45=[exp1_45(:,2);exp2_45(:,2);exp3_45(:,2);exp4_45(:,2)];
curv1_0 =[ exp1_0(:,3); exp2_0(:,3); exp3_0(:,3); exp4_0(:,3)];
curv1_15=[exp1_15(:,3);exp2_15(:,3);exp3_15(:,3);exp4_15(:,3)];
curv1_30=[exp1_30(:,3);exp2_30(:,3);exp3_30(:,3);exp4_30(:,3)];
curv1_45=[exp1_45(:,3);exp2_45(:,3);exp3_45(:,3);exp4_45(:,3)];
%% plot curvatures
h1=figure(1);
subplot(2,2,1); plot(density_0,curv0_0,'r^','MarkerFaceColor','r','MarkerSize',7);
hold on; plot(density_0,curv1_0,'dk','MarkerFaceColor','k','MarkerSize',7); hold off
axis([90 600 -.05 2.1]); set(gca,'FontSize',11,'Box','on','LineWidth',1.5)
ylabel('Curvature [mm^-^1]','FontSize',11); %xlabel('Density [cells/mm^2]','FontSize',11);
subplot(2,2,2); plot(density_15,curv0_15,'r^','MarkerFaceColor','r','MarkerSize',7);
hold on; plot(density_15,curv1_15,'dk','MarkerFaceColor','k','MarkerSize',7); hold off
axis([90 600 -.05 2.1]); set(gca,'FontSize',11,'Box','on','LineWidth',1.5)
% xlabel('Density [cells/mm^2]','FontSize',11);ylabel('Curvature [mm^-^1]','FontSize',11);
subplot(2,2,3); plot(density_30,curv0_30,'r^','MarkerFaceColor','r','MarkerSize',7);
hold on; plot(density_30,curv1_30,'dk','MarkerFaceColor','k','MarkerSize',7); hold off
axis([90 600 -.05 2.1]); set(gca,'FontSize',11,'Box','on','LineWidth',1.5)
xlabel('Density [cells/mm^2]','FontSize',11);ylabel('Curvature [mm^-^1]','FontSize',11);
subplot(2,2,4); plot(density_45,curv0_45,'r^','MarkerFaceColor','r','MarkerSize',7);
hold on; plot(density_45,curv1_45,'dk','MarkerFaceColor','k','MarkerSize',7); hold off
axis([90 600 -.05 2.1]); set(gca,'FontSize',11,'Box','on','LineWidth',1.5)
xlabel('Density [cells/mm^2]','FontSize',11);%ylabel('Curvature [mm^-^1]','FontSize',11);
%% get sigf data
sigf0_0 =[ exp1_0(:,4); exp2_0(:,4); exp3_0(:,4); exp4_0(:,4)];
sigf0_15=[exp1_15(:,4);exp2_15(:,4);exp3_15(:,4);exp4_15(:,4)];
sigf0_30=[exp1_30(:,4);exp2_30(:,4);exp3_30(:,4);exp4_30(:,4)];
sigf0_45=[exp1_45(:,4);exp2_45(:,4);exp3_45(:,4);exp4_45(:,4)];
sigf0_0(sigf0_0==0)=NaN; sigf0_15(sigf0_15==0)=NaN; sigf0_30(sigf0_30==0)=NaN; sigf0_45(sigf0_45==0)=NaN;
sigf1_0 =[ exp1_0(:,5); exp2_0(:,5); exp3_0(:,5); exp4_0(:,5)];
sigf1_15=[exp1_15(:,5);exp2_15(:,5);exp3_15(:,5);exp4_15(:,5)];
sigf1_30=[exp1_30(:,5);exp2_30(:,5);exp3_30(:,5);exp4_30(:,5)];
sigf1_45=[exp1_45(:,5);exp2_45(:,5);exp3_45(:,5);exp4_45(:,5)];
sigf1_0(sigf1_0==0)=NaN; sigf1_15(sigf1_15==0)=NaN; sigf1_30(sigf1_30==0)=NaN; sigf1_45(sigf1_45==0)=NaN;
%% plot sigf
h5=figure(5); plot(density_0,sigf0_0,'r^','MarkerFaceColor','r','MarkerSize',8);
hold on; plot(density_0,sigf1_0,'dk','MarkerFaceColor','k','MarkerSize',8); hold off
axis([100 600 -.2 9]); set(gca,'FontSize',14,'Box','on','LineWidth',2)
xlabel('Density [cells/mm^2]','FontSize',14);ylabel('\sigma_f [kPa]','FontSize',14);
h6=figure(6); plot(density_15,sigf0_15,'r^','MarkerFaceColor','r','MarkerSize',8);
hold on; plot(density_15,sigf1_15,'dk','MarkerFaceColor','k','MarkerSize',8); hold off
axis([100 600 -.2 9]); set(gca,'FontSize',14,'Box','on','LineWidth',2)
xlabel('Density [cells/mm^2]','FontSize',14);ylabel('\sigma_f [kPa]','FontSize',14);
h7=figure(7); plot(density_30,sigf0_30,'r^','MarkerFaceColor','r','MarkerSize',8);
hold on; plot(density_30,sigf1_30,'dk','MarkerFaceColor','k','MarkerSize',8); hold off
axis([100 600 -.2 9]); set(gca,'FontSize',14,'Box','on','LineWidth',2)
xlabel('Density [cells/mm^2]','FontSize',14);ylabel('\sigma_f [kPa]','FontSize',14);
h8=figure(8); plot(density_45,sigf0_45,'r^','MarkerFaceColor','r','MarkerSize',8);
hold on; plot(density_45,sigf1_45,'dk','MarkerFaceColor','k','MarkerSize',8); hold off
axis([100 600 -.2 9]); set(gca,'FontSize',14,'Box','on','LineWidth',2)
xlabel('Density [cells/mm^2]','FontSize',14);ylabel('\sigma_f [kPa]','FontSize',14);
%% get sigmax data
sigmax0_0 =[ exp1_0(:,6); exp2_0(:,6); exp3_0(:,6); exp4_0(:,6)];
sigmax0_15=[exp1_15(:,6);exp2_15(:,6);exp3_15(:,6);exp4_15(:,6)];
sigmax0_30=[exp1_30(:,6);exp2_30(:,6);exp3_30(:,6);exp4_30(:,6)];
sigmax0_45=[exp1_45(:,6);exp2_45(:,6);exp3_45(:,6);exp4_45(:,6)];
sigmax0_0(sigmax0_0==0)=NaN; sigmax0_15(sigmax0_15==0)=NaN; sigmax0_30(sigmax0_30==0)=NaN; sigmax0_45(sigmax0_45==0)=NaN;
sigmax1_0 =[ exp1_0(:,7); exp2_0(:,7); exp3_0(:,7); exp4_0(:,7)];
sigmax1_15=[exp1_15(:,7);exp2_15(:,7);exp3_15(:,7);exp4_15(:,7)];
sigmax1_30=[exp1_30(:,7);exp2_30(:,7);exp3_30(:,7);exp4_30(:,7)];
sigmax1_45=[exp1_45(:,7);exp2_45(:,7);exp3_45(:,7);exp4_45(:,7)];
sigmax1_0(sigmax1_0==0)=NaN; sigmax1_15(sigmax1_15==0)=NaN; sigmax1_30(sigmax1_30==0)=NaN; sigmax1_45(sigmax1_45==0)=NaN;
%% plot sigmax
h9=figure(9); plot(density_0,sigmax0_0,'r^','MarkerFaceColor','r','MarkerSize',8);
hold on; plot(density_0,sigmax1_0,'dk','MarkerFaceColor','k','MarkerSize',8); hold off
axis([100 600 -.3 11]); set(gca,'FontSize',14,'Box','on','LineWidth',2)
xlabel('Density [cells/mm^2]','FontSize',14);ylabel('\sigma_m_a_x [kPa]','FontSize',14);
h10=figure(10); plot(density_15,sigmax0_15,'r^','MarkerFaceColor','r','MarkerSize',8);
hold on; plot(density_15,sigmax1_15,'dk','MarkerFaceColor','k','MarkerSize',8); hold off
axis([100 600 -.3 11]); set(gca,'FontSize',14,'Box','on','LineWidth',2)
xlabel('Density [cells/mm^2]','FontSize',14);ylabel('\sigma_m_a_x [kPa]','FontSize',14);
h11=figure(11); plot(density_30,sigmax0_30,'r^','MarkerFaceColor','r','MarkerSize',8);
hold on; plot(density_30,sigmax1_30,'dk','MarkerFaceColor','k','MarkerSize',8); hold off
axis([100 600 -.3 11]); set(gca,'FontSize',14,'Box','on','LineWidth',2)
xlabel('Density [cells/mm^2]','FontSize',14);ylabel('\sigma_m_a_x [kPa]','FontSize',14);
h12=figure(12); plot(density_45,sigmax0_45,'r^','MarkerFaceColor','r','MarkerSize',8);
hold on; plot(density_45,sigmax1_45,'dk','MarkerFaceColor','k','MarkerSize',8); hold off
axis([100 600 -.3 11]); set(gca,'FontSize',14,'Box','on','LineWidth',2)
xlabel('Density [cells/mm^2]','FontSize',14);ylabel('\sigma_m_a_x [kPa]','FontSize',14);
%%
h13=figure;
subplot(2,4,1); plot(density_0,sigf0_0,'r^','MarkerFaceColor','r','MarkerSize',6);
hold on; plot(density_0,sigf1_0,'dk','MarkerFaceColor','k','MarkerSize',6); hold off
axis([90 450 -.5 9]); set(gca,'FontSize',10,'Box','on','LineWidth',1.5)
ylabel('\sigma_f [kPa]','FontSize',10);
subplot(2,4,2); plot(density_15,sigf0_15,'r^','MarkerFaceColor','r','MarkerSize',6);
hold on; plot(density_15,sigf1_15,'dk','MarkerFaceColor','k','MarkerSize',6); hold off
axis([90 450 -.5 9]); set(gca,'FontSize',10,'Box','on','LineWidth',1.5)
% xlabel('Density [cells/mm^2]','FontSize',10);%ylabel('\sigma_f [kPa]','FontSize',10);
subplot(2,4,3); plot(density_30,sigf0_30,'r^','MarkerFaceColor','r','MarkerSize',6);
hold on; plot(density_30,sigf1_30,'dk','MarkerFaceColor','k','MarkerSize',6); hold off
axis([90 450 -.5 9]); set(gca,'FontSize',10,'Box','on','LineWidth',1.5)
% xlabel('Density [cells/mm^2]','FontSize',10);%ylabel('\sigma_f [kPa]','FontSize',10);
subplot(2,4,4); plot(density_45,sigf0_45,'r^','MarkerFaceColor','r','MarkerSize',6);
hold on; plot(density_45,sigf1_45,'dk','MarkerFaceColor','k','MarkerSize',6); hold off
axis([90 450 -.5 9]); set(gca,'FontSize',10,'Box','on','LineWidth',1.5)
% xlabel('Density [cells/mm^2]','FontSize',10);%ylabel('\sigma_f [kPa]','FontSize',10);
subplot(2,4,5); plot(density_0,sigmax0_0,'r^','MarkerFaceColor','r','MarkerSize',6);
hold on; plot(density_0,sigmax1_0,'dk','MarkerFaceColor','k','MarkerSize',6); hold off
axis([90 450 -.5 11]); set(gca,'FontSize',10,'Box','on','LineWidth',1.5)
xlabel('Density [cells/mm^2]','FontSize',10);ylabel('\sigma_m_a_x [kPa]','FontSize',10);
subplot(2,4,6); plot(density_15,sigmax0_15,'r^','MarkerFaceColor','r','MarkerSize',6);
hold on; plot(density_15,sigmax1_15,'dk','MarkerFaceColor','k','MarkerSize',6); hold off
axis([90 450 -.5 11]); set(gca,'FontSize',10,'Box','on','LineWidth',1.5)
xlabel('Density [cells/mm^2]','FontSize',10);%ylabel('\sigma_m_a_x [kPa]','FontSize',10);
subplot(2,4,7); plot(density_30,sigmax0_30,'r^','MarkerFaceColor','r','MarkerSize',6);
hold on; plot(density_30,sigmax1_30,'dk','MarkerFaceColor','k','MarkerSize',6); hold off
axis([90 450 -.5 11]); set(gca,'FontSize',10,'Box','on','LineWidth',1.5)
xlabel('Density [cells/mm^2]','FontSize',10);%ylabel('\sigma_m_a_x [kPa]','FontSize',10);
subplot(2,4,8); plot(density_45,sigmax0_45,'r^','MarkerFaceColor','r','MarkerSize',6);
hold on; plot(density_45,sigmax1_45,'dk','MarkerFaceColor','k','MarkerSize',6); hold off
axis([90 450 -.5 11]); set(gca,'FontSize',10,'Box','on','LineWidth',1.5)
xlabel('Density [cells/mm^2]','FontSize',10);%ylabel('\sigma_m_a_x [kPa]','FontSize',10);
%% sigmax/density
sigmaxCell0_0 = sigmax0_0./ density_0; sigmaxCell1_0 = sigmax1_0./ density_0;
sigmaxCell0_15 = sigmax0_15./density_15; sigmaxCell1_15 = sigmax1_15./density_15;
sigmaxCell0_30 = sigmax0_30./density_30; sigmaxCell1_30 = sigmax1_30./density_30;
sigmaxCell0_45 = sigmax0_45./density_45; sigmaxCell1_45 = sigmax1_45./density_45;
sigmaxCellTot=[sigmaxCell0_0;sigmaxCell1_0;sigmaxCell0_15;sigmaxCell1_15;sigmaxCell0_30;sigmaxCell1_30;sigmaxCell0_45;sigmaxCell1_45];
%% boxplot labels
label0time=[repmat({'0 h'},size(sigmaxCell0_0,1),1);repmat({'1 h'},size(sigmaxCell1_0,1),1)];
label15time=[repmat({'0 h'},size(sigmaxCell0_15,1),1);repmat({'1 h'},size(sigmaxCell1_15,1),1)];
label30time=[repmat({'0 h'},size(sigmaxCell0_30,1),1);repmat({'1 h'},size(sigmaxCell1_30,1),1)];
label45time=[repmat({'0 h'},size(sigmaxCell0_45,1),1);repmat({'1 h'},size(sigmaxCell1_45,1),1)];
deg=sprintf(char(176));
label0angle =repmat({['0',deg]},2*size(sigmaxCell0_0,1),1);
label15angle=repmat({['15',deg]},2*size(sigmaxCell0_15,1),1);
label30angle=repmat({['30',deg]},2*size(sigmaxCell0_30,1),1);
label45angle=repmat({['45',deg]},2*size(sigmaxCell0_45,1),1);
g1=[label0time;label15time;label30time;label45time];
g2=[label0angle;label15angle;label30angle;label45angle];
%% boxplot
h14=figure(14);
bp=boxplot(sigmaxCellTot*1000,{g2,g1},'factorgap',[10,0],'symbol','o','colors','rk');
set(bp,'LineWidth',2);
ylabel('\sigma_c_e_l_l [Pa]','FontSize',12)
%% sigf/density
sigfCell0_0 = sigf0_0./ density_0; sigfCell1_0 = sigf1_0./ density_0;
sigfCell0_15 = sigf0_15./density_15; sigfCell1_15 = sigf1_15./density_15;
sigfCell0_30 = sigf0_30./density_30; sigfCell1_30 = sigf1_30./density_30;
sigfCell0_45 = sigf0_45./density_45; sigfCell1_45 = sigf1_45./density_45;
sigfCellTot=[sigfCell0_0;sigfCell1_0;sigfCell0_15;sigfCell1_15;sigfCell0_30;sigfCell1_30;sigfCell0_45;sigfCell1_45];
%% boxplot
h15=figure(15);
bp=boxplot(sigfCellTot*1000,{g2,g1},'factorgap',[10,0],'symbol','o','colors','rk');
set(bp,'LineWidth',2);
ylabel('\sigma_f [Pa]','FontSize',12)
clear all; close all
%% import nuclear orientation on films
angles=xlsread('compleet overzicht experiment 1 - 4.xlsx','orientation0','C1:CN1');
deg0=xlsread('compleet overzicht experiment 1 - 4.xlsx','orientation0','C2:CN33');
deg15=xlsread('compleet overzicht experiment 1 - 4.xlsx','orientation15','C2:CN32');
deg30=xlsread('compleet overzicht experiment 1 - 4.xlsx','orientation30','C2:CN30');
deg45=xlsread('compleet overzicht experiment 1 - 4.xlsx','orientation45','C2:CN27');
meanD0=mean(deg0,1); semD0=std(deg0,0,1)/sqrt(size(deg0,1)); stdD0=std(deg0,0,1);