直方图均衡化及OpenCV实现

一、直方图均衡化原理：

出自Sonka “Image processing, analysis and machine vision”, 3rd ed, 2008.

直方图均衡化是想把左边的直方图变换为右边的直方图，使各个灰度值上的像素个数均匀分布，这样图像整体的亮度分布较均匀，不会过亮或过暗，同时可以增加图像的对比度。

假设变换前后图像的灰度等级都为k，则变换过程为：

G(q i) k

i=0=H(p i)

k

i=0

其中，p i为变换前的灰度值，q i为变换后的灰度值，H(p)为变换前的直方图，G(p)为变换后的直方图。假设图像大小为N*N，那么G q i应为常量值：N2

q k?q0

那么用精确的连续函数表达，有：

其中，称为累积直方图(cumulative histogram).

转换为离散形式，得到：

二、算法实现

1) 声明数组H，数组长度为图像的灰度等级G(一般为256)，图像大小为N*M;

2) 扫描图像，建立灰度直方图： H(pi) = H(pi) + 1;

3) 建立累积直方图Hc(pi)：

Hc(p0) = H(p0);

Hc(pi) = Hc(pi-1) + H(pi);

4) 扫描图像，计算新像素值qi，也可以建立一个查找表：

qi=(G-1)Hc(pi)/(N*M)

三、代码实现：

Mat& MyHistEqual(Mat& I)

{

CV_Assert(I.data);

// accept only char type matrices

CV_Assert(I.depth() != sizeof(uchar));

const int channels = I.channels();

switch(channels)

{

case 1:

{

int Hist[256];

memset(Hist, 0, sizeof(Hist));

// form image histogramp

MatIterator_ it, end;

for( it = I.begin(), end = I.end(); it != end; it++ )

Hist[(*it)] ++;

// form cumulative histograme

for( int i = 1; i < 256; i++ )

Hist[i] += Hist[i-1];

// histogram equalization

int nRows = I.rows;

int nCols = I.cols;

for( it = I.begin(), end = I.end(); it != end; it++ )

*it = (uchar)( 255 * Hist[(*it)] / (nRows * nCols) );

break;

}

case 3:

{

MatIterator_ it, end;

for( it = I.begin(), end = I.end(); it != end; it++ )

{

(*it)[0] = 255 - (*it)[0];

(*it)[1] = 255 - (*it)[1];

(*it)[2] = 255 - (*it)[2];

}

}

}

return I;

}

四、OpenCV2.4.3中的源代码：

CV_IMPL void cvEqualizeHist( const CvArr* srcarr, CvArr* dstarr )

{

CvMat sstub, *src = cvGetMat(srcarr, &sstub);

CvMat dstub, *dst = cvGetMat(dstarr, &dstub);

CV_Assert( CV_ARE_SIZES_EQ(src, dst) && CV_ARE_TYPES_EQ(src, dst) && CV_MAT_TYPE(src->type) == CV_8UC1 );

CvSize size = cvGetMatSize(src);

if( CV_IS_MAT_CONT(src->type & dst->type) )

{

size.width *= size.height;

size.height = 1;

}

int x, y;

const int hist_sz = 256;

int hist[hist_sz];

memset(hist, 0, sizeof(hist));

for( y = 0; y < size.height; y++ )

{

const uchar* sptr = src->data.ptr + src->step*y;

for( x = 0; x < size.width; x++ )

hist[sptr[x]]++;

}

float scale = 255.f/(size.width*size.height);

int sum = 0;

uchar lut[hist_sz+1];

for( int i = 0; i < hist_sz; i++ )

{

sum += hist[i];

int val = cvRound(sum*scale);

lut[i] = CV_CAST_8U(val);

}

lut[0] = 0;

for( y = 0; y < size.height; y++ )

{

const uchar* sptr = src->data.ptr + src->step*y; uchar* dptr = dst->data.ptr + dst->step*y;

for( x = 0; x < size.width; x++ )

dptr[x] = lut[sptr[x]];

}

}