function sober
% This function builds on top of ourRankOpt with emphasis on the following
% points: 1) no default setting of evaluation bias if a predicate is not
% observed

fail = load('failtrace'); % positive samples
pass = load('passtrace'); % negative samples

disp('done with loading data files');

% [m, col] = size(fail); % get number of columns from failing matrix, i.e., the smaller one, to expedite computation
[n, col] = size(pass);
N = col / 3; 

i = 0: (N-1);

% getting model from correct runs
pass_true = pass(:, 3 * i + 1);
pass_sum = pass(:, 3 * i + 3);

pass_bias = pass_true ./ pass_sum; % this step can be further optimized if we inline pass_true and pass_sum
nan_index = isnan(pass_bias); % index for units that are effective, i.e., not NaN 
pass_bias(nan_index) = 0.5;
% miu_p = zeros(1, N);
% sigma_p = zeros(1, N);

miu_p = mean(pass_bias);
sigma_p = std(pass_bias);


fail_true = fail(:, 3 * i + 1);
fail_sum = fail(:, 3 * i + 3);

fail_bias = fail_true ./ fail_sum;
nan_index = isnan(fail_bias);
effect_index_fail = ~nan_index; 
fail_bias(nan_index) = 0.5;
effect_counts = sum(effect_index_fail);
prune1 = (effect_counts == 0); % pruning rule #1: predicate p is of no interests if it is never observed


scores = zeros(1, N);
c = zeros(1, N);

[m, dummy] = size(fail_bias);

c = (sum(fail_bias, 1) - m * miu_p) / sqrt(m);% c is defined the same as in paper when proving the incoporation of software invariants i.e., c = \frac{\sum_{i=1}^m x_i - m \miu_p}{\sqrt{m}}.
scores = log(sigma_p) - 0.5 * log(m) + 0.5 * (c ./ sigma_p).^2; % Originally, score should be calculated in this way, for the jth predicate, scores(j) = log(sigma_p(j)) - 0.5 * log(mj) + 0.5 * (c(j)/sigma_p(j))^2; However, since m is uniform for all predicates, the above calculation does not take m into account. 


sigma_p_index = (sigma_p == 0);
c_index_eqzero = (c == 0);
c_index_notzero = not(c_index_eqzero);

prune2 = sigma_p_index & c_index_eqzero;

% scores(sigma_p_index & c_index_notzero) = +inf;
scores(prune1 | prune2) = -inf; % prune it
% till now, scores should have no NaN any more.

sort_score = -scores; % for descending order sort
[dummy, index] = sort(sort_score);

values = scores(index);

fid = fopen('rankResultWe', 'w');
% fid2 = fopen('ourWithValue', 'w');
for j = 0: (N/2 - 1) % due to complementaryness
   i = 2 * j + 1;
   fprintf(fid, '%d\n', index(i));
  % fprintf(fid2, '%d\t%f\n', index(i), values(i));
end
fclose(fid);

return;