% msf4.m
%clear

load fnnwb4
%load plantv
t=0.05;
e=-2:0.05:2;
T0=0.005:0.005:20;
in=0.6104*sin(0.1*pi*T0);  %20000/32768=0.6104 is plant can track
T=in(1:100:4000);  

[n,N]=size(e);
%k3=50;
ec=e;
P=[e;ec];
NWB=wbg4;

i=1;
 W1=[NWB(i,1) 0;NWB(i,2) 0;NWB(i,3) 0;NWB(i,4) 0;NWB(i,5) 0;0 NWB(i,6);0 NWB(i,7);0 NWB(i,8)];
  B1=[NWB(i,9);NWB(i,10);NWB(i,11);NWB(i,12);NWB(i,13);NWB(i,14);NWB(i,15);NWB(i,16)];
  W2=[NWB(i,17) NWB(i,18) NWB(i,19) NWB(i,20) NWB(i,21)...
       NWB(i,22) NWB(i,23) NWB(i,24) NWB(i,25)...
       NWB(i,26) NWB(i,27) NWB(i,28) NWB(i,29) NWB(i,30)...
       NWB(i,31) NWB(i,32)];

%Calculate the fitness values

% NETWORK PARAMETERS
 %W1 = [W10(:,1) zeros(5,1); zeros(5,1) W10(:,2)] ;

lfun=gause(W1*P,B1);
xp=multil(lfun);
cm=purelin(W2*xp);


plot(cm),ylabel('control value'),
pause,

%plot([e'*327.68 ec'*327.68]),
%ylabel('error'),

%plot([e(5:30)'*655.36 ec(5:30)'*655.36]),
%ylabel('error'),pause
plot(e,lfun(1,:)),
hold on
plot(e,lfun(2,:)),
plot(e,lfun(3,:)),
plot(e,lfun(4,:)),

hold off,
pause,

plot(e,lfun(5,:)),
hold on
plot(e,lfun(6,:)),
plot(e,lfun(7,:)),
plot(e,lfun(8,:)),

hold off
pause

m=1;
for e=-2:0.2:2
    n=1;
    for ec=-2:0.2:2
        lfun=gause(W1*[e;ec],B1);
	xp=multil(lfun);
	UR(m,n)=purelin(W2*xp);      
        n=n+1;
    end
    m=m+1;
end



%E=[-2 -1.5 -1 -0.5 0 -0.5 1.0 1.5 2.];
%EC=E;
e=-2:0.2:2;
ec=e;

x=ec;
y=e;
z=UR;

mesh(x,y,z);
pause
surface(x,y,z);  % good
pause,
meshz(x,y,z);	% the best