OpenCV--021：直方图规定化

直方图规定化

---

就是通过一个灰度映像函数，将原灰度直方图改造成所希望的直方图。所以直方图修正的关键就是灰度映像函数。


直方图规定化是用于产生处理后有特殊直方图的图像方法。

直方图均衡化能自动增强图像的整体对比度，但是往往结果难以受到控制。实际中常常需要增强某个特定灰度值范围内的对比度或使图像灰度值的分布满足特定需求。这个时候使用直方图规定化会有较好的结果。

直方图规定化就是要调整原始图像的直方图去逼近规定的目标直方图。M为原始图的灰度级数，N为目标图的灰度级数，且M>N。

直方图规定化的基本原理

---

根据直方图均衡化理论，首先对原始图像进行直方图均衡化处理.即求变换函数

$s=T(r)={\int }_0^r{P}_r(x)dx$

现假定直方图规定化的目标图像已经实现，因此，对于目标图像也采用同样的方法进行均衡化处理，因而:

$v=G(z)={\int }_0^z{P}_z(x)dx$

上式的逆变换为

$z=G^{-1}(v)$

上式表明，可通过均衡化后的灰度级v求出目标函数的灰度级z。由于对目标图像和原始图像都进行了均衡化处理，因此具有相同的分布密度，即

$P_s(s)=P_v(v)$

因而可以用原始图像均衡化以后的灰度级s代表v，即

$z=G^{-1}(v)=G^{-1}(s)$

所以可以依据原始图像均衡化后的图像的灰度值得到目标图像的灰度级z。

直方图规定化方法

---

直方图规定化处理的一般步骤是：

1. 根据直方图均衡化原理，对原图像的直方图进行灰度均衡化处理。
2. 按照目标图像的概率密度函数$P_z(z)$,求解目标图像进行均衡化处理的变换函数$G(z)$;
3. 用原始图像均衡化中得到的灰度级s代替v,求解逆变换$z=G^{-1}(s)$.

经过上述处理得到的目标图像的灰度级将具有事先规定的概率密度$P_z(z)$。上述变换过程中所包含的两个变换函数$T（r）$和$G^{-1}(s)$可形成复合函数，即可表示为

$z=G^{-1}(s)=G^{-1}[T(r)]$

由此可知，无需进行直方图均衡化运算就可以直接实现直方图规定化处理，通过复合函数关系有效的简化了直方图规定化处理过程，求出$T(r)$和$G^{-1}(s)$之间的复合函数关系就可以直接对原始图像进行变换。

举例

---

假设一幅大小为64像素（MN=4096）的3比特图像（L=8）的灰度分布如表3.1所示，其中灰度级范围[0,L-1]=[0,7]中的整数。

表3.1 大小为64*64像素的3比特数字图像的灰度分布和直方图值

$r_k$	$n_k$	$p_r(r_k)=n_k/MN$
$r_0$	790	0.19
$r_1$	1023	0.25
$r_2$	850	0.21
$r_3$	656	0.16
$r_4$	329	0.08
$r_5$	245	0.06
$r_6$	122	0.03
$r_7$	81	0.02

假设图像的直方图如图3.19(a)所示。直方图累计函数分布:
$$s_0=T(r_0)=7\sum _{j=0}^0{p}_r(r_j)=7p_r(r_0)=1.33$$$$s_1=T(r_1)=7\sum _{j=0}^1{p}_r(r_j)=7p_r(r_0)+7p_r(r_1)=7\ast 0.19+7\ast 0.25=3.08$$$$s_2=T(r_2)=7\sum _{j=0}^2{p}_r(r_j)=7p_r(r_0)+7p_r(r_1)+7p_r(r_2)=7\ast 0.19+7\ast 0.25+7\ast 0.21=4.55$$
······
类似可得出：
$s_3=5.67,s_4=6.23,s_5=6.65,s_6=6.86,s_7=7.00$。
该变换函数的形状为阶梯状态,如图3.19(b)。


由于s是由概率和求得，所以它一直是一个小数，因此要把它们改成最接近的整数：

$s_k$	$G(z_k)$
$s_0$	1
$s_1$	3
$s_2$	5
$s_3$	6
$s_4$	6
$s_5$	7
$s_6$	7
$s_7$	7

这些是均衡后的直方图的值。很明显，只有5个不同的灰度级。因为$r_0=0$被映射为$s_0=1$,在均衡后的图像中有790个像素具有该值(见表3.1)。另外，在图像中有1023个像素取$s_1=3$这个值，有850个像素取$s_2=5$这个值。由于$r_3$和$r_4$都被映射为同一个值6，所以均衡后的图像中有（656+329）=985个像素取这个值。类似的，在均衡后的图像中有（245+122+81）=448个像素取7这个值。
均衡后直方图取值：

表3.11

$s_k$	$n_k$	$p_k(s_k)=n_k/MN$
$s_0=1$	790	0.19
$s_1=3$	1023	0.25
$s_2=5$	850	0.21
$s_{3,4}=6$	985	0.24
$s_{5,6,7}=7$	448	0.11

再次考虑上面64*64的假设图像，其直方图重复显示在图3.22(a)中。希望变换该直方图，以便使其具有表3.2中第二列规定的值。图3.22(b)显示了该直方图的大概形状。

表3.2 规定直方图和实际直方图

$z_q$	规定的$p_z(Z_q)$
$z_0=0$	0.00
$z_1=1$	0.00
$z_2=2$	0.00
$z_3=3$	0.15
$z_4=4$	0.20
$z_5=5$	0.30
$z_6=6$	0.20
$z_7=7$	0.15

计算直方图均衡后的累计分布值：

$$G(z_0)=7\sum _{j=0}^0{p}_z(z_j)=7\ast 0.00=0.00$$
$$G(z_1)=7\sum _{j=0}^1{p}_z(z_j)=7\ast 0.00+7\ast 0.00=0.00$$
$$G(z_2)=7\sum _{j=0}^2{p}_z(z_j)=0$$
$$G(z_3)=7\sum _{j=0}^3{p}_z(z_j)=7\ast 0.15=1.05$$
$$G(z_4)=7\sum _{j=0}^4{p}_z(z_j)=7\ast 0.15+7\ast 0.20=2.45$$
$$G(z_5)=7\sum _{j=0}^5{p}_z(z_j)=7\ast 0.15+7\ast 0.20+7\ast 0.3=4.55$$
$$G(z_6)=7\sum _{j=0}^6{p}_z(z_j)=4.55+7\ast 0.20=5.95$$
$$G(z_7)=7\sum _{j=0}^7{p}_z(z_j)=5.95+7\ast 0.15=7.00$$

把这些分数值转换为有效区间[0,7]内的整数。结果是：

结果形成表3.3，则可得变换函数3.22©。

表3.3

$z_q$	$G(z_q)$
$z_0$	0
$z_1$	0
$z_2$	0
$z_3$	1
$z_4$	2
$z_5$	5
$z_6$	6
$z_7$	7

第三步，找到$z_q$最小值。所以值$G(Z_q)$最接近$s_k$。我们对每一个$s_k$值这样做，以产生从s到z所需要的映射。

• $s_0=1$,我们看到$G(z_3)=1$,在这种情况下，这是最完美的匹配，因此有对应$s_0=>z_3$。也就是说，直方图均衡后的图像中的每一个值为1的像素映射为直方图规定化后的图像中（相应位置上）的值为3的像素。
• $s_1=3$,$G(z_q)$中不存在对应的3，则找最接近3的值，则为2或者4。则有$s_1=>z_4$.
• $s_2=5$,$G(z_5)=5$与之完美对应，则有$s_2=>z_5$.
• $s_3=6$,$G(z_6)=6$与之完美对应，则有$s_3=>z_6$.

由此可以得到映射关系：

表3.4

$s_k$	—>	$z_q$
1	—>	3
$3$	—>	4
$5$	—>	5
$6$	—>	6
$7$	—>	7

结合3.11和3.4中的映射把直方图均衡后的图像中的每个像素映射为新创建的直方图规定化图像中相应的像素。结果如下：

$z_q$	规定的$p_z(Z_q)$	实际的$p_z(z_k)$
$z_0=0$	0.00	0.00
$z_1=1$	0.00	0.00
$z_2=2$	0.00	0.00
$z_3=3$	0.15	0.19
$z_4=4$	0.20	0.25
$z_5=5$	0.30	0.21
$z_6=6$	0.20	0.24
$z_7=7$	0.15	0.11

代码实现

直方图规定化中要注意两点：

---

实际操作中不会进行两次均衡化。在推导中发现，假如sk 规定化后的对应灰度是zm的话，需要满足的条件是sk的累积概率和zm的累积概率是最接近的。所以可以根据计算累计密度的差值来进行映射。

	//1.传入原图和待处理的目标图
	Mat src1, src2;
	src1 = imread("D:/test/src1.jpg");
	src2 = imread("D:/test/src2.jpg");

	//2. 将原图和目标图转化为灰度图像
	Mat graySrc1, graySrc2;
	cvtColor(src1, graySrc1, COLOR_BGR2GRAY);
	cvtColor(src2, graySrc2, COLOR_RGB2GRAY);

	//3. 将图像进行均衡化处理
	equalizeHist(graySrc1, graySrc1);
	equalizeHist(graySrc2, graySrc2);

	//4. 计算图像的直方图
	/*
	calcHist( const Mat* images, int nimages,
                          const int* channels, 
						  InputArray mask, OutputArray hist, 
						  int dims, const int* histSize,
                          const float** ranges, 
						  bool uniform = true, bool accumulate = false );
	*/
	int channels = 0;
	Mat dst(src2.cols, src2.rows, CV_8UC3);
	//计算均衡化后的图像的直方图
	int dims = 1;
	int histSize = 256;
	float range[] = { 0,255 };
	const float *histRange[] = { range };
	Mat src1Hist, src2Hist;
	calcHist(&graySrc1, 1, &channels, Mat(), src1Hist, dims, &histSize, histRange);
	calcHist(&graySrc2, 1, &channels, Mat(), src2Hist, dims, &histSize, histRange);

	//计算图像的累计直方图
	float src1_cdf[256] = { 0 };
	float src2_cdf[256] = { 0 };
	for (int i = 0; i < 256; i++) {
		if (i == 0) {
			src1_cdf[i] = src1Hist.at<float>(i);
			src2_cdf[i] = src2Hist.at<float>(i);
		}
		else
		{
			src1_cdf[i] = src1_cdf[i - 1] + src1Hist.at<float>(i);
			src2_cdf[i] = src2_cdf[i - 1] + src2Hist.at<float>(i);
		}
	}

	//对图像目标进行规定化处理
	// 1. 计算累计概率的差值
	float diff_cdf[256][256];
	for (int i = 0;i < 256;i++) {
		for (int j = 0; j < 256; j++)
		{
			diff_cdf[i][j] = fabs(src1_cdf[i] - src2_cdf[j]);
		}
	}
	//2. 构建灰度级映射表
	Mat lut(1, 256, CV_8U);
	for (int i = 0; i < 256; i++)
	{
		//查找源灰度级为i的映射灰度和i的累积概率差最小的规定化灰度
		float min = diff_cdf[i][0];
		int index = 0;
		for (int j = 0; j < 256; j++) {
			if (min > diff_cdf[i][j]) {
				min = diff_cdf[i][j];
				index = j;
			}
		}
		lut.at<uchar>(i) = static_cast<uchar>(index);
	}
	LUT(graySrc2, lut, dst);

	waitKey(0);
	return 0;

学习：
OpenCV 直方图处理：直方图均衡和规定化（匹配）