冷軋帶鋼缺陷表面形狀恢復(fù)新算法

楊永敏，李 戈，趙 杰

(哈爾濱工業(yè)大學(xué)機器人技術(shù)與系統(tǒng)國家重點實驗室，150080哈爾濱)

目前，大多數(shù)表面缺陷檢測系統(tǒng)只能檢測二維信息，缺乏對三維信息的檢測.三維信息不僅能提供更加直觀的缺陷表達，同時更利于提高缺陷分類識別性能，因此對提取表面缺陷三維信息進行研究是十分必要的.三維重構(gòu)是提取三維信息的有效手段，其中從陰影恢復(fù)形狀［1］(Shape-from-Shading，SFS)的方法只需要單幅圖像就能夠完成三維重構(gòu)，該方法已在很多領(lǐng)域得到了應(yīng)用［2-3］.

SFS方法本身存在不確定性，解決該問題需從兩方面考慮:正確光照模型的建立以及有效SFS數(shù)值算法的設(shè)計.光照模型主要有漫反射模型、鏡面反射模型以及混合表面反射模型，其中朗伯體模型是理想的漫反射模型［4］，絕大多數(shù)SFS方法是基于該光照模型推導(dǎo)的，但是該模型與實際情況相差較大，因此提出了各種改進光照模型，如Oren-Nayar漫反射模型［5］，Phong鏡面反射模型［6］、Torrance-Sparrow模型［7］以及Healey-Binford模型［8］等.SFS數(shù)值經(jīng)典算法大致分為三類:優(yōu)化迭代方法［9］、局部方法［10］和傳播方法［11］等.近年來又出現(xiàn)了一些新的算法［12-14］，建立了更接近實際情況的光照模型，并引入新數(shù)學(xué)工具消除了SFS方法的不確定性.

本文根據(jù)冷軋帶鋼表面缺陷檢測系統(tǒng)的自身特點，建立了改進的Oren-Nayar漫反射型光照模型，推導(dǎo)出了透視投影下相應(yīng)的數(shù)學(xué)模型，利用高階Lax-Friedrichs漢密爾頓函數(shù)與牛頓迭代相結(jié)合的快速掃描算法完成了方程的求解，最終通過實驗驗證了SFS方法的正確性和有效性.

1光照模型的建立

1.1 冷軋帶鋼表面缺陷檢測系統(tǒng)分析

冷軋帶鋼表面缺陷檢測系統(tǒng)放置在酸洗工藝之后，此時表面已經(jīng)十分光滑潔凈，其表面粗糙度大概在1.2 μm左右，其鏡面反射率在60%左右.該檢測系統(tǒng)采用中心波長為700 nm的高亮度LED線光源，光源發(fā)出平行光在帶鋼表面形成一條照度均勻的光帶.圖像采集器采用的是敏感波長為700 nm的線陣CCD，同時光學(xué)鏡頭中安裝有690～710 nm的帶通濾光片，從而確保不受外界環(huán)境光的影響.光路設(shè)計中使得入射角等于反射角.這種照明方式下，三維表面缺陷會產(chǎn)生漫反射，且隨著表面形狀的變化反射到采集傳感器上的照度也會隨之變化.照度的變化會使采集到的圖像灰度發(fā)生變化，從而檢測出三維缺陷.由上述分析可知，缺陷圖像灰度變化是由漫反射變化引起的，因此利用SFS方法來恢復(fù)缺陷三維信息時采用漫反射光照模型.

1.2 改進光照模型

根據(jù)對冷軋帶鋼表面缺陷檢測系統(tǒng)的分析，選擇更接近實際情況的Oren-Nayar漫反射模型為原型進行改進，建立更加適合本系統(tǒng)的光照模型.

Oren-Nayar漫反射模型在建模過程中將觀察方向也考慮在內(nèi)，令光源方向為(θi，φi)，觀察方向為(θr，φr)，并定義α=max［θr，θi］，β=min［θr，θi］，則在計算精度要求不高的情況下，當(dāng)入射光強度為E0時，漫反射率為ρ的物體表面的輻照度為

式中:A=1.0-0.5σ2/(σ2+0.33)，B=0.45σ2/(σ2+0.09)，σ為高斯分布的標(biāo)準(zhǔn)偏差，其大小由表面粗糙度決定.

在冷軋帶鋼表面缺陷檢測系統(tǒng)中，光源為平行光入射，對于帶鋼表面的微面元來說，可將光源視為與其距離為r的點光源.圖像采集傳感器的光軸與光源位于同一平面之內(nèi)，光源和圖像采集傳感器之間的夾角是一個固定的較小值，可以近似認為光源為前端照明，即認為光源位于圖像采集傳感器的光學(xué)中心.此時光源方向與圖像采集傳感器的光軸方向一致，應(yīng)該有θi=θr=α=β=γ，同時考慮距離r所起到的衰減作用，則更適合本系統(tǒng)的光照模型具體表達式為

式中:η為常數(shù)，與光源的輻照度、帶鋼表面粗糙度和反射率參數(shù)、以及光學(xué)成像系統(tǒng)的各種參數(shù)有關(guān)，→L為光源方向，→N為表面法向.

2 數(shù)學(xué)模型的建立與求解

2.1 透視投影數(shù)學(xué)模型

在光照模型建立的過程中，假設(shè)圖像采集傳感器為透視投影，光源位于圖像采集傳感器的光學(xué)中心O，如圖1所示，圖中圖像采集傳感器由光學(xué)中心和成像平面組成，其焦距為f＞0.

圖1 透視投影成像模型

假設(shè)物體成像后由表面S來表達，并定義Ω為代表圖像定義域的實數(shù)集R2的一個開集，為一個矩形區(qū)域，其大小代表了圖像的尺寸，即:［0，X］×［0，Y］.用從集合Ω的閉包Ω到R3的映射來表示表面坐標(biāo)x到物體表面S(x)的映射，即表面S的形狀信息表示為

其中:u(x)為物體表面的連續(xù)函數(shù)，是未知量.對于這樣的表面，任意一點S(x)的法向量可以表示為

此點的照明方向單位向量可以表示為

用Ⅰ(x)替代Ei/η，光源到物體表面的距離為r=fu(x)，帶入到光照模型表達式中得到反射圖方程:

式中的(→L·→N)為光源方向矢量和表面法向矢量之間的內(nèi)積，將其數(shù)學(xué)表達式代入可得

式中:

將式(1)帶入反射圖方程中，并利用變量代換方法進行簡化，令v=lnu(x)，并通過整理可得

式(2)為靜態(tài)H-J方程，該方程通常情況下不存在經(jīng)典確定解.根據(jù)黏性解的概念，為了確保該方程具有唯一解，需要給定Dirichlet邊界條件［15］，得到如下的方程組:

該方程組中的H-J方程是非凸的，而Lax-Friedrichs快速掃描方法不用考慮漢密爾頓函數(shù)的凸凹性，且計算簡單，速度快，因此用該方法對方程組進行求解.

2.2 高階Lax-Friedrichs快速掃描方法

Lax-Friedrichs快速掃描策略是Kao等［16］提出的基于Lax-Friedrichs單調(diào)數(shù)值漢密爾頓函數(shù)的快速掃描方法，可用來近似任何空間維度中的任意靜態(tài)H-J方程的黏性解.高階快速掃描的本質(zhì)就是利用高階加權(quán)本質(zhì)無振蕩格式來進行偏導(dǎo)數(shù)的近似，同時結(jié)合Lax-Friedrichs單調(diào)數(shù)值漢密爾頓函數(shù)以及Guass-Seidel迭代來完成4個方向上的快速掃描.

考慮二維空間中如下的靜態(tài)H-J方程組:

其中:Ω為R中的計算域，Γ是Ω的邊界子集，

p=?φ/?x，q=?φ/?y，漢密爾頓函數(shù)H為一個非線性的Lipschitz連續(xù)函數(shù).應(yīng)用到該靜態(tài)方程組的Lax-Friedrichs單調(diào)數(shù)值漢密爾頓函數(shù)為

其中σx和σy為人工黏性因子，且σx≥max｜?H/?p｜，σy≥max｜?H/?q｜，偏導(dǎo)數(shù)p和q用高階加權(quán)本質(zhì)無振蕩格式來近似.下面以三階加權(quán)本質(zhì)無振蕩格式來進行分析，其插值模板如圖2所示.

圖2 三階加權(quán)本質(zhì)無振蕩策略模板

網(wǎng)格點(i，j)處前向和后向差分p±的三階加權(quán)本質(zhì)無振蕩近似表示為［17］

其中ε為防止分母為零的正常數(shù)，后續(xù)的計算中選擇10-5.類似可以推導(dǎo)出網(wǎng)格點(i，j)處q±的三階加權(quán)本質(zhì)無振蕩的近似表示.

將式(3)、(4)帶入Lax-Friedrichs單調(diào)數(shù)值漢密爾頓函數(shù)，經(jīng)過整理得靜態(tài)H-J方程組的數(shù)值粘性解的迭代計算公式為

其中，φ是φi，j的更新函數(shù)值為同一網(wǎng)格點上當(dāng)前的函數(shù)值.式(5)右側(cè)的偏導(dǎo)數(shù)近似值未使用上標(biāo)，其上標(biāo)由Gauss-Seidel迭代的交替掃描方向決定，確保式中的，總是使用函數(shù)φ在網(wǎng)格點上的最新當(dāng)前值來計算.

2.3 基于快速掃描的SFS算法

根據(jù)改進的光照模型推導(dǎo)出的數(shù)學(xué)模型為靜態(tài)H-J方程，為了應(yīng)用高階Lax-Friedrichs快速掃描方法進行求解，將其改寫成如下的形式:

式中:

首先對圖像區(qū)域Ω進行離散化處理，假設(shè)圖像的大小為M×N，可以將圖像區(qū)域Ω離散化為大小為M×N的矩形網(wǎng)格區(qū)域，網(wǎng)格上的點即為像素點{(xi，yj)，i=1，2，…，M;j=1，2，…，N}，均一化網(wǎng)格的大小等同于像素的大小(Δx，Δy)，實際計算中認為Δx=Δy.每個像素點的高度值是未知的，離散化的目的就是要獲取每個像素點上高度值的近似解Vi，j=v(iΔx，jΔy).

根據(jù)高階Lax-Friedrichs掃描策略，可以得到近似解的迭代更新公式如下:

對于Vi，j來說，式(6)是非線性的，需要使用牛頓迭代法對其進行求解.

令t=Vi，j，a=-1/(σx/Δx+σy/Δy)，且將式(6)的左側(cè)用b表示，則上式轉(zhuǎn)化為

函數(shù)f的一階導(dǎo)數(shù)為

令初始值t0=Vni，j，利用如下的公式來計算t的更新值

較少迭代步驟后，就能夠完成迭代公式的更新計算.

迭代更新公式中的黏性因子σx和σy按照參考文獻［18］中所提到的非線性黏性因子的方法進行確定:

綜上所述，整個算法具體實現(xiàn)步驟如下.

1)初始化.給定圖像Ⅰ(x，y)，并設(shè)置圖像中邊界點的高度值，根據(jù)邊界條件外推計算附近點的高度值，這些高度值在后續(xù)的迭代運算過程中保持不變.為其他像素點的高度值分配一個大的正值C，使其大于所有像素點對應(yīng)的高度值的最大值.

2)交替掃描.在第n+1次迭代中，除了高度數(shù)值固定的像素點外，其他所有像素點均按照迭代更新公式來計算高度值Vni，+j1.當(dāng)且僅當(dāng)?shù)玫降闹敌∮谏弦淮蔚嬎阒礦n i，j時，將像素點的高度值Vi，j更新為Vni，+j1.每個迭代過程分別從4個不同的掃描方向上進行:(1)i=1∶M，j=1∶N;(2)i=M∶1，j=1∶N;(3)i=1∶M，j=N∶1;(4)i=M∶1，j=N∶1.

3)更新邊界條件.每次掃描之后，更新計算域外部虛擬點的高度數(shù)值的計算公式如下［16］:

對于表面缺陷圖像，其計算域邊界往往屬于正常部分，像素點高度值是固定不變的，可以認為其計算域外部虛擬像素點的高度也是固定的，因此在實際運算中邊界條件的更新可以省略.

4)計算結(jié)束.如果‖Vn+1i，j-Vni，j‖≤δ(δ為給定的收斂閾值，表示連續(xù)兩次函數(shù)值計算結(jié)果的差異)，或者迭代次數(shù)到達給定的次數(shù)，則算法結(jié)束并輸出像素點高度函數(shù)的數(shù)值.

3 實驗分析

下面對人工合成圖像以及實際表面缺陷圖像分別進行實驗，驗證改進光照模型的正確性，以及所采用算法的有效性.實驗運行的硬件條件為Intel Core i5 760 2.8 GHz，RAM為2 GB.下面對基于SFS方法的三維重構(gòu)實驗進行詳細分析.

3.1 人工合成圖像實驗

選擇最通用半球和花瓶合成圖像來進行實驗.合成圖像采用本文提出的基于透視投影的改進光照模型，圖像按照如下的參數(shù)生成:光源位于光心位置(0，0，0)T，圖像采集傳感器的焦距為100像素，取參數(shù)σ=0.195.

半球人工合成圖像的三維形狀生成公式為

花瓶人工合成圖像的三維形狀生成公式為

按式(7)、(8)并結(jié)合反射圖方程，生成分辨率為120×120的合成圖像.由于基于透視投影的反射圖方程具有伸縮不變性，因此每次迭代運算后對得到的高度函數(shù)進行歸一化處理，迭代30次后得到的半球和花瓶的三維形狀重構(gòu)結(jié)果如圖3、4所示.

圖3 合成半球的三維重構(gòu)結(jié)果

圖4 合成花瓶的三維重構(gòu)結(jié)果

同時引入高度平均誤差(ME)、高度均方根誤差(RMSE)以及CPU運行時間等標(biāo)準(zhǔn)對結(jié)果進行分析和評價.圖5、6分別為合成半球和花瓶重構(gòu)的三維形狀高度平均誤差和均方根誤差隨迭代次數(shù)的變化情況.表1給出了迭代30次時重構(gòu)結(jié)果.

圖5 合成半球的三維重構(gòu)結(jié)果

表1 歸一化高度平均誤差、均方根誤差及運行時間

由以上圖標(biāo)分析可知，該方法能夠獲得較高精度的重構(gòu)結(jié)果.實驗過程中發(fā)現(xiàn)，迭代次數(shù)超過20次，高度的平均誤差和均方根誤差的下降不再明顯，基本保持不變，這應(yīng)該是數(shù)值計算誤差引起的，包括圖像離散化以及差分近似計算引起的誤差.

3.2 實際表面缺陷圖像實驗

選擇孔洞、夾雜兩種表面缺陷進行實驗，系統(tǒng)的參數(shù)設(shè)置如下:光源位于光心位置(0，0，0)T，圖像采集傳感器焦距為20 mm，取參數(shù)σ=0.428，設(shè)置迭代次數(shù)為20作為結(jié)束條件.對表面缺陷分別在朗伯體反射模型以及本文提出的改進光照模型下進行實驗，圖7為表面缺陷灰度圖像，實驗結(jié)果如圖8和9所示.為更清晰地表達三維重構(gòu)后的結(jié)果，孔洞缺陷本來是凹陷的，將其進行翻轉(zhuǎn)以凸起的形式來表達.

圖7 表面缺陷圖像

圖8 孔洞缺陷三維重構(gòu)結(jié)果

圖9 夾雜缺陷三維重構(gòu)結(jié)果

獲取冷軋帶鋼表面缺陷的確切三維高度信息在目前的實驗條件下是比較困難的，在無法進行定量分析的情況下，只能從視覺效果上對它們進行對比分析.由圖9中可以明顯看到，基于改進光照模型恢復(fù)的三維形狀要比基于朗伯體光照模型恢復(fù)的三維形狀更加接近于缺陷的實際形狀，表明改進的光照模型能夠消除凸凹二義性，而通用的朗伯體模型無法做到這一點.總之，基于改進光照模型和高階Lax-Friedrichs和牛頓迭代相結(jié)合的快速掃描SFS算法是正確且有效的，能夠完成實際表面缺陷的三維重構(gòu).

4 結(jié)論

1)在對冷軋帶鋼表面缺陷檢測系統(tǒng)的光路布置以及Oray-Nayar漫反射模型深入分析的基礎(chǔ)上，根據(jù)光路布置對該漫反射模型進行了合理簡化，同時考慮到光源到物體表面距離所起到的衰減作用，提出了一種適用于表面缺陷檢測系統(tǒng)的改進光照模型，并在實驗中得到了驗證.

2)將透視投影模型與改進光照模型相結(jié)合，推導(dǎo)出了透視投影下的反射圖方程，完成了數(shù)學(xué)模型的建立，針對該反射圖方程，提出了一種高階Lax-Friedrichs漢密爾頓函數(shù)和牛頓迭代相結(jié)合的快速掃描算法來實現(xiàn)該方程的求解.

3)利用人工合成圖像和表面缺陷圖像進行了實驗驗證，人工合成圖像實驗結(jié)果表明本文所設(shè)計的三維重構(gòu)算法能夠獲得較高的重構(gòu)精度，同時表面缺陷圖像實驗結(jié)果進一步驗證了改進光照模型以及改進高階Lax-Friedrichs快速掃描算法的正確性與有效性，充分說明了利用基于SFS原理的三維重構(gòu)算法能夠有效的實現(xiàn)表面缺陷的三維重構(gòu).

［1］ZHANG Ruo，TSAI P S，EDWIN J，et al.Shape from shading:a survey［J］.IEEE Transactions on Pattern Analysis and Machine Intelligence，1999，21(8):690-706.

［2］CARTER C N，PUSATERI R J，CHEN Dongqing，et al.Shape from shading for hybrid surfaces as applied to tooth reconstruction［C］//Proceedings of 2010 IEEE 17thInternational Conference on Image Processing.Hong Kong:IEEE Computer Society，2010:4049-4052.

［3］趙輝煌，周德儉，黃春躍，等.基于SFS原理的SMT焊點表面三維重構(gòu)技術(shù)研究［J］.計算機科學(xué)，2009，36(7):288-291.

［4］HORN B K.Height and gradient from shading［J］.International Journal of Computer Vision，1990，5(1):37-75.

［5］OREN M，NAYAR S K.Generalization of lambert’s reflectance model［C］//Proceedings of the 21stAnnual Conference on Computer Graphics and Interactive Techniques.Orlando，USA:ACM，1994:239-246.

［6］PHONG B T.Illumination for computer generated pictures［J］.Communications of the ACM，1975，18(6):311-317.

［7］TORRANCE K E，SPARROW E M.Theory for off-specular refelction from roughed surfaces［J］.Journal of the Optical Society of America，1967，57(9):1105-1114.

［8］HEALEY G，BINFORD T O.Local shape from specularity［J］.Computer Vision，Graphics and Image Processing，1988，42(1):62-86.

［9］FRANKOT R T，CHELLAPPA R.A method for enforcing integrability in shape from shading algorithms［J］.IEEE Transactions on PAMI，1988，10(4):439-451.

［10］LEE C H，ROSENFELD A.Improved methods of estimating shape from shading using the light source coordinate system［J］.Artificial Intelligence，1985，26(2):125-143.

［11］DUPUIS P，OLIENSIS J.Direct method for reconstructing shape from shading［C］//Proceedings of the image Understanding Workshop.San Diego，CA，USA:Morgan Kaufmann Publ Inc，1992:453-458.

［12］ZHANG Li，YIP A M，TAN C L.Shape from shading based on lax-friedrichs sweeping and regularization technique with application to document image restoration［C］//Proceedings of 2007 IEEE Computer Society Conference on CVPR.Minneapolis，MN，USA:IEEE Computer Society，2007:1180-1187.

［13］YUEN S Y，TSUI Y Y，Chow C K.A fast marching formulation of perspective shape from shading under frontal illumination［J］.Pattern Recognition Letters，2007，28(7):806-824.

［14］PARDOS E，CAMILLI F，F(xiàn)AUGERAS O.A unifying and rigorous shape from shading method adapted to realistic data and applications［J］.Journal of Mathematical Imaging and Vision，2006，259(3):307-328.

［15］BARDI M，CAPUZZO-DOLCETTA I.Optimal control and viscosity solutions of hamilton-jacobi-bellman equations［M］.［S.l.］:Birkh?user Boston，1997:1-23.

［16］KAO C Y，OSHER S，QIAN J.Lax-friedrichs sweeping scheme for static hamilton-jacobi equations［J］.Journal of Computational Physics，2004，196(1):367-391.

［17］ZHANG Yongtao，ZHAO Hongkai，QIAN Jianliang.High order fast sweeping methods for static hamilton-jacobi equations［J］.Journal of Scientific Computing，2006，29(1):25-56.

［18］劉瑞玲，韓九強.透視投影下的鏡面反射表面形狀恢復(fù)新算法［J］.西安交通大學(xué)學(xué)報，2009，43(2):6-10.