劉綱,李孟珠,蔣偉,張維慶
(1.山地城鎮(zhèn)建設與新技術教育部重點實驗室(重慶大學),重慶 400045;2.重慶大學 土木工程學院,重慶 400045)
中國古建筑具有悠久的歷史傳統(tǒng),是當時科學技術水平和人民社會生活狀況的真實見證,具有重要的歷史、藝術和科學價值[1].在以木結構為主的古建筑中,受木材蠕變、腐蝕、蟲蛀和長期荷載作用[2-3],木質構件變形過大問題較為普遍[4],木結構變形監(jiān)測可用于判斷變形發(fā)展趨勢、提供結構受力計算所需變形數(shù)據(jù),是古建筑木結構預防性保護的重要手段.
傳統(tǒng)上古建筑變形測量以接觸式方法為主,例如在結構構件上布置觀測點或粘貼光纖傳感器[4-5],這些方法對古建筑構件有一定程度的破壞,不利于古建筑保護.雖然近年來無需在古建筑上布設測試點的免棱鏡技術得到快速發(fā)展,但是其測試精度受較多因素限制,在部分條件下較難滿足工程測試要求.近年來,非接觸式變形測量方法得到快速發(fā)展,其中,基于數(shù)字圖像相關(Digital image correlation,DIC)原理[6]的測試方法隨著攝像技術、計算能力的提升得到了大力發(fā)展,已在航天、機械和土木工程領域得到實際應用.DIC 法最早可追溯到20 世紀80 年代,美國南卡羅萊納大學Peters 等[7]針對材料表面均勻拉伸應變測量提出基于DIC 進行位移測量的思路;Sutton 等[8]針對全場平面內變形的整像素位移搜索,提出粗-細搜索法,從而提高整像素級圖像匹配的計算速度.Gencturk 等[9]將DIC 法應用于預應力混凝土位移測量,結果表明該方法可準確測量預應力混凝土變形信息.奧村運明等[10]應用DIC 法進行古建筑墻體裂縫觀測,結果顯示該方法能滿足工程上裂縫測量的要求.
DIC 法通過在被測結構表面噴涂具有一定特征的散斑圖,然后利用結構位移前后散斑圖的變化來識別結構變形[11],而中國古建筑木結構的表面往往裝飾有彩繪圖案[12],故利用彩繪圖也可實現(xiàn)位移測量.但應指出的是,與特制散斑圖相比,彩繪圖灰度梯度分布不均且存在大塊灰度近似區(qū),故直接采用現(xiàn)有DIC 位移解算方法將無法保證位移測試精度.針對這一問題,在現(xiàn)有自適應十字模式搜索方法(Adaptive Rood Pattern Search,ARPS)的基礎上,提出局部窮舉-自適應十字模式搜索方法(Improve-APRS,IARPS),提高整像素位移搜索精度,從而增強古建筑彩繪圖位移辨識精度,并提升位移辨識計算的穩(wěn)定性和效率.通過實驗驗證所提方法的適用性,并與常用整像素位移解算方法進行對比,為古建筑木結構非接觸式位移測量提供新的解算方法.
基于DIC 原理的位移測試技術在被測構件上噴涂特制散斑圖,借用數(shù)字攝影技術將散斑圖在相機感光組件上成像,然后通過匹配變形前后數(shù)字圖像像素的變化,從而獲得被測物體表面各點的位移信息,該方法的位移辨識原理如圖1 所示.
圖1 二維DIC 位移辨識原理Fig.1 Displacement identification principle of two-dimensional DIC
在圖1 中,一個小方格代表數(shù)字圖像上的一個像素,變形前的圖像稱為參考圖像,變形后的圖像稱為目標圖像.針對參考圖像、目標圖像建立統(tǒng)一的坐標系x-y,參考圖像中(x,y)處像素點的灰度值為(fx,y);目標圖像中(x′,y′)處的像素點灰度值為g(x′,y′).設參考圖像中P 點坐標為(xp,yp),目標圖像中與P 點對應的P′點坐標為(x′p,y′p),DIC 方法計算P′點到P 的距離可分為以下幾個步驟:
1)在參考圖像中,選取整數(shù)M,然后選擇一個計算點P(xp,yp),以該點為中心選?。?M+1)×(2M+1)像素大小的區(qū)域為計算子集,其中2M+1 稱為子集半徑.
2)在目標圖像中,任選一點并作為搜索點,以該點為中心,框選(2M+1)×(2M+1)的計算子集,再計算相關系數(shù)CZNSSD:
式中:fm表示參考圖像子集的灰度平均值;gm表示目標圖像子集的灰度平均值.然后,通過一定的搜索方法計算多個點的相關系數(shù),并將相關系數(shù)最值點對應的搜索點作為P′點的匹配點.此時,可確定P 點整像素x 向位移up、y 向位移vp.
3)根據(jù)求得的匹配點的坐標,在(xp+up,yp+vp)像素點相鄰1 個像素范圍內通過一定方法進行亞像素位移搜索,得到P 點分別在x、y 向亞像素位移Δup、Δvp.
4)根據(jù)求得的整、亞像素位移,最終解算得到P點x 向位移u、y 向位移v 分別為:
從以上步驟可知,DIC 法先執(zhí)行整像素位移搜索,并將搜索結果作為亞像素搜索的初始值,故整像素搜索算法不但決定整像素搜索精度,也將影響亞像素搜索精度.因此,整像素搜索算法對位移辨識精度的影響較大.
需指出的是,第一次在目標圖像中任選一點時,若該點距圖像邊緣的距離小于2M+1,則無法以該點為中心組成計算子集,稱為邊緣問題.為避免這一問題,可在初次選擇計算點時進行預判:若計算點不存在邊緣問題,則進行以上步驟的搜索;若不滿足,則重新選取計算點,直至找到不存在邊緣問題的計算點.
目前,DIC 中最常用的整像素位移搜索算法有粗-細搜索法(Coarse-Fine Search,CFS)[13]、三步搜索法[14]及菱形搜索法[15]等,其中CFS 算法的精度和計算效率相對較高,其位移搜索的基本思路為:將整個目標圖像視為搜索區(qū)域,首先采用較大的計算步長計算相關系數(shù),以相關系數(shù)極小值對應的搜索點為本輪搜索的匹配點;然后以此匹配點為中心,縮小搜索區(qū)域及計算步長選取搜索點進行下一輪搜索;最后,當計算步長縮小為1 像素時,可得P′點的匹配點,進而計算得到該點的整像素位移.
設參考圖像和目標圖像的像素均為I×K,則兩幅圖像的均方誤差MES 定義為:
兩幅圖像的峰值信噪比(Peak-Signal-to-Noise-Ratio,PSNR)定義為:
式中:MAXI表示圖像灰度的最大值,對于采樣點用8 進制表示的灰度圖像取255.峰值信噪比可評價兩幅圖像子集的近似性,近似程度越高,PSNR 系數(shù)越高,即搜索算法的精度越高.
CFS 算法假定每一輪最佳匹配點處的相關函數(shù)值CZNSSD單調減小,即由相關系數(shù)組成的二維曲面是單峰的.通常情況下,由散斑圖得到的二維曲面可能存在幾個局部極值現(xiàn)象,如圖2(a)所示.此時,只要計算步長選擇合適,CFS 仍可搜索到正確結果.但對于彩繪圖像,由于其灰度分布不均,所得相關函數(shù)曲面往往存在多峰現(xiàn)象,在CFS 搜索過程中的單調性不復存在,如圖2(b)所示.此時,初始計算步長對CFS 算法精度甚至正誤影響很大,而初始計算步長通常由人工主觀選取,具有很大的隨機性.
為克服計算步長選取的隨機性,Nie 等[16]提出了自適應十字模式搜索法(Adaptive Rood Pattern Search,ARPS).其基本思路是,假設物體的變形是連續(xù)的,在進行整像素位移搜索時,采用相鄰已知位移像素點的位移作為搜索初值,其具體搜索步驟為:
圖2 單峰與多峰曲面示意圖Fig.2 Schematic figure of unimodal and multimodal surfaces
1)對參考圖像中的某計算點P,設其在x-y 坐標系中的坐標為(u1,v1);記P 相鄰像素點已知整像素位移為(u,v),記R=max(|u|,|v|),R=(u,v).在目標圖像中,以坐標(u1,v1)所在的像素點為中心點O,若是第一個計算點,則在該點沿x-y 軸向相距2 像素的上下左右各選1 個點,形成十字形搜索模式;若不是第一個計算點,則在該點沿x-y 軸向相距R 像素的上下左右各選1 個點,形成十字形搜索模式,同時選取R 為一個搜索點,共計5 個搜索點,如圖3(a)中黑色小圓點所示.
2)按式(1)計算圖3(a)中5 個搜索點的相關系數(shù),選取相關系數(shù)最小值所在點為初始匹配點.以初始匹配點為中心,選取該點及上下左右4 個相鄰像素點形成單位十字搜索模式,并將這5 個點作為新的搜索點,如圖3(b)所示.
3)按式(1)計算單位十字搜索模式中5 個搜索點的相關系數(shù),若最佳匹配位置為單位十字搜索模式的中心點,則結束搜索,該點即為計算點P 的最佳匹配點P′(u2,v2);否則,將單位十字中心移至新的相關系數(shù)最小點,重復進行單位十字模式搜索,直至相關系數(shù)最小點為十字搜索模式的中心為止,此時十字模式的中心點即為計算點的最佳匹配點P′.
4)根據(jù)求得的最佳匹配點P′的坐標(u2,v2),可解算P 點的整像素位移(u2-u1,v2-v1).
當計算參考圖像中的第一點時,其相鄰點位移為未知,故在第1)步中僅采用十字形搜索模式的4個點進行搜索,其余步驟與以后各計算點的步驟相同.
ARPS 法根據(jù)相鄰像素已知位移自適應調整步長,有利于解決相關系數(shù)曲面存在多個極值的問題,但對第一個位移搜索點的位移初值采取零位移假設(僅采用4 個搜索點進行計算),可能帶來一定的誤差.為解決該問題,可先預估第一個位移搜索點的位移最大值,并將其作為初始搜索半徑,進行局部窮舉搜索,從而為后續(xù)各搜索點的位移計算提供準確的初值.局部窮舉搜索具體步驟為:
圖3 整像素搜索法示意圖Fig.3 Schematic figure of integer pixel search
1)對于參考圖形中的某搜索點P,已知其在x-y坐標系中的坐標為(u1,v1);在目標圖像中,以坐標(u1,v1)所在的點為中心點,選擇一個初始搜索半徑,如圖3(c)所示,在該方形區(qū)域內對每個點進行局部窮舉搜索.
2)針對方形區(qū)域內的每一個點,計算所有點的相關系數(shù).
3)選取所有搜索點中相關系數(shù)最小點為第一個初始計算點的最佳匹配點,然后按照ARPS 法的第2)~4)步計算該點位移.
4)按照ARPS 方法的第1)~4)步計算剩余各計算點的位移.
采用窮舉搜索可確保搜索到精確的位移初始向量,避免陷入圖2(b)所示的錯誤局部最優(yōu)解,將增加局部窮舉搜索后的ARPS 法稱為修正的自適應十字模式搜索法(Improved-ARPS,IARPS).
分別通過模擬散斑圖發(fā)生剛體位移、實驗室木梁彩繪圖變形實驗對比IARPS 算法與CFS、ARPS 算法的區(qū)別.
3.1.1 散斑圖位移模擬
采用Zhou 等[17]提出的散斑圖生成算法得到256×256 像素的散斑圖作為參考圖像,其中高斯光斑數(shù)目為1 000 個,高斯光斑大小為4 像素,高斯光斑的中心光強為0.7.將參考圖像在y 向移動6 個像素作為目標圖像,即目標圖像各整像素點的位移向量均為(0,6).為模擬測試噪聲,在目標圖像中加入均值為0,標準差為2 的高斯白噪聲.為節(jié)約篇幅,將參考圖像、目標圖像分左右兩幅給出,如圖4 所示.
圖4 模擬散斑圖像Fig.4 Simulated speckle image
在散斑圖中,僅計算圖4 中虛線所在像素點發(fā)生的位移,該虛線由90 個像素點組成,其像素坐標自上而下定義為1~90 pixel.
3.1.2 彩繪圖實驗
在實驗室中采用長寬高分別為2 000 mm×40 mm×25 mm 的木梁進行實驗,將規(guī)則花紋粘貼于梁中部側面.采用重物分兩次在跨中加載,使梁產(chǎn)生撓度,分別記為工況1 和工況2,并在跨中安裝千分表對辨識位移的精度進行驗證,具體的加載裝置如圖5所示.實驗中采用的圖像采集系統(tǒng)為VIC-2D,相機型號為GZL-CL-41C6M-C,其分辨率為2 048×2 048像素.將相機放置在彩繪梁正前方約2 m 處,并使其光軸對準梁的形心且與彩繪圖案面垂直,將相機與計算機相連并用測試系統(tǒng)的軟件進行圖像采集.
考慮到窮舉搜索法(Exhaustive Search,ES)將計算搜索指定范圍內所有可能的像素點,其計算精度在所有的搜索算法中最好,故采用該算法作為基準驗證其他整像素解算方法的精度.
圖5 實驗加載裝置Fig.5 Experimental loading device
為方便與千分表測試結果對比,取千分表正上方虛線所對應像素進行位移辨識,共包括90 個像素點,從上向下依次定義為1~90 pixel.一個像素對應的實際長度為25 mm/90=0.278 mm.計算子集中M取為15,則一個子集包含31×31 像素.
CFS 及IARPS 法計算的PSNR 結果如圖6 所示.由圖6(a)可知,對于散斑圖的整像素位移解算,CFS、IARPS 法均可達到與ES 法相同的精度,這主要是因為散斑圖的灰度分布較均勻,其相關函數(shù)組成的二維曲面是單峰的,相關系數(shù)往往單調減小,故所有方法均可達到一樣的精度;但對于ARPS 算法而言,在像素坐標1 處的精度較低,這是由于假定的初始假設位移為(0,0),不能精確解算第一個搜素點的位移.從圖6(b)可看出,在不同像素坐標處,CFS 法對彩繪圖位移的解算精度會出現(xiàn)降低現(xiàn)象,例如從像素坐標13 到20;APRS 法在像素坐標1 處的精度也較低,而IARPS 法對所有像素坐標的解算精度均與ES 法相同,這主要是因為IARPS 法可精確計算第一個搜索點的整像素位移,且繼承了ARPS 高效的自適應性,因此,可有效克服傳統(tǒng)算法局部計算不精確的缺陷,比APRS、CFS 更適合灰度不均勻彩繪梁位移測量.
圖6 不同整像素識別算法的PSNRFig.6 PSNR of different integer pixel recognition algorithms
采用算法對相關函數(shù)的調用次數(shù)作為計算效率的評價指標,調用次數(shù)越多,計算效率越低.4 種方法的計算效率比較如圖7 所示,由于ES 方法是對于全局搜索的,而且ES 方法只起對比作用,因此為使圖像美觀,ES 方法也選擇跟其他算法一樣的8 像素的搜索半徑.
從圖7 可知,窮舉法在每個像素點處均需289次相關函數(shù)運算;CFS 法在不同的搜索點處大約需43 次相關函數(shù)運算;而IARPS 法除在第一個像素點處需局部窮舉計算,需298 次相關函數(shù)運算之外(除局部的窮舉搜索外還進行了自適應的搜索,因此搜索次數(shù)大于ES 法搜索次數(shù)),其余各像素點則需9次相關函數(shù)運算.針對90 個整像素點的位移解算,IARPS 法較CFS 法的計算效率提高71.6%,因此IARPS 法總體上具有較高的計算效率.
圖7 計算效率對比Fig.7 Comparison about calculation efficiency
在整像素搜索的基礎上,采用常用梯度亞像素搜索法進行亞像素位移辨識,得到的位移識別結果與千分表的測試結果見表1.當僅采用整像素的識別結果時,相對誤差高達6.05%,這主要是由于1 個像素對應的位移精度僅為0.278 mm,而整像素搜索結果誤差必然大于該精度,故誤差較大.采用亞像素搜索后,相對誤差可降低到0.064%,小于千分之一.因此,所提IARPS 算法與常規(guī)的亞像素識別算法結合后可用于古建筑構件的非接觸式位移測試.
表1 位移識別結果對比Tab.1 Comparison of displacement recognition results mm
3.3.1 初始搜索半徑大小
對于彩繪圖實驗,采用90 個像素點PSNR 系數(shù)之和——SUM-PSNR 來衡量不同初始搜索半徑下不同算法的解算精度,結果如圖8 所示.
由圖8 可知,CFS 法的初始搜索半徑要達到6像素以上才能有較高的SUM-PSNR 系數(shù),即初始搜索半徑大于像素點實際位移時才能較精確地識別整像素位移.取不同搜索半徑時,ARPS 法的識別精度較IARPS 低,這主要是因為初始計算點假設選取的緣故,IARPS 法的搜索半徑僅為4 像素時就能達到與ES 相同的精度,這是由于IARPS 法中盡管第一步局部窮舉方法識別結果存在少許誤差,但是之后的自適應搜索部分依然可以進行快速糾正,并保證計算結果的精確性.
圖8 初始搜索半徑對識別結果的影響Fig.8 Effect of initial search radius on the recognition results
3.3.2 圖像噪聲分析
采用初始搜索半徑為8 像素,給模擬散斑圖、彩繪圖分別添加均值為0,標準差為1~10 的高斯白噪聲,以測試圖像噪聲對整像素位移識別的影響,結果如圖9 所示.
由圖9 可看出,不管是散斑圖還是彩繪圖,隨著噪聲標準差的逐漸增大,SUM-PSNR 系數(shù)逐漸降低,表明噪聲將降低目標圖像與參考圖像之間的相關性.對散斑圖而言,在噪聲影響下,CFS、IARPS 法的識別精度與ES 法相當,這主要是因為散斑圖灰度分布較均勻,CFS 法可獲得較好的解算精度,但ARPS的精度略低,這是由于該方法在計算第一個搜索點位移時有較大誤差所致;對彩繪圖而言,僅IARPS 法與ES 法的識別精度相當,且CFS 法的精度比ARPS法更差,這主要是因為相同噪聲水平下,灰度分布不均勻導致相關系數(shù)曲面出現(xiàn)多峰情況,導致CFS 精度大幅下降所致.
圖9 測試噪聲對位移識別的影響Fig.9 Effect of test noise on displacement recognition
基于DIC 方法基本原理,結合古建筑彩繪梁圖案灰度分布不均的特點,引入自適應十字模式搜索法,并針對第一個位移搜索點精度較低的缺點,提出修正的自適應十字模式搜索法(IARPS).模擬散斑圖和彩繪梁變形實驗表明:
1)在整像素位移搜索方面,IARPS 法具有較高的計算精度,在不同的像素坐標以及不同計算步長條件下均能得到正確解算結果,有效克服了CFS 法無法適應灰度不均、ARPS 法在第一個位移搜索點精度較低的缺陷.
2)IARPS 法繼承了ARPS 法高效計算的優(yōu)點,在解算90 個整像素點位移時,計算效率較CFS 法提高71.6%.
3)在不同搜索半徑、不同噪聲水平下,IARPS 法的位移解算能力優(yōu)于ARPS 法和CFS 法,且IARPS法的計算精度與精確的完全窮舉法相當.