直方圖均衡化原理與實現

    直方圖均衡化(Histogram Equalization) 又稱直方圖平坦化,實質上是對圖像進行非線性拉伸,重新分配圖像象元值,使一定灰度范圍內象元值的數量大致相等。這樣,原來直方圖中間的峰頂部分對比度得到增強,而兩側的谷底部分對比度降低,輸出圖像的直方圖是一個較平的分段直方圖:如果輸出數據分段值較小的話,會產生粗略分類的視覺效果。

    直方圖是表示數字圖像中每一灰度出現頻率的統計關系。直方圖能給出圖像灰度范圍、每個灰度的頻度和灰度的分布、整幅圖像的平均明暗和對比度等概貌性描述。灰度直方圖是灰度級的函數, 反映的是圖像中具有該灰度級像素的個數, 其橫坐標是灰度級r, 縱坐標是該灰度級出現的頻率( 即像素的個數) pr( r) , 整個坐標系描述的是圖像灰度級的分布情況, 由此可以看出圖像的灰度分布特性, 即若大部分像素集中在低灰度區域, 圖像呈現暗的特性; 若像素集中在高灰度區域, 圖像呈現亮的特性。

    圖1所示就是直方圖均衡化, 即將隨機分布的圖像直方圖修改成均勻分布的直方圖。基本思想是對原始圖像的像素灰度做某種映射變換, 使變換后圖像灰

度的概率密度呈均勻分布。這就意味着圖像灰度的動態范圍得到了增加, 提高了圖像的對比度。

 

圖1 直方圖均衡化

通過這種技術可以清晰地在直方圖上看到圖像亮度的分布情況, 並可按照需要對圖像亮度調整。另外,這種方法是可逆的, 如果已知均衡化函數, 就可以恢復原始直方圖。

設變量r 代表圖像中像素灰度級。對灰度級進行歸一化處理, 則0≤r≤1, 其中r= 0表示黑, r= 1表示白。對於一幅給定的圖像來說, 每個像素值在[ 0,1] 的灰度級是隨機的。用概率密度函數來表示圖像灰度級的分布。

為了有利於數字圖像處理, 引入離散形式。在離散形式下, 用 代表離散灰度級, 用 代表 , 並且下式成立:

其中, 0≤≤1, k=0, 1, 2, …, n-1。式中 為圖像中出現這種灰度的像素數, n是圖像中的像素總數, 而就是概率論中的頻數。圖像進行直方圖均衡化的函數表達式為:

式中, k為灰度級數。相應的反變換為:

四、算法實現及結果分析

4.1核心算法

#define HDIM 256
#define SRC 0
#define DST 1
int main(int argc, char** argv)
{
    IplImage *src = 0, *dst = 0;
    int n[] = {HDIM,HDIM,HDIM};
    int r[256] = {0}, g[256] = {0}, b[256] = {0};
 
    if(argc!=2 || (src = cvLoadImage(argv[1],3))== NULL)    return -1;
 
    cvNamedWindow("source",1);
    cvNamedWindow("result",1);
 
    int width = src->width;        
    int height = src->height;
    int sum = width * height;       //圖像中的像素點綜合
    int i,j;
 
    //分別統計直方圖的RGB分布
    for(i=0; i<height; i++)
        for(j=0; j<width; j++)
        {
            b[((uchar*)(src->imageData+i*src->width))[j*src->nChannels+0]]++;
            g[((uchar*)(src->imageData+i*src->width))[j*src->nChannels+1]]++;
            r[((uchar*)(src->imageData+i*src->width))[j*src->nChannels+2]]++;
        }
 
    ////構建直方圖的累計分布方程，用於對直方圖進行均衡化
    double val[3] = {0};
    for(i=0; i<HDIM; i++)
    {
        val[0] += b[i];
        val[1] += g[i];
        val[2] += r[i];
        b[i] = val[0]*255/sum;
        g[i] = val[1]*255/sum;
        r[i] = val[2]*255/sum;
    }
 
    dst = cvCreateImage(cvSize(width,height),8,3);
    //歸一化直方圖
    for(i=0; i<height; i++)
        for(j=0; j<width; j++)
        {
        ((uchar*)(dst->imageData+i*dst->widthStep))[j*dst->nChannels+0]=b[((uchar*)(src->imageData+i*src->widthStep))[j*src->nChannels+0]];
        ((uchar*)(dst->imageData+i*dst->widthStep))[j*dst->nChannels+1]=g[((uchar*)(src->imageData+i*src->widthStep))[j*src->nChannels+1]];
        ((uchar*)(dst->imageData+i*dst->widthStep))[j*dst->nChannels+2]=r[((uchar*)(src->imageData+i*src->widthStep))[j*src->nChannels+2]];
        }
    cvShowImage("source",src);
    cvShowImage("result",dst);
    cvSaveImage("out.jpg",dst);
    cvWaitKey(0);
 
    cvDestroyWindow("source");
    cvDestroyWindow("result");
    cvReleaseImage(&src);
    cvReleaseImage(&dst);
    cvReleaseHist(&hist);
 
    return 0;
}

　　

4.2結果展示

 

對於第一幅圖像，首先我們在RGB空間對其進行直方圖均衡化得如下結果：

可見右上角由較大的光影，我們對其進行同台濾波，如下圖所示，但是效果也不太好。

如下圖所示為直方圖均衡前后圖像RGB分量的直方圖，可以看到收到了較好的效果，由於直方圖均衡第二個作業已經做過，我們在此不再贅述。

對於HSI空間的直方圖均衡化，首先我們要進行RGB和HIS空間顏色分量的轉化，代碼如下：

/**************************************************************

函數功能：對圖像進行由RGB空間到HSI空間的轉化

輸入參數：源圖像src；目標圖像des；圖像參數width,height,nChannels；

輸出參數：目標圖像

**************************************************************/

void rgb_hsi(unsigned char* des, const unsigned char* src, int width, int height, int nChannels)

{

for(int y=0; y<height; y++)

{

for (int x=0; x<width; x++)

{

double B= src[y * width * nChannels + x * nChannels ] ;

double G= src[y * width * nChannels + x * nChannels + 1] ;

double R= src[y * width * nChannels + x * nChannels + 2] ;

double H,S,I=0;//H色調、S飽和度（純度）、I強度

double mx,mi;

mx=max(max(R,G),B);

mi=min(min(R,G),B);

if (mx==mi) //如果RGB相等

{

k=k+1;

H=0; //H分量為0

S=0; //S分量為0

I=mi;

}

else

{

if (B<=G)

{

H=acos((0.5*((R-B)+(R-G)))/(sqrt(1.0*((R-G)*(R-G)+(R-B)*(G-B)))));

}

else

{

H=360-acos((0.5*((R-B)+(R-G)))/(sqrt(1.0*((R-G)*(R-G)+(R-B)*(G-B)))));

}

S=(3*mi)/(R+B+G);

S=1-S;

I=(R+B+G)/3;

}

des[y * width * nChannels + x * nChannels + 0]= int(H);

des[y * width * nChannels + x * nChannels + 1]= int (S*255);

des[y * width * nChannels + x * nChannels + 2]=int(I);

}

}

}

/**************************************************************

函數功能：對圖像進行由HSI空間到RGB空間的轉化

輸入參數：源圖像src；目標圖像des；圖像參數width,height,nChannels；

輸出參數：目標圖像

**************************************************************/

void hsi_rgb(unsigned char* des, const unsigned char* src, int width, int height, int nChannels)

{

for(int y=0; y<height; y++)

{

for (int x=0; x<width; x++)

{

double H= src[y * width * nChannels + x * nChannels ] ;

//printf("H%d",H);

double S= src[y * width * nChannels + x * nChannels + 1]/255 ;

//printf("S%d",S);

double I= src[y * width * nChannels + x * nChannels + 2] ;

//printf("I%d",I);

double R,G,B;//H色調、S飽和度（純度）、I強度

if((H>=0)&&(H<120))

{

B=I*(1-S);

R=I*(1+S*cos(H)/cos(60-H));

G=3*I-(R+B);

}

else if((H>=120)&&(H<240))

{

H=H-120;

R=I*(1-S);

G=I*(1+S*cos(H)/cos(60-H));

B=3*I-(R+G);

}

else

{

H=H-240;

G=I*(1-S);

B=I*(1+S*cos(H)/cos(60-H));

R=3*I-(B+G);

}

des[y * width * nChannels + x * nChannels + 0]= int(B);

des[y * width * nChannels + x * nChannels + 1]= int(G);

des[y * width * nChannels + x * nChannels + 2]= int(R);

for(int n=0;n<nChannels;n++)

{

int val=des[y * width * nChannels + x * nChannels + n];

if(val<0)

des[y * width * nChannels + x * nChannels + n]=0;

else if(val>255)

des[y * width * nChannels + x * nChannels + n]=255;

else

des[y * width * nChannels + x * nChannels + n]=des[y * width * nChannels + x * nChannels + n];

}

}

}

}

　　

但是效果並不是很好，結果如下圖所示，可見變換后的圖像雖然比以前增強了，但是幾乎變成了灰度圖像。