clear
clc

deltaT = 0.1;
h_set = 0.5;
hinit = 0.5;
uinit = 0.5;
qininit = 0.5;

qin_set1 = 0.5;
qin_set2 = 0.7;
qin_set3 = 0.9;
qin_set4 = 0.95;

load('MPC5_newsim.mat')
mpc_case5_h_ISE=trapz( 0:deltaT:200, (h_set-[hinit; sim.x]).^2);
mpc_case5_h_IAE=trapz( 0:deltaT:200, abs(h_set-[hinit; sim.x]))
mpc_case5_qout_IAE=trapz(     0:deltaT:50,  abs(qin_set1-[uinit sim.u(1:500)]))+...
                   trapz(  50.1:deltaT:100, abs(qin_set2-sim.u(501:1000)))+...
                   trapz( 100.1:deltaT:150, abs(qin_set3-sim.u(1001:1500)))+...
                   trapz( 150.1:deltaT:200, abs(qin_set4-sim.u(1501:2000)))
mpc_case5_qin_IAE=trapz(     0:deltaT:50,  abs(qin_set1-[qininit qin(1:500)]))+...
                  trapz(  50.1:deltaT:100, abs(qin_set2-qin(501:1000)))+...
                  trapz( 100.1:deltaT:150, abs(qin_set3-qin(1001:1500)))+...
                  trapz( 150.1:deltaT:200, abs(qin_set4-qin(1501:2000)))

%load('MPC1.mat')
%load('MPC1_newtunings.mat')
%load('MPC1_newsim.mat')
load('MPC1_newsim1000.mat')
mpc1_h_ISE=trapz( 0:deltaT:200, (h_set-[hinit; sim.x]).^2);
mpc1_h_IAE=trapz( 0:deltaT:200, abs(h_set-[hinit; sim.x]))
% mpc1_qout_IAE=trapz(     0:deltaT:50,  abs(qin_set1-[uinit sim.u(1:100)]))+...
%               trapz(  50.5:deltaT:100, abs(qin_set2-sim.u(101:200)))+...
%               trapz( 100.5:deltaT:150, abs(qin_set3-sim.u(201:300)))+...
%               trapz( 150.5:deltaT:200, abs(qin_set4-sim.u(301:400)))+...
%               trapz( 200.5:deltaT:400, abs(qin_set5-sim.u(401:800)))
% mpc1_qin_IAE=trapz(     0:deltaT:50,  abs(qin_set1-[qininit qin(1:100)]))+...
%              trapz(  50.5:deltaT:100, abs(qin_set2-qin(101:200)))+...
%              trapz( 100.5:deltaT:150, abs(qin_set3-qin(201:300)))+...
%              trapz( 150.5:deltaT:200, abs(qin_set4-qin(301:400)))+...
%              trapz( 200.5:deltaT:400, abs(qin_set5-qin(401:800)))
mpc1_qout_IAE=trapz(     0:deltaT:50,  abs(qin_set1-[uinit sim.u(1:500)]))+...
              trapz(  50.1:deltaT:100, abs(qin_set2-sim.u(501:1000)))+...
              trapz( 100.1:deltaT:150, abs(qin_set3-sim.u(1001:1500)))+...
              trapz( 150.1:deltaT:200, abs(qin_set4-sim.u(1501:2000)))
mpc1_qin_IAE=trapz(     0:deltaT:50,  abs(qin_set1-[qininit qin(1:500)]))+...
             trapz(  50.1:deltaT:100, abs(qin_set2-qin(501:1000)))+...
             trapz( 100.1:deltaT:150, abs(qin_set3-qin(1001:1500)))+...
             trapz( 150.1:deltaT:200, abs(qin_set4-qin(1501:2000)))

clear sim.x;
clear sim.u;
clear qin;

%load('MPC2.mat')
% load('MPC2_newtunings.mat')
% load('MPC2_newsim.mat')
load('MPC2_newsim1000.mat')
mpc2_h_ISE=trapz( 0:deltaT:200, (h_set-[hinit; sim.x]).^2);
mpc2_h_IAE=trapz( 0:deltaT:200, abs(h_set-[hinit; sim.x]))
mpc2_qout_IAE=trapz(     0:deltaT:50,  abs(qin_set1-[uinit sim.u(1:500)]))+...
              trapz(  50.1:deltaT:100, abs(qin_set2-sim.u(501:1000)))+...
              trapz( 100.1:deltaT:150, abs(qin_set3-sim.u(1001:1500)))+...
              trapz( 150.1:deltaT:200, abs(qin_set4-sim.u(1501:2000)))
mpc2_qin_IAE=trapz(     0:deltaT:50,  abs(qin_set1-[qininit qin(1:500)]))+...
             trapz(  50.1:deltaT:100, abs(qin_set2-qin(501:1000)))+...
             trapz( 100.1:deltaT:150, abs(qin_set3-qin(1001:1500)))+...
             trapz( 150.1:deltaT:200, abs(qin_set4-qin(1501:2000)))      
clear sim.x;
clear sim.u;
clear qin;

%load('MPC3.mat')
% load('MPC3_newsim.mat')
load('MPC3_newsim1000.mat')
mpc3_h_ISE=trapz( 0:deltaT:200, (h_set-[hinit; sim.x]).^2);
mpc3_h_IAE=trapz( 0:deltaT:200, abs(h_set-[hinit; sim.x]))
mpc3_qout_IAE=trapz(     0:deltaT:50,  abs(qin_set1-[uinit sim.u(1:500)]))+...
              trapz(  50.1:deltaT:100, abs(qin_set2-sim.u(501:1000)))+...
              trapz( 100.1:deltaT:150, abs(qin_set3-sim.u(1001:1500)))+...
              trapz( 150.1:deltaT:200, abs(qin_set4-sim.u(1501:2000)))
mpc3_qin_IAE=trapz(     0:deltaT:50,  abs(qin_set1-[qininit qin(1:500)]))+...
             trapz(  50.1:deltaT:100, abs(qin_set2-qin(501:1000)))+...
             trapz( 100.1:deltaT:150, abs(qin_set3-qin(1001:1500)))+...
             trapz( 150.1:deltaT:200, abs(qin_set4-qin(1501:2000)))
clear sim.x;
clear sim.u;
clear qin;

% load('MPC4.mat')
% load('MPC4_newsim.mat')
load('MPC4_newsim1000.mat')
mpc4_h_ISE=trapz( 0:deltaT:200, (h_set-[hinit; sim.x]).^2);
mpc4_h_IAE=trapz( 0:deltaT:200, abs(h_set-[hinit; sim.x]))
mpc4_qout_IAE=trapz(     0:deltaT:50,  abs(qin_set1-[uinit sim.u(1:500)]))+...
              trapz(  50.1:deltaT:100, abs(qin_set2-sim.u(501:1000)))+...
              trapz( 100.1:deltaT:150, abs(qin_set3-sim.u(1001:1500)))+...
              trapz( 150.1:deltaT:200, abs(qin_set4-sim.u(1501:2000)))
mpc4_qin_IAE=trapz(     0:deltaT:50,  abs(qin_set1-[qininit qin(1:500)]))+...
             trapz(  50.1:deltaT:100, abs(qin_set2-qin(501:1000)))+...
             trapz( 100.1:deltaT:150, abs(qin_set3-qin(1001:1500)))+...
             trapz( 150.1:deltaT:200, abs(qin_set4-qin(1501:2000)))
clear sim.x;
clear sim.u;
clear qin;


