clear
format compact
nntwarn off

w = -1 : 0.1 : -0.5;
g = 0.5 : 0.1 : 1;

[x,y] = meshgrid(w,g);
P = [x(:) y(:)]';
T = P(1,:)./P(2,:);

[R, Q] = size(P);
S1 = 20;
S2 = 10;
[S3, Q] = size(T);

[W1, B1] = nwtan(S1, R);
[W2, B2] = nwtan(S2, S1);
[W3, B3] = rands(S3, S2);

[dW1, dB1] = nwtan(S1, R);
[dW2, dB2] = nwtan(S2, S1);
[dW3, dB3] = rands(S3, S2);

A1 = tansig(W1*P, B1);
A2 = tansig(W2*A1, B2);
A3 = purelin(W3*A2, B3);

E = T - A3;

lr = 0.1;
mc = 0.95;
max_epoch = 50000;
disp_freq = 100;
err_goal = 0.1;

SSE = sumsqr(E);

for epoch = 1:max_epoch
    if SSE < err_goal
        epoch = epoch - 1; break
    end

    D3 = deltalin(A3, E);
    D2 = deltatan(A2, D3, W3);
    D1 = deltatan(A1, D2, W2);

    [dW1, dB1] = learnbpm(P, D1, lr, mc, dW1, dB1);
    [dW2, dB2] = learnbpm(A1, D2, lr, mc, dW2, dB2);
    [dW3, dB3] = learnbpm(A2, D3, lr, mc, dW3, dB3);

    W1 = W1 + dW1;
    B1 = B1 + dB1;

    W2 = W2 + dW2;
    B2 = B2 + dB2;

    W3 = W3 + dW3;
    B3 = B3 + dB3;

    A1 = tansig(W1*P, B1);
    A2 = tansig(W2*A1, B2);
    A3 = purelin(W3*A2, B3);

    E = T - A3;

    SSE = sumsqr(E);

    if(rem(epoch, disp_freq) == 0)
        fprintf('\n %d epoch, SSE = %g', epoch, SSE)
        mesh(x,y, reshape(T, [length(x), length(y)]), 'edgecolor', 'g')
        hold on
        mesh(x,y,reshape(A3, [length(x), length(y)]), 'edgecolor', 'b')
        hold on
        mesh(x,y,reshape(A2, [length(x), length(y)]), 'edgecolor', 'y')
        hold off
        pause(1e-10)
    end
end
fprintf('\n %d epoch, SSE = %g', epoch, SSE)