曹鵬輝，呂書(shū)強(qiáng)，侯妙樂(lè)，趙林毅，汪萬(wàn)福

(1. 北京建筑大學(xué)測(cè)繪與城市空間信息學(xué)院，北京 100044； 2. 北京市建筑遺產(chǎn)精細(xì)重構(gòu)與健康監(jiān)測(cè)重點(diǎn)實(shí)驗(yàn)室，北京 100044；3. 敦煌研究院保護(hù)研究所，甘肅酒泉 736200； 4. 國(guó)家古代壁畫(huà)與土遺址保護(hù)工程技術(shù)研究中心，甘肅酒泉 736200)

0 引 言

中國(guó)古代壁畫(huà)存世數(shù)量巨大、色彩瑰麗、內(nèi)容豐富，具有極高的藝術(shù)和研究?jī)r(jià)值，是中國(guó)文化遺產(chǎn)的重要組成部分。然而，由于存世年代久遠(yuǎn)，自然環(huán)境以及人為等因素導(dǎo)致壁畫(huà)出現(xiàn)各種病害。由于壁畫(huà)的不可再生性，所以如何對(duì)壁畫(huà)進(jìn)行無(wú)損留存與保護(hù)成為現(xiàn)在研究的熱點(diǎn)之一。而壁畫(huà)的線(xiàn)狀特征決定著畫(huà)面的整體空間的劃分、畫(huà)面主題定位及物體間的位置關(guān)系等各元素。它們構(gòu)成了每一幅畫(huà)面的雛形，也是將畫(huà)面的生命力和感染力傾注在畫(huà)面上的奠基者[1]，在壁畫(huà)的保護(hù)與修復(fù)中占據(jù)著重要地位。

近年來(lái)高光譜技術(shù)快速發(fā)展[2-4]，因其分析手段無(wú)損、波段數(shù)目多、波段覆蓋范圍廣、波段分布連續(xù)、光譜分辨率高，以及“圖譜合一”等特點(diǎn)，可以更好地表達(dá)彩繪文物的完整信息，所以被廣泛應(yīng)用在彩繪文物的數(shù)字化保護(hù)中。例如在信息挖掘方面，郭新蕾等[5]利用主成分分析技術(shù)提取了古畫(huà)頭冠周?chē)磕ê圹E。Herens等[6]利用高光譜技術(shù)對(duì)凱斯·凡·東根的一幅畫(huà)作進(jìn)行顏料分析，發(fā)現(xiàn)了一個(gè)隱藏起來(lái)的女性肖像。吳太夏等[7]通過(guò)光譜分析，成功地提取了中國(guó)黑墨水的特征，在墓葬的兩個(gè)重疊木片的隱藏面上發(fā)現(xiàn)了60多個(gè)漢字。史寧昌等[8]利用短波紅外波段，對(duì)部分故宮館藏書(shū)畫(huà)文物進(jìn)行分析，在增強(qiáng)印章的文字信息、發(fā)現(xiàn)書(shū)畫(huà)的涂改涂抹痕跡和提取底稿線(xiàn)等方面得到了一定進(jìn)展。在顏料識(shí)別方面，王功明等[9]采用光譜表示顏料的顏色，然后借助獨(dú)立成分分析，設(shè)計(jì)了一種分離混合顏料的方法，采用蒙賽爾色卡的光譜信息進(jìn)行模擬實(shí)驗(yàn)，恢復(fù)出原始顏料的種類(lèi)及比例。Rohani等[10]利用高光譜稀疏模型對(duì)埃及發(fā)掘的一幅古畫(huà)進(jìn)行了顏料解混，并識(shí)別了不同顏料的種類(lèi)。目前，常用于增強(qiáng)線(xiàn)狀特征方法有通過(guò)對(duì)高光譜數(shù)據(jù)進(jìn)行主成分分析(Principle Component Analysis，PCA)或最小噪聲分離(Minimum Noise Fraction，MNF)變換后找尋線(xiàn)狀特征清晰的特征波段[11-12]，利用模糊聚類(lèi)算法對(duì)數(shù)字影像的模糊邊緣進(jìn)行增強(qiáng)[13]，利用Laplacian算子進(jìn)行邊緣增強(qiáng)[14]等方法。但是，利用MNF及PCA方法在選擇特征波段的時(shí)候需要人為干預(yù)，影響增強(qiáng)結(jié)果。其他方法增強(qiáng)的對(duì)象限于真彩色數(shù)字影像，其波段有限，無(wú)法完全留存壁畫(huà)信息，增強(qiáng)結(jié)果會(huì)有一定的信息缺失。

小波變換也被頻繁應(yīng)用于圖像增強(qiáng)、圖像融合、人臉識(shí)別等領(lǐng)域[15-19]，但小波變換處理的圖像大多為可見(jiàn)光范圍內(nèi)的數(shù)字影像，無(wú)法有效利用近紅外波段線(xiàn)狀特征信息豐富的優(yōu)點(diǎn)。

利用高光譜影像覆蓋波長(zhǎng)寬、信息豐富的優(yōu)點(diǎn)，結(jié)合MNF變換與Haar小波變換，提出了一種壁畫(huà)線(xiàn)狀特征的增強(qiáng)方法，以青海省瞿曇寺的壁畫(huà)局部為研究實(shí)例進(jìn)行了線(xiàn)狀特征增強(qiáng)，并與PCA線(xiàn)狀特征增強(qiáng)進(jìn)行了效果對(duì)比。

1 研究對(duì)象與數(shù)據(jù)預(yù)處理

1.1 研究對(duì)象

瞿曇寺位于青海省海東市樂(lè)都區(qū)南21 km處，始建于明洪武二十五年(1392年)，1982年被評(píng)為第二批全國(guó)重點(diǎn)文物保護(hù)單位。該寺現(xiàn)存明清兩代壁畫(huà)1 338 m2[20]，為明清兩代宮廷畫(huà)師所作，是瞿曇寺的藝術(shù)瑰寶，具有極高的藝術(shù)價(jià)值。目前，由于存世久遠(yuǎn)，部分壁畫(huà)存在明顯的空鼓、起甲、顏料層脫落、褪色等病害，對(duì)其保護(hù)修復(fù)亟待進(jìn)行。因此，利用高光譜MNF變換與Haar小波變換結(jié)合，選擇瞿曇寺的東西回廊部分壁畫(huà)進(jìn)行線(xiàn)狀特征的增強(qiáng)處理，期待為壁畫(huà)的保護(hù)修復(fù)提供更豐富更直觀(guān)的參考信息。

1.2 數(shù)據(jù)獲取與預(yù)處理

1.2.1數(shù)據(jù)獲取 所用儀器為VNIR400H型高光譜成像儀，光譜范圍為400～1 000 nm，光譜分辨率為2.8 nm，相機(jī)為內(nèi)置掃描方式，波段數(shù)為1 040個(gè)，成像畫(huà)幅為1 392×1 000個(gè)像元。拍攝時(shí)間為2018年8月28日下午，現(xiàn)場(chǎng)無(wú)直射自然光，環(huán)境較暗，采用人工光源。儀器距壁畫(huà)109 cm，光圈為4.0，曝光時(shí)間為80 ms。主要實(shí)驗(yàn)區(qū)域的正射影像如圖1所示。紅色區(qū)域?yàn)檠芯繀^(qū)域。

圖1 西回廊十五區(qū)正射影像圖Fig.1 Orthophoto of the 15th district of west corridor

1.2.2預(yù)處理 高光譜成像儀采集的原始數(shù)據(jù)是目標(biāo)物反射的輻射亮度，其中包含壁畫(huà)本身的輻射信息和噪聲。另外考慮到地面高光譜成像的物距只有1 m左右，大氣輻射傳輸對(duì)輻射亮度的影響可以忽略。因此，只需要對(duì)高光譜原始數(shù)據(jù)進(jìn)行反射率重建，如式(1)。

(1)

式中，Rref為反射率重建后數(shù)據(jù)；Rdate為高光譜原始影像數(shù)據(jù)；Rwhite為同等環(huán)境下白板原始影像數(shù)據(jù)；Rdark為暗電流噪聲數(shù)據(jù)。

白板使用反射率為99%的標(biāo)準(zhǔn)白板，為減少燈光和儀器拍攝角度對(duì)白板數(shù)據(jù)的影響，采集白板的環(huán)境與采集壁畫(huà)時(shí)的環(huán)境和儀器參數(shù)保持一致。暗電流數(shù)據(jù)在關(guān)閉光源，蓋上鏡頭蓋時(shí)進(jìn)行數(shù)據(jù)獲取，反映了儀器在沒(méi)有輻射能量輸入時(shí)儀器的噪聲數(shù)據(jù)。

2 研究方法

針對(duì)壁畫(huà)線(xiàn)狀特征經(jīng)常會(huì)褪色而難以辨認(rèn)的問(wèn)題，綜合利用MNF變換能分離噪聲、小波變換能分離低頻和高頻信號(hào)的特點(diǎn)，設(shè)計(jì)了壁畫(huà)線(xiàn)狀特征增強(qiáng)的流程(圖2)。對(duì)反射率重建后的高光譜影像進(jìn)行MNF變換，計(jì)算MNF變換后的前m個(gè)波段的平均梯度(Average Gradient，AG)，選擇AG最高的波段作為最優(yōu)波段，該波段線(xiàn)狀特征最豐富。同時(shí)對(duì)MNF逆變換的影像進(jìn)行RGB組合并轉(zhuǎn)化為灰度圖像。利用Haar小波對(duì)最優(yōu)波段和MNF逆變換合成的灰度圖像進(jìn)行小波分解，將兩者獲得的低頻信號(hào)相融合，并與MNF逆變換的灰度影像分解的高頻信號(hào)進(jìn)行Haar小波逆變換，得到線(xiàn)狀特征明顯的增強(qiáng)影像。

圖2 線(xiàn)狀特征增強(qiáng)流程圖Fig.2 Flow chart of the algorithm proposed

2.1 MNF變換

MNF變換[21]是一種高光譜數(shù)據(jù)常用的降維方法，本質(zhì)是進(jìn)行兩次主成分變換。第一次變換分離出原始數(shù)據(jù)中的噪聲，如式(2)，使得變換后的噪聲數(shù)據(jù)具有最小的方差且各波段不相關(guān)。第二次變換是對(duì)噪聲數(shù)據(jù)進(jìn)行標(biāo)準(zhǔn)主成分變換，如式(3)。

X=S+N

(2)

式中，X為反射率重建后數(shù)據(jù)；S為無(wú)噪聲數(shù)據(jù)；N為噪聲數(shù)據(jù)。

Y=U×X

(3)

式中，Y為經(jīng)過(guò)MNF變換后數(shù)據(jù)；U為N的協(xié)方差的逆矩陣與X的協(xié)方差矩陣相乘的特征向量。

原始高光譜數(shù)據(jù)通過(guò)MNF變換以后，獲得各波段不相關(guān)且信噪比由高到低排列的數(shù)據(jù)。

2.2 平均梯度法

平均梯度的是描述圖像在垂直和水平方向上灰度的平均變化率[22]。由于線(xiàn)狀特征的灰度值與其他物質(zhì)的灰度值差距較大，所以若影像中線(xiàn)狀特征越豐富、清晰，該影像的平均梯度就越大。利用這種方法對(duì)降維后的高光譜數(shù)據(jù)進(jìn)行最優(yōu)波段選擇，選擇出線(xiàn)狀特征豐富并且清晰的波段，如式(4)。

2.3 Haar小波變換

Haar小波是小波分析中具有緊支撐的正交小波函數(shù)。Haar小波函數(shù)定義如式(5)，其尺度函數(shù)如式(6)。

(5)

(6)

利用Haar小波對(duì)圖像進(jìn)行一級(jí)分解可以得到一個(gè)低頻信號(hào)(LL)和3個(gè)不同方向的高頻信號(hào)，分別為水平方向分量(HL)，豎直方向分量(LH)和對(duì)角線(xiàn)分量(HH)。

圖像分解后，主要信息集中在低頻信號(hào)上，而噪聲和細(xì)節(jié)信息集中在高頻信號(hào)，所以只需要對(duì)低頻信號(hào)進(jìn)行增強(qiáng)，對(duì)高頻信號(hào)進(jìn)行弱化即可以達(dá)到圖像增強(qiáng)的目的。

2.4 線(xiàn)狀特征增強(qiáng)方法

使用經(jīng)過(guò)反射率重建的數(shù)據(jù)進(jìn)行MNF變換，選擇信噪比較高的前m個(gè)波段，分別利用式(4)計(jì)算各波段的平均梯度，選擇平均梯度值最大的波段為最優(yōu)波段。利用前m個(gè)波段信噪比高的特點(diǎn)進(jìn)行MNF逆變換獲得低噪聲的高光譜數(shù)據(jù)。對(duì)其進(jìn)行RGB組合，利用式(7)將RGB影像轉(zhuǎn)化為灰度影像。

I=R×0.299+G×0.587+B×0.114

(7)

式中，I為合成灰度影像；R為紅波段影像；G為綠波段影像；B為藍(lán)波段影像。

由于最優(yōu)波段影像與MNF逆變換后合成的灰度影像的灰度區(qū)間不同，所以對(duì)最優(yōu)波段影像和MNF逆變換后合成的灰度影像進(jìn)行歸一化，保證其灰度區(qū)間在[0，1]內(nèi)，如式(8)。

(8)

式中，sn為歸一化后的影像；s為待歸一化影像；smin為待歸一化影像s中灰度的最小值；smax為待歸一化影像s中灰度的最大值。

利用Haar小波對(duì)進(jìn)行歸一化處理后的灰度影像和最優(yōu)波段影像進(jìn)行分解。雖然最優(yōu)波段影像的低頻信號(hào)可以清晰地反映出絕大多數(shù)線(xiàn)狀特征信息并且沒(méi)有噪聲的影響，但單一波段無(wú)法完全保留壁畫(huà)的線(xiàn)狀特征，需要與MNF逆變換的灰度影像的低頻信號(hào)相融合，從而確保線(xiàn)狀特征的完整度。MNF逆變換影像比原始影像的噪聲更低。在經(jīng)過(guò)小波分解后，高頻信號(hào)中的噪聲信息降低，突出了圖像細(xì)節(jié)信息。這不但達(dá)到削弱高頻信號(hào)的目的，也保留了圖像的細(xì)節(jié)信息。利用增強(qiáng)后的低頻信號(hào)與削弱后的高頻信號(hào)后進(jìn)行重構(gòu)，獲得增強(qiáng)后影像。在低頻信號(hào)融合的過(guò)程中，應(yīng)按式(9)對(duì)兩幅影像的低頻信號(hào)進(jìn)行加權(quán)。具體流程圖如圖3。

圖3 融合示意圖Fig.3 Fusion flow chart

L=λL1+(1-λ)L2

(9)

式中，L1、L2分別為最優(yōu)波段影像與灰度影像利用Haar小波分解的低頻部分；L為低頻部分融合結(jié)果；λ為權(quán)重系數(shù)。

3 結(jié)果與分析

以瞿曇寺西回廊15區(qū)的部分壁畫(huà)為例，進(jìn)行實(shí)驗(yàn)并分析結(jié)果。在ENVI 5.1上對(duì)原始影像進(jìn)行MNF變換、MNF逆變換、RGB組合與灰度圖像轉(zhuǎn)換。其余所有處理均使用Matlab 2013 a實(shí)現(xiàn)。

3.1 實(shí)驗(yàn)過(guò)程

3.1.1MNF變換與MNF逆變換 對(duì)輻射校正后的影像進(jìn)行MNF變換獲得信噪比由高到低排列的各波段不相關(guān)的數(shù)據(jù)，其特征值隨波段增加而減小，選擇特征值高的前10個(gè)波段作為MNF變換的結(jié)果。利用前10波段進(jìn)行MNF逆變換，去除高光譜圖像中的噪聲。

3.1.2基于平均梯度的最優(yōu)波段選擇 對(duì)MNF變換后的10個(gè)波段進(jìn)行歸一化后分別計(jì)算平均梯度，結(jié)果如圖4。MNF變換后的第二波段的平均梯度明顯高于其他波段，利用目視法進(jìn)行檢驗(yàn)，證明該波段為線(xiàn)狀特征最為豐富、清晰的波段，所以選擇該波段為最優(yōu)波段。

圖4 MNF變換前10波段的平均梯度Fig.4 Average gradient of the first 10 bands

3.1.3基于Haar小波變換的線(xiàn)狀特征增強(qiáng) 選擇MNF逆變換后的band 383(640 nm)作為紅光波段，band 241(550 nm)為綠色波段，band 94(460 nm)為藍(lán)色波段進(jìn)行RGB組合，并利用式(7)將其轉(zhuǎn)換為灰度影像(圖5a)。

利用Haar小波對(duì)歸一化后最優(yōu)波段和MNF逆變換后的灰度圖像進(jìn)行分解、重構(gòu)，得到線(xiàn)狀特征的增強(qiáng)影像(圖5b)。

圖5 結(jié)果對(duì)比Fig.5 Results of comparison

為同時(shí)利用MNF逆變換后合成的灰度影像與最優(yōu)波段影像的低頻信號(hào)，在低頻信號(hào)融合時(shí)權(quán)重均為0.5。

3.2 結(jié)果分析

如圖5，對(duì)比增強(qiáng)后灰度影像與合成灰度影像發(fā)現(xiàn)區(qū)域1和區(qū)域2中的線(xiàn)狀特征明顯增強(qiáng)，區(qū)域3中被顏料覆蓋的線(xiàn)狀特征通過(guò)本研究方法增強(qiáng)后也清晰可見(jiàn)。將數(shù)據(jù)歸一化后，截取區(qū)域1～3利用Sobel算子進(jìn)行邊緣提取后，計(jì)算其平均梯度，結(jié)果如表1?？梢?jiàn)增強(qiáng)后影像平均梯度均大于原始灰度影像，進(jìn)一步證明線(xiàn)狀特征獲得增強(qiáng)。

表1 平均梯度對(duì)比結(jié)果Table 1 Results of average gradient

為了證明該方法的普適性，選擇瞿曇寺東回廊三區(qū)和西回廊十五區(qū)的另一部分壁畫(huà)進(jìn)行實(shí)驗(yàn)，實(shí)驗(yàn)結(jié)果如圖6～7所示。

圖6 東回廊三區(qū)實(shí)驗(yàn)結(jié)果Fig.6 Experimental results of the third district of east corridor

圖7 西回廊十五區(qū)另一部分實(shí)驗(yàn)結(jié)果Fig.7 Experimental results of another area in the 15th district of west corridor

通過(guò)對(duì)比真彩色影像與增強(qiáng)后灰度影像發(fā)現(xiàn)東回廊三區(qū)人物的腰帶和衣服細(xì)節(jié)等區(qū)域的線(xiàn)狀特征明顯增強(qiáng)。西回廊十五區(qū)的另一部分除增強(qiáng)線(xiàn)狀特征外，還使被顏料覆蓋的文字信息清晰可見(jiàn)。證明該方法在壁畫(huà)的線(xiàn)狀特征增強(qiáng)和隱含信息的挖掘的方面具有良好的普適性。

3.3 線(xiàn)狀特征增強(qiáng)對(duì)比

文獻(xiàn)[11]中，提出基于PCA變換的線(xiàn)狀特征增強(qiáng)方法。利用該方法對(duì)圖1所示區(qū)域高光譜數(shù)據(jù)進(jìn)行PCA，選擇最優(yōu)波段，然后與圖5中所得結(jié)果進(jìn)行對(duì)比，如圖8。另外，計(jì)算局部區(qū)域的平均梯度與結(jié)構(gòu)相似性(Structural similarity index measurement，SSIM)進(jìn)行對(duì)比。對(duì)數(shù)據(jù)進(jìn)行歸一化后截取區(qū)域1～3，分別計(jì)算基于PCA增強(qiáng)、MNF-Haar小波結(jié)合增強(qiáng)兩種線(xiàn)狀特征增強(qiáng)方法與降噪后合成灰度影像的結(jié)構(gòu)相似性，結(jié)果如表2。利用Sobel算子對(duì)歸一化后的區(qū)域1～3進(jìn)行邊緣檢測(cè)后，分別計(jì)算兩種方法不同區(qū)域的平均梯度，結(jié)果如表2。區(qū)域1、區(qū)域2由于裂縫的影響，使用文獻(xiàn)[11]中方法增強(qiáng)后裂縫過(guò)于明顯，所以導(dǎo)致計(jì)算平均梯度值高于MNF-Haar小波結(jié)合增強(qiáng)方法。對(duì)于無(wú)裂縫的區(qū)域3，MNF-Haar小波結(jié)合增強(qiáng)方法的平均梯度較好，證明該方法可以更好地實(shí)現(xiàn)線(xiàn)狀特征的增強(qiáng)，且增強(qiáng)后圖像與原圖具有較高的結(jié)構(gòu)相似性。

圖8 方法對(duì)比Fig.8 Methods contrast

表2 不同方法的結(jié)構(gòu)相似性與平均梯度Table 2 Structural similarity index measurement and average gradient of different methods

4 結(jié) 論

利用高光譜MNF變換與Haar小波結(jié)合，提出了一種壁畫(huà)線(xiàn)狀特征的增強(qiáng)方法。首先，對(duì)高光譜影像進(jìn)行MNF變換，選擇信息量集中的前幾個(gè)波段進(jìn)行MNF逆變換實(shí)現(xiàn)降噪處理。其次，對(duì)重構(gòu)后的影像選擇真彩色波段變換為灰度圖像，對(duì)其進(jìn)行Haar小波分解。然后，對(duì)MNF變換后的影像，利用最大平均梯度法進(jìn)行最優(yōu)波段選擇，挑選出線(xiàn)狀特征信息豐富的波段進(jìn)行Haar小波分解，利用分解后的低頻信號(hào)與灰度圖像分解后的低頻信號(hào)相融合，高頻信號(hào)使用重構(gòu)后的降噪圖像代替。最后對(duì)優(yōu)化組合的低頻和高頻信號(hào)進(jìn)行Haar小波逆變換得到結(jié)果圖像。以青海省瞿曇寺壁畫(huà)局部為研究實(shí)例，利用提出的方法進(jìn)行了線(xiàn)狀特征增強(qiáng)，經(jīng)過(guò)與原始灰度影像、主成分分析線(xiàn)狀特征增強(qiáng)方法進(jìn)行對(duì)比，結(jié)果表明MNF變換與Haar小波結(jié)合的壁畫(huà)線(xiàn)狀特征增強(qiáng)方法具有較好的效果，驗(yàn)證了提出方法的有效性。研究表明，利用高光譜與小波變換相結(jié)合的方法可以有效增強(qiáng)壁畫(huà)模糊的線(xiàn)狀特征，實(shí)現(xiàn)隱含信息提取，能夠?yàn)楸诋?huà)修復(fù)提供線(xiàn)狀特征參考，使得修復(fù)有據(jù)可依。