% eda11_12
%
% Use a neural net to approximate a 3D step function

% ---Step function definitions---

% Nc: the number of centers
Nc = 1;

% Xc: the center of the step function
Xc = [1.5];
Yc = [1.5];

% Dx: the width of the step function
Dx = [1];
Dy = [1];

% h: the height of the step function
h = [1];

% slope: slope of the step function
slope = 999/4;

% ---Neural Net definitions---

% Nmax: number of layers;
Lmax = 4;

% Lmax: maximum number of neurons in any layer
Nmax = Nc * 4;

[N,w,b,a] = eda_net_init(Lmax,Nmax);

% z: a(1:N(i),i) is a column-vector that gives the
% argument of the activitites of the N(i) neurons on the i-th layer
z = zeros( Nmax, Lmax );

N = [2, Nc*4, Nc, 1]';

for ncvalue = 1:Nc

    w12 = 4*slope;
    w22 = 4*slope;
    w32 = 4*slope;
    w42 = 4*slope;
    w( ((ncvalue-1)*4)+1:((ncvalue-1)*4)+2, 1 ,2) = [w12 ; w22];
    w( ((ncvalue-1)*4)+3:((ncvalue-1)*4)+4, 2 ,2) = [w32 ; w42];
    
    w13 = 4*slope;
    w23 = -4*slope;
    w33 = 4*slope;
    w43 = -4*slope;
    w(ncvalue, ((ncvalue-1)*4)+1:((ncvalue-1)*4)+4,3) = [w13;w23;w33;w43];
    
    w14 = h(ncvalue);
    w(1, ncvalue, 4) = w14;

    b12 = -w12 * ( Xc(ncvalue) - Dx(ncvalue)/2);
    b22 = -w22 * ( Xc(ncvalue) + Dx(ncvalue)/2);
    b32 = -w32 * ( Yc(ncvalue) - Dy(ncvalue)/2);
    b42 = -w42 * ( Yc(ncvalue) + Dy(ncvalue)/2);
    b( ((ncvalue-1)*4)+1:((ncvalue-1)*4)+4 ,2) = [b12;b22;b32;b42];
    
    b13 = -(3*(4*slope))/2;
    b( ncvalue,3) = b13;
    
end 

% ---Picking the points to plot---
figure(5);
clf;
eda_net_view(N,w,b);
% Nx: the number of points we want to plot
Nx = 50;
Ny = 50;
% xmax: the maximum value of the points we want to plot
xmax = 3.1;
ymax = 3.1;
% xmin: the minimum value of the points we want to plot
xmin = 0;
ymin = 0;
% x: the points we are plotting
x = xmin + (xmax-xmin)* [0:Nx-1]'/(Nx-1);
y = ymin + (ymax-ymin)* [0:Ny-1]'/(Ny-1);
m= zeros(Nx,Ny);
zprev = zeros(Nx,Ny);

% ---Calculating the Activities---

for kx = [1:Nx]
    for ky = [1:Ny]
    
        a(1:N(1),1) = [x(kx);x(ky)]; 

        for i=[2:Lmax]
            z(1:N(i),i) = w(1:N(i),1:N(i-1),i) * a(1:N(i-1),i-1) + b(1:N(i),i);
            if i < Lmax
                a(1:N(i),i) = 1.0./(1.0+exp(-z(1:N(i),i)));
            else
                a(1:N(i),i) = z(1:N(i),i);
            end
        end
        zprev(kx,ky) = z(1,3);
        m(kx,ky) = a(1:N(3),4);
    end
end

hold on;

eda_draw(zprev,'caption z(1,3)','   ',m,'caption a(1,4)'); %mtrue,'caption True','   ',



