[Matlab] Histogram Equalization

% BI_histeq.m
% Histogram Equalization operation
% Written by Barak Itkin, http://barak-itkin.blogspot.com

%%%%%%%%%%%%%%%%%%%% Image commands "cheat sheet" %%%%%%%%%%%%%%%%%%%%%%%%%

% im = imread (path);           % read an image file (try 'cameraman.tif'
%                               % and 'peppers.png')
% info = imfinfo (path);        % Get information about the image (use the
%                               % same path!). Try printing this ;)
%
% if info.ColorType == 'truecolor' % Check if the image is RGB type
%    img = rgb2gray(im);        % If so, convert to gray 0-255
% else                          % Otherwise (if gray), do nothing
%    img = im;                  % Nothing :)
% end                           % End of if statement
%
% imgd = double(img) / 255;     % convert from range 0-255 to 0-1
% imshow (imgd);                % Display some image in rgb/gray/0-1/0-255

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%% BI_histeq (img, outputLevels) %%%
% img           = Input image, one gray channel with values in range 0-1
% outputLevels  = Maximal amount of output levels in the equalized
%                 histogram, 256 if no value is explicitly specified
% imeq (return) = an image with equalized histogram

% Example useage - imshow (BI_histeq (imgd)
function imeq = BI_histeq (img, outputLevels)

% Check if number of arguments was less than 2. If so, then no value was
% passed for the outputLevels parameter - so use 256
if nargin < 2
   outputLevels = 256;
end

% Create a result image, initialize it to zeros
imeq = zeros (size(img));

% Total number of pixels (in source image)
total = numel (img);

% Create a sorted list of the output levels in ascending order
outstep = 1/(outputLevels-1);
b = 0:outstep:1;

% Create a sorted list of the intput levels in ascending order
a = unique(img);
% count the amount of input pixels
jlim = numel(a);

% Initialize loop parameters
% An index over the input levels
j = 1;
% A counter of how many pixels we already converted (range [0-total])
sigma = 0;


% Pass over all the output levels
for bi=b

   % Pass over input levels we haven't covered yet
   % Note that we are adding one step, to make sure that for output level
   % 0, we accumulate pixels until we reached the brightness of the next
   % step. This way, we won't have 0 pixels in output level 0
   while sigma < (bi + outstep) * total && j <= jlim

       % Create a matrix with 0's and 1's, indicating where in the source
       % image we have the current input level
       temp = (img == a(j));

       % Add to the result image, areas with the current input level
       % colored with the current output level
       imeq = imeq + temp * bi;

       % Add the amount of pixels we converted in this step, to the
       % total amount of pixels we converted
       % MATLAB: M(:) makes a vector containing all the elements in the
       %         matrix M. Remember also that temp had 0's and 1's
       %         indicating where we matched the input level (match=1)
       sigma = sigma + sum (temp(:));

       % Progress to the next input level
       j = j + 1;

   % End the iteration over input levels
   end

% End the iteration over output levels

end