%% General plotting options
fs=10;
scrsz = get(0,'ScreenSize');
figure('Units','centimeters','Position',[8 2 27 27])
set(gcf,'papersize',[27,25])
fig = gcf;
style = hgexport('factorystyle');
style.Bounds = 'tight';
%hgexport(fig,'-clipboard',style,'applystyle', true);
drawnow;
%% Run initialisation script
Initialize;
%% Run simulation for power profile 1
P_D = y_exp_P1;
fault_B1_enable=[0 0];
fault_B2_enable=[0 0];
noise_B1_enable=[1 1];
noise_B2_enable=[1 1];
disp("Running simulation for power profile 1\n");
P1=sim("Modular_Controller.slx","TimeOut",t_tot);
%% Power profile 1 plots
% Power profile vs achieved comparison [Profile 1]
figure
plot(P1.Demand_PD.Time/60,P1.Demand_PD.Data,'b-');
hold on
plot(P1.Power_Split_1.Time/60,eta_T*(P1.Power_Split_1.Data(:,1)+P1.Power_Split_1.Data(:,2)),'m--');
hold on
plot(P1.Demand_Paux.Time/60,P1.Demand_Paux.Data,'k-');
hold on
plot(P1.Power_Split_1.Time/60,P1.Power_Split_1.Data(:,3)+P1.Power_Split_1.Data(:,4)+P1.Power_Split_1.Data(:,5)-P1.Power_Split_1.Data(:,2)/(eta_IM*eta_PC),'c--');
grid('minor');
legend("P_D^{ref}","P_D","P_{aux}^{ref}","P_{aux}");
xlabel('Time [min]','FontName','Times New Roman','FontSize',24,'interpreter','latex')
ylabel('Power [kW]','FontName','Times New Roman','FontSize',24)
% Propulsion power split - loading of ice and induction motor [Profile 1]
figure
plot(P1.Power_Split_1.Time/60,P1.Power_Split_1.Data(:,1),'k-');
hold on
plot(P1.Power_Split_1.Time/60,P1.Power_Split_1.Data(:,2),'m-');
legend("ICE","Induction Motor")
grid('minor');
ylim([0,900]);
xlabel('Time [min]','FontName','Times New Roman','FontSize',24,'interpreter','latex')
ylabel('Power [kW]','FontName','Times New Roman','FontSize',24)
% Generator power split - loading of gensets [Profile 1]
figure
plot(P1.Power_Split_1.Time/60,P1.Power_Split_1.Data(:,3),'k-');
hold on
plot(P1.Power_Split_1.Time/60,P1.Power_Split_1.Data(:,4),'m-');
hold on
plot(P1.Power_Split_1.Time/60,P1.Power_Split_1.Data(:,5),'b-');
legend("Genset 1","Genset 2","Battery 1")
grid('minor');
ylim([-200,600]);
xlabel('Time [min]','FontName','Times New Roman','FontSize',24,'interpreter','latex');
ylabel('Power [kW]','FontName','Times New Roman','FontSize',24);
% Battery power split - loading of batteries and SOC estimation on the right [Profile 1]
figure
plot(P1.Measurements_B1.Time/60,P1.Measurements_B1.Data(:,3),'r-');
hold on;
plot(P1.Measurements_B1.Time/60,100*ones(size(P1.Measurements_B1.Data(:,3))),'k--');
hold on;
plot(P1.Measurements_B1.Time/60,10*ones(size(P1.Measurements_B1.Data(:,3))),'k--');
ylabel('SOC [%]','FontName','Times New Roman','FontSize',24)
legend("$\hat{x}^{(5)}$","SOC_{min}","SOC_{max}",'intepreter','latex');
xlabel("Time [min]",'FontName','Times New Roman','FontSize',24,'interpreter','latex');
grid('minor');
ylim([0 110])
%% Run simulation for power profile 2
P_D = y_exp_P2;
fault_B1_enable=[1 0];
fault_B2_enable=[0 1];
noise_B1_enable=[1 1];
noise_B2_enable=[1 1];
time1=[20*60 50];
time2=[20 60*60];
t_fault2= 100*60;
disp("Running simulation for power profile 2");
%%
P2=sim("Modular_Controller.slx","TimeOut",t_tot);
%% Power profile 2 plots
%figure
plot(P2.Demand_PD.Time/60,P2.Demand_PD.Data,'b-');
hold on
ach_prop_P2=P2.Power_Split_2.Data(P2.Power_Split_2.Time<=t_fault2,1)+P2.Power_Split_2.Data(P2.Power_Split_2.Time<=t_fault2,2);
ach_prop_P1=P2.Power_Split_1.Data(P2.Power_Split_2.Time>t_fault2,1)+P2.Power_Split_1.Data(P2.Power_Split_2.Time>t_fault2,2);
plot(P2.Power_Split_2.Time/60,eta_T*[ach_prop_P2' ach_prop_P1'],'m--');
hold on
plot(P2.Demand_Paux.Time/60,P2.Demand_Paux.Data,'k-');
hold on
plot(P2.Power_Split_2.Time/60,P2.Power_Split_2.Data(:,3)+P2.Power_Split_2.Data(:,4)+P2.Power_Split_2.Data(:,5)+P2.Power_Split_2.Data(:,6)-P2.Power_Split_2.Data(:,2)/(eta_IM*eta_PC),'r-.');
grid('minor');
leg4=legend("$P_D^{ref}$","$P_D$","$P_{aux}^{ref}$","$P_{aux}$");
xlabel('Time [min]','FontName','Times New Roman','FontSize',24,'interpreter','latex')
ylabel('Power [kW]','FontName','Times New Roman','FontSize',24,'interpreter','latex')
set(leg4,'Interpreter','latex');
set(leg4,'FontSize',17);
set(gca,'FontSize',20)
%%
% Propulsion power split - loading of ice and induction motor [Profile 1]
%figure
plot(P2.Power_Split_2.Time/60,[P2.Power_Split_2.Data(P2.Power_Split_2.Time<=t_fault2,1)' P2.Power_Split_1.Data(P2.Power_Split_2.Time>t_fault2,1)'],'k-.');
hold on
plot(P2.Power_Split_2.Time/60,[P2.Power_Split_2.Data(P2.Power_Split_2.Time<=t_fault2,2)' P2.Power_Split_1.Data(P2.Power_Split_2.Time>t_fault2,2)'],'m-');
leg3=legend("ICE","Induction Motor")
grid('minor');
ylim([0,900]);
yticks([0:200:900])
xlabel('Time [min]','FontName','Times New Roman','FontSize',24,'interpreter','latex')
ylabel('Power [kW]','FontName','Times New Roman','FontSize',24,'interpreter','latex')
set(leg3,'Interpreter','latex');
set(leg3,'FontSize',17);
set(gca,'FontSize',20)
%%
% Generator power split - loading of gensets [Profile 1]
plot(P2.Power_Split_2.Time/60,[P2.Power_Split_2.Data(P2.Power_Split_2.Time<=t_fault2,3)' P2.Power_Split_1.Data(P2.Power_Split_2.Time>t_fault2,3)'],'k-');
hold on
plot(P2.Power_Split_2.Time/60,[P2.Power_Split_2.Data(P2.Power_Split_2.Time<=t_fault2,4)' P2.Power_Split_1.Data(P2.Power_Split_2.Time>t_fault2,4)'],'m-.');
hold on
plot(P2.Power_Split_2.Time/60,[P2.Power_Split_2.Data(P2.Power_Split_2.Time<=t_fault2,5)' P2.Power_Split_1.Data(P2.Power_Split_2.Time>t_fault2,5)'],'b--');
hold on
plot(P2.Power_Split_2.Time/60,[P2.Power_Split_2.Data(P2.Power_Split_2.Time<=t_fault2,6)' zeros(1,length(P2.Power_Split_2.Time)-length(P2.Power_Split_2.Data(P2.Power_Split_2.Time<=t_fault2,6)))],'r:','LineWidth',2);
leg2=legend("Genset 1","Genset 2","Battery 1","Battery 2")
grid('minor');
ylim([-200,800]);
yticks([-200:200:800])
xlabel('Time [min]','FontName','Times New Roman','FontSize',24,'interpreter','latex');
ylabel('Power [kW]','FontName','Times New Roman','FontSize',24,'interpreter','latex');
set(leg2,'Interpreter','latex');
set(leg2,'FontSize',17);
set(gca,'FontSize',20)
% Battery power split - loading of batteries and SOC estimation on the right [Profile 1]
figure
plot(P2.Measurements_B1.Time/60,P2.Measurements_B1.Data(:,3),'m-');
hold on;
plot(P2.Measurements_B2.Time/60,P2.Measurements_B2.Data(:,3),'r-');
hold on;
plot(P2.Measurements_B1.Time/60,100*ones(size(P2.Measurements_B1.Data(:,3))),'k--');
hold on;
plot(P2.Measurements_B1.Time/60,10*ones(size(P2.Measurements_B1.Data(:,3))),'k--');
ylabel('SOC [%]')
leg1=legend("$\hat{x}^{(5)}$","$\hat{x}^{(6)}$","SOC_{min}","SOC_{max}");
set(leg1,'Interpreter','latex');
set(leg1,'FontSize',17);
grid('minor')
xlabel("Time [min]",'FontName','Times New Roman','FontSize',24,'interpreter','latex');
ylim([0 110])
%%
% Switching signals for battery sensor faults [Profile 2]
subplot(2,3,1);
yyaxis left
plot(P2.kI_bat1.Time/60,P2.kI_bat1.Data(:,1),'k-');
ylim([0 3])
xlim([0 140])
xticks([0:20:140])
yticks([0:1:3])
grid('minor');
ylabel("$\sigma$(1)",'FontName','Times New Roman','FontSize',24,'interpreter','latex');
hold on
yyaxis right
plot(P2.kI_bat1.Time/60,[zeros(size(P2.kI_bat1.Time(P2.kI_bat1.Data<=1)))' ones(size(P2.kI_bat1.Time(P2.kI_bat1.Data>1)))'],'r-.')
ylim([0 1.05])
yticks([0:1:1.05])
xlabel("Time [min]",'FontName','Times New Roman','FontSize',24,'interpreter','latex');
ylabel("$D^{(5,1)}$",'FontName','Times New Roman','FontSize',24,'interpreter','latex');
title('(a)','FontName','Times New Roman','FontSize',24)
set(gca,'FontSize',20)
subplot(2,3,2);
title('(b)','FontName','Times New Roman','FontSize',24)
yyaxis left
plot(P2.kV_bat1.Time/60,P2.kV_bat1.Data(:,1),'k-');
ylim([0 3])
yticks([0:1:3])
xlim([0 140])
xticks([0:20:140])
grid('minor');
ylabel("$\sigma$(2)",'FontName','Times New Roman','FontSize',24,'interpreter','latex');
hold on
yyaxis right
plot(P2.kV_bat1.Time/60,[zeros(size(P2.kV_bat1.Time(P2.kV_bat1.Data<=1)))' ones(size(P2.kV_bat1.Time(P2.kV_bat1.Data>1)))'],'r-.')
ylim([0 1.05])
yticks([0:1:1.05])
xlabel("Time [min]",'FontName','Times New Roman','FontSize',24,'interpreter','latex');
ylabel("$D^{(5,2)}$",'FontName','Times New Roman','FontSize',24,'interpreter','latex');
set(gca,'FontSize',20)
subplot(2,3,4);
yyaxis left
plot(P2.kI_bat2.Time/60,ones(size(P2.kI_bat2.Time/60)),'k-');%P2.kI_bat2.Data(:,1)
ylim([0 3])
yticks([0:1:3])
grid('minor');
ylabel("$\sigma$(3)",'FontName','Times New Roman','FontSize',24,'interpreter','latex');
hold on
yyaxis right
plot(P2.kI_bat2.Time/60,0*[zeros(size(P2.kI_bat2.Time(P2.kI_bat2.Data<=1)))' ones(size(P2.kI_bat2.Time(P2.kI_bat2.Data>1)))'],'r-.')
ylim([0 1.05])
yticks([0:1:1.05])
xlim([0 100])
xticks([0:20:100])
grid('minor');
grid('minor');
xlabel("Time [min]",'FontName','Times New Roman','FontSize',24,'interpreter','latex');
ylabel("$D^{(6,1)}$",'FontName','Times New Roman','FontSize',24,'interpreter','latex');
title('(c)','FontName','Times New Roman','FontSize',24)
set(gca,'FontSize',20)
subplot(2,3,5);
yyaxis left
plot(P2.kV_bat2.Time/60,P2.kV_bat2.Data(:,1),'k-');
ylim([0 3])
yticks([0:1:3])
grid('minor');
ylabel("$\sigma$(4)",'FontName','Times New Roman','FontSize',24,'interpreter','latex');
hold on
yyaxis right
plot(P2.kV_bat2.Time/60,[zeros(size(P2.kV_bat2.Time(P2.kV_bat2.Data<=1)))' ones(size(P2.kV_bat2.Time(P2.kV_bat2.Data>1)))'],'r-.')
ylim([0 1.05])
yticks([0:1:1.05])
xlim([0 100])
xticks([0:20:100])
grid('minor');
grid('minor');
xlabel("Time [min]",'FontName','Times New Roman','FontSize',24,'interpreter','latex');
ylabel("$D^{(6,2)}$",'FontName','Times New Roman','FontSize',24,'interpreter','latex');
title('(d)','FontName','Times New Roman','FontSize',24)
set(gca,'FontSize',20)
subplot(2,3,6);
Controller_change=[2*ones(size(P2.Power_Split_2.Time(P2.Power_Split_2.Time<t_fault2)))' ones(size(P2.Power_Split_2.Time(P2.Power_Split_2.Time>=t_fault2)))'];
yyaxis left
plot(P2.Power_Split_2.Time/60,Controller_change,'k-')
xlabel("Time [min]",'FontName','Times New Roman','FontSize',24,'interpreter','latex');
ylabel("$\sigma$(5)",'FontName','Times New Roman','FontSize',24,'interpreter','latex');
ylim([0 3])
yticks([0:1:3])
xlim([0 140])
xticks([0:20:140])
grid('minor');
hold on
change_signal=[ones(size(P2.Power_Split_2.Time(P2.Power_Split_2.Time<t_fault2)))' 0 ones(size(P2.Power_Split_2.Time(P2.Power_Split_2.Time>t_fault2)))']
yyaxis right
plot(P2.Power_Split_2.Time/60,change_signal,'r-.')
yticks([0:1:1.05])
ylim([0 1.05])
ylabel("$\Sigma\times\mathcal{F}\mapsto\mathcal{I}$","FontSize",24,"Interpreter",'latex')
title('(e)','FontName','Times New Roman','FontSize',24)
set(gca,'FontSize',20)
%figure
% Switching signals for battery DoS events [Profile 2]
%figure
%%
plot(t_exp/60,y_exp_P1,'b-','LineWidth',2)
hold on
plot(t_exp/60,y_exp_P2,'m--','LineWidth',2)
hold on
plot(t_exp/60,P_aux,'r-.','LineWidth',2)
leg1=legend("$P_D$ (Mission 1)","$P_D$ (Mission 2)","$P_{aux}$")
grid('minor')
xlabel("Time[min]",'FontName','Times New Roman','FontSize',24,'interpreter','latex')
ylabel('Power [kW]','FontName','Times New Roman','FontSize',24,'interpreter','latex')
set(gca,'FontSize',20)
set(leg1,'Interpreter','latex');
set(leg1,'FontSize',17);