function filename = SVM(fname)
%
% This function runs the support vector machine trainer--meaning it
% calculates the parameters for the support vector machine and stores them
% in the file specified by FNAME.
%
% See also CLASSIFYSVM
%

%% Check type
if (nargin ~= 1), error('Too few (or too many) input parameters. Need a model filename.'); end
if (~ischar(fname)), error('Filename ''fname'' must be a valid filename.'); end

%% Check to make sure SVM package is "installed"
whereIsSVM = which('svm_learn');
if(~exist(whereIsSVM))
    error('Cannot find SVM_light package.');
end

[path func ext] = fileparts(whereIsSVM);

%% Set SVM Params
params = '''ExecPath'', path, ''Kernel'', 0, ''C'', 1'; %''Verbosity'',2

%% Provide data
[xe xs xh y] = prepareData('SVM');

norm_xe = [];
norm_xs = [];
norm_xh = [];

for M = 1:size(xe,1),
    norm_xe = [norm_xe; norm(xe(M,:))];
    norm_xs = [norm_xs; norm(xs(M,:))];
    norm_xh = [norm_xh; norm(xh(M,:))];
end

xe = xe./repmat(norm_xe,1,size(xe,2));
xs = xs./repmat(norm_xs,1,size(xs,2));
xh = xh./repmat(norm_xh,1,size(xh,2));

x={xe, xs, xh};
type = {'CIELAB' 'sCIELAB' 'hasler'};
tmp_file = 'temp.tmp';

for K = 1:length(type),
    svmlwrite(tmp_file, x{K}, y);
    [Y, X] = svmlread(tmp_file);

    % call the SVM routine
    eval(['net = svml([pwd ''\'' fname ''.' type{K} '''], ' params ');']);

    net = svmltrain(net, X, Y);
    eval(['SVMmodels.net_' type{K} ' = net;']);
    delete(tmp_file);
end
%delete *.transduction;

save SVMmodels.mat SVMmodels;
filename = 'SVMmodels.mat';

%% Real Code Ends Here
%% Make Pretty Plots
% modify the assignment to temp_x to display different vectors--either x
% (CIELAB error) or xs (sCIELAB error)

% index = find(y==0);
% temp_x = x(index,:);
% plot(temp_x(:,4), temp_x(:,1),'g^');
% hold on
data = x;
labels = y;
for K = 1:length(type),
    y = labels;
    x=data{K};
    figure()
    if (~isempty(find(y==0)))
        fid = fopen([fname '.' type{K} '.transduction']);
        update = textscan(fid,'%n%n%n%n','delimiter',':');
        fclose(fid);
        y_new = update{2};
        y(find(labels==0)) = y_new;
    end

    index = find(y==-1);
    temp_x = x(index,:);
    % plot the data; change the index into temp_x to see comparisons among
    % different dimensions; temp_x(:,1) is constant so it's useless to plot
    % against it (just ignore it)
    plot(temp_x(:,4), temp_x(:,1), 'bo');
    hold on

    index = find(y==1);
    temp_x = x(index,:);
    plot(temp_x(:,4), temp_x(:,1),'rx');
    hold on

    % index = find(y==3);
    % temp_x = x(index,:);
    % plot(temp_x(:,5), temp_x(:,2),'rx');

    grid on

    v = axis;
    x1 = v(3):.01:v(4);

    [Y X]= svmlread([fname '.' type{K}]);
    b = Y(9); % almost always the 9th entry is the "b" value
    alpha = Y(10:end);
    v = X(10:end,:);
    w = alpha'*v;
    

    %[x1 x2]= meshgrid(v(3):v(4), v(5):v(6));
    %x1 = [v(5):.01:v(6)];
    %x2 = [v(3):.01:v(4)];
    % Since we're looking for h_t = 0.5, that is equivalent to (t*x')' = 0. So
    % we want t(1) + t(2)*x1 + t(3)*x2 + t(4)*x3 + t(5)*x4 = 0. Ignoring the x2 
    % and x3 dimensions, and using the x4 vector above, we get the x1 vector
    % below:

    x4 = (b - w(1)*x1)/w(4);
    x5 = (b - w(1)*x1-1)/w(4);
    x6 = (b - w(1)*x1+1)/w(4);

    plot(x4,x1,x5,x1,'r--',x6,x1,'g--')
    hold off
    title(type{K});
        xlabel('Metric');
    ylabel('Distance (cm)');
end