%% Art3_v6_Fig5and6.m
% Roeland C. van de Vijsel
% NIOZ Yerseke, The Netherlands

clear all
close all
clc

% V1 - 20/04/2021: first, complete version
% V2 - 29/09/2021: final version for JGR publication: y-axes starting at 0
% and exclude labels in the graph

%%
lw = 2;
fs = 17;
varia = 4;
%%
if varia == 1;
%% RAI4
VAR = 'RAI4_ORIGINAL';
% summary = csvread(['Figure5/rawdata/output/',VAR,'_summary.csv'],1,2);
summary = csvread(['Figure5_R1/rawdata/output/',VAR,'_summary.csv'],1,3);

% Convert total content (micro-gram) into gram per m2
T0N_mea = summary(1,2);
T0Y_mea = summary(2,2);
T1N_mea = summary(3,2);
T1Y_mea = summary(4,2);
T2N_mea = summary(5,2);
T2Y_mea = summary(6,2);
T3N_mea = summary(7,2);
T3Y_mea = summary(8,2);

T0N_sem = summary(1,4);
T0Y_sem = summary(2,4);
T1N_sem = summary(3,4);
T1Y_sem = summary(4,4);
T2N_sem = summary(5,4);
T2Y_sem = summary(6,4);
T3N_sem = summary(7,4);
T3Y_sem = summary(8,4);

X1 = [0.5 1.5];
Y1 = [T3Y_mea, T3N_mea];
Z1 = [T3Y_sem, T3N_sem];

X2 = [0,2,4,10];
Y2 = [T0Y_mea,T1Y_mea,T2Y_mea,T3Y_mea];
Z2 = [T0Y_sem,T1Y_sem,T2Y_sem,T3Y_sem];
X3 = [0,2,4,10];
Y3 = [T0N_mea,T1N_mea,T2N_mea,T3N_mea];
Z3 = [T0N_sem,T1N_sem,T2N_sem,T3N_sem];

%% Plot RAI4 at end (Fig5)
close(figure(1));figure(1)
bar(X1(1),Y1(1),'FaceColor',[0 0.8 0],'EdgeColor',[0 0 0]); hold on
bar(X1(2),Y1(2),'FaceColor',[0.6 0.6 0.6],'EdgeColor',[0 0 0]);
errorbar(X1,Y1,Z1,'.k','LineWidth',lw)
set(gca,'xtick',X1,'xticklabel',{'With algae','Control'},'FontSize',fs)
% axis([0 2 1 4.8])
axis([0 2 0 4.8])
% box off
pbaspect([4. 1 1])
set(gca,'XGrid', 'off','YGrid', 'on')
% ylabel('Sediment deformation ratio [-]')
% title('(a) Sediment deformation ratio [-]', 'Units', 'normalized', 'Position', [0.3, 0.8, 0])
% saveas(gcf,'Art3_v6_Fig5_RAI4_ORGINAL.pdf')
saveas(gcf,'Art3_v11_Fig5_RAI4_ORGINAL.pdf')

%% Plot RAI4 through time (Fig6)
close(figure(1));figure(1)
ms = 6
errorbar(X2,Y2,Z2,'o-','Color',[0 0.8 0],'LineWidth',lw,'MarkerSize',ms); hold on
errorbar(X3,Y3,Z3,'o-','Color',[0.6 0.6 0.6],'LineWidth',lw,'MarkerSize',ms);
% axis([-.5 10.5 1 7.5])
axis([-.5 10.5 0 7.5])
set(gca,'xtick',[0:2:10])
set(gca,'FontSize',fs+3)
box off
set(gca,'XGrid', 'on','YGrid', 'on')
% legend('With algae','Control')
% legend boxoff
ylabel('Sediment deformation ratio [-]')
xlabel('Experiment time [weeks]')
% title('(a) Sediment deformability', 'Units', 'normalized', 'Position', [0.5, 0.9, 0])
% saveas(gcf,'Art3_v6_Fig6_RAI4_ORGINAL.pdf')
saveas(gcf,'Art3_v11_Fig6_RAI4_ORGINAL.pdf')

elseif varia == 2
%% CHLA
% VAR = 'CHLA_ORIGINAL';
VAR = 'CHLA'; % For CHLA, the data is not transformed, so this is the same as "ORIGINAL"
% summary = csvread(['Figure6/rawdata/output/',VAR,'_summary.csv'],1,2);
summary = csvread(['Figure6_R1/rawdata/output/',VAR,'_summary.csv'],1,3);

% Convert total content (micro-gram) into gram per m2
T0N_mea = summary(1,2);
T0Y_mea = summary(1,2);
T1N_mea = summary(2,2);
T1Y_mea = summary(3,2);
T2N_mea = summary(4,2);
T2Y_mea = summary(5,2);
T3N_mea = summary(6,2);
T3Y_mea = summary(7,2);

T0N_sem = summary(1,4);
T0Y_sem = summary(1,4);
T1N_sem = summary(2,4);
T1Y_sem = summary(3,4);
T2N_sem = summary(4,4);
T2Y_sem = summary(5,4);
T3N_sem = summary(6,4);
T3Y_sem = summary(7,4);

X1 = [0.5 1.5];
Y1 = [T3Y_mea, T3N_mea];
Z1 = [T3Y_sem, T3N_sem];

X2 = [0,2,4,10];
Y2 = [T0Y_mea,T1Y_mea,T2Y_mea,T3Y_mea];
Z2 = [T0Y_sem,T1Y_sem,T2Y_sem,T3Y_sem];
X3 = [0,2,4,10];
Y3 = [T0N_mea,T1N_mea,T2N_mea,T3N_mea];
Z3 = [T0N_sem,T1N_sem,T2N_sem,T3N_sem];

%% Plot CHLA
close(figure(1));figure(1)
bar(X1(1),Y1(1),'FaceColor',[0 0.8 0],'EdgeColor',[0 0 0]); hold on
bar(X1(2),Y1(2),'FaceColor',[0.6 0.6 0.6],'EdgeColor',[0 0 0]);
errorbar(X1,Y1,Z1,'.k','LineWidth',lw)
set(gca,'xtick',X1,'xticklabel',{'With algae','Control'},'FontSize',fs)
axis([0 2 0 700])
pbaspect([3 1 1])
set(gca,'XGrid', 'off','YGrid', 'on')
% box off
% ylabel('Chlorophyll-a concentration [\mug / g]')
% title('(b) Chlorophyll-a concentration [\mug / g]', 'Units', 'normalized', 'Position', [0.65, 0.8, 0])
% saveas(gcf,'Art3_v6_Fig5_CHLA_ORGINAL.pdf')
saveas(gcf,'Art3_v11_Fig5_CHLA_ORGINAL.pdf')

%% Plot CHLA through time (Fig6)
close(figure(1));figure(1)
ms = 6
errorbar(X2,Y2,Z2,'o-','Color',[0 0.8 0],'LineWidth',lw,'MarkerSize',ms); hold on
errorbar(X3,Y3,Z3,'o-','Color',[0.6 0.6 0.6],'LineWidth',lw,'MarkerSize',ms);
axis([-.5 10.5 0 700])
set(gca,'xtick',[0:2:10])
set(gca,'FontSize',fs+3)
box off
legend('With algae','Control','Location','West')
legend boxoff
set(gca,'XGrid', 'on','YGrid', 'on')
ylabel('Chlorophyll-a concentration [\mug / g]')
xlabel('Experiment time [weeks]')
% title('(b) Chlorophyll-a', 'Units', 'normalized', 'Position', [0.5, 0.9, 0])
% saveas(gcf,'Art3_v6_Fig6_CHLA_ORGINAL.pdf')
saveas(gcf,'Art3_v11_Fig6_CHLA_ORGINAL.pdf')

elseif varia == 3
%% ElevChange

% 3mL = 4.4mm
facML = 4.4/3;
VAR = 'ElevChange_ORIGINAL';
% summary = csvread(['Figure6/rawdata/output/',VAR,'_summary.csv'],1,2);
summary = csvread(['Figure6_R1/rawdata/output/',VAR,'_summary.csv'],1,3);

% Convert total content (micro-gram) into gram per m2
T0N_mea = 0;
T0Y_mea = 0;
T1N_mea = summary(1,2)*facML;
T1Y_mea = summary(2,2)*facML;
T2N_mea = summary(3,2)*facML;
T2Y_mea = summary(4,2)*facML;
T3N_mea = summary(5,2)*facML;
T3Y_mea = summary(6,2)*facML;

T0N_sem = 0;
T0Y_sem = 0;
T1N_sem = summary(1,4)*facML;
T1Y_sem = summary(2,4)*facML;
T2N_sem = summary(3,4)*facML;
T2Y_sem = summary(4,4)*facML;
T3N_sem = summary(5,4)*facML;
T3Y_sem = summary(6,4)*facML;

X1 = [0.5 1.5];
Y1 = [T3Y_mea, T3N_mea];
Z1 = [T3Y_sem, T3N_sem];

X2 = [0,2,4,10];
Y2 = [T0Y_mea,T1Y_mea,T2Y_mea,T3Y_mea];
Z2 = [T0Y_sem,T1Y_sem,T2Y_sem,T3Y_sem];
X3 = [0,2,4,10];
Y3 = [T0N_mea,T1N_mea,T2N_mea,T3N_mea];
Z3 = [T0N_sem,T1N_sem,T2N_sem,T3N_sem];

%% Plot ElevChange
close(figure(1));figure(1)
bar(X1(1),Y1(1),'FaceColor',[0 0.8 0],'EdgeColor',[0 0 0]); hold on
bar(X1(2),Y1(2),'FaceColor',[0.6 0.6 0.6],'EdgeColor',[0 0 0]);
errorbar(X1,Y1,Z1,'.k','LineWidth',lw)
set(gca,'xtick',X1,'xticklabel',{'With algae','Control'},'FontSize',fs)
axis([0 2 -5 0])
pbaspect([3 1 1])
set(gca,'XGrid', 'off','YGrid', 'on')
% box off
% ylabel('Elevation change relative to t_0 [mm]')
% title('(c) Substrate elevation change', 'Units', 'normalized', 'Position', [0.35, 0.1, 0])
% saveas(gcf,'Art3_v6_Fig5_ElevChange_ORGINAL.pdf')
saveas(gcf,'Art3_v11_Fig5_ElevChange_ORGINAL.pdf')

%% Plot ElevChange through time (Fig6)
close(figure(1));figure(1)
ms = 6
errorbar(X2,Y2,Z2,'o-','Color',[0 0.8 0],'LineWidth',lw,'MarkerSize',ms); hold on
errorbar(X3,Y3,Z3,'o-','Color',[0.6 0.6 0.6],'LineWidth',lw,'MarkerSize',ms);
axis([-.5 10.5 -5 0])
set(gca,'xtick',[0:2:10])
set(gca,'FontSize',fs+3)
box off
set(gca,'XGrid', 'on','YGrid', 'on')
% legend('With algae','Control')
ylabel('Elevation change relative to t_0 [mm]')
xlabel('Experiment time [weeks]')
% title('(c) Substrate elevation change', 'Units', 'normalized', 'Position', [0.5, 0.9, 0])
% saveas(gcf,'Art3_v6_Fig6_ElevChange_ORGINAL.pdf')
saveas(gcf,'Art3_v11_Fig6_ElevChange_ORGINAL.pdf')

elseif varia == 4
%% WC
VAR = 'WC_ORIGINAL';
% summary = csvread(['Figure6/rawdata/output/',VAR,'_summary.csv'],1,2);
summary = csvread(['Figure6_R1/rawdata/output/',VAR,'_summary.csv'],1,3);

% Convert total content (micro-gram) into gram per m2
T0N_mea = summary(1,2);
T0Y_mea = summary(1,2);
T1N_mea = summary(2,2);
T1Y_mea = summary(3,2);
T2N_mea = summary(4,2);
T2Y_mea = summary(5,2);
T3N_mea = summary(6,2);
T3Y_mea = summary(7,2);

T0N_sem = summary(1,4);
T0Y_sem = summary(1,4);
T1N_sem = summary(2,4);
T1Y_sem = summary(3,4);
T2N_sem = summary(4,4);
T2Y_sem = summary(5,4);
T3N_sem = summary(6,4);
T3Y_sem = summary(7,4);

X1 = [0.5 1.5];
Y1 = [T3Y_mea, T3N_mea];
Z1 = [T3Y_sem, T3N_sem];

X2 = [0,2,4,10];
Y2 = [T0Y_mea,T1Y_mea,T2Y_mea,T3Y_mea];
Z2 = [T0Y_sem,T1Y_sem,T2Y_sem,T3Y_sem];
X3 = [0,2,4,10];
Y3 = [T0N_mea,T1N_mea,T2N_mea,T3N_mea];
Z3 = [T0N_sem,T1N_sem,T2N_sem,T3N_sem];

%% Plot WC
close(figure(1));figure(1)
bar(X1(1),Y1(1),'FaceColor',[0 0.8 0],'EdgeColor',[0 0 0]); hold on
bar(X1(2),Y1(2),'FaceColor',[0.6 0.6 0.6],'EdgeColor',[0 0 0]);
errorbar(X1,Y1,Z1,'.k','LineWidth',lw)
set(gca,'xtick',X1,'xticklabel',{'With algae','Control'},'FontSize',fs)
% axis([0 2 50 80])
axis([0 2 0 80])
pbaspect([3 1 1])
set(gca,'XGrid', 'off','YGrid', 'on')
% ylabel('Sediment water content [%]')
% title('(d) Water content', 'Units', 'normalized', 'Position', [0.75, 0.8, 0])
% saveas(gcf,'Art3_v6_Fig5_WC_ORGINAL.pdf')
saveas(gcf,'Art3_v11_Fig5_WC_ORGINAL.pdf')

%% Plot CHLA through time (Fig6)
close(figure(1));figure(1)
ms = 6
errorbar(X2,Y2,Z2,'o-','Color',[0 0.8 0],'LineWidth',lw,'MarkerSize',ms); hold on
errorbar(X3,Y3,Z3,'o-','Color',[0.6 0.6 0.6],'LineWidth',lw,'MarkerSize',ms);
% axis([-.5 10.5 45 80])
axis([-.5 10.5 0 80])
set(gca,'xtick',[0:2:10])
set(gca,'FontSize',fs+3)
box off
set(gca,'XGrid', 'on','YGrid', 'on')
% legend('With algae','Control')
ylabel('Sediment water content [%]')
xlabel('Experiment time [weeks]')
% title('(d) Water content', 'Units', 'normalized', 'Position', [0.5, 0.9, 0])
% saveas(gcf,'Art3_v6_Fig6_WC_ORGINAL.pdf')
saveas(gcf,'Art3_v11_Fig6_WC_ORGINAL.pdf')

end

