【VS开发】【图像处理】直方图均衡与平台直方图

直方图均衡化(Histogram Equalization

直方图均衡化处理的“中心思想”是把原始图像的灰度直方图从比较集中的某个灰度区间变成在全部灰度范围内的均匀分布。直方图均衡化就是对图像进行非线性拉伸,重新分配图像像素值,使一定灰度范围内的像素数量大致相同。直方图均衡化就是把给定图像的直方图分布改变成“均匀”分布直方图分布。

直方图均衡化的主要过程

  • 统计每一个灰度值的像素点数和所占的百分比
  • 累加每一个灰度值的百分比
  • 利用新的灰度值的百分比,重新计算每一个灰度值的变化
  • 将新的的对应的关系映射到新的图像中

一个简单的例子


首先需要说明的是,如果你说的是一道完整的题目,则这道题目没有唯一解,因为题目中没有说明原始图像的灰度级数(比如原始图像是16个灰度级的,或者是32个灰度级的,等等)。为了给你提供一个解题思路,现在人为假设原始图像是16个灰度级的,其它灰度级的解法类似。
1、图像的灰度直方图求法为:
(1)先计算图像中各个灰度级的出现频率,用h(i)表示灰度级i的出现频率,其值等于灰度级出现次数/图像像素个数:
h(0)=2/16
h(1)=1/16
h(2)=3/16
h(3)=2/16
h(4)=0/16
h(5)=1/16
h(6)=4/16
h(7)=1/16
h(8)=1/16
h(9)=1/16
h(10)=h(11)=h(12)=h(13)=h(14)=h(15)=0/16。
然后以灰度级i为横轴,出现频率h(i)为纵轴即可绘制出图像对应的直方图。
(2)图像进行直方图均衡化处理的过程为:
先计算累积分布,用r(i)表示灰度级i的累积分布:
r(0)=h(0)=2/16
r(1)=r(0)+h(1)=2/16+1/16=3/16
r(2)=r(1)+h(2)=3/16+3/16=6/16
r(3)=r(2)+h(3)=6/16+2/16=8/16
r(4)=r(3)+h(4)=8/16+0/16=8/16
r(5)=r(4)+h(5)=8/16+1/16=9/16
r(6)=r(5)+h(6)=9/16+4/16=13/16
r(7)=r(6)+h(7)=13/16+1/16=14/16
r(8)=r(7)+h(8)=14/16+1/16=15/16
r(9)=r(8)+h(9)=15/16+1/16=16/16=1
r(10)=r(11)=r(12)=r(13)=r(14)=r(15)=1
将累积分布进行量化(量化时需要用到原始图像的灰度级数,这也是为什么前面需要说明的原因),量化后的灰度级用rq(i)表示,量化公式为rq(i)=ROUND(r(i)*15),(说明:量化公式中的15等于原始图像灰度级数减1),可得:
rq(0)=ROUND(r(0)*15)=2
rq(1)=ROUND(r(1)*15)=3
rq(2)=ROUND(r(2)*15)=6
rq(3)=ROUND(r(3)*15)=8
rq(4)=ROUND(r(4)*15)=8
rq(5)=ROUND(r(5)*15)=8
rq(6)=ROUND(r(6)*15)=12
rq(7)=ROUND(r(7)*15)=13
rq(8)=ROUND(r(8)*15)=14
rq(9)=ROUND(r(9)*15)=15
rq(10)=ROUND(r(10)*15)=15
rq(11)=ROUND(r(11)*15)=15
rq(12)=ROUND(r(12)*15)=15
rq(13)=ROUND(r(13)*15)=15
rq(14)=ROUND(r(14)*15)=15
rq(15)=ROUND(r(15)*15)=15
因此,原始图像中的灰度级和均化后图像中的灰度级之间的对应关系为:
0->2
1->3
2->6
3->8
4->8
5->8
6->12
7->13
8->14
9->15
10->15
11->15
12->15
13->15
14->15
15->15
将原始图像中对应的灰度值安装上述对应关系替换成相应的灰度值,即可得到均化图像,结果如下:
3  8  13  8
6  12  2  12
14  6  12  8
15  6  12  2

关键的代码实现

[cpp] view plain copy
 在CODE上查看代码片派生到我的代码片
  1. </pre><pre name="code" class="cpp"><span style="white-space:pre"> </span>Mat src,dest;  
  2. <span style="white-space:pre">    </span>Src.copyTo(src);  
  3.   
  4.     if (src.channels() > 1)  
  5.     {  
  6.         cvtColor(src, src, CV_BGR2GRAY);  
  7.     }  
  8.   
  9.     MatND  hist;  
  10.     const int histSize = 256;  
  11.     float range[] = { 0, 255 };  
  12.     const float *ranges[] = { range };  
  13.     const int channels = 0;  
  14.     cv::calcHist(&src, 1, &channels, Mat(), hist, 1, &histSize, ranges);  
  15.     float total = src.size().width* src.size().height;  
  16.       
  17.     float bins[histSize] = { 0 };  
  18.     float binsAcc[histSize] = { 0 };  
  19.     Mat lut(1, 256, CV_8U);  
  20.     vector<float> vectorBins;  
  21.     vector<float> maxBins;  
  22.     float sumBins = 0.0;  
  23.     int countMax = 0;  
  24.     float TValue = 0;  
  25.       
  26.     // Find the mapping table  
  27.     for (int i = 0; i<histSize; i++)  
  28.     {  
  29.         float bin_val = hist.at<float>(i); // 第i灰度级上的数  
  30.         bins[i] = bin_val / total;  
  31.   
  32.         if (bins[i] > 0)  
  33.         {  
  34.             vectorBins.push_back(bins[i]);  
  35.         }  
  36.           
  37.         if (i>0)  
  38.         {  
  39.             binsAcc[i] =binsAcc[i-1] + bins[i];  
  40.         }  
  41.         else  
  42.         {  
  43.             binsAcc[0] = bins[0];  
  44.         }  
  45.         lut.at<uchar>(i) = static_cast<uchar>(cvRound(binsAcc[i] * 255));  
  46.     }  
利用openCV中均衡化的代码只要一句话:
[cpp] view plain copy
 在CODE上查看代码片派生到我的代码片
  1. Mat dest2;  
  2. equalizeHist(src, dest2);  
  3. imshow("equlization2", dest2);  



平台直方图及均衡化

平台直方图的概念

这种方法是基于直方图均衡算法的改进算法,它是通过在直方图分布中设计一个阈值T 对原来直方图进行改造。
这种方法的思路是:如果原直方图中某些灰度级对应的值大于阈值T,就将该处的值设置为T。要是该处的值不大于T时,那么该处的值保持不变。



其中 Pt(k)是平台直方图, Pr(k)是原来的直方图
和直方图均衡类似,对平台直方图进行累加计算,得到累积函数
然后通过均衡化把原图灰度值改为新的灰度值
在这种方法中T的确定是最主要的任务,它的确定方法介绍如下

平台阈值的确定




关键代码实现

[cpp] view plain copy
 在CODE上查看代码片派生到我的代码片
  1.        Mat src,dest;  
  2. Src.copyTo(src);  
  3.   
  4. if (src.channels() > 1)  
  5. {  
  6.     cvtColor(src, src, CV_BGR2GRAY);  
  7. }  
  8.   
  9.   
  10.   
  11. MatND  hist;  
  12. const int histSize = 256;  
  13. float range[] = { 0, 255 };  
  14. const float *ranges[] = { range };  
  15. const int channels = 0;  
  16. cv::calcHist(&src, 1, &channels, Mat(), hist, 1, &histSize, ranges);  
  17. float total = src.size().width* src.size().height;  
  18.   
  19. float bins[histSize] = { 0 };  
  20. float binsAcc[histSize] = { 0 };  
  21. Mat lut(1, 256, CV_8U);  
  22. vector<float> vectorBins;  
  23. vector<float> maxBins;  
  24. float sumBins = 0.0;  
  25. int countMax = 0;  
  26. float TValue = 0;  
  27.   
  28. // Find the mapping table  
  29. for (int i = 0; i<histSize; i++)  
  30. {  
  31.     float bin_val = hist.at<float>(i); // 第i灰度级上的数  
  32.     bins[i] = bin_val / total;  
  33.   
  34.     if (bins[i] > 0)  
  35.     {  
  36.         vectorBins.push_back(bins[i]);  
  37.     }  
  38.       
  39.   
  40. }  
  41.   
  42. // Calculate the Meadin value by 3 sapce  
  43. for (int i = 1; i < vectorBins.size() - 1; i++)  
  44. {  
  45.     if (vectorBins[i] < vectorBins[i - 1] && vectorBins[i - 1] < vectorBins[i + 1] || vectorBins[i] > vectorBins[i - 1] && vectorBins[i - 1]  > vectorBins[i + 1])  
  46.     {  
  47.         vectorBins[i] = vectorBins[i - 1];  
  48.     }  
  49.     else if (vectorBins[i] < vectorBins[i + 1] && vectorBins[i + 1] < vectorBins[i - 1] || vectorBins[i] > vectorBins[i + 1] && vectorBins[i + 1] > vectorBins[i - 1])  
  50.     {  
  51.         vectorBins[i] = vectorBins[i + 1];  
  52.     }  
  53. }  
  54.   
  55. // Calculate the max peak value  
  56. for (int i = 1; i < vectorBins.size() - 1; i++)  
  57. {  
  58.     if (vectorBins[i] - vectorBins[i - 1] >= 0 && vectorBins[i+1] - vectorBins[i] <= 0)  
  59.     {  
  60.         maxBins.push_back(vectorBins[i]);  
  61.         sumBins += vectorBins[i];   
  62.         countMax++;  
  63.     }  
  64. }  
  65.   
  66. TValue = sumBins / countMax;  
  67.   
  68.   
  69. // Find the mapping table  
  70. for (int i = 0; i<histSize; i++)  
  71. {  
  72.   
  73.     if (bins[i] > TValue)  
  74.     {  
  75.         bins[i] = TValue;  
  76.     }  
  77.   
  78.     if (i>0)  
  79.     {  
  80.         binsAcc[i] = binsAcc[i - 1] + bins[i];  
  81.     }  
  82.     else  
  83.     {  
  84.         binsAcc[0] = bins[0];  
  85.     }  
  86. }  
  87.   
  88. for (int i = 0; i < histSize; i++)  
  89. {  
  90.     lut.at<uchar>(i) = static_cast<uchar>(cvRound(binsAcc[i] * 255 / binsAcc[255]));  
  91. }  
  92.   
  93. LUT(src, lut, dest);  
  94. imshow("src", src);  
  95. imshow("equlization", dest);  

欢迎转载,请注明本文的出处:http://blog.csdn.net/fioletfly/article/details/51011399  谢谢!

posted @ 2016-05-10 10:56  ZhangPYi  阅读(280)  评论(0编辑  收藏  举报