OpenCV图像处理专栏二十 | 复现2003年一篇古老的低照度图像增强论文

1. 前言

今天分享一篇2003年的低照度图像增强论文《Adaptive Logarithmic Mapping For Displaying High Contrast Scenes》，论文地址为：https://domino.mpi-inf.mpg.de/intranet/ag4/ag4publ.nsf/0/53A4B81D590A3EEAC1256CFD003CE441/$file/logmap.pdf。

2. 原理

首先论文第一个重要的公式是：

其中代表每个像素的显示亮度，代表原图像亮度，是原图像亮度的最大值。

第二个重要的公式是：

这就是高中学过的换底公式了。

第三个重要公式是：

这个公式的来源为：经过实验表明人眼对亮度的适应比较符合对数曲线，为了使得对数变换变得"平滑"，使用了上述的变换。这个变换就是将一个数值做一个指数变换，来达到调节的目的，当时，即是，当取时，得到的调整函数最接近的伽马矫正结果，论文还尝试，，，的不同恢复结果，最后在代码实现部分选择了，这个值看起来是最优秀的。

第4个重要的公式为：

其中是设定的一个比例因子，根据不同的显示器需要进行调整，CRT显示器可以取，这部分的细节和的取值建议去看原论文，最后作者在实现过程中使用的是将原始的rgb转成ciexyz色彩空间，然后在ciexyz空间进行变换后再转回rgb颜色空间。

按照这几个公式以及自己的理解去复现了一把，和论文展示的结果是比较类似的，应该不会有大问题。

3. 代码复现

代码语言：javascript复制

double Transform(double x)
{
	if (x <= 0.05)return x * 2.64;
	return 1.099*pow(x, 0.9 / 2.2) - 0.099;
}
struct zxy {
	double x, y, z;
}s[2500][2500];

int work(cv::Mat input_img, cv::Mat out_img) {
	int rows = input_img.rows;
	int cols = input_img.cols;
	double r, g, b;
	double lwmax = -1.0, base = 0.75;
	for (int i = 0; i < rows; i  ) {
		for (int j = 0; j < cols; j  ) {
			b = (double)input_img.at<Vec3b>(i, j)[0] / 255.0;
			g = (double)input_img.at<Vec3b>(i, j)[1] / 255.0;
			r = (double)input_img.at<Vec3b>(i, j)[2] / 255.0;
			s[i][j].x = (0.4124*r   0.3576*g   0.1805*b);
			s[i][j].y = (0.2126*r   0.7152*g   0.0722*b);
			s[i][j].z = (0.0193*r   0.1192*g   0.9505*b);
			lwmax = max(lwmax, s[i][j].y);
		}
	}
	for (int i = 0; i < rows; i  ) {
		for (int j = 0; j < cols; j  ) {
			double xx = s[i][j].x / (s[i][j].x   s[i][j].y   s[i][j].z);
			double yy = s[i][j].y / (s[i][j].x   s[i][j].y   s[i][j].z);
			double tp = s[i][j].y;
			//修改CIE:X,Y,Z
			s[i][j].y = 1.0 * log(s[i][j].y   1) / log(2   8.0*pow((s[i][j].y / lwmax), log(base) / log(0.5))) / log10(lwmax   1);
			double x = s[i][j].y / yy*xx;
			double y = s[i][j].y;
			double z = s[i][j].y / yy*(1 - xx - yy);

			//转化为用RGB表示
			r = 3.2410*x - 1.5374*y - 0.4986*z;
			g = -0.9692*x   1.8760*y   0.0416*z;
			b = 0.0556*x - 0.2040*y   1.0570*z;

			if (r < 0)r = 0; if (r>1)r = 1;
			if (g < 0)g = 0; if (g>1)g = 1;
			if (b < 0)b = 0; if (b>1)b = 1;

			//修正补偿
			r = Transform(r), g = Transform(g), b = Transform(b);
			out_img.at<Vec3b>(i, j)[0] = int(b * 255);
			out_img.at<Vec3b>(i, j)[1] = int(g * 255);
			out_img.at<Vec3b>(i, j)[2] = int(r * 255);
		}
	}
	return 0;
}