function [alpha,theta,solution,t,alpha1,alpha2]=ekozinec2(X,I,epsilon,tmax)% EKOZINEC2 finds epsilon-optimal decision hyperplane, dichotomy.%  [alpha,theta,solution,t,alpha1,alpha2]=ekozinec2(X,I,epsilon,tmax)%% EKOZINEC2 is faster version of EKOZINEC (see help ekozinec).%%  Input:%   X [DxM] contains M training points in D-dimensional the feature %      space. X=[x1,x2,..XM] where xi is i-th column vectors.%   I [1xM] contains class labels of the points in X. A class label %      has  to be 1 for the first class and 2 for the second class.%   epsilon [1x1] desired quality of the solution.%   tmax [1x1] maximal number of iterations.%%  Output:%   alpha [Dx1] found normal vector of the hyperplane or, if tmax==-1,%      badly classified point.%   theta [1x1] found threshold of the hyperplane.%   solution [1x1], 1 ... solution is found,%                   0 ... solution is not found.%   t [1x1] is number of steps the algorithm performed.%% See also EKOZINEC, KOZINEC, PERCEPTR, LINSVM.%% Statistical Pattern Recognition Toolbox, Vojtech Franc, Vaclav Hlavac% (c) Czech Technical University Prague, http://cmp.felk.cvut.cz% Written Vojtech Franc (diploma thesis) 02.11.1999, 13.4.2000% Modifications% 21-apr-2001, V.Franc, created% Process input argumentsif nargin < 5,   t=0;   alpha1=0;   alpha2=0;endif nargin < 4,   tmax=inf;endX1=X(:,find(I == 1));X2=X(:,find(I == 2));n1 = size(X1,2);n2 = size(X2,2);alpha1=X1(:,1);alpha2=X2(:,2);solution=0;t = 0;while solution == 0 & tmax > t,  t = t + 1;  solution = 1;  da1 = (alpha1 - alpha2);  nda1 = norm( da1 );    XA2 = X1 - repmat( alpha2, 1, n1 );    [minXDA1, inx] = min( XA2' * da1 /nda1 );    if minXDA1 <= (nda1 - epsilon/2 ),    xt = X1(:, inx );    k = min(1, (da1'*(alpha1-xt)) / ((alpha1-xt)'*(alpha1-xt)));        alpha1 = alpha1 * (1 - k) + xt * k;            solution = 0;  else    da2 = (alpha2 - alpha1 );    nda2 = norm( da2 );      XA1 = X2 - repmat( alpha1, 1, n2 );      [minXDA2, inx ] = min( XA1' * da2  / nda2 );      if minXDA2 <= (nda2 - epsilon/2 ),      xt = X2(:, inx );      k = min(1, (da2'*(alpha2-xt)) / ((alpha2-xt)'*(alpha2-xt)));      alpha2 = alpha2 * (1 - k) + xt * k;            solution = 0;    end  end    endtheta = 0.5 * ( alpha1'*alpha1 - alpha2'*alpha2 );alpha = alpha1 - alpha2;return;