屈乾龍，朱慶偉，艾衛(wèi)濤，溫 波，馬肇祥

(1.西安科技大學(xué) 測(cè)繪科學(xué)與技術(shù)學(xué)院，陜西 西安 710054；2.河南城建學(xué)院 測(cè)繪與城市空間信息學(xué)院，河南 平頂山 467000)

0 引 言

隨著時(shí)間的推移，古代高塔如今都是“十塔九斜”，塔體受基礎(chǔ)和結(jié)構(gòu)本身的影響及諸多外界因素而產(chǎn)生彎曲變形，綜合表現(xiàn)為整個(gè)建筑體傾斜變形[1-3]。傾斜度是衡量施工技術(shù)質(zhì)量和后期建筑物安全維護(hù)的重要指標(biāo)，如何快速準(zhǔn)確測(cè)量塔類(lèi)建筑物傾斜度一直都是變形觀測(cè)領(lǐng)域的研究熱點(diǎn)。三維激光掃描技術(shù)以其無(wú)需直接接觸即可快速獲取物體表面高精度三維點(diǎn)云數(shù)據(jù)的特點(diǎn)被廣泛應(yīng)用于塔類(lèi)建筑物的變形監(jiān)測(cè)[4-7]。國(guó)內(nèi)對(duì)塔形建筑物的變形監(jiān)測(cè)研究從未停止過(guò)，梁華等基于重慶某鐵塔點(diǎn)云數(shù)據(jù)，對(duì)鐵塔結(jié)構(gòu)底部特征進(jìn)行圓柱擬合獲取鐵塔中心坐標(biāo)，再分別提取關(guān)鍵部位點(diǎn)云，利用特征擬合獲取幾何參數(shù)，計(jì)算鐵塔上、中、下端中心到鐵塔軸線偏差距離[8]。辛星等基于三維激光掃描技術(shù)對(duì)塔形構(gòu)筑物進(jìn)行傾斜監(jiān)測(cè)，利用Imageware軟件對(duì)點(diǎn)云數(shù)據(jù)提取構(gòu)筑物中軸線，計(jì)算出構(gòu)筑物的傾斜值。最后與全站儀前方交會(huì)精度實(shí)驗(yàn)對(duì)比驗(yàn)證三維激光掃描技術(shù)在塔型構(gòu)筑物傾斜監(jiān)測(cè)中應(yīng)用的可靠性[9]。王莉等以銅川市延昌塔為研究對(duì)象，采用三維激光掃描技術(shù)獲取塔體點(diǎn)云數(shù)據(jù)，利用Geomagic軟件進(jìn)行建模，然后利用MATLAB軟件對(duì)點(diǎn)云數(shù)據(jù)擬合計(jì)算，得到塔體中軸線，繼而分析判斷塔體的傾斜變化狀況[10]。楊永林等以西安萬(wàn)壽寺塔為例，基于三維激光掃描技術(shù)對(duì)塔體進(jìn)行變形監(jiān)測(cè)，提出最小二乘塔體傾斜度計(jì)算方法。該方法通過(guò)提取每層塔體的輪廓點(diǎn)云，根據(jù)最小二乘原理擬合出每層塔體的圓心，對(duì)每層圓心進(jìn)行線性擬合提取塔體中軸線，繼而計(jì)算出塔體的傾斜量[11]。

以上學(xué)者在采用三維激光掃描技術(shù)做傾斜監(jiān)測(cè)時(shí)，無(wú)論是橫斷面提取中心軸還是提取輪廓點(diǎn)擬合中心軸，均采用擬合中心軸的方法，數(shù)據(jù)處理方案較為單一，且均是將塔體看作一個(gè)剛體計(jì)算整體傾斜量，不能逐層計(jì)算。筆者提出一種點(diǎn)云輪廓點(diǎn)Z坐標(biāo)傾斜監(jiān)測(cè)數(shù)據(jù)處理方案，與提取中心軸法共同對(duì)西安大雁塔點(diǎn)云數(shù)據(jù)進(jìn)行傾斜量計(jì)算，最后將計(jì)算的傾斜量進(jìn)行對(duì)比。計(jì)算結(jié)果顯示，點(diǎn)云數(shù)據(jù)處理方法對(duì)塔體建筑物傾斜監(jiān)測(cè)具有可行性和工程實(shí)際意義。

1 研究對(duì)象概括

大雁塔位于中國(guó)陜西省西安市南郊大慈恩寺內(nèi)，大雁塔是磚仿木結(jié)構(gòu)的四方形樓閣式磚塔，由塔基、塔身、塔剎組成，現(xiàn)通高為64.517 m。塔基高4.200 m，南北約48.700 m，東西45.700 m；塔體呈方錐形，平面呈正方形，底邊長(zhǎng)為25.500 m，塔身高59.900 m，塔剎高4.870 m[12]。據(jù)記載，陜西測(cè)繪部門(mén)自1985年開(kāi)始對(duì)其進(jìn)行測(cè)量，其時(shí)大雁塔的傾斜速度處于加快過(guò)程中，到1996年左右塔傾斜達(dá)到最大程度，傾斜度達(dá)到1.011 m。如今，具有1 300多年歷史的大雁塔整體屬于動(dòng)態(tài)平衡。大雁塔歷史悠久且傾斜明顯，故以大雁塔為例驗(yàn)證文中方法的可行性，對(duì)大雁塔的安全運(yùn)營(yíng)及維修加固具有一定的參考意義。

2 點(diǎn)云數(shù)據(jù)獲取與處理

2.1 數(shù)據(jù)來(lái)源

三維激光掃描數(shù)據(jù)采集的誤差來(lái)源大致分為3種：與目標(biāo)物反射面有關(guān)的誤差、儀器誤差、外界環(huán)境條件。王俊杰等通過(guò)以Trimble GX三維激光掃描儀為例研究數(shù)據(jù)采集誤差對(duì)數(shù)據(jù)精度的影響[13]。最終實(shí)驗(yàn)得知測(cè)距誤差對(duì)點(diǎn)位精度影響微乎其微；由于掃描角的誤差，點(diǎn)位誤差隨測(cè)量距離的增大而增大，可通過(guò)點(diǎn)位精度估計(jì)公式來(lái)計(jì)算出掃描最大點(diǎn)位誤差，根據(jù)實(shí)際工作要求的點(diǎn)位精度來(lái)選取適當(dāng)?shù)膾呙杈嚯x。文中數(shù)據(jù)是于2016年所測(cè)大雁塔點(diǎn)云數(shù)據(jù)。該數(shù)據(jù)是由Scan Station C10掃描儀采集的點(diǎn)云數(shù)據(jù)。設(shè)站情況如圖1所示，站點(diǎn)分布在塔的四周。近距離站點(diǎn)共8站，主要是獲取塔的中下部點(diǎn)云數(shù)據(jù)；遠(yuǎn)距離的站點(diǎn)共4站，主要是獲取塔的頂部點(diǎn)云數(shù)據(jù)[12]。數(shù)據(jù)采集時(shí)，X軸在橫向掃描面內(nèi)，Y軸在橫向掃描面內(nèi)與X軸垂直，Z軸與橫向掃描面垂直。文中所涉及的所有Z值均是指Z軸坐標(biāo)值。

圖1 設(shè)站分布

2.2 數(shù)據(jù)預(yù)處理

在獲得原始點(diǎn)云數(shù)據(jù)后，第一步要對(duì)多站掃描數(shù)據(jù)進(jìn)行拼接配準(zhǔn)處理。由于三維激光掃描過(guò)程中難免會(huì)掃描到塔體周?chē)鞣N地物而產(chǎn)生巨量的噪聲點(diǎn)，為了防止拼接過(guò)程中因?yàn)楹Ａ奎c(diǎn)云數(shù)據(jù)而導(dǎo)致軟件崩潰以及提高各站點(diǎn)云數(shù)據(jù)拼接速度，在原始數(shù)據(jù)導(dǎo)入Leica Cyclone軟件后，對(duì)每站除待測(cè)塔體和標(biāo)靶球區(qū)域外大部分噪聲點(diǎn)云，直接在Cyclone軟件中進(jìn)行區(qū)域裁剪，快速除去大量噪聲點(diǎn)云，這樣后續(xù)點(diǎn)云數(shù)據(jù)拼接和融合速度大大提高[14-16]。然后逐站標(biāo)記同名點(diǎn)和標(biāo)靶球，通過(guò)標(biāo)靶球和同名點(diǎn)混合拼接方法進(jìn)行數(shù)據(jù)拼接。為了減少個(gè)別標(biāo)靶球或同名點(diǎn)的誤差傳播，采用每2站數(shù)據(jù)拼接一次，然后再把拼接好的結(jié)果進(jìn)行拼接[17]，嚴(yán)格把控各站拼接精度，對(duì)于拼接誤差較大的測(cè)站，將該站標(biāo)靶球和同名點(diǎn)刪除后重新標(biāo)記同名點(diǎn)再進(jìn)行拼接，直至拼接誤差在0.006 m左右。最終將各站基于儀器內(nèi)部獨(dú)立坐標(biāo)系統(tǒng)的點(diǎn)云數(shù)據(jù)轉(zhuǎn)化到同一坐標(biāo)系中，形成一個(gè)整體。

拼接成一個(gè)整體后，仍存在大量的重復(fù)點(diǎn)云，所以需要通過(guò)設(shè)置平均采樣間隔或其它參數(shù)設(shè)置來(lái)進(jìn)行點(diǎn)云數(shù)據(jù)融合處理，去除大量重復(fù)冗余點(diǎn)云[18-20]。

3 傾斜監(jiān)測(cè)方法

3.1 已有三維激光掃描傾斜監(jiān)測(cè)方法

文獻(xiàn)[9]提出基于中軸線上節(jié)點(diǎn)坐標(biāo)偏移的方法，其思路是提取中軸線上節(jié)點(diǎn)的坐標(biāo)[21-22]?；玖鞒倘鐖D2所示。

圖2 切片中心軸傾斜監(jiān)測(cè)基本流程

文獻(xiàn)[10]利用Geomagic軟件進(jìn)行建模，構(gòu)建OBJ三維格網(wǎng)模型；再通過(guò)軟件提取每層塔體特征點(diǎn)坐標(biāo)，計(jì)算出該層中心點(diǎn)坐標(biāo)，利用 MATLAB軟件擬合計(jì)算塔的中心軸線，最終計(jì)算分析塔體傾斜變化狀況，基本原理如圖3所示。

圖3 擬合中心軸計(jì)算原理

文獻(xiàn)[11]提出基于最小二乘的塔體中軸線提取方法，其思路是首先提取出每層塔體外輪廓點(diǎn)云的圓心坐標(biāo)，然后基于最小二乘原理對(duì)所有層圓心坐標(biāo)進(jìn)行線性擬合獲取塔體中軸線的方程，最后計(jì)算2期塔體的傾斜度，實(shí)現(xiàn)塔體的傾斜監(jiān)測(cè)。

以上方法無(wú)論是對(duì)塔體進(jìn)行區(qū)域分割還是基于塔體塔檐提取特征輪廓點(diǎn)，最后都是基于點(diǎn)云數(shù)據(jù)提取塔體中心軸的方法來(lái)確定塔體傾斜度，文中所提出的方法與該類(lèi)方法有所不同，不用提取中心點(diǎn)即可計(jì)算出塔體傾斜度。

3.2 輪廓點(diǎn)Z坐標(biāo)傾斜監(jiān)測(cè)方法

筆者提出基于輪廓點(diǎn)云Z坐標(biāo)提取計(jì)算法進(jìn)行傾斜監(jiān)測(cè)?；舅悸肥菙M合提取塔體固定輪廓特征點(diǎn)坐標(biāo)，比如每一層的塔檐。如圖4所示，提取出每層塔檐外輪廓點(diǎn)云的坐標(biāo)，A，B，C，D是塔檐的4個(gè)角，我們只需在提取出的外輪廓點(diǎn)云坐標(biāo)中篩選出其Z坐標(biāo)值最大點(diǎn)D和Z坐標(biāo)值最小點(diǎn)B，此時(shí)LDO為D、B這2點(diǎn)Z坐標(biāo)之差ΔZ，LDB可根據(jù)2點(diǎn)間距離公式求出，最后利用正弦定理公式(1)和反三角函數(shù)公式(2)即可算出該2點(diǎn)連線LDB與水平面夾角∠OBD，該角在數(shù)值上等于塔體傾斜角。

圖4 基于輪廓點(diǎn)Z坐標(biāo)傾斜監(jiān)測(cè)法原理

(1)

又∠DOB=90°，故公式(1)可化成

(2)

如圖5所示，這種是塔體傾斜特殊情況，塔體塔檐連線嚴(yán)格向另一側(cè)傾斜，這時(shí)LAD，LBC線上各點(diǎn)的Z坐標(biāo)總體相同，此時(shí)可直接選取同側(cè)塔角點(diǎn)A，B，或者點(diǎn)D，C為Z坐標(biāo)最大值和最小值點(diǎn)計(jì)算∠OEF。

圖5 塔體傾斜特殊情況

本方法同樣適用于圓形或六邊形等規(guī)則塔體建筑，原理相同，在此不再贅述。

4 實(shí)驗(yàn)數(shù)據(jù)計(jì)算對(duì)比分析

4.1 擬合提取中心軸線傾斜監(jiān)測(cè)法

擬合提取中心軸線傾斜監(jiān)測(cè)法具體操作在此不再詳細(xì)介紹，基本操作流程如圖6所示。用Leica Cyclone軟件將點(diǎn)云數(shù)據(jù)預(yù)處理后，將拼接后的點(diǎn)云數(shù)據(jù)導(dǎo)入Imageware軟件中進(jìn)行區(qū)域等距離切割點(diǎn)云，并將切割后的點(diǎn)云輪廓數(shù)據(jù)提取出來(lái)[23-25]，為了提高精度，準(zhǔn)確計(jì)算出塔體各中心坐標(biāo)，因此首先將輪廓點(diǎn)云擬合成等間距圓，如圖7所示，共提取7條輪廓線。

圖6 中心軸法實(shí)驗(yàn)流程

圖7 切割輪廓點(diǎn)提取擬合圓

再計(jì)算出各層中心點(diǎn)圓心坐標(biāo)并將計(jì)算得出的坐標(biāo)和輪廓點(diǎn)坐標(biāo)數(shù)據(jù)導(dǎo)入軟件具體效果如圖8所示。

圖8 頂點(diǎn)及各層中心點(diǎn)效果

計(jì)算出各層中心點(diǎn)坐標(biāo)導(dǎo)入MATLAB軟件顯示后明顯發(fā)現(xiàn)各中心點(diǎn)與塔頂中心點(diǎn)并不在同一直線上，這是因?yàn)樗w建筑傾斜形變并不是線性形變的，因此中心軸法需要將塔體作為一個(gè)剛體整體分析。將中心點(diǎn)坐標(biāo)和塔頂點(diǎn)坐標(biāo)通過(guò)MATLAB編寫(xiě)程序，基于最小二乘擬合原理提取塔體中軸線，結(jié)果如圖9所示。

圖9 塔體中心圓點(diǎn)擬合直線(單位：m)

最終得到擬合后頂點(diǎn)和各層中心點(diǎn)坐標(biāo)，見(jiàn)表1。

表1 擬合后頂點(diǎn)和各層中心點(diǎn)坐標(biāo)

最終通過(guò)擬合中心坐標(biāo)計(jì)算，該年大雁塔塔體大約整體傾斜0.960 m。

4.2 基于輪廓點(diǎn)Z坐標(biāo)傾斜監(jiān)測(cè)法

流程如圖10所示。在點(diǎn)云拼接、去噪后，為了保證輪廓點(diǎn)提取精度，采用區(qū)域分割分層提取檐輪點(diǎn)云塔廓線，區(qū)域分割時(shí)，只需保證塔檐整體部分在分割區(qū)域內(nèi)即可，如圖11所示。

圖10 輪廓點(diǎn)Z坐標(biāo)提取法實(shí)驗(yàn)流程

圖11 區(qū)域分割

因?yàn)樵跀?shù)據(jù)融合時(shí)為增加后期數(shù)據(jù)處理速度，故設(shè)置采樣間隔進(jìn)行點(diǎn)云抽稀，為防止頂點(diǎn)和塔檐四角點(diǎn)云誤刪而導(dǎo)致Z坐標(biāo)值最大最小點(diǎn)提取不準(zhǔn)確，故在Cloudcompare軟件中選取各層塔檐四角和頂點(diǎn)坐標(biāo)加入輪廓線提取坐標(biāo)群，如圖12所示。

圖12 選取的塔角四點(diǎn)及塔頂

根據(jù)該方法進(jìn)行輪廓點(diǎn)提取，最終各層塔檐共提取點(diǎn)云數(shù)量見(jiàn)表2。

表2 各層塔檐輪廓提取點(diǎn)云數(shù)量

提取各層輪廓點(diǎn)點(diǎn)云坐標(biāo)后，導(dǎo)入表格從大到小排序，將Z坐標(biāo)值最大最小點(diǎn)坐標(biāo)單獨(dú)導(dǎo)出，根據(jù)公式(1)、(2)，求得各層傾斜角度以及基于該層傾斜度乘以塔的通高所求的塔整體傾斜量xi，詳細(xì)數(shù)據(jù)見(jiàn)表3。

表3 各層傾斜度及基于該層傾斜度的整體傾斜量

在求塔體整體傾斜量時(shí)，對(duì)表3中基于各層傾斜度整體傾斜量直接取算術(shù)平均值不能滿足測(cè)量中的精度要求，為了提高塔體整體傾斜計(jì)算精度，引入權(quán)要素，對(duì)表3中整體傾斜量做加權(quán)處理。對(duì)于塔形建筑物而言，底層較頂層相對(duì)傾斜量較小，且對(duì)整體傾斜影響最大，由于激光掃描點(diǎn)精度相同，且輪廓點(diǎn)提取點(diǎn)云密度直接影響能否提取出真正Z坐標(biāo)值最大最小點(diǎn)，每一個(gè)輪廓點(diǎn)都相當(dāng)于同精度觀測(cè)一次，因此用各層塔檐輪廓提取點(diǎn)數(shù)作為Ni次同精度觀測(cè)值參與定權(quán)計(jì)算。7次相同精度觀測(cè)塔體求解塔體傾斜量的過(guò)程，故可求同精度傾斜值算術(shù)平均值的權(quán)，具體公式見(jiàn)下式

(3)

(4)

(5)

式中Pi為第i層傾斜量的權(quán)；Ni為第i層輪廓點(diǎn)提取點(diǎn)云數(shù)，數(shù)據(jù)見(jiàn)表2；σ點(diǎn)為每一個(gè)點(diǎn)云的中誤差，每層測(cè)量精度相同，其中誤差為σ0，σi為第i層傾斜量中誤差。

將公式(4)(5)導(dǎo)入公式(3)中得

(6)

式中C為單位權(quán)輪廓點(diǎn)的點(diǎn)數(shù)，個(gè)。

取第一層輪廓點(diǎn)提取數(shù)做單位權(quán)觀測(cè)值，即C=2054。則根據(jù)公式(6)，可計(jì)算出每層傾斜量的權(quán)Pi。

(7)

由公式(7)可得傾斜量加權(quán)平均值即塔體整體傾斜量，結(jié)果為0.962 m。

4.3 數(shù)據(jù)分析與對(duì)比

通過(guò)提取的輪廓點(diǎn)最大值最小值點(diǎn)坐標(biāo)數(shù)據(jù)對(duì)比后發(fā)現(xiàn)，文中方法所提取的Z坐標(biāo)值最大、最小值點(diǎn)均為塔檐對(duì)角點(diǎn)，可根據(jù)文中方法原理可大致判斷出塔體傾斜方向?yàn)闁|南向西北方向傾斜?；谖闹蟹椒ǎ谀芘袛嗨w整體大致傾斜走向情況下，可直接從塔檐各角坐標(biāo)中提取Z坐標(biāo)值最大、最小點(diǎn)計(jì)算，省去提取輪廓點(diǎn)坐標(biāo)這一步驟，減少計(jì)算量，可大大提高數(shù)據(jù)處理效率。為了驗(yàn)證該規(guī)律的可靠性，在明確塔體大致傾斜方向后，隨機(jī)提取某層塔柱沿該傾斜方向?qū)屈c(diǎn)來(lái)進(jìn)行驗(yàn)證。最終，對(duì)塔體第3層和第4層塔柱的兩對(duì)角點(diǎn)進(jìn)行提取坐標(biāo)計(jì)算塔體的傾斜量以此來(lái)驗(yàn)證該規(guī)律的可靠性，具體數(shù)據(jù)見(jiàn)表4。同樣，和提取中心軸法相比文中方法原理易懂，數(shù)據(jù)計(jì)算處理更為簡(jiǎn)便。

根據(jù)表5各層傾斜度數(shù)據(jù)可以發(fā)現(xiàn)，每層傾斜角度并不相同，塔體底層傾斜量較小，隨著塔體升高各層傾斜量不一，這也是提取中心軸法不便將塔體分層計(jì)算傾斜量原因，塔各層的變形并不是線性的，采用不同層數(shù)據(jù)進(jìn)行計(jì)算整體傾斜量會(huì)有不同的結(jié)果。因此現(xiàn)有基于點(diǎn)云數(shù)據(jù)提取中心軸法進(jìn)行塔形建筑物傾斜監(jiān)測(cè)分析時(shí)，均是將塔體假設(shè)為一個(gè)剛體進(jìn)行分析，利用塔體中心點(diǎn)數(shù)據(jù)擬合出塔體的中軸線數(shù)據(jù)再計(jì)算塔體傾斜量。而文中提出的方法是通過(guò)分層計(jì)算塔體傾斜角度，再基于該傾斜角度進(jìn)行整體偏移量的加權(quán)計(jì)算，既可計(jì)算塔整體傾斜量，亦可計(jì)算各層傾斜量。通過(guò)表4可以看出，無(wú)論是從塔檐提取輪廓點(diǎn)還是從塔柱提取輪廓點(diǎn)，都可以滿足文中輪廓點(diǎn)Z坐標(biāo)傾斜監(jiān)測(cè)方法。結(jié)合表5看出，2種方法計(jì)算塔體傾斜量盡管存在一些差異，所得出的大雁塔整體傾斜量相差0.002 m，約0.26%，但是在測(cè)量誤差范圍之內(nèi)?？傮w來(lái)說(shuō)，通過(guò)與擬合中心軸法整體傾斜量計(jì)算對(duì)比表明，文中所提的基于輪廓點(diǎn)Z坐標(biāo)的塔形建筑物傾斜監(jiān)測(cè)方法切實(shí)可行。

表4 數(shù)據(jù)特征驗(yàn)證

表5 兩方法傾斜量計(jì)算結(jié)果對(duì)比

5 結(jié) 論

1)對(duì)大雁塔三維激光掃描數(shù)據(jù)進(jìn)行研究計(jì)算，驗(yàn)證文中所提出的基于輪廓點(diǎn)Z坐標(biāo)傾斜監(jiān)測(cè)法的有效性，該方法可以對(duì)基于三維激光掃描塔形建筑物傾斜監(jiān)測(cè)提供一種新的數(shù)據(jù)處理方案。

2)采用基于塔體各層傾斜角度加權(quán)計(jì)算塔形建筑物整體傾斜量的方法，因此可以明確獲得塔體各層傾斜度，為塔體某層具體修繕維護(hù)提供數(shù)據(jù)參考。

3)在明確塔體基本傾斜方向后，可以直接提取該傾斜方向各層塔檐對(duì)角點(diǎn)云坐標(biāo)數(shù)據(jù)，作為Z坐標(biāo)最大最小值點(diǎn)參與傾斜量計(jì)算，減少點(diǎn)云計(jì)算量，可大大提高數(shù)據(jù)處理效率。

4)該方法適用性廣泛，圓形、六邊形等同類(lèi)項(xiàng)目，可基于文中方法進(jìn)行傾斜監(jiān)測(cè)，對(duì)無(wú)中心軸不規(guī)則建筑物傾斜監(jiān)測(cè)項(xiàng)目也具有一定的參考價(jià)值。