Local histogram equalization manually

I am implementing local histogram equalization manually but the result is no satisfactory What i want and what i am getting is in picture ...

here is code:

I=rgb2gray(imread('peppers.png'));

if (isa(I,'uint8'))
  I=double(I)/255;
end
if (size(I,3)==3)
  I=(I(:,:,1)+I(:,:,2)+I(:,:,3))/3; % average the RGB channels
end

windowsize=17;

% Create an empty array
Ieq = zeros(size(I,1),size(I,2));
% Apply this over a NxN region around each pixel (N is odd)
n = floor(windowsize/2); % <-- N is the half size of the region ie. [-N:N,-N:N]

for r=1+n:size(I,1)-n
  for c=1+n:size(I,2)-n
    % -- INSERT YOUR CODE BELOW ------------------------------------------
    % NOTE: For pixels near the boundary, ensure the NxN neighbourhood is still
    % inside the image (this means for pixels near the boundary the pixel may
    % not be at the centre of the NxN neighbourhood).

      if r-n <=1
      fromrow=1;
      torow=r+n;
      else
          fromrow=abs(r-n);
          if n+r >= size(I,1)
              torow=size(I,1);
          else
              torow=r+n;
          end

      end

      if c-n <= 1
          fromcol=1;
          tocol=c+n;
      else
          fromcol=abs(c-n);
          if c+n > size(I,2);
              tocol=size(I,2);
          else
              tocol=c+n;
          end
      end

      neighbour = I(fromrow:torow,fromcol:tocol);

      lessoreq=neighbour(neighbour<=I(r,c));

      sumofval=sum(lessoreq);

      pixval=sumofval/(size(neighbour,1)*size(neighbour,2));

      Ieq(r,c)=pixval;

      % -- INSERT YOUR CODE ABOVE ------------------------------------------

    end
  end
  imshow(Ieq);

Answer 1

[Image Analyst]

I'm assuming this isn't your homework since you didn't label it as homework, so this will work, however you can do it much faster with adapthisteq() or by using intlut() to make up a look up table rather than using a double for loop.

clc;    % Clear the command window.
close all;  % Close all figures (except those of imtool.)
clear;  % Erase all existing variables. Or clearvars if you want.
workspace;  % Make sure the workspace panel is showing.
format long g;
format compact;
fontSize = 20;

grayImage = rgb2gray(imread('peppers.png'));
subplot(1,2,1);
imshow(grayImage);
axis on;
drawnow;

if (isa(grayImage,'uint8'))
  grayImage=double(grayImage)/255;
end
if (size(grayImage,3)==3)
  grayImage=(grayImage(:,:,1)+grayImage(:,:,2)+grayImage(:,:,3))/3; % average the RGB channels
end

windowsize = 17;

% Create an empty array
Ieq = zeros(size(grayImage,1),size(grayImage,2));
% Apply this over a NxN region around each pixel (N is odd)
n = floor(windowsize/2); % <-- N is the half size of the region ie. [-N:N,-N:N]

for r=1+n:size(grayImage,1)-n
  fprintf('Processing row #%d\n', r);
  for c=1+n:size(grayImage,2)-n
    % -- INSERT YOUR CODE BELOW ------------------------------------------
    % NOTE: For pixels near the boundary, ensure the NxN neighbourhood is still
    % inside the image (this means for pixels near the boundary the pixel may
    % not be at the centre of the NxN neighbourhood).
    
    if r-n <=1
      fromrow=1;
      torow=r+n;
    else
      fromrow=abs(r-n);
      if n+r >= size(grayImage,1)
        torow=size(grayImage,1);
      else
        torow=r+n;
      end
      
    end
    
    if c-n <= 1
      fromcol=1;
      tocol=c+n;
    else
      fromcol=abs(c-n);
      if c+n > size(grayImage,2);
        tocol=size(grayImage,2);
      else
        tocol=c+n;
      end
    end
    neighbors = grayImage(fromrow:torow, fromcol:tocol);
    counts = imhist(neighbors, 256);
    cdf = cumsum(counts);
    cdf = cdf ./ cdf(end); % Normalize.
    % Get the original value of the central pixel that we need to replace.
    centerGrayLevel = grayImage(r, c);
    % Now see what it gets mapped to.
    index = ceil(256 * centerGrayLevel);
    newPixelValue = cdf(index);
    % Do the assignment
    Ieq(r,c) = newPixelValue;
    
    % -- INSERT YOUR CODE ABOVE ------------------------------------------
    
  end
end
subplot(1,2,2);
imshow(Ieq, []);
axis on;
msgbox('Done');

[Rashid Hussain]

Thanks a lot :) but on problem description these things we need to follow :

"That is, for each local neighbourhood we simply find the number of pixels of lesser or equal value to the pixel I(r,c) and use this to compute the new value for I(r,c). There is no need to compute any histograms or local mapping table m()."

[Image Analyst]

Then try this:

neighbors = grayImage(fromrow:torow, fromcol:tocol);
thisGreyLevel = grayImage(r, c);
theSum = sum(neighbors(:) <= thisGreyLevel);
newPixelValue = theSum / numel(neighbors);
% Do the assignment
Ieq(r,c) = newPixelValue;

[John BG]

and the question reads:

.. local histogram equalization manually ..

manually seems to be the key word here

[Rashid Hussain]

Image Analyst : It worked :) yaayy many many thanks

i didn't get what was the problem even procedure was same

[Image Analyst]

John, yeah, "manually" is always an ill-defined/ambiguous definition. I assume it meant he couldn't use the built-in histogram equalization functions histeq() or adapthisteq().

Answer 2

[Image Analyst]

Try using []:

imshow(Ieq, []);

also try using different windowSize's until you get the look you want.

[Rashid Hussain]

I have tried using this : imshow(Ieq, []); and also used different window but still no luck

[Image Analyst]

The problem is this:

lessoreq=neighbour(neighbour<=grayImage(r,c));

All that does is to set your new value equal to the average of the values in the neighborhood that are less than the central pixel's value. You need to get the histogram and then use cumsum() on the counts to get the cdf, then do an inverse look up to get the new value.