曹鹏辉,吕书强,侯妙乐,赵林毅,汪万福
(1. 北京建筑大学测绘与城市空间信息学院,北京 100044; 2. 北京市建筑遗产精细重构与健康监测重点实验室,北京 100044;3. 敦煌研究院保护研究所,甘肃酒泉 736200; 4. 国家古代壁画与土遗址保护工程技术研究中心,甘肃酒泉 736200)
中国古代壁画存世数量巨大、色彩瑰丽、内容丰富,具有极高的艺术和研究价值,是中国文化遗产的重要组成部分。然而,由于存世年代久远,自然环境以及人为等因素导致壁画出现各种病害。由于壁画的不可再生性,所以如何对壁画进行无损留存与保护成为现在研究的热点之一。而壁画的线状特征决定着画面的整体空间的划分、画面主题定位及物体间的位置关系等各元素。它们构成了每一幅画面的雏形,也是将画面的生命力和感染力倾注在画面上的奠基者[1],在壁画的保护与修复中占据着重要地位。
近年来高光谱技术快速发展[2-4],因其分析手段无损、波段数目多、波段覆盖范围广、波段分布连续、光谱分辨率高,以及“图谱合一”等特点,可以更好地表达彩绘文物的完整信息,所以被广泛应用在彩绘文物的数字化保护中。例如在信息挖掘方面,郭新蕾等[5]利用主成分分析技术提取了古画头冠周围涂抹痕迹。Herens等[6]利用高光谱技术对凯斯·凡·东根的一幅画作进行颜料分析,发现了一个隐藏起来的女性肖像。吴太夏等[7]通过光谱分析,成功地提取了中国黑墨水的特征,在墓葬的两个重叠木片的隐藏面上发现了60多个汉字。史宁昌等[8]利用短波红外波段,对部分故宫馆藏书画文物进行分析,在增强印章的文字信息、发现书画的涂改涂抹痕迹和提取底稿线等方面得到了一定进展。在颜料识别方面,王功明等[9]采用光谱表示颜料的颜色,然后借助独立成分分析,设计了一种分离混合颜料的方法,采用蒙赛尔色卡的光谱信息进行模拟实验,恢复出原始颜料的种类及比例。Rohani等[10]利用高光谱稀疏模型对埃及发掘的一幅古画进行了颜料解混,并识别了不同颜料的种类。目前,常用于增强线状特征方法有通过对高光谱数据进行主成分分析(Principle Component Analysis,PCA)或最小噪声分离(Minimum Noise Fraction,MNF)变换后找寻线状特征清晰的特征波段[11-12],利用模糊聚类算法对数字影像的模糊边缘进行增强[13],利用Laplacian算子进行边缘增强[14]等方法。但是,利用MNF及PCA方法在选择特征波段的时候需要人为干预,影响增强结果。其他方法增强的对象限于真彩色数字影像,其波段有限,无法完全留存壁画信息,增强结果会有一定的信息缺失。
小波变换也被频繁应用于图像增强、图像融合、人脸识别等领域[15-19],但小波变换处理的图像大多为可见光范围内的数字影像,无法有效利用近红外波段线状特征信息丰富的优点。
利用高光谱影像覆盖波长宽、信息丰富的优点,结合MNF变换与Haar小波变换,提出了一种壁画线状特征的增强方法,以青海省瞿昙寺的壁画局部为研究实例进行了线状特征增强,并与PCA线状特征增强进行了效果对比。
瞿昙寺位于青海省海东市乐都区南21 km处,始建于明洪武二十五年(1392年),1982年被评为第二批全国重点文物保护单位。该寺现存明清两代壁画1 338 m2[20],为明清两代宫廷画师所作,是瞿昙寺的艺术瑰宝,具有极高的艺术价值。目前,由于存世久远,部分壁画存在明显的空鼓、起甲、颜料层脱落、褪色等病害,对其保护修复亟待进行。因此,利用高光谱MNF变换与Haar小波变换结合,选择瞿昙寺的东西回廊部分壁画进行线状特征的增强处理,期待为壁画的保护修复提供更丰富更直观的参考信息。
1.2.1数据获取 所用仪器为VNIR400H型高光谱成像仪,光谱范围为400~1 000 nm,光谱分辨率为2.8 nm,相机为内置扫描方式,波段数为1 040个,成像画幅为1 392×1 000个像元。拍摄时间为2018年8月28日下午,现场无直射自然光,环境较暗,采用人工光源。仪器距壁画109 cm,光圈为4.0,曝光时间为80 ms。主要实验区域的正射影像如图1所示。红色区域为研究区域。
图1 西回廊十五区正射影像图Fig.1 Orthophoto of the 15th district of west corridor
1.2.2预处理 高光谱成像仪采集的原始数据是目标物反射的辐射亮度,其中包含壁画本身的辐射信息和噪声。另外考虑到地面高光谱成像的物距只有1 m左右,大气辐射传输对辐射亮度的影响可以忽略。因此,只需要对高光谱原始数据进行反射率重建,如式(1)。
(1)
式中,Rref为反射率重建后数据;Rdate为高光谱原始影像数据;Rwhite为同等环境下白板原始影像数据;Rdark为暗电流噪声数据。
白板使用反射率为99%的标准白板,为减少灯光和仪器拍摄角度对白板数据的影响,采集白板的环境与采集壁画时的环境和仪器参数保持一致。暗电流数据在关闭光源,盖上镜头盖时进行数据获取,反映了仪器在没有辐射能量输入时仪器的噪声数据。
针对壁画线状特征经常会褪色而难以辨认的问题,综合利用MNF变换能分离噪声、小波变换能分离低频和高频信号的特点,设计了壁画线状特征增强的流程(图2)。对反射率重建后的高光谱影像进行MNF变换,计算MNF变换后的前m个波段的平均梯度(Average Gradient,AG),选择AG最高的波段作为最优波段,该波段线状特征最丰富。同时对MNF逆变换的影像进行RGB组合并转化为灰度图像。利用Haar小波对最优波段和MNF逆变换合成的灰度图像进行小波分解,将两者获得的低频信号相融合,并与MNF逆变换的灰度影像分解的高频信号进行Haar小波逆变换,得到线状特征明显的增强影像。
图2 线状特征增强流程图Fig.2 Flow chart of the algorithm proposed
MNF变换[21]是一种高光谱数据常用的降维方法,本质是进行两次主成分变换。第一次变换分离出原始数据中的噪声,如式(2),使得变换后的噪声数据具有最小的方差且各波段不相关。第二次变换是对噪声数据进行标准主成分变换,如式(3)。
X=S+N
(2)
式中,X为反射率重建后数据;S为无噪声数据;N为噪声数据。
Y=U×X
(3)
式中,Y为经过MNF变换后数据;U为N的协方差的逆矩阵与X的协方差矩阵相乘的特征向量。
原始高光谱数据通过MNF变换以后,获得各波段不相关且信噪比由高到低排列的数据。
平均梯度的是描述图像在垂直和水平方向上灰度的平均变化率[22]。由于线状特征的灰度值与其他物质的灰度值差距较大,所以若影像中线状特征越丰富、清晰,该影像的平均梯度就越大。利用这种方法对降维后的高光谱数据进行最优波段选择,选择出线状特征丰富并且清晰的波段,如式(4)。
Haar小波是小波分析中具有紧支撑的正交小波函数。Haar小波函数定义如式(5),其尺度函数如式(6)。
(5)
(6)
利用Haar小波对图像进行一级分解可以得到一个低频信号(LL)和3个不同方向的高频信号,分别为水平方向分量(HL),竖直方向分量(LH)和对角线分量(HH)。
图像分解后,主要信息集中在低频信号上,而噪声和细节信息集中在高频信号,所以只需要对低频信号进行增强,对高频信号进行弱化即可以达到图像增强的目的。
使用经过反射率重建的数据进行MNF变换,选择信噪比较高的前m个波段,分别利用式(4)计算各波段的平均梯度,选择平均梯度值最大的波段为最优波段。利用前m个波段信噪比高的特点进行MNF逆变换获得低噪声的高光谱数据。对其进行RGB组合,利用式(7)将RGB影像转化为灰度影像。
I=R×0.299+G×0.587+B×0.114
(7)
式中,I为合成灰度影像;R为红波段影像;G为绿波段影像;B为蓝波段影像。
由于最优波段影像与MNF逆变换后合成的灰度影像的灰度区间不同,所以对最优波段影像和MNF逆变换后合成的灰度影像进行归一化,保证其灰度区间在[0,1]内,如式(8)。
(8)
式中,sn为归一化后的影像;s为待归一化影像;smin为待归一化影像s中灰度的最小值;smax为待归一化影像s中灰度的最大值。
利用Haar小波对进行归一化处理后的灰度影像和最优波段影像进行分解。虽然最优波段影像的低频信号可以清晰地反映出绝大多数线状特征信息并且没有噪声的影响,但单一波段无法完全保留壁画的线状特征,需要与MNF逆变换的灰度影像的低频信号相融合,从而确保线状特征的完整度。MNF逆变换影像比原始影像的噪声更低。在经过小波分解后,高频信号中的噪声信息降低,突出了图像细节信息。这不但达到削弱高频信号的目的,也保留了图像的细节信息。利用增强后的低频信号与削弱后的高频信号后进行重构,获得增强后影像。在低频信号融合的过程中,应按式(9)对两幅影像的低频信号进行加权。具体流程图如图3。
图3 融合示意图Fig.3 Fusion flow chart
L=λL1+(1-λ)L2
(9)
式中,L1、L2分别为最优波段影像与灰度影像利用Haar小波分解的低频部分;L为低频部分融合结果;λ为权重系数。
以瞿昙寺西回廊15区的部分壁画为例,进行实验并分析结果。在ENVI 5.1上对原始影像进行MNF变换、MNF逆变换、RGB组合与灰度图像转换。其余所有处理均使用Matlab 2013 a实现。
3.1.1MNF变换与MNF逆变换 对辐射校正后的影像进行MNF变换获得信噪比由高到低排列的各波段不相关的数据,其特征值随波段增加而减小,选择特征值高的前10个波段作为MNF变换的结果。利用前10波段进行MNF逆变换,去除高光谱图像中的噪声。
3.1.2基于平均梯度的最优波段选择 对MNF变换后的10个波段进行归一化后分别计算平均梯度,结果如图4。MNF变换后的第二波段的平均梯度明显高于其他波段,利用目视法进行检验,证明该波段为线状特征最为丰富、清晰的波段,所以选择该波段为最优波段。
图4 MNF变换前10波段的平均梯度Fig.4 Average gradient of the first 10 bands
3.1.3基于Haar小波变换的线状特征增强 选择MNF逆变换后的band 383(640 nm)作为红光波段,band 241(550 nm)为绿色波段,band 94(460 nm)为蓝色波段进行RGB组合,并利用式(7)将其转换为灰度影像(图5a)。
利用Haar小波对归一化后最优波段和MNF逆变换后的灰度图像进行分解、重构,得到线状特征的增强影像(图5b)。
图5 结果对比Fig.5 Results of comparison
为同时利用MNF逆变换后合成的灰度影像与最优波段影像的低频信号,在低频信号融合时权重均为0.5。
如图5,对比增强后灰度影像与合成灰度影像发现区域1和区域2中的线状特征明显增强,区域3中被颜料覆盖的线状特征通过本研究方法增强后也清晰可见。将数据归一化后,截取区域1~3利用Sobel算子进行边缘提取后,计算其平均梯度,结果如表1。可见增强后影像平均梯度均大于原始灰度影像,进一步证明线状特征获得增强。
表1 平均梯度对比结果Table 1 Results of average gradient
为了证明该方法的普适性,选择瞿昙寺东回廊三区和西回廊十五区的另一部分壁画进行实验,实验结果如图6~7所示。
图6 东回廊三区实验结果Fig.6 Experimental results of the third district of east corridor
图7 西回廊十五区另一部分实验结果Fig.7 Experimental results of another area in the 15th district of west corridor
通过对比真彩色影像与增强后灰度影像发现东回廊三区人物的腰带和衣服细节等区域的线状特征明显增强。西回廊十五区的另一部分除增强线状特征外,还使被颜料覆盖的文字信息清晰可见。证明该方法在壁画的线状特征增强和隐含信息的挖掘的方面具有良好的普适性。
文献[11]中,提出基于PCA变换的线状特征增强方法。利用该方法对图1所示区域高光谱数据进行PCA,选择最优波段,然后与图5中所得结果进行对比,如图8。另外,计算局部区域的平均梯度与结构相似性(Structural similarity index measurement,SSIM)进行对比。对数据进行归一化后截取区域1~3,分别计算基于PCA增强、MNF-Haar小波结合增强两种线状特征增强方法与降噪后合成灰度影像的结构相似性,结果如表2。利用Sobel算子对归一化后的区域1~3进行边缘检测后,分别计算两种方法不同区域的平均梯度,结果如表2。区域1、区域2由于裂缝的影响,使用文献[11]中方法增强后裂缝过于明显,所以导致计算平均梯度值高于MNF-Haar小波结合增强方法。对于无裂缝的区域3,MNF-Haar小波结合增强方法的平均梯度较好,证明该方法可以更好地实现线状特征的增强,且增强后图像与原图具有较高的结构相似性。
图8 方法对比Fig.8 Methods contrast
表2 不同方法的结构相似性与平均梯度Table 2 Structural similarity index measurement and average gradient of different methods
利用高光谱MNF变换与Haar小波结合,提出了一种壁画线状特征的增强方法。首先,对高光谱影像进行MNF变换,选择信息量集中的前几个波段进行MNF逆变换实现降噪处理。其次,对重构后的影像选择真彩色波段变换为灰度图像,对其进行Haar小波分解。然后,对MNF变换后的影像,利用最大平均梯度法进行最优波段选择,挑选出线状特征信息丰富的波段进行Haar小波分解,利用分解后的低频信号与灰度图像分解后的低频信号相融合,高频信号使用重构后的降噪图像代替。最后对优化组合的低频和高频信号进行Haar小波逆变换得到结果图像。以青海省瞿昙寺壁画局部为研究实例,利用提出的方法进行了线状特征增强,经过与原始灰度影像、主成分分析线状特征增强方法进行对比,结果表明MNF变换与Haar小波结合的壁画线状特征增强方法具有较好的效果,验证了提出方法的有效性。研究表明,利用高光谱与小波变换相结合的方法可以有效增强壁画模糊的线状特征,实现隐含信息提取,能够为壁画修复提供线状特征参考,使得修复有据可依。