% 2006/1/12
% noisy cancellation
%
function [E,NetDef,W1,W2]=noisecan(x,y,M)
%---------- Generate training and test set ----------
%PHI=0:0.25:6;
%Y=sin(PHI);
%PHI1 = PHI(:,1:2:25);
%Y1   = Y(:,1:2:25);
%PHI2 = PHI(:,2:2:24);
%Y2   = Y(:,2:2:24);

PHI=x';
[d,N]=size(PHI);
Y=y';
PHI1 = PHI(:,1:1:N);
Y1   = Y(:,1:1:N);
%PHI2 = PHI(:,2:2:N-1);
%Y2   = Y(:,2:2:N-1);
%figure(1);
%plot(PHI1,Y1,'+')
%title('Data set for training')
%drawnow


% ----------- Training the network ----------
k=5;
while k~=7,

  %k = menu('Choose one of the training algorithms below:',...
  %          'Back Propagation (Batch)',...
  %          'Back Propagation (Recursive)',...
  %          'Recursive Prediction Error Method (forgetting factor)',...
  %          'Recursive Prediction Error Method (Constant Trace)',...
  %          'Recursive Prediction Error Method (EFRA)',...
  %          'Levenberg-Marquardt method',...
  %          'quit');
         
  %chvec = get(0,'Children'); % Close all figures except for figure(1)
  %chvec = chvec(find(chvec~=1));       
  %close(chvec);
  %---------- Define network structure and initialize weights ----------
  rand('seed',0);
  W1         = rand(M,d+1);  % Weights to hidden layer 
  W2         = rand(1,M+1);  % Weights to output
  NetDef1=[];
  NetDef2=[];
  for m=1:M
      NetDef1=[NetDef1 'H'];
      if m==1
          NetDef2=[NetDef2 'L'];
      else
          NetDef2=[NetDef2 '-'];
      end
  end
  NetDef = [NetDef1;NetDef2];
%  NetDef = ['HHHH'       
%            'L---'];       % The top row represents the structure of 
                           % the hidden layer. The bottom row represents
                           % the output structure.
  trparms = settrain;      % Set training parameters to default values

  % ----- Back propagation (Batch) -----
  if k==1,
    trparms=settrain(trparms,'maxiter',1000,'eta',0.01);
    [W1,W2,PI_vector,iter]=batbp(NetDef,W1,W2,PHI1,Y1,trparms);

  
  % ----- Back propagation (Recursive) -----
  elseif k==2,
    trparms=settrain(trparms,'maxiter',1000,'eta',0.01);
    [W1,W2,PI_vector,iter]=incbp(NetDef,W1,W2,PHI1,Y1,trparms);
    

  % ----- RPE algorithm (Forgetting factor) -----
  elseif k==3,
    trparms=settrain(trparms,'maxiter',200,'p0',200,'lambda',0.98);
    [W1,W2,PI_vector,iter]=rpe(NetDef,W1,W2,PHI1,Y1,trparms);

    
  % ----- RPE algorithm (Constant Trace) -----
  elseif k==4,
    trparms=settrain(trparms,'method','ct');
    trparms=settrain(trparms,'maxiter',200,'alpha_max',200,'alpha_min',0.001);
    [W1,W2,PI_vector,iter]=rpe(NetDef,W1,W2,PHI1,Y1,trparms);

  
  % ----- RPE algorithm (EFRA) -----
  elseif k==5,  
    trparms=settrain(trparms,'method','efra','maxiter',200,'alpha',1);
    [W1,W2,PI_vector,iter]=rpe(NetDef,W1,W2,PHI1,Y1,trparms);

  
  % ----- Marquardt algorithm -----
  elseif k==6,
    trparms=settrain(trparms,'maxiter',200);
    [W1,W2,PI_vector,iter,lambda]=marq(NetDef,W1,W2,PHI1,Y1,trparms);
  end

  if k~=7,
    
    % -----------  Validate Network  -----------
    [Y_sim,E,PI] = nneval(NetDef,W1,W2,PHI1,Y1,1);
    fprintf('MSE for training data %f\n',mean(E.^2));
    fprintf('ME for training data %f',mean(abs(E)));
    % -----------  Plot Cost function  -----------
    %figure
    %semilogy(PI_vector)
    %title('Criterion evaluated after each iteration')
    %xlabel('Iteration (epoch)')
    %ylabel('Criterion')
    %grid
  end
  k=7;
end
