function [bin] = classifyNB(type, data)
% This function uses the Naive Bayes classifier to take in the features
% from prepareData() and outputs the bin it belongs in. There are three
% bins: bin 1 holds the more non-colored points, bin 2 holds points that
% may be seen as colored or non-colored, bin 3 holds points that are more
% colored points.
%
% Inputs: type = 'CIELAB', 'sCIELAB', or 'HASLER', data is a matrix where 
% each row is a different point, and columns contain: distance in cm, 
% frequency/fov, angle in degrees, cielab or scielab error value
%
% Outputs: vector of bin numbers in which each point has been classified


%% Check type and etc inputs
if (nargin ~= 2), error('Too few input parameters. Need a type and a data matrix.'); end
if (~ischar(type)), error('Data ''type'' must be a string: ''sCIELAB'' or ''CIELAB'''); end
if (strcmp(upper(type),'SCIELAB') == 0 & strcmp(upper(type),'CIELAB') == 0 & strcmp(upper(type),'HASLER') == 0)
    error('You must specify the metric to use: ''sCIELAB'', ''CIELAB'', or ''HASLER''');
end

% load Naive Bayes training data
load NBmodels.mat

datasize = size(data);

for i = 1:datasize(1),
    
    d_in = data(i,1);
    f_in = data(i,2);
    a_in = data(i,3);
    err_in = data(i,4);
    
    % discretize input values
%     if(d_in < 9*2.54)
%         d_in = 6*2.54;
%     else if(d_in < 18*2.54)
%             d_in = 12*2.54;
%         else if(d_in < 30*2.54)
%                 d_in = 24*2.54;
%             else d_in = 36*2.54;
%             end
%         end
%     end
    
    if(d_in < 9*2.54)
        d_in = 1;
    else if(d_in < 18*2.54)
            d_in = 2;
        else if(d_in < 30*2.54)
                d_in = 3;
            else d_in = 4;
            end
        end
    end
    
    if(f_in < 3) f_in = 2;
    else if(f_in < 5) f_in = 4;
        else if(f_in < 7) f_in = 6;
            else if(f_in < 9) f_in = 8;
                else if(f_in < 11) f_in = 10;
                    else if(f_in < 13) f_in = 12;
                        else if(f_in < 15) f_in = 14;
                            else if(f_in < 17) f_in = 16;
                                else if(f_in < 19) f_in = 18;
                                    else f_in = 20;
                                    end
                                end
                            end
                        end
                    end
                end
            end
        end
    end

    if(a_in < 5) a_in = 0;
    else if(a_in < 15) a_in = 10;
        else if(a_in < 25) a_in = 20;
            else if(a_in < 35) a_in = 30;
                else if(a_in < 45) a_in = 40;
                    else if(a_in < 55) a_in = 50;
                        else if(a_in < 65) a_in = 60;
                            else if(a_in < 75) a_in = 70;
                                else if(a_in < 85) a_in = 80;
                                    else a_in = 90;
                                    end
                                end
                            end
                        end
                    end
                end
            end
        end
    end

    for c = 1:NUMBINS,
%        I_d = find(dist_cm == d_in);    % index of dist input
        I_d = d_in;
        I_f = find(freq == f_in);       % index of freq input
        I_a = find(ang == a_in);        % index of angle input
    
        if(strcmp(upper(type),'CIELAB') == 1) % for cielab error
            % discretize error
            if(err_in <= xe2)
                I_e = 1;
            else if(err_in <= xe3)
                    I_e = 2;
                else if(err_in <= xe4)
                        I_e = 3;
                    else if(err_in <= xe5)
                            I_e = 4;
                        else if(err_in <= xe6)
                                I_e = 5;
                            else
                                I_e = 6;
                            end
                        end
                    end
                end
            end
            eval([ 'out(c) = p_c(c)*p_dist_' num2str(c) '(I_d)*p_freq_' num2str(c) '(I_f)*p_ang_' num2str(c) '(I_a)*p_err_' num2str(c) '(I_e);' ]);
        elseif(strcmp(upper(type),'SCIELAB') == 1) % for scielab error
            % discretize error
            if(err_in <= xs2)
                I_e = 1;
            else if(err_in <= xs3)
                    I_e = 2;
                else if(err_in <= xs4)
                        I_e = 3;
                    else if(err_in <= xs5)
                            I_e = 4;
                        else if(err_in <= xs6)
                                I_e = 5;
                            else
                                I_e = 6;
                            end
                        end
                    end
                end
            end
            eval([ 'out(c) = p_c(c)*p_dist_' num2str(c) '_s(I_d)*p_freq_' num2str(c) '_s(I_f)*p_ang_' num2str(c) '_s(I_a)*p_err_' num2str(c) '_s(I_e);' ]);
        else    % for Hasler metric
            % discretize value
            if(err_in <= xh2)
                I_e = 1;
            else if(err_in <= xh3)
                    I_e = 2;
                else if(err_in <= xh4)
                        I_e = 3;
                    else if(err_in <= xh5)
                           I_e = 4;
                        else if(err_in <= xh6)
                                I_e = 5;
                            else
                                I_e = 6;
                            end
                        end
                    end
                end
            end
            eval([ 'out(c) = p_c(c)*p_dist_' num2str(c) '_h(I_d)*p_freq_' num2str(c) '_h(I_f)*p_ang_' num2str(c) '_h(I_a)*p_err_' num2str(c) '_h(I_e);' ]);
        end
        
    [Y,bin_i] = max(out);
    bin(i,:) = bin_i - 1;
    end
end

% % Plot bin results
% 
%     figure()
%     index = find(bin==0);
%     temp_x = data(index,:);
% 
%     % plot the data; change the index into temp_x to see comparisons among
%     % different dimensions; 
%     plot(temp_x(:,4), temp_x(:,1), 'bo');
%     hold on
% 
%     index = find(bin==1);
%     temp_x = data(index,:);
%     plot(temp_x(:,4), temp_x(:,1),'g^');
%     hold on
% 
%     index = find(bin==2);
%     temp_x = data(index,:);
%     plot(temp_x(:,4), temp_x(:,1),'rx');
%     hold on
% 
%     grid on
% 
%     hold off
%     title(['Classification Results']);
%     xlabel('Metric');
%     ylabel('Distance (cm)');
