采用主成分分析与梯度金字塔的高动态范围图像生成方法

张淑芳, 丁文鑫, 韩泽欣, 刘孟娅, 郭志鹏

(天津大学电气自动化与信息工程学院, 300072, 天津)

高动态范围(HDR)图像能够记录具有很高亮度动态范围的景象,能够呈现细腻的色彩信息以及清晰的亮度层次。因此,HDR图像广泛应用于摄影、游戏、医疗和计算机视觉等对图像质量要求较高的领域。

目前,获取HDR图像主要有两种方法。第一种是通过采集多曝光低动态范围(LDR)图像来融合生成HDR图像。Mann和Picard通过设置不同曝光时间来获取多曝光LDR图像,进而生成HDR图像[1];Goshtasby等人提出基于最优块的HDR图像融合算法[2],但上述算法对采集设备要求高且不适用于动态场景。第二种是基于单幅LDR图像生成HDR图像,主要分为多曝光LDR图像获取和HDR图像融合两个过程。针对多曝光LDR图像获取过程中的映射问题,Tumblin等人借助Stevens定律[3]从全局的角度压缩亮度动态范围[4],但只适用于信息简单且有条理的灰度图像;方华猛等人从主成分分析(Principal Component Analysis,PCA)和引导滤波角度进行色调映射[5],该算法对图像的色彩和边缘保持较好,但变换较为繁杂。Mertens等人根据对比度、饱和度及曝光度计算图像标量权重图,借助金字塔变换将对应分辨率的LDR图像和权重图融合得到HDR图像[6],该算法更加切合实际,但当多曝光LDR图像的亮度范围变化较大时,会出现光晕;Shen等人将局部权重、全局权重和基于恰可察觉差异的显著权重组合生成新的曝光权重图,使用曝光权重对图像的细节和基础信号进行增强,并结合拉普拉斯金字塔进行图像融合[7],融合后的图像很好地保留了颜色和纹理结构信息,但只适用于静态图像;Qin等人考虑颜色信息,能量优化过程中尺度、平移和旋转等几何信息,构造新的能量成本函数,并利用新的迭代方法,同时进行匹配和融合过程,使用随机游走算法保留输入曝光序列的颜色和纹理细节[8]。该算法可用于捕捉现实、HDR图像和视频处理。

本文研究基于单幅LDR图像融合生成HDR图像的方法。首先,利用PCA变换提取输入LDR图像的亮度图,对该亮度图进行改进的S曲线全局映射和Retinex局部色调映射处理,利用PCA逆变换获取多张曝光值不同的LDR图像;然后,利用多曝光LDR图像的对比度、饱和度与最佳曝光度和亮度信息,构造多曝光LDR图像的标量权重图;最后,利用梯度金字塔模型生成HDR图像。本文方法适用于动态场景,生成的HDR图像能够较好地展现真实场景,有较好的鲁棒性。

1　单幅LDR图像生成HDR图像方法

与直接获取多曝光LDR图像相比,基于单幅LDR图像生成多幅曝光值不同的LDR图像对采集设备要求低,且适用于动态场景。图1为该类方法的示意图。

本文借鉴了Wang等人提出的基于局部区域调整的伪多曝光融合方法[9],该方法用改进的S曲线对单幅LDR图像进行映射获取多曝光LDR图像,并采用亮度信息权重计算函数分别计算多曝光LDR图像对应的权重图,最后将多曝光LDR图像与权重图进行融合生成HDR图像。

图1　基于单幅LDR图像生成HDR图像示意图

1.1　改进的S曲线

由于S曲线具有增强高光、阴影和重要的中间亮度区域对比度的效果,符合视觉系统早期阶段的全局自适应特性,因此Wang等人将改进的S曲线作为全局映射曲线[9]。改进的S曲线亮度值为

(1)

式中:Ls,k(i,j)表示经改进的S曲线映射的第k幅LDR图像(i,j)位置处的亮度值;Pk是LDR图像的调节因子,10-Pk用于控制相邻LDR图像之间的亮度差,Pk越大,亮度差越大;La,k表示第k幅LDR图像的平均亮度值,满足La,k=1+exp(λEk),其中Ek表示第k幅LDR图像的曝光值,Ek越大,LDR图像的平均亮度越大,λ是可调节的常数,取λ=0.85;Lr(i,j)表示参考LDR图像中(i,j)位置处的亮度值,亮度范围为[0,1];η是个亮度常量,实验中取值为382.5;Lmax,k表示第k幅LDR图像的最大亮度,由参数Pk和Ek共同控制。实验表明,将Pk设为(1.6,1.3,1,0.85,0.75)、Ek设为(-1,-0.5,0,0.5,1)时生成的LDR图像效果最好。

1.2　亮度信息权重计算函数

为了增强明亮图像的暗区细节,以及黑暗图像中明亮区域的细节,Wang等人为不同曝光度的图像赋予不同的权重[9],其权重函数为

(2)

式中:wk(i,j)表示第k张图像(i,j)位置处的亮度权重值;Ls,k(i,j)表示经过改进的S曲线映射后的第k张图像(i,j)位置处的亮度值,由式(1)得到;Lmid,k是第k张图像的中间亮度值,由第k幅图像的最大和最小亮度值Lmax,k和Lmin,k求得,满足Lmid,k=(Lmax,k-Lmin,k)/2;αk表示用于调节的常数;αk和Lmid,k都是与k相关的。为了提高较暗图像中较亮区域的亮度,选择较亮区域的中间亮度值作为Lmid,k。相似地,为了降低较亮图像中较暗区域的亮度,选择较暗区域的中间亮度值作为Lmid,k。图像亮暗区域划分用下面的直方图法。

(1)多曝光图像选取。取5幅曝光度不同的LDR亮度图像组成多曝光序列,分别记为L1,…,L5。L1和L2为不同程度的欠曝光图像,L3为曝光度正常的图像,L4和L5为过曝光的图像。

(2)局部区域划分。局部区域划分是基于曝光良好的图像的直方图,即L3的直方图。首先,使用直方图的中间亮度Lmid将L3划分为两个区域。对于由Lmid分成的两部分,再使用Llow和Lhigh分别将上述两个区域再分成两部分。这样,就将整幅图像划分成4个部分,记为b1、b2、b3和b4。

(3)各区域中间亮度值。对于L3的上述4个区域,利用公式Lmid=(Lmax-Lmin)/2,求出各个区域的中间亮度值,记为Lmid,bi,其中i=1,2,3,4。

(4)分配中间亮度值。对于最暗的图像L1,将中间亮度值Lmid,1设为Lmid,b4。对于最亮的图像L5,将中间亮度值Lmid,5设为Lmid,b1。以此类推,Lmid,2=Lmid,b3,Lmid,4=Lmid,b2,Lmid,3=(Lmax,3-Lmin,3)/2。

(5)确定调整常数αk。αk用来调整第k张图像的亮度测量因子高斯曲线的斜率。当αk的值较小时,权重值变化较快;当αk的值较大时,权重值变化较慢。为了减弱最亮和最暗图像对融合图像的影响,应为最亮图像和最暗图像设定较大的αk。经过实验,αk取(1,0.1,0.1,0.1,1)时效果最好。

2　基于亮色分离和梯度金字塔的HRD图像生成方法

2.1　基于PCA变换的多曝光LDR图像获取

本文提出用PCA变换提取图像亮度成分,结合改进的S曲线全局映射和改进的Retinex色调映射算法获取多曝光LDR图像,获取过程如图2所示。

图2　基于PCA变换的多曝光LDR图像获取过程示意图

(1)PCA变换提取亮度图像。PCA变换能够完美地去除自然图像亮度和色度之间的相关性,降低由亮度变动引起的色差,进而防止光晕和泛灰现象。本文采用PCA变换算法[10],对输入LDR彩色图像I进行PCA变换,将变换后的第一主成分作为亮度图像。

(2)改进的S曲线映射。常用于调整亮度的色调映射曲线包括线性曲线、γ曲线和S曲线。其中,S曲线更符合视觉系统早期阶段的全局自适应特性。利用1.1节改进的S曲线对PCA变换后的亮度图像进行全局映射,得到(i,j)位置处的亮度值L(i,j)。

(3)基于Retinex的局部自适应滤波。经过全局处理后的图像细节及清晰度并不是很高。本文使用基于Retinex的色调映射算法,通过局部自适应地改变滤波器的形状使之与高对比度边缘匹配,从而消除光晕现象[10]。基于Retinex的局部自适应滤波后的亮度值为

Lad(i,j)=log(L(i,j))-β(i,j)log(mask(i,j))

(3)

式中:Lad(i,j)表示(i,j)位置处经过Retinex局部自适应滤波后的亮度值;L(i,j)表示位置在(i,j)处的经过改进的S曲线映射后的亮度值,为浮点型数据,范围为[0,1]。为了进行对数运算,先将L(i,j)范围扩展到[0.01,100],然后再将对数编码后的L(i,j)缩小到[0,1],计算公式如下

(4)

本文采用(i,j)位置处像素点来衡量掩模的权重,其基于sigmoid函数的计算公式为

(5)

式(5)可以让中间亮度值不受约束地变化,对于高强度像素,(i,j)处像素点的掩模会被一个接近0的值进行加权。低强度像素会被一个接近1的值进行加权。

对于每一个像素来说,滤波器是不同的,所以不能使用卷积的方法。本文通过自适应地改变滤波器的形状使明亮区域对相邻暗淡区域的影响减小,从而消除光晕。掩模采用下面公式计算

(6)

式中:mask(i,j)为(i,j)处像素点的掩模;θ表示径向角度;d表示与中心像素的距离;σθ,d的定义如下

(7)

mask(i,j)以径向方式连续地选择像素。首先,第一个像素即为(i,j)处的中心像素;然后,沿径向方向的所有像素进行累加,将空间常数为σ0的高斯函数作为加权函数,如果一个边缘沿着径向方向跨过,σθ,d会取值为相对较小的σ1,并会保持这一值直至d=dmax,其中dmax=3σ0;不断对每个方向的像素进行加权和,直到环绕完成。对于新的径向方向,σθ,d重置为初始值σ0。像素的权重会被权重和归一化,即使该像素被边缘环绕,也能保证每一个像素的掩膜是同等效果。σ0和σ1的值由图像的尺寸确定。由于篇幅限制,具体参数确定可参考文献[10]。

(4)PCA逆变换。将上一步得到的亮度图插入到输入图像I中,参考Wang方法[9],本文生成了5张曝光值不同的LDR图像,记为I1,I2,I3,I4,I5。通常设置Pk=(2.4,2.1,1.8,1.4,1.2)可以规避噪声。

2.2　权重图的获取

目前,广泛使用的融合算法将对比度、饱和度与最佳曝光度作为测量因子。在此基础上,本文将上述3个测量因子与亮度信息组合,构造多曝光LDR图像的标量权重图。

对2.1节中得到的多曝光LDR图像对应的灰度图进行拉普拉斯滤波,将(i,j)位置处的滤波器响应的绝对值作为对比度测量因子[6],记为C(i,j)。计算过程如下

C(i,j)=|imfilter(g(i,j))|

(8)

式中:imfilter表示拉普拉斯滤波,g(i,j)表示LDR图像(i,j)位置处的灰度值。

将多曝光LDR图像(i,j)位置处R、G、B这3个色彩通道的色彩值的标准差作为饱和度测量因子[6],记为S(i,j)。饱和度的计算过程如下

(9)

式中:ΔR(i,j)、ΔG(i,j)和ΔB(i,j)分别表示红、绿和蓝分量值相对于均值μ(i,j)的偏差;R(i,j)、G(i,j)和B(i,j)分别表示3个颜色通道的值;μ(i,j)表示R(i,j)、G(i,j)和B(i,j)分量的平均值。

使用高斯曲线来衡量多曝光LDR图像的R、G和B通道(i,j)位置处像素值距离0.5的远近程度(靠近0表示曝光不足,靠近1表示曝光过度)[6],然后将3个通道的结果相乘得到曝光度因子Ex(i,j),计算公式为

(10)

式中:Lexp,R(i,j)、Lexp,G(i,j)和Lexp,B(i,j)分别表示第k张LDR图像的R(i,j)、G(i,j)和B(i,j)分量经过高斯曲线处理后的曝光值;σ为标准差,根据文献[6],σ取值为0.2。

亮度测量因子的计算过程为:首先,得到多曝光LDR图像对应的亮度图像;然后,利用1.2节的亮度信息权重函数来获取第k幅图像(i,j)位置处的实际亮度信息Lad(i,j)。

由4个测量因子得到LDR图像标量权重的计算公式为

W(i,j)=

(C(i,j))ωC(S(i,j))ωS(Ex(i,j))ωE(Lad(i,j))ωL

(11)

式中:W(i,j)表示(i,j)位置处的标量权重值;ωC、ωS、ωE和ωL分别用来控制对比度测量因子C(i,j)、饱和度测量因子S(i,j)、曝光度测量因子Ex(i,j)以及亮度测量因子Lad(i,j)对标量权重图W(i,j)的影响程度。如果指数ω等于0,则没有考虑相应的度量因子。参考Mertens等人提出的算法[6],文中使用的图片都是由JPEG编码生成的,伽马校正和相机的响应曲线都是未知的。因此,本方法采用同样的加权质量度量ωC=ωS=ωE=1。亮度信息Lad(i,j)的取值范围为[0,1]。通过实验发现,ωL的取值越小,融合后图像的亮度范围越小,不能很好地实现HDR图像显示更高动态范围亮度的目标。因此,为保证融合生成的HDR图像有较高的亮度范围,本文将ωL设为1。

对式(11)进行归一化,得到第k张多曝光LDR图像在(i,j)位置处的权重值为

(12)

2.3　基于梯度金字塔的HDR图像融合过程

(13)

3　实验结果与分析

3.1　基于PCA变换的多曝光LDR图像获取实验

本实验使用的CCD数码相机型号为佳能PowerShot Pro90 IS,该款相机能够还原较高的亮度动态范围。通过调整光圈和快门,获取曝光良好的LDR图像。用基于S型曲线的色调映射算法和本文提出的基于PCA变换的多曝光LDR图像获取方法分别对同一幅LDR图像进行处理,生成多幅曝光值不同的LDR图像,如图3所示。

(a)基于S型曲线的色调映射算法[9]

(b)本文提出的色调映射算法图3　两种色调映射算法对同一幅LDR图像的处理结果对比

从图3a可以看出,基于S型曲线的色调映射处理在高亮区域容易出现光晕现象,生成的LDR图像较为模糊,如图3a第5张图像的窗户边缘模糊;低对比度区域容易出现泛灰,如图3a第1张图像,树边缘很模糊,天空颜色泛灰。本文提出的方法,利用PCA变换去除亮度和色度之间的相关性,降低由于亮度变动引起的色差;使用改进S曲线实现亮度图像动态范围的首次压缩;用改进的Retinex色调映射算法来消除光晕。从图3b可以看出,第1张图像的树边缘清晰,天空颜色也较饱和;第5张图像保留了窗户细节,且动态范围较高。因而,本方法获取的多曝光LDR图像不仅能获得较高的动态范围,而且能够有效规避光晕和泛灰。

3.2　改进的梯度金字塔融合生成HDR实验

本实验使用图4所示的曝光度良好的LDR图像,参考HDR图像是由Photomatix生成的,如图5所示,Photomatix是采用多曝光融合技术生成HDR图像的工具软件。使用的电脑为64 bit操作系统,处理器为Intel(R) Core(TM) i5-6500 CPU @ 3.20 GHz,内存为8.00 GB,显示器是分辨率为1920×1080像素的BENQ GW2270-T;使用软件为MATLAB R2015b。

(a)建筑 (b)教堂 (c)厨房

(d)窗外 (e)书房 (f)书桌图4　LDR图像

(a)建筑 (b)教堂 (c)厨房

(d)窗外 (e)书房 (f)书桌图5　Photomatix生成的HDR图像

(a)拉普拉斯算法　(b)梯度算法　(c)改进的梯度算法图6　3种基于金字塔的融合算法实验结果

用2.1节方法处理图4的LDR图像,获得多幅曝光值不同的LDR图像,然后分别用拉普拉斯金字塔融合算法、梯度金字塔融合算法及本文提出的改进梯度金字塔融合算法来生成HDR图像。实验结果如图6所示。与原LDR图像对比,上述3种算法都较好地展现了原始场景。如图6a所示,拉普拉斯金字塔融合算法的多分辨率和多尺度特性,生成的HDR图像能够有效地规避光晕,保留图像细节。图6b中,由梯度金字塔生成的HDR图像与图6a对应的HDR图像相比,更加充分地展示了图像的纹理和边缘信息。本文改进的梯度金字塔融合算法,引入亮度高斯函数,生成的HDR图像,如图6c所示,相比另外两种算法,进一步提高了图像对比度,亮暗区域更加分明且边缘细节更加清晰。因此,本文提出的改进梯度金字塔融合算法,对不同场景下的图像都能取得较理想的效果,具有较好的鲁棒性和实用性。

采用高动态范围图像可见差异预测2.2方法(High Dynamic Range Visible Difference Predictor,HDR-VDP-2.2)[12]和Deng等人提出的相对熵模型[13],对本文提出的融合方法、拉普拉斯和梯度金字塔融合方法进行客观质量评价。HDR-VDP-2.2的统计评价是基于测试HDR图像和参考HDR图像,当视觉差异大于75%时其图像像素分布概率结果如表1所示。用相对熵模型评价图像的局部对比度,平均相对熵与噪声程度比值(用符号ξ表示)越大,则该方法合成的图像质量越高[13],结果如表2所示。

表1　视觉差异大于75%时图像像素分布概率

从表1的数据来看:在视觉差异概率大于75%时,由拉普拉斯金字塔、梯度金字塔和本文提出的改进梯度金字塔算法融合生成的HDR图像与参考HDR图像的平均像素分布概率分别为16.94%、36.85%和6.1%;对于窗外和书房两个场景,由于涉及屋内和屋外环境,亮度范围差异较大,从而由单幅LDR图像生成的多幅LDR曝光图像差异较大,因此利用3种融合方法生成的HDR图像效果较好,并且由于原始LDR图像具有门窗等结构度较强的信息,所以无梯度信息的拉普拉斯方法融合效果较好,但与本文所提算法相比,效果差异很小;3种算法中,本文所提融合算法对建筑图像的合成效果最好,和梯度算法相比,其图像像素分布概率降低了55.47%,降低率为98.96%。

表2　3种算法的平均相对熵与噪声程度比值

从表2可看出,本文算法的ξ与拉普拉斯金字塔和梯度金字塔算法相比,ξ分别提高了0.542 1和0.508 9,说明本文方法合成的HDR图像信息损失更小,合成效果较好。

单次运行时间能够体现算法的复杂度。使用2.1节的方法获取多曝光LDR图像,对3种金字塔融合算法的运行时间进行测量,结果如表3所示。

表3　3种融合算法的运行时间对比

从表3中看出,对于不同尺寸的图像,本文提出的改进梯度金字塔,因为引入了亮度测量因子,导致运行时间相对较长。下一步工作应研究降低运算复杂度,规范代码,减少算法运行时间。

4　结　论

本文提出一种基于单曝光LDR图像的HDR图像生成方法。首先利用改进的S曲线全局映射算法和基于Retinex的自适应滤波算法对PCA变换后的亮度图像进行处理,并结合PCA逆变换获取多曝光LDR图像;其次利用多曝光LDR图像的对比度、饱和度、曝光度和亮度信息构造标量权重图;最后利用梯度金字塔生成HDR图像。实验结果表明,本文方法生成的多曝光LDR图像能有效地规避光晕和泛灰,融合生成的HDR图像亮暗区域分明且细节丰富。总体来讲,本文提出的方法能够较好地保留真实场景的亮度动态范围特性,具有一定的参考价值。

参考文献:

[1]MANN S, PICARD R W. On being ‘undigital’ with digital cameras: extending dynamic range by combining differently exposed pictures [C]∥Proceedings of the 1995 IS & I’s 48th Annual Conference. Springfield, VA, USA: Society for Imaging Science and Technology, 1995: 442-448.

[2]GOSHTASBY A. Fusion of multi-exposure images [J]. Image & Vision Computing, 2005, 23(6): 611-618.

[3]ZWISLOCKI J, HELLMAN R P. On the “psychophysical law” [J]. The Journal of the Acoustical Society of America, 1960, 32(7): 924-924.

[4]TUMBLIN J, RUSHMEIER H. Tone reproduction for realistic images [J]. IEEE Computer Graphics and Applications, 1993, 13(6): 42-48.

[5]方华猛, 易本顺, 赵继勇. 基于主成分分析和引导滤波的色调映射算法 [J]. 光电子·激光, 2014, 25(12): 2423-2429.

FANG Huameng, YI Benshun, ZHAO Jiyong. A tone mapping algorithm based on PCA and guided filter [J]. Journal of Optoelectronics·Laser, 2014, 25(12): 2423-2429.

[6]MERTENS T, KAUTZ J, VAN REETH F. Exposure fusion: a simple and practical alternative to high dynamic range photography [J]. Computer Graphics Forum, 2009, 28(1): 161-171.

[7]SHEN Jianbin, ZHAO Ying, YAN Shuicheng, et al. Exposure fusion using boosting Laplacian pyramid [J]. IEEE Transactions on Cybernetics, 2014, 44(9): 1579-1590.

[8]QIN Xiameng, SHEN Jianbin, MAO Xiaoyang, et al. Robust match fusion using optimization [J]. IEEE Transactions on Cybernetics, 2015, 45(8): 1549-1560.

[9]WANG T H, CHIU C W, WU W C, et al. Pseudo multiple exposure-based tone fusion with local region adjustment [J]. IEEE Transactions on Multimedia, 2015, 17(4): 470-484.

[10] MEYLAN L, SUSSTRUNK S. High dynamic range image rendering with a Retinex-based adaptive filter [J]. IEEE Transactions on Image Processing, 2006, 15(9): 2820-2830.

[11] 李建林, 俞建成, 孙胜利. 基于梯度金字塔图像融合的研究 [J]. 科学技术与工程, 2007, 7(22): 5818-5822.

LI Jianlin, YU Jiancheng, SUN Shengli. The image fusion based on gradient pyramid [J]. Science Technology and Engineering, 2007, 7(22): 5818-5822.

[12] NARWARIA M, MANTIUK R, CALLET P L. HDR-VDP-2.2: a calibrated method for objective quality prediction of high-dynamic range and standard images [J]. Journal of Electronic Imaging, 2015, 24(1): 010501.

[13] DENG G. An entropy interpretation of the logarithmic image processing model with application to contrast enhancement [J]. IEEE Transactions on Image Processing, 2009, 18(5): 1135-1140.