多尺度輪廓波分解的群稀疏壁畫修復算法

陳 永，趙夢雪，陶美風

(1.蘭州交通大學 電子與信息工程學院，甘肅 蘭州 730070；2.甘肅省人工智能與圖形圖像處理工程研究中心，甘肅 蘭州 730070)

敦煌莫高窟壁畫是享譽世界的文化瑰寶，具有極高的研究價值。然而，由于自然氣候腐蝕、人為破壞等因素導致大量壁畫出現(xiàn)空鼓、起甲、酥堿等亟待解決的病害問題。數(shù)字化修復可以克服人工物理修復風險大且不可逆的問題，是目前的研究熱點[1]。

目前，傳統(tǒng)圖像修復方法可以分為：基于序列修復[2]和基于稀疏表示修復[3]。其中，基于序列修復方法又包括擴散法[4-5]和樣本合成法[6-7]。擴散法將完好信息向破損區(qū)域動力學擴散，適用于細小區(qū)域的修復。樣本合成法則是根據(jù)塊匹配原則，在完好區(qū)域選擇相似樣本塊對破損區(qū)域進行復制填充，但易出現(xiàn)樣本塊填充錯誤等問題。基于稀疏表示修復方法則是在過完備字典上計算破損塊的稀疏編碼，從而將稠密樣本轉化為模型復雜度較低的稀疏表達形式，借助稀疏編碼和過完備字典對受損區(qū)域進行重構[8-9]。ZHANG等[10]提出了組稀疏圖像修復方法，通過局部自相似性選擇匹配塊來構建字典，但該方法未考慮全局圖像特征，易出現(xiàn)結構斷裂問題。CHEN等[11]將自適應分組稀疏表示引入到圖像修復中，但該方法采用歐式距離計算樣本塊之間的分組差異，未考慮圖像紋理和結構信息的差異性。ZHA等[12]提出樣本塊聯(lián)合組稀疏重構破損圖像的方法，但該方法建立聯(lián)合字典的過程也僅在局部區(qū)域進行字典學習，導致修復結果存在紋理模糊的問題。馮象初等[13]利用博弈理論和稀疏先驗性提出了一種圖像修復模型，但該方法在完好區(qū)域與破損區(qū)域修復博弈時，未考慮樣本塊局部自相似性特點。

綜上所述，現(xiàn)有稀疏表示在壁畫圖像修復時，未考慮壁畫結構信息與紋理信息的差異性，導致修復結果易出現(xiàn)紋理模糊和結構線條斷裂等問題。針對上述問題，筆者提出了一種基于多尺度輪廓波分解的群稀疏壁畫修復算法。主要工作有：① 通過非下采樣輪廓波變換對待修復壁畫圖像進行多尺度分解，將其分解為紋理重復低頻分量和結構輪廓高頻分量，克服了現(xiàn)有稀疏表示壁畫修復時，未考慮壁畫結構和紋理信息差異性的不足；② 對低頻分量采用提出的改進群稀疏算法，完成對壁畫紋理分量的修復；③ 采用三次立方卷積插值算法對高頻分量壁畫結構信息進行插值修復，完成壁畫結構信息的修復；④ 采用非下采樣輪廓波逆變換實現(xiàn)對修復后各尺度分量的融合重構，從而完成壁畫修復。通過對敦煌壁畫修復實驗，結果證明所提方法與對比方法相比，在主觀及客觀定量指標評價中均獲得了更好的修復效果。

1 文中算法

敦煌壁畫中不僅包含灰度值過渡平緩、較規(guī)則的紋理信息，還包含豐富的線條輪廓等像素值變化劇烈的結構信息[14]。為了充分利用壁畫紋理和結構信息，筆者提出了一種基于多尺度輪廓波分解的群稀疏壁畫修復算法，將壁畫分解為紋理重復變化緩慢的低頻分量和結構輪廓變化劇烈的高頻分量，從而將壁畫修復視點從像素域轉移至頻域，利用頻率差異[15]對輸入圖像進行多尺度分解，繼而對各分解分量分別進行修復，并最終融合重構，從而恢復圖像原始頻率分量。

1.1 基于非下采樣輪廓波變換的多尺度壁畫分解

非下采樣輪廓波變換(NonSubsampled Contourlet Transform，NSCT)使用非下采樣操作，使分解圖像具備位移不變、多尺度、多方向等特性，消除了因缺乏位移不變性出現(xiàn)在奇異點附近的偽吉布斯現(xiàn)象，其分解系數(shù)對曲線、奇異點等細節(jié)紋理信息具有更稀疏的表達[16]。首先，通過由雙通道二維濾波器組和采樣矩陣構成的非下采樣金字塔濾波器(NonSubsampled Pyramid，NSP)對尺寸大小為M×M的壁畫圖像進行a層多分辨率特征提取，將待修復壁畫圖像分解為高頻分量和低頻分量。對得到的高頻分量繼續(xù)進行方向濾波，而低頻分量則通過金字塔NSP濾波器繼續(xù)進行下一尺度迭代分解，直至達到設定分解層數(shù)。以2尺度2方向分解為例，其分解結構框架如圖1所示。

圖1 壁畫分解結構框架示意圖

與其他圖像單一分解方法不同，為了獲得豐富的結構信息，筆者對圖1中分解得到的高頻分量，繼續(xù)采用由雙通道扇形濾波器組和梅花矩陣構成的非下采樣方向濾波器組(NonSubsampled Directional Filter Bank，NSDFB)進行兩級方向濾波疊加，其四通道濾波結構圖以及濾波頻域分布如圖2所示。第1級為扇形濾波器組，首先對高頻分量進行垂直方向U0(z)和水平方向U1(z)濾波，然后分別將濾波結果送入第2級象限濾波器中，與各象限濾波方向的重疊部分即為該方向子圖的頻域分布范圍，輸出重疊區(qū)域的方向子圖yk，完成方向濾波過程。經(jīng)過濾波各個尺度的高通分量分別被分解為2bi(i=1，2，…，a)個方向子帶圖像，i為分解層數(shù)，bi為當前層高頻子帶濾波方向數(shù)，各方向頻域分布如圖2(b)所示。

(a) 四通道濾波結構圖

該濾波結構輸出等價公式可表示為

(1)

NSP和NSDFB對圖像的疊加過程中均無下采樣環(huán)節(jié)，使其獲得平移不變性，以圖3敦煌莫高窟第57窟被譽為唐代美人的“菩薩圖”的壁畫分解結果為例，可以看出，分解后的低頻分量圖3(b)中分布著大量像素變化緩慢的紋理特征；而高頻分量圖3(c)～圖3(h)則攜帶著明顯的圖像變化劇烈的結構特征，并包含了壁畫圖像的主要輪廓結構信息，其中圖3(c)和圖3(d)分別對應第一尺度扇形濾波器組水平方向和垂直方向的濾波結果，圖3(e)～圖3(h)為二尺度四個方向的結構特征，較細尺度分解子圖則依據(jù)濾波劃分頻帶提取四個方向的結構特征，將分布在同方向的奇異點合成NSCT分解系數(shù)。從圖3可以發(fā)現(xiàn)，經(jīng)過雙通道扇形濾波器組和梅花矩陣對其進行聯(lián)合作用后，可以在多個方向的高頻分量上提取到充分的結構特征，可以在重構壁畫時具有足夠完備的結構特征信息選擇。

圖3 壁畫分解結果示意圖

1.2 改進群稀疏的壁畫低頻分量修復

在采用NSCT對壁畫分解為低頻分量和高頻分量的基礎上，其低頻分量為大量重復和規(guī)則變化的紋理信息，根據(jù)其局部自相似性和全局不同的特點，采用改進群稀疏方法對其進行重構修復。組稀疏算法以圖像結構組作為稀疏表示的原子，在求解模型時首先考慮樣本塊之間的相關性，即構造相似圖像塊集合，相似集合由結構紋理相似并且尺寸一致的匹配塊矩陣構成[3，10]。但傳統(tǒng)組稀疏算法構建結構組時使用歐式距離作為衡量匹配塊相似度的基準，衡量標準單一，未考慮圖像紋理和結構信息的差異性。

針對上述不足，根據(jù)分解后壁畫低頻分量為大量重復和局部相似的特征，以及紋理信息像素之間存在強相關性的特點，在構造相似集合時，提出將結構相似性作為匹配原則，將圖像結構信息以及對比度、亮度等更多衡量因素考慮進去，充分利用先驗信息來提高字典的構建準確性。如圖4所示，對低頻分量壁畫圖像群稀疏修復時，首先將壁畫的低頻分量劃分為n個N×N大小的重疊圖像塊Xi；然后將Xi表示為Bi×1維列矩陣，形成待匹配樣本塊候選集；接著計算匹配塊之間的結構相似性，選擇匹配度最高的c個樣本塊組成集合SXi}。其結構相似性的計算公式為

(2)

再將集合{SXi}內的圖像塊橫向展開即轉化為大小為Bi×c的相似群XGi，其中XGi為第Xi個圖像塊對應的相似群。

在構造完字典后，進行字典學習與系數(shù)更新。對于已經(jīng)構建的群稀疏模型，每個樣本塊都由相似群域中的原子結合稀疏系數(shù)近似重構，通過rGi奇異值分解估計各相似群字典，即

(3)

(4)

圖4 壁畫結構組構造示意圖

圖5 低頻分量與訓練字典

其中，DGi為相似群XGi的對應字典；dGi×k為第i個群字典中的第k個原子列，k=1，2，…，n；uGi、vGi分別為UGi、VGi的列向量。通過提出的多尺度輪廓波分解變換，將大量的紋理信息保留在低頻分量上，而將結構信息分離在高頻分量上，從而避免了傳統(tǒng)組稀疏表示方法在字典構造時，未考慮結構和紋理信息差異的問題。對于分解后的低頻分量，采用式(4)構建得到的訓練字典，如圖5所示。

對每個樣本塊相似群XGi的訓練字典DGi，采用分裂伯格曼算法求解其對應的稀疏向量αGi=[αGi?1，αGi?2，…，αGi?c]，其數(shù)學定義為

(5)

(6)

1.3 壁畫高頻分量修復及融合重構

對分解后的高頻結構信息，采用三次立方卷積插值算法進行插值修復。利用鄰域已知點像素值進行卷積操作，增強其本身攜帶的點奇異特征信息，提高各尺度高頻分量的分辨率。三次立方卷積插值算法利用浮點坐標(i+u，j+v)鄰域內16個坐標像素值，通過下式求得待修復像素值：

f(i+u，j+v)=[A][B][C] ，

(7)

其中，f(i+u，j+v)為待修復坐標像素值，[A]=[S(u+1)S(u)S(u-1)S(u-2)]，[B]為目標點16個鄰點像素值構成的4×4維矩陣，[C]=[S(v+1)S(v)S(v-1)S(v-2)]T，S是對sin(xπ)/(xπ)的逼近函數(shù)。

最后將低頻分量的修復結果圖XL與各尺度各方向高頻分量的插值修復結果圖XH進行逆變換得到最終修復圖像

X=XL+XH。

(8)

2 實驗結果與分析

為了驗證所提方法的有效性，進行不同修復方法對比實驗。對比實驗均在相同的軟硬件環(huán)境下進行：MATLAB R2016a編程實現(xiàn)，硬件環(huán)境為Intel(R)Core(TM)i5-1155G7 @ 2.50 GHz，16.0 GB RAM，NVIDIA GeForce MX450。采用評價準則峰值信噪比PSNR和結構相似性SSIM分別進行人為添加破損和真實敦煌壁畫修復實驗對比及客觀性量化評價，并與曲率驅動擴散(Curvature Driven Diffusions，CDD)算法、Criminisi圖像修復算法、文獻[10]組稀疏和文獻[12]聯(lián)合樣本塊稀疏表示方法進行對比分析。

2.1 人為添加破損壁畫修復

選擇4幅完好壁畫圖像經(jīng)人為添加破損后進行修復比較，結果如圖6所示。圖6(a)為完好壁畫，圖6(b)為添加破損圖像。圖6(c)為CDD算法修復結果，因該算法采用流體擴散修復，其擴散能力有限，出現(xiàn)了修復不徹底的現(xiàn)象，如壁畫1右上菩薩頭光部分未完成修復。圖6(d)為 Criminisi算法結果，該算法采用塊匹配復制填充的方法，易出現(xiàn)填充錯誤的現(xiàn)象，如壁畫1左側矩形框內結構線條出現(xiàn)了錯誤匹配，壁畫4矩形框也出現(xiàn)了紋理紊亂的現(xiàn)象。圖6(e)為文獻[10]組稀疏修復結果，因該方法未考慮壁畫紋理和結構的差異性，壁畫中線條信息修復后存在斷裂及不連貫的問題，如圖6(e)壁畫3“狩獵·野牛圖”的身體部分出現(xiàn)線條結構斷裂。圖6(f)為文獻[12]聯(lián)合樣本塊稀疏修復結果，在壁畫1和壁畫4中存在結構斷裂和塊狀模糊問題。圖6(g)為文中算法結果，與其他算法相比，修復后圖像線條連貫更自然，獲得了更好的視覺效果。

為了進一步定量比較分析，采用PSNR和SSIM對圖6的修復結果定量比較，結果見表1。其中PSNR和SSIM分別表示圖像失真程度和結構相似性對比，其值越大，說明失真越小，修復效果越好。從表1可以看出，文中方法在 PSNR 和 SSIM 上均優(yōu)于對比算法，表明筆者所提方法修復質量更好。

表1 人為添加破損壁畫修復客觀比較

2.2 真實破損壁畫修復

下面繼續(xù)對4組真實破損壁畫圖像進行修復對比實驗，修復結果如圖7所示。其中，圖7(c)為CDD算法修復結果，對于壁畫2和壁畫4，該算法均無法完成擴散修復。圖7(d)為Criminisi算法修復結果，其中壁畫1和壁畫3矩形框內出現(xiàn)修復線條錯位以及塊匹配出錯現(xiàn)象，導致修復結果平滑度欠缺，不滿足人眼視覺效果。圖7(e)為文獻[10]組稀疏修復結果，該方法對于細小裂紋有較好的修復效果，但對于較大區(qū)域其修復結果中會出現(xiàn)結構模糊的問題，如壁畫4矩形框內存在模糊和暈染現(xiàn)象。圖7(f)為文獻[12]修復結果，該方法存在線條不連續(xù)的和塊效應問題，如壁畫3和壁畫4矩形框內出現(xiàn)了結構斷裂和紋理模糊的現(xiàn)象。圖7(g)為文中算法修復結果，因采用結構和紋理分解修復，其對真實壁畫的修復效果相較對比方法，取得了較好的視覺修復效果。

圖7 真實破損壁畫修復實驗對比圖

對于真實破損壁畫修復結果的評價，因缺少相應的標準參考圖像，故采用無參考評價方法進行評價。首先，采用無參考主觀評價指標(Mean Opinion Scores，MOS)進行評價[17]。MOS值代表主觀評分值，其值越大，表明修復效果越滿足人眼視覺效果，其與破損量對應關系見表2。為了對主觀評價進行量化分析，進一步通過計算皮爾森線性相關系數(shù)(Pearson Linear Correlation Coefficient，PLCC)進行量化評價。通過PLCC系數(shù)計算MOS值與信息熵的一致程度，該值越大，代表修復效果越好，比較結果見表3。從表3中可以看出，筆者所提方法修復效果優(yōu)于對比算法。

表2 真實破損壁畫修復主觀評價及量化表

表3 真實破損壁畫修復結果主觀評價對比

在完成真實壁畫主觀評價后，下面進一步采用無參考客觀評價指標(Sum of Modulus of gray Difference，SMD)和信息熵對圖7修復結果進行量化比較，結果見表4。其中SMD和信息熵分別代表灰度方差函數(shù)和圖像信息豐富程度，二者值越大，表明修復后信息越豐富，修復結果更清晰。從表4看出，文中方法的SMD和信息熵評價值均優(yōu)于對比算法，進一步驗證了文中方法對于真實破損壁畫修復的有效性。

表4 真實破損壁畫修復客觀比較

2.3 算法泛化性能分析

通過平均準確率對不同方法進行泛化能力對比。為了便于比較，采用自建數(shù)據(jù)集中有參照人為添加破損壁畫進行量化比較。在平均準確率度量時，使用均方誤差 (Mean Square Error，MSE)和平均峰值信噪比(Mean Peak Signal Noise Ratio，MPSNR)進行評價，二者均表示修復后壁畫圖像的失真程度，其中MSE值越小，代表修復誤差越小，MPSNR值越大，表明修復效果更好，泛化能力更強。各算法比較結果如表5 所示。可以看出，在5種方法中，文中方法均方根誤差值最小，且MPSNR值最大，說明文中方法修復后，壁畫失真最少，其泛化性能優(yōu)于其他比較算法。

表5 不同算法MSE和MPSNR比較

3 結束語

針對稀疏表示未考慮壁畫結構信息與紋理信息的差異性，導致修復結果易出現(xiàn)紋理模糊和結構線條斷裂等問題，提出了一種基于多尺度輪廓波分解的群稀疏壁畫修復算法。利用非下采樣輪廓波變換將待修復壁畫分解為紋理低頻分量和結構高頻分量。由于低頻分量包含著大量待修復塊的可匹配塊，對于低頻分量采用改進群稀疏算法進行修復，以結構相似性作為塊匹配準則，同時曲線細節(jié)特征被提取分離至高頻分量上，因此避免了稀疏算法破壞壁畫輪廓特征導致線條斷裂的現(xiàn)象，并有效克服了稀疏算法修復結果中出現(xiàn)的塊效應問題。采用三次立方卷積插值算法加強紋理特征，獲得高分辨率的輪廓線條。通過NSCT逆變換將各分量融合重構得到最終修復效果。通過對比，筆者所提方法的修復結果在主客觀評價上均優(yōu)于比較算法，更加符合視覺一致性。