
%% Testing polyfitn and wpolyfitn

close all; clear all; clf; clc

    u1 = 5*rand(20,1).*(-1).^(round(10*rand(20,1)));
    u2 = 5*sort(rand(20,1)).*(-1).^(round(10*rand(20,1)));
    
            figure(1)
            plot(u1,u2, 'b*')
            vector_position = [ 53         379        1067         569];
            hold on;
            for j=1:size(u1,1)
                if mod(j,2)==0
                    text(u1(j,1),u2(j,1),[' \leftarrow',num2str(j)],...
                        'VerticalAlignment','middle',...
                        'HorizontalAlignment','left',...
                        'FontSize',8,'Color','r')
                else
                    text(u1(j,1),u2(j,1),[num2str(j),' \rightarrow'],...
                        'VerticalAlignment','middle',...
                        'HorizontalAlignment','right',...
                        'FontSize',8,'Color','r')
                end
                
            end
            grid on;
            xlabel('u1')
            ylabel('u2')
            title('Distribution of u1 and u2')
            
print -dpng values_u1_and_u2

%% PLOT TRUE_VALUE FUNCTION

val_origin = 2;
for i=1:size(u1,1)
    
%     true_val(i) = u1(i).^2-u1(i)*u2(i)+2*u1(i)-u2(i).^2+3*u2(i)+3
%            z(i) = u1(i).^2-u1(i)*u2(i)+2*u1(i)-u2(i).^2+3*u2(i)+3; % the "fake" function values
      true_val(i) = exp(u1(i))-u2(i)+1;
             z(i) = exp(u1(i))-u2(i)+1 + 10*randn(1);
end


xfor_plot = [-6: 0.1: 6]';
yfor_plot = [-6: 0.1: 6]';
[mesh_x, mesh_y] = meshgrid(xfor_plot,yfor_plot);

for i=1:size(mesh_x,1)
    
%true_val_for_plot(:,i) =  mesh_x(:,i).^2-mesh_x(:,i).*mesh_y(:,i)+2*mesh_x(:,i)-mesh_y(:,i).^2+3*mesh_y(:,i)+3;
 true_val_for_plot(:,i) =  exp(mesh_x(:,i))-mesh_y(:,i)+1;
end

f2= figure(2)

mesh(xfor_plot,yfor_plot,true_val_for_plot);
xlabel('u1')
ylabel('u2')
zlabel('Z')

grid on
axis square
set(f2,'paperpositionmode','auto','Visible','on');
set(f2,'Colormap',[0 0 1;0.0322580635547638 0.0322580635547638 1;0.0645161271095276 0.0645161271095276 1;0.0967741906642914 0.0967741906642914 1;0.129032254219055 0.129032254219055 1;0.161290317773819 0.161290317773819 1;0.193548381328583 0.193548381328583 1;0.225806444883347 0.225806444883347 1;0.25806450843811 0.25806450843811 1;0.290322571992874 0.290322571992874 1;0.322580635547638 0.322580635547638 1;0.354838699102402 0.354838699102402 1;0.387096762657166 0.387096762657166 1;0.419354826211929 0.419354826211929 1;0.451612889766693 0.451612889766693 1;0.483870953321457 0.483870953321457 1;0.516129016876221 0.516129016876221 1;0.548387110233307 0.548387110233307 1;0.580645143985748 0.580645143985748 1;0.612903237342834 0.612903237342834 1;0.645161271095276 0.645161271095276 1;0.677419364452362 0.677419364452362 1;0.709677398204803 0.709677398204803 1;0.74193549156189 0.74193549156189 1;0.774193525314331 0.774193525314331 1;0.806451618671417 0.806451618671417 1;0.838709652423859 0.838709652423859 1;0.870967745780945 0.870967745780945 1;0.903225779533386 0.903225779533386 1;0.935483872890472 0.935483872890472 1;0.967741906642914 0.967741906642914 1;1 1 1;1 0.96875 0.96875;1 0.9375 0.9375;1 0.90625 0.90625;1 0.875 0.875;1 0.84375 0.84375;1 0.8125 0.8125;1 0.78125 0.78125;1 0.75 0.75;1 0.71875 0.71875;1 0.6875 0.6875;1 0.65625 0.65625;1 0.625 0.625;1 0.59375 0.59375;1 0.5625 0.5625;1 0.53125 0.53125;1 0.5 0.5;1 0.46875 0.46875;1 0.4375 0.4375;1 0.40625 0.40625;1 0.375 0.375;1 0.34375 0.34375;1 0.3125 0.3125;1 0.28125 0.28125;1 0.25 0.25;1 0.21875 0.21875;1 0.1875 0.1875;1 0.15625 0.15625;1 0.125 0.125;1 0.09375 0.09375;1 0.0625 0.0625;1 0.03125 0.03125;1 0 0]);
clim([-8 8])
title('Plot true function')

print -dpng polyfitn_true_function



%% GET WEIGHTS


z_normalized=z'./val_origin;
w = get_weights(z_normalized)


f3=figure(3)

plot(z_normalized,w,'k*')
ylim ([0, max(w)+min(w)])
xlabel ('Znorm value')
ylabel('Weights')
title('Z-value dependent weights')

print -dpng distribution_of_weights_23


%% PLOT UN_WEIGHTED FUNCTION
z=z';

[polymodel] = polyfitn([u1,u2],z,2);

f = polyvaln(polymodel,[u1,u2]);

f4 = figure(4);

        plot3(u1,u2,true_val,'*g', 'MarkerSize',15);
        hold on;

        for j=1:size(w,1)
            if w(j)==1
            plot3(u1(j),u2(j),z(j),'*y', 'MarkerSize',15);
            hold on;
            else
            plot3(u1(j),u2(j),z(j),'*k', 'MarkerSize',15);
            hold on;
            end
        hold on;
        end
        xlabel('u1')
        ylabel('u2')
        zlabel('Z')
        grid on
        axis square
        hold on;

        for i=1:size(mesh_x,1)
            ffor_plot(:,i) = polyvaln(polymodel,[mesh_x(:,i), mesh_y(:,i)]);
        end
        mesh(xfor_plot,yfor_plot,ffor_plot)

        set(f4,'paperpositionmode','auto','Visible','on');
        set(f4,'Colormap',[0 0 1;0.0322580635547638 0.0322580635547638 1;0.0645161271095276 0.0645161271095276 1;0.0967741906642914 0.0967741906642914 1;0.129032254219055 0.129032254219055 1;0.161290317773819 0.161290317773819 1;0.193548381328583 0.193548381328583 1;0.225806444883347 0.225806444883347 1;0.25806450843811 0.25806450843811 1;0.290322571992874 0.290322571992874 1;0.322580635547638 0.322580635547638 1;0.354838699102402 0.354838699102402 1;0.387096762657166 0.387096762657166 1;0.419354826211929 0.419354826211929 1;0.451612889766693 0.451612889766693 1;0.483870953321457 0.483870953321457 1;0.516129016876221 0.516129016876221 1;0.548387110233307 0.548387110233307 1;0.580645143985748 0.580645143985748 1;0.612903237342834 0.612903237342834 1;0.645161271095276 0.645161271095276 1;0.677419364452362 0.677419364452362 1;0.709677398204803 0.709677398204803 1;0.74193549156189 0.74193549156189 1;0.774193525314331 0.774193525314331 1;0.806451618671417 0.806451618671417 1;0.838709652423859 0.838709652423859 1;0.870967745780945 0.870967745780945 1;0.903225779533386 0.903225779533386 1;0.935483872890472 0.935483872890472 1;0.967741906642914 0.967741906642914 1;1 1 1;1 0.96875 0.96875;1 0.9375 0.9375;1 0.90625 0.90625;1 0.875 0.875;1 0.84375 0.84375;1 0.8125 0.8125;1 0.78125 0.78125;1 0.75 0.75;1 0.71875 0.71875;1 0.6875 0.6875;1 0.65625 0.65625;1 0.625 0.625;1 0.59375 0.59375;1 0.5625 0.5625;1 0.53125 0.53125;1 0.5 0.5;1 0.46875 0.46875;1 0.4375 0.4375;1 0.40625 0.40625;1 0.375 0.375;1 0.34375 0.34375;1 0.3125 0.3125;1 0.28125 0.28125;1 0.25 0.25;1 0.21875 0.21875;1 0.1875 0.1875;1 0.15625 0.15625;1 0.125 0.125;1 0.09375 0.09375;1 0.0625 0.0625;1 0.03125 0.03125;1 0 0]);
        clim([-8 8])

        title('Plot unweighted function')

print -dpng polyfitn_unweighted_function_23

%% PLOT WEIGHTED FUNCTION

[polymodel_weighted] = wpolyfitn([u1,u2],z,w,2);

f_weighted = polyvaln(polymodel_weighted,[u1,u2]);

for i=1:size(mesh_x,1)
ffor_plot_weighted(:,i) = polyvaln(polymodel_weighted,[mesh_x(:,i), mesh_y(:,i)]);
end

f5=figure(5);

        plot3(u1,u2,true_val,'*g', 'MarkerSize',15);
        hold on;

        for j=1:size(w,1)
            if w(j)==1
            plot3(u1(j),u2(j),z(j),'*y', 'MarkerSize',15);
            hold on;
            else
            plot3(u1(j),u2(j),z(j),'*k', 'MarkerSize',15);
            hold on;
            end
        hold on;
        end
        xlabel('u1')
        ylabel('u2')
        zlabel('Z')
        grid on
        axis square
        hold on;

        mesh(xfor_plot,yfor_plot,ffor_plot_weighted)
        set(f5,'paperpositionmode','auto','Visible','on');
        set(f5,'Colormap',[0 0 1;0.0322580635547638 0.0322580635547638 1;0.0645161271095276 0.0645161271095276 1;0.0967741906642914 0.0967741906642914 1;0.129032254219055 0.129032254219055 1;0.161290317773819 0.161290317773819 1;0.193548381328583 0.193548381328583 1;0.225806444883347 0.225806444883347 1;0.25806450843811 0.25806450843811 1;0.290322571992874 0.290322571992874 1;0.322580635547638 0.322580635547638 1;0.354838699102402 0.354838699102402 1;0.387096762657166 0.387096762657166 1;0.419354826211929 0.419354826211929 1;0.451612889766693 0.451612889766693 1;0.483870953321457 0.483870953321457 1;0.516129016876221 0.516129016876221 1;0.548387110233307 0.548387110233307 1;0.580645143985748 0.580645143985748 1;0.612903237342834 0.612903237342834 1;0.645161271095276 0.645161271095276 1;0.677419364452362 0.677419364452362 1;0.709677398204803 0.709677398204803 1;0.74193549156189 0.74193549156189 1;0.774193525314331 0.774193525314331 1;0.806451618671417 0.806451618671417 1;0.838709652423859 0.838709652423859 1;0.870967745780945 0.870967745780945 1;0.903225779533386 0.903225779533386 1;0.935483872890472 0.935483872890472 1;0.967741906642914 0.967741906642914 1;1 1 1;1 0.96875 0.96875;1 0.9375 0.9375;1 0.90625 0.90625;1 0.875 0.875;1 0.84375 0.84375;1 0.8125 0.8125;1 0.78125 0.78125;1 0.75 0.75;1 0.71875 0.71875;1 0.6875 0.6875;1 0.65625 0.65625;1 0.625 0.625;1 0.59375 0.59375;1 0.5625 0.5625;1 0.53125 0.53125;1 0.5 0.5;1 0.46875 0.46875;1 0.4375 0.4375;1 0.40625 0.40625;1 0.375 0.375;1 0.34375 0.34375;1 0.3125 0.3125;1 0.28125 0.28125;1 0.25 0.25;1 0.21875 0.21875;1 0.1875 0.1875;1 0.15625 0.15625;1 0.125 0.125;1 0.09375 0.09375;1 0.0625 0.0625;1 0.03125 0.03125;1 0 0]);
        clim([-8 8])
        title('Plot weighted function')

print -dpng polyfitn_weighted_function_23

f6=figure(6)

    subplot(3,1,1)
    bar([f-true_val',f_weighted-true_val']);
    legend ('Non-weighted errors','Weighted errors')
    xtck = 1:size(w,1);
    set(gca,'Xtick',xtck)
    xlabel('Points')
    ylabel('Polynomial value - true value')
    hold on;
    
    subplot(3,1,2)
    bar(w,0.4,'m');
    legend ('Weights')
    xtck = 1:size(w,1);
    set(gca,'Xtick',xtck)
    xlabel('Points')
    ylabel('Weights')
    
    subplot(3,1,3)
    bar(z_normalized,0.4,'b');
    legend ('Znorm. value')
    xtck = 1:size(w,1);
    set(gca,'Xtick',xtck)
    xlabel('Points')
    ylabel('Normalized Z value')
    
set(f6,'paperpositionmode','auto','Visible','on');
set(f6,'position',[280         172        1166         736]);
  
print -dpng illustration_errors_versus_weights_23

 Coef_init = polymodel.Coefficients
 Coef_weighted = polymodel_weighted.Coefficients
    
 RMSE_init = polymodel.RMSE
 RMSE_weighted = polymodel_weighted.RMSE
 
 R2_init = polymodel.R2
 R2_weighted = polymodel_weighted.R2
