
clc
close all
clear fig
clear norm_comb_mean_1

% Calculation of several values
Np.a = 25;
Np.c = 5;
opt.Np = 5e3;
Np.tot = Np.a:Np.c:opt.Np;
n.PLS  = (opt.Np-Np.a)/Np.c+1;
% Np.tot = Np.a:Np.c:itg.t;
% n.PLS  = (itg.t-Np.a)/Np.c+1;
opt.normtresh = 0.01*Np.c;
opt.comptresh = 0.02;
n.norm = 25./Np.c;
n.inlet = 197;
n.comb = (1000-Np.a)/Np.c+1;

Ver.Ind = 'MolFrac';

% Selection of the chosen data points
switch Ver.Ind
case 'MolFlow'
    U.tra.sel = U.tra.samp;
case 'MolFrac'
    U.tra.sel = U.tra.Molfrac;
case 'PartPress'
    U.tra.sel = U.tra.PartPress;
end

% Calculation of the norms
if ~exist('norm_comb_mean_1','var') || isempty(norm_comb_mean_1{2,1})
    
% Predefinition of the matrices
if true
    xload      = cell(n.PLS,n.red);
    yload      = xload;
    xloads     = xload;
    load_used  = xload;
    xscores    = xload;
    yscores    = xload;
    stats      = xload;
    beta       = xload;
    pctvar     = xload;
    sumpct     = xload;
    diff_loads = xload;
    diff_mat   = cell(n.in,n.red);
    
    % 1-norm
    norm_comb_1           = xload;
    norm_comb_mean_1      = xloads;
    
    % 2-norm
    norm_comb_2           = xload;
    norm_comb_mean_2      = xloads;
    
    % inf-norm
    norm_comb_i           = xload;
    norm_comb_mean_i      = xloads;
    
    % frob-norm
    norm_comb_f           = xload;
    norm_comb_mean_f      = xloads;
    
    norm_loads = cell(n.PLS,n.red,n.in);
    norm_loads_mean = norm_loads;
    optpls = statset('UseParallel',true);
end
tic
for k = 1:n.PLS
    
    for j = 1:n.red     
        
        % Scaling of the independent variables
        [U.tra.scaled,Var.U.mu,Var.U.std] = zscore(U.tra.sel(1:Np.tot(k),:));

        % Scaling of the dependent variables
        [x.tra.scaled,Var.x.mu,Var.x.std]  = zscore(x.tra.red(1:Np.tot(k),j));
        
        % Calculation of the Partial least square regression
        [xload{k,j},yload{k,j},xscores{k,j},yscores{k,j},...
                 beta{k,j},pctvar{k,j},~,stats{k,j}] = plsregress(U.tra.scaled,x.tra.scaled,n.in,'options',optpls);
         
        load_used{k,j} = stats{k,j}.W;xload{k,j};%/(xload{k,j}'*xload{k,j});  
        
        
        % Rearrangement and normalization of the loads
%         xloads{k,j} = load_used{k,j}./(ones(n.in,1)*max(abs(load_used{k,j})));
        xloads{k,j} = load_used{k,j}./(ones(n.in,1)*sqrt(diag(stats{k,j}.W'*stats{k,j}.W)'));
        sumpct{k,j} = sum(pctvar{k,j}(2,:) > opt.comptresh);
                
        % Definition of the index of the maximum absolute and
        % adjustment of the loads
        [~,in.m] = max(abs(xloads{j}));
        in.m = in.m+(0:n.in:n.in*(n.in-1));
        in.c = sign(xloads{j}(in.m)) ~= sign(xloads{max(j-Np.c,1)}(in.m));
        xloads{j}(:,in.c) = -xloads{j}(:,in.c);
        
        if k > 1
            diff_loads{k,j} = (xloads{k,j}-xloads{k-1,j});
            
            % 1-norm
            norm_comb_1{k,j}   = norm(diff_loads{k,j}(:,1:min(sumpct{k,j},sumpct{k-1,j})),1);
            norm_comb_mean_1{k,j} = mean(cell2mat(norm_comb_1(max(k-n.norm,1)+1:k,j)));
            
            % 2-norm
            norm_comb_2{k,j}   = norm(diff_loads{k,j}(:,1:min(sumpct{k,j},sumpct{k-1,j})),2);
            norm_comb_mean_2{k,j} = mean(cell2mat(norm_comb_2(max(k-n.norm,1)+1:k,j)));
            
            % inf-norm
            norm_comb_i{k,j}   = norm(diff_loads{k,j}(:,1:min(sumpct{k,j},sumpct{k-1,j})),inf);
            norm_comb_mean_i{k,j} = mean(cell2mat(norm_comb_i(max(k-n.norm,1)+1:k,j)));
            
            % inf-norm
            norm_comb_f{k,j}   = norm(diff_loads{k,j}(:,1:min(sumpct{k,j},sumpct{k-1,j})),'fro');
            norm_comb_mean_f{k,j} = mean(cell2mat(norm_comb_f(max(k-n.norm,1)+1:k,j)));
            for l = 1:n.in
                norm_loads{k,j,l} = norm(diff_loads{k,j}(:,l),'fro');
                norm_loads_mean{k,j,l} = mean(cell2mat(norm_loads(max(k-n.norm,1)+1:k,j,l)));
                diff_mat{l,j}(:,k) = diff_loads{k,j}(:,l);
            end
        end
        test2{j}(:,k) = xloads{k,j}(:,1);
        test3{j}(:,k) = mean(test2{j}(:,max(k-n.norm,1)+1:k),2);
    end
    
    if mod(Np.tot(k),100) == 0
        fprintf('%i points are sampled.\n',Np.tot(k))
    end
    
end
toc
end

% Definition of integer values
itg.fig = 1;
itg.lab = 1;
itg.fda = 1;
itg.axe = 1;
itg.ann = 1;

% Definition of figure and values
fig.dimw12 = 12.3;  % 12.5 for ESCAPE
fig.dimh12 = 5.5;
fig.dimw21 = 8.4;
fig.dimh21 = 12;
fig.dimw = 8.4;
fig.dimh = 6;
fig.Font = [10 8];
fig.label = [];
fig.ann   = [];
fig.dat   = [];
%%
% Definition of the box values
box_x = [0 Np.tot(n.inlet)];
box_y = 10.^[-2 1];
box_x_v = box_x([1 2 2 1 1]);
box_y_v = box_y([1 1 2 2 1]);

% Definition of a colour string
rk1 = {'k','k','k','r','r','r','r'};%[0.93 0.69 0.13]}
ls1 = {'-','--','-.','-','--','-.',':'};
rk2 = {'k','k','r','r'};%[0.93 0.69 0.13]}
ls2 = {'-','--','-','--'};

% Plot of the differences for the individual weights
fig.in(itg.fig) = figure(itg.fig);
set(fig.in(itg.fig),'units','centimeters','position', [5+itg.fig,5+itg.fig, fig.dimw, fig.dimh]);
fig.axes(itg.axe) = gca;
fig.axes(itg.axe).Position(3) = 0.68;
fig.axes(itg.axe).Position(2) = 0.15;
fig.axes(itg.axe).Position(1) = 0.15;
hold on
for l = 1:n.in %sumpct{end,j}
    fig.dat(itg.fda) = plot(Np.tot(2:n.PLS),cell2mat(norm_loads_mean(1:n.PLS,1,l)),...
                            'Color',rk1{l},'LineStyle',ls1{l}); itg.fda = itg.fda+1;
end
[~,fig.leg1] = legend('$i = 1$','$i = 2$','$i = 3$','$i = 4$','$i = 5$','$i = 6$','$i = 7$');
set(fig.axes(itg.axe).Legend,'Location','best','FontSize',fig.Font(2),'Interpreter','latex','Box','off')
fig.label(itg.lab) = xlabel('Number of sampling points $n_p$'); itg.lab = itg.lab+1;
fig.label(itg.lab) = ylabel('Norm $||\Delta\mathbf{w}_{i}^k||_{F}^{av}$'); itg.lab = itg.lab+1;
fig.axes(itg.axe).YScale = 'log';
fig.axes(itg.axe).YTick = 10.^(-4:1:1);
fig.axes(itg.axe).YLim = [10^-4 10];
itg.axe = itg.axe+1;
itg.fig = itg.fig+1;

% Plot of the differences for the insignificant weights
fig.dat(itg.fda) = plot(box_x_v,box_y_v,'k'); itg.fda = itg.fda+1;
fig.in(itg.fig) = figure(itg.fig);
set(fig.in(itg.fig),'units','centimeters','position', [5+itg.fig,5+itg.fig, fig.dimw, fig.dimh]);
fig.axes(itg.axe) = gca;
hold on
for l = 1:n.in-sumpct{end}
    fig.dat(itg.fda) = plot(Np.tot(2:n.inlet),cell2mat(norm_loads_mean(1:n.inlet,1,l+sumpct{end})),...
                            'Color',rk1{l+sumpct{end}},'LineStyle',ls1{l+sumpct{end}}); itg.fda = itg.fda+1;
%                             'Color',rk2{l},'LineStyle',ls2{l}); itg.fda = itg.fda+1;
end
[~,fig.leg2] = legend('$i = 4$','$i = 5$','$i = 6$','$i = 7$');
set(fig.axes(itg.axe).Legend,'Location','best','FontSize',fig.Font(2),'Interpreter','latex','Box','off')
fig.label(itg.lab) = xlabel('Number of sampling points $n_p$'); itg.lab = itg.lab+1;
fig.label(itg.lab) = ylabel('Norm $||\Delta\mathbf{w}_{i}^k||_{F}^{av}$'); itg.lab = itg.lab+1;
fig.axes(itg.axe).YScale = 'log';
% fig.axes(itg.axe).YTick = 10.^(-4:1:1);
% fig.axes(itg.axe).YLim = [10^-4 10];
itg.axe = itg.axe+1;
itg.fig = itg.fig+1;

% Plot of the overall difference
fig.in(itg.fig) = figure(itg.fig);
set(fig.in(itg.fig),'units','centimeters','position', [5+itg.fig,5+itg.fig, fig.dimw, fig.dimh]);
fig.axes(itg.axe) = gca;
fig.dat(itg.fda) = plot(Np.tot(2:n.comb),cell2mat(norm_comb_mean_f(2:n.comb)),'r'); itg.fda = itg.fda+1;
hold on
fig.dat(itg.fda) = plot([0 Np.tot(n.comb)],[1 1]*opt.normtresh,'k--'); itg.fda = itg.fda+1;
fig.ann(itg.ann) = annotation('textbox',[.80 .24 .05 .05],'String','$\gamma$'); itg.ann = itg.ann+1;
set(fig.axes(itg.axe).Legend,'Location','Best','FontSize',fig.Font(1),'Interpreter','latex','Box','off')
fig.label(itg.lab) = xlabel('Number of sampling points $n_p$'); itg.lab = itg.lab+1;
fig.label(itg.lab) = ylabel('Norm $||\Delta\mathbf{W}^k||_{F}^{av}$'); itg.lab = itg.lab+1;
itg.axe = itg.axe+1;
itg.fig = itg.fig+1;
    
for k = n.in+1:2:3*n.in
    fig.leg1(k).XData = [fig.leg1(k).XData(1) 0.3];
end
fig.axes(1).Legend.Position = [0.85 0.4 0.2 0.5];


set(fig.ann, 'Interpreter', 'Latex');
set(fig.ann, 'Fontsize', fig.Font(1));
set(fig.ann, 'LineStyle', 'none');
set(fig.label,'Fontsize', fig.Font(1));
set(fig.label,'Interpreter', 'Latex');
set(fig.axes,'FontSize', fig.Font(2));
set(fig.axes,'Box','on');
set(fig.axes,'TickLabelInterpreter', 'Latex');
set(fig.axes,'LineWidth',.75);
set(fig.axes,'XMinorTick', 'on');
set(fig.axes,'YMinorTick', 'on');
set(fig.dat,'LineWidth',.75)
% fig.axes(2).FontSize = 6;

for k = 1:n.in
    fig.leg1(k).Position = [.35 fig.leg1(k).Position(2:3)];
end
%%
for k = 1:n.in
    fig.leg1(k).Interpreter = 'latex';
    fig.leg1(k).FontSize = fig.Font(2);
end

for k = 1:4
    fig.leg2(k).Interpreter = 'latex';
    fig.leg2(k).FontSize = fig.Font(2);
end
