Saturday, January 12, 2013

Histogram equalization using C++: Image Processing

Theory

The histogram equalization is an approach to enhance a given image. The approach is to design a transformation T such that the gray values in the output is uniformly distributed in [0, 1].

Algorithm

Compute a scaling factor, α= 255 / number of pixels
Calculate histogram of the image
Create a look up table LUT with
    LUT[0] =  α * histogram[0]
for all remaining grey levels, i, do
    LUT[i] = LUT[i-1] + α * histogram[i]
end for
for all pixel coordinates, x and  y, do
    g(x, y) = LUT[f(x, y)]
end for

Source Code : C++

#include <iostream>
#include <opencv2/highgui/highgui.hpp>
#include <opencv2/imgproc/imgproc.hpp>

using std::cout;
using std::cin;
using std::endl;

using namespace cv;

void imhist(Mat image, int histogram[])
{

    // initialize all intensity values to 0
    for(int i = 0; i < 256; i++)
    {
        histogram[i] = 0;
    }

    // calculate the no of pixels for each intensity values
    for(int y = 0; y < image.rows; y++)
        for(int x = 0; x < image.cols; x++)
            histogram[(int)image.at<uchar>(y,x)]++;

}

void cumhist(int histogram[], int cumhistogram[])
{
    cumhistogram[0] = histogram[0];

    for(int i = 1; i < 256; i++)
    {
        cumhistogram[i] = histogram[i] + cumhistogram[i-1];
    }
}

void histDisplay(int histogram[], const char* name)
{
    int hist[256];
    for(int i = 0; i < 256; i++)
    {
        hist[i]=histogram[i];
    }
    // draw the histograms
    int hist_w = 512; int hist_h = 400;
    int bin_w = cvRound((double) hist_w/256);

    Mat histImage(hist_h, hist_w, CV_8UC1, Scalar(255, 255, 255));

    // find the maximum intensity element from histogram
    int max = hist[0];
    for(int i = 1; i < 256; i++){
        if(max < hist[i]){
            max = hist[i];
        }
    }

    // normalize the histogram between 0 and histImage.rows

    for(int i = 0; i < 256; i++){
        hist[i] = ((double)hist[i]/max)*histImage.rows;
    }


    // draw the intensity line for histogram
    for(int i = 0; i < 256; i++)
    {
        line(histImage, Point(bin_w*(i), hist_h),
                              Point(bin_w*(i), hist_h - hist[i]),
             Scalar(0,0,0), 1, 8, 0);
    }

    // display histogram
    namedWindow(name, CV_WINDOW_AUTOSIZE);
    imshow(name, histImage);
}



int main()
{
    // Load the image
    Mat image = imread("scene.jpg", CV_LOAD_IMAGE_GRAYSCALE);

    // Generate the histogram
    int histogram[256];
    imhist(image, histogram);

    // Caluculate the size of image
    int size = image.rows * image.cols;
    float alpha = 255.0/size;

    // Calculate the probability of each intensity
    float PrRk[256];
    for(int i = 0; i < 256; i++)
    {
        PrRk[i] = (double)histogram[i] / size;
    }

    // Generate cumulative frequency histogram
    int cumhistogram[256];
    cumhist(histogram,cumhistogram );

    // Scale the histogram
    int Sk[256];
    for(int i = 0; i < 256; i++)
    {
        Sk[i] = cvRound((double)cumhistogram[i] * alpha);
    }


    // Generate the equlized histogram
    float PsSk[256];
    for(int i = 0; i < 256; i++)
    {
        PsSk[i] = 0;
    }

    for(int i = 0; i < 256; i++)
    {
        PsSk[Sk[i]] += PrRk[i];
    }

    int final[256];
    for(int i = 0; i < 256; i++)
        final[i] = cvRound(PsSk[i]*255);


    // Generate the equlized image
    Mat new_image = image.clone();

    for(int y = 0; y < image.rows; y++)
        for(int x = 0; x < image.cols; x++)
            new_image.at<uchar>(y,x) = saturate_cast<uchar>(Sk[image.at<uchar>(y,x)]);

   // Display the original Image
    namedWindow("Original Image");
    imshow("Original Image", image);

    // Display the original Histogram
    histDisplay(histogram, "Original Histogram");

    // Display equilized image
    namedWindow("Equilized Image");
    imshow("Equilized Image",new_image);

    // Display the equilzed histogram
    histDisplay(final, "Equilized Histogram");

    waitKey();
    return 0;
}
 
Note: OpenCV is used for read and display image only.

Output
Original Image and Histogram



Equalized Image and Histogram
















10 comments:

  1. i have some problems with visual studio 2013 , please help , anybody ?

    ReplyDelete
    Replies
    1. i have a question. In my case equalize image is only about one third of the mine, and the remaining two thirds of the original image appears I don't know what is wrong.

      Delete
  2. This comment has been removed by the author.

    ReplyDelete
  3. actually you use opencv too in:
    histogram[(int)image.at(y,x)]++;

    ReplyDelete
  4. Is there a mistak at the line 126
    I think it should be:
    final[i] = cvRound(PsSk[i] * 255.0 * 255.0);

    ReplyDelete
  5. can i say this is the answer for, implementation of algorithm of histogram equalization which works on grayscale image and single channel of color image.
    if not,i would like to get an answer for this please

    ReplyDelete
  6. can i call this as , the algorithm of histogram equalization ,which is working on grayscale image and single channel of colour image?
    if not,I need a coding for this particular thing

    ReplyDelete