基于多層中軸點(diǎn)擬合的古塔變形檢測(cè)

劉世杰， 郭成成， 王穗輝， 童小華

(1.同濟(jì)大學(xué) 測(cè)繪與地理信息學(xué)院，上海 200092;2.同濟(jì)大學(xué) 土木工程防災(zāi)國家重點(diǎn)實(shí)驗(yàn)室， 上海 200092)

古塔是具有歷史價(jià)值、藝術(shù)價(jià)值和人文價(jià)值的國家文物，是國家不可再生的文化資源.現(xiàn)存的古塔對(duì)研究我國古代建筑史、宗教史和地方史具有重要的意義.但由于歷史、地質(zhì)、戰(zhàn)爭等原因，使得古塔傾斜、損壞的現(xiàn)象十分普遍，所謂“十塔九歪”.因此只有對(duì)古塔質(zhì)量建立科學(xué)的評(píng)價(jià)體系，進(jìn)行及時(shí)的保護(hù)和修繕，才能延長古塔壽命，繼續(xù)發(fā)揮其獨(dú)特的價(jià)值[1].

在古塔保護(hù)、維修的過程中，首先應(yīng)該對(duì)古塔的傾斜變形進(jìn)行評(píng)估.傳統(tǒng)的古塔傾斜測(cè)量方法首先在古塔周圍建立局部平面坐標(biāo)系，在塔底進(jìn)行測(cè)量得到塔底多邊形，采用多邊形擬合或?qū)蔷€相交的方法確定底部幾何中心，再將測(cè)量得到的塔頂三維坐標(biāo)投影到底部，計(jì)算出投影點(diǎn)與底部中心的距離，根據(jù)此距離與塔高的比值計(jì)算傾斜程度[2].傳統(tǒng)方法僅用塔底幾何中心和塔頂?shù)倪B線來表征古塔姿態(tài)，并未反映出古塔的分層變形情況.劉琦等[3]提出將古塔每層的8個(gè)觀測(cè)點(diǎn)投影到水平面上，然后應(yīng)用最小二乘技術(shù)進(jìn)行圓擬合獲取每層中心，但在實(shí)際的測(cè)量中，塔身及塔體棱角多會(huì)受到破壞而模糊，難以對(duì)棱角點(diǎn)進(jìn)行準(zhǔn)確觀測(cè).周千[4]提出把到各層觀測(cè)點(diǎn)距離平方和最小的點(diǎn)作為各層中心點(diǎn)，進(jìn)而對(duì)中心點(diǎn)擬合空間軸線，其方法對(duì)觀測(cè)點(diǎn)的質(zhì)量有較高要求且難以確定實(shí)際塔體測(cè)量點(diǎn).Tang等[5]通過測(cè)量點(diǎn)擬合表征古塔的正交橫截面，將橫截面交線作為塔體空間軸線，進(jìn)而計(jì)算整體傾斜程度.

本文提出了一種新的古塔變形觀測(cè)與計(jì)算方法，對(duì)古塔的姿態(tài)進(jìn)行全面客觀的描述.使用全站儀對(duì)古塔墻面進(jìn)行點(diǎn)位觀測(cè)，獲取墻面觀測(cè)點(diǎn)位的三維坐標(biāo)數(shù)據(jù)，通過多層中軸點(diǎn)擬合計(jì)算古塔傾斜和扭曲變形.基于本文提出的方法，以安徽省安慶市的振風(fēng)塔為例，進(jìn)行了數(shù)據(jù)采集和變形分析，驗(yàn)證了方法的有效性.

1 實(shí)驗(yàn)對(duì)象和測(cè)量數(shù)據(jù)

實(shí)驗(yàn)對(duì)象為安徽省安慶市的振風(fēng)塔，坐落在長江之畔，建于公元1570年，至今已有400多年歷史，享有“萬里長江第一塔”、“過了安慶不說塔”的美譽(yù)，是七層八角樓閣式建筑，如圖1所示.

圖1 振風(fēng)塔Fig.1 Zhenfeng Pagoda

首先在振風(fēng)塔四周建立了包含15條測(cè)邊的閉合導(dǎo)線網(wǎng)進(jìn)行控制測(cè)量，建立獨(dú)立坐標(biāo)系，選取某一地面控制點(diǎn)為原點(diǎn)，以南北方向的兩個(gè)控制點(diǎn)所在直線向北為x軸方向，與這兩點(diǎn)所在直線垂線往東向?yàn)閥軸方向.然后在布設(shè)的導(dǎo)線控制網(wǎng)中，將全站儀架設(shè)在觀測(cè)條件良好的控制點(diǎn)上對(duì)古塔各層的每一墻面觀測(cè)均勻分布的5～6個(gè)點(diǎn)位以進(jìn)行平面擬合，并測(cè)量了塔尖三維坐標(biāo)，墻面平整度較差時(shí)需增加觀測(cè)點(diǎn)數(shù)量以減小表面方程表征墻面的誤差.部分觀測(cè)數(shù)據(jù)如表1所示.為了便于點(diǎn)位觀測(cè)數(shù)據(jù)管理，對(duì)觀測(cè)點(diǎn)位進(jìn)行了編號(hào)，其中觀測(cè)點(diǎn)編號(hào)的首位表示對(duì)應(yīng)的古塔墻面序號(hào)，中間位表示塔層，末位表示相應(yīng)塔層墻面的觀測(cè)點(diǎn)號(hào).如表1中872號(hào)觀測(cè)點(diǎn)，表示古塔第7層第8面的2號(hào)觀測(cè)點(diǎn).

表1 全站儀部分觀測(cè)數(shù)據(jù)Tab.1 Part of the observed data of total station

振風(fēng)塔各層墻面測(cè)量點(diǎn)位分布如圖2所示.

圖2 振風(fēng)塔墻面觀測(cè)點(diǎn)分布Fig.2 Distribution of measurement pointson Zhenfeng Pagoda

2 多層中軸點(diǎn)擬合變形檢測(cè)方法

本文提出的基于多層中軸點(diǎn)擬合的古塔變形檢測(cè)流程如圖3所示，包括平面擬合與棱線確定、各層棱角點(diǎn)和中軸點(diǎn)坐標(biāo)計(jì)算、變形參數(shù)計(jì)算等步驟.

圖3 基于多層中軸點(diǎn)擬合的古塔變形檢測(cè)流程圖Fig.3 Flowchart of pagoda deformation estimationbased on multilayer central points’ fitting

2.1 平面擬合與棱線確定

整體最小二乘平差(TLS)方法在參數(shù)求解中同時(shí)考慮系數(shù)矩陣和觀測(cè)向量存在的誤差.TLS是以系數(shù)矩陣和觀測(cè)向量的殘差改正數(shù)的F范數(shù)極小為約束準(zhǔn)則，求解未知參數(shù).考慮到系數(shù)矩陣包含隨機(jī)誤差，整體最小二乘平差方法在擬合精度上要優(yōu)于普通最小二乘法[6-8].其估計(jì)準(zhǔn)則如下：

(1)

對(duì)全站儀測(cè)量數(shù)據(jù)進(jìn)行處理.分別對(duì)每一層每一面的平面方程進(jìn)行擬合，由于系數(shù)矩陣由觀測(cè)值構(gòu)成，存在測(cè)量誤差，故在求解平面方程系數(shù)的最優(yōu)估值時(shí)采用整體最小二乘平差方法.得到各面的平面方程后，將相鄰墻面平面方程聯(lián)立得到對(duì)應(yīng)塔層的棱線方程，如式2所示.

(2)

式中：(a1,b1,c1)、(a2,b2,c2)分別為兩平面的法向量.則相應(yīng)的棱線方程為

(3)

通過棱線方程確定的不同高度處的棱角點(diǎn)坐標(biāo)決定著后續(xù)傾斜計(jì)算的精度.這里對(duì)x關(guān)于z的棱線方程進(jìn)行微分，分析平面系數(shù)誤差與棱角點(diǎn)坐標(biāo)誤差的關(guān)系.

x=(b1-b2)g+(c1b2-b1c2)zg

計(jì)算x關(guān)于平面方程系數(shù)的導(dǎo)數(shù)，結(jié)果如下：

(4)

(5)

同理可得x、y關(guān)于a2、b2、c2的導(dǎo)數(shù)，則有如下關(guān)系式：

(6)

將古塔墻面擬合系數(shù)代入式(6)計(jì)算導(dǎo)數(shù)值，可發(fā)現(xiàn)在高程值確定的條件下，棱角點(diǎn)坐標(biāo)對(duì)平面擬合系數(shù)的變化較為敏感，因此須保證墻面觀測(cè)點(diǎn)的精度.通過平面擬合殘差檢測(cè)粗差點(diǎn)(偏離擬合平面三倍中誤差)并剔除，以提高古塔墻面擬合的精度，從而保證棱線和棱角點(diǎn)精度.

2.2 各層棱角點(diǎn)和中軸點(diǎn)坐標(biāo)計(jì)算

得到各層相鄰墻面棱線方程后，給定高程值即可得到對(duì)應(yīng)高程面的棱角點(diǎn)坐標(biāo).實(shí)驗(yàn)中選取了各層所有觀測(cè)點(diǎn)的高程均值作為相應(yīng)塔層棱角點(diǎn)計(jì)算的高程值，然后對(duì)每層的8個(gè)棱角點(diǎn)坐標(biāo)進(jìn)行圓或橢圓擬合，計(jì)算各層中軸點(diǎn).

圓曲線方程如下：

(7)

式中：(xi,yi)為棱角點(diǎn)坐標(biāo)；(x,y)為待求中心點(diǎn)坐標(biāo)；R為半徑.

古塔發(fā)生傾斜變形后，同一高程上的各個(gè)棱角點(diǎn)會(huì)因變形而構(gòu)成一個(gè)橢圓，橢圓曲線方程如下：

Ax2+Bxy+Cy2+Dx+Ey+1=0

(8)

式中：A、B、C、D、E為系數(shù).

橢圓中心坐標(biāo)為

(9)

圓或橢圓平差擬合過程中系數(shù)矩陣包含了坐標(biāo)測(cè)量誤差項(xiàng)和常數(shù)項(xiàng)，故采用混合最小二乘法進(jìn)行參數(shù)求解，其求解原則是系數(shù)矩陣的誤差項(xiàng)和整體殘差的F范數(shù)最小,如下式所示：

(10)

式中：A2為設(shè)計(jì)矩陣誤差項(xiàng)；y為觀測(cè)量.

2.3 古塔變形參數(shù)計(jì)算

得到古塔各層中軸點(diǎn)坐標(biāo)后，對(duì)表征古塔整體傾斜度的中心軸線進(jìn)行整體最小二乘擬合[9-10]，并計(jì)算古塔變形參數(shù)來評(píng)估古塔變形狀態(tài).古塔變形參數(shù)包括古塔傾斜度、偏移量、各層傾斜度以及扭曲度等.傾斜度用古塔各層中心擬合軸線與豎直方向的夾角來表示.偏移量定義為塔尖偏離塔底中心的距離，可以分解為東西和南北兩個(gè)方向偏移分量.為了細(xì)致描述古塔變形，可以計(jì)算古塔的各層傾斜度，即每一層的幾何中心與底部中心連線向量與豎直向量的夾角.古塔的各層為多邊形結(jié)構(gòu)，通過計(jì)算各層棱角點(diǎn)構(gòu)成的多邊形之間的旋轉(zhuǎn)角度得到古塔相鄰層之間的扭曲度，進(jìn)而對(duì)其扭曲狀態(tài)進(jìn)行詳細(xì)評(píng)估.

3 結(jié)果與分析

基于古塔各層墻面觀測(cè)點(diǎn)擬合空間平面，在平面擬合計(jì)算過程中，檢測(cè)到第7層第8面中存在一個(gè)粗差觀測(cè)值，導(dǎo)致該平面擬合誤差較大，最大距離達(dá)200 mm，經(jīng)剔除后，平面擬合誤差在5 mm以內(nèi).所有有效觀測(cè)點(diǎn)到相應(yīng)擬合平面的距離分布，如圖4所示.由圖4可以看出，大部分點(diǎn)到平面的距離在1 mm以內(nèi)，3 mm以內(nèi)的觀測(cè)點(diǎn)占90%，所有有效觀測(cè)點(diǎn)到相應(yīng)擬合平面的最大距離不超過6 mm.

圖4 墻面觀測(cè)點(diǎn)到擬合平面距離分布圖Fig.4 Distribution of distance betweenobservation points and fitting plane

在求解各層相鄰墻面棱線方程的基礎(chǔ)上，將高程值帶入解算棱角點(diǎn)坐標(biāo)，進(jìn)而計(jì)算得到古塔各層中軸點(diǎn)坐標(biāo)，如表2所示.表2列出了基于圓擬合和橢圓擬合的結(jié)果及均方根誤差，可以看出，兩種擬合方法得到的各層中心坐標(biāo)幾乎一致.均方根誤差顯示橢圓擬合精度高于圓擬合精度，這也符合古塔在發(fā)生傾斜時(shí)，等高程的各角點(diǎn)的分布更接近橢圓曲線.

表2 各層中心坐標(biāo)及擬合誤差Tab.2 Center coordinates of each layer and standard deviation of fitting

利用中軸點(diǎn)及塔尖點(diǎn)坐標(biāo)進(jìn)行整體最小二乘擬合得到古塔中心軸線為

各層中軸點(diǎn)、塔尖點(diǎn)及擬合空間軸線如圖5所示.

擬合中心空間軸線向量為(0.000 6，0.000 4,1.000 0)，計(jì)算其傾斜角為2′28.74″，傾斜率約為1/1 400，此為振風(fēng)塔的整體傾斜度.由于軸線擬合的平差特性，可能會(huì)使得計(jì)算出的傾斜度偏小.利用傳統(tǒng)的塔頂、塔底中心的水平投影偏移與塔高的比值來計(jì)算整體傾斜度，其結(jié)果為7′30.2″，傾斜率約為1/460.通過第1層中軸點(diǎn)和塔尖坐標(biāo)計(jì)算古塔偏移量，南北方向?yàn)?.037 m，東西方向?yàn)?.121 m，整體向東北方向傾斜.由于塔尖單點(diǎn)測(cè)量的不準(zhǔn)確性和塔尖結(jié)構(gòu)的晃動(dòng)，塔尖所反映的變形往往誤差較大.古塔各層傾斜度及扭曲度計(jì)算結(jié)果如表3所示.

變形計(jì)算結(jié)果顯示，古塔傾斜量很小，說明基礎(chǔ)及塔體依然穩(wěn)固.古塔各層傾斜程度并不完全一致，其中第1層到第3層的傾斜度相對(duì)較大，達(dá)16′，傾斜率約為1/210.其原因有待后續(xù)進(jìn)一步觀測(cè)確定.如果后續(xù)觀測(cè)呈現(xiàn)一樣的傾斜特征，則可能該傾斜特征在該塔建造之時(shí)即存在，即在建造完第3層發(fā)現(xiàn)相對(duì)傾斜后對(duì)上層予以糾正，從而呈現(xiàn)計(jì)算得出的傾斜分布.此外，扭曲度結(jié)果顯示，相鄰塔層多邊形結(jié)構(gòu)的扭曲程度不大，均在8′以內(nèi).目前對(duì)古塔傾斜的允許值還未建立相應(yīng)的標(biāo)準(zhǔn)，應(yīng)定期對(duì)古塔進(jìn)行傾斜測(cè)量，在古塔周圍環(huán)境發(fā)生較大變化時(shí)應(yīng)當(dāng)增加觀測(cè)頻率，通過古塔觀測(cè)的序列數(shù)據(jù)來分析變形趨勢(shì)及制訂相應(yīng)的保護(hù)措施.

圖5 多層中軸點(diǎn)分布及空間軸線擬合示意圖Fig.5 Distribution of multilayer centralpoints and spatial axis fitting表3 古塔各層傾斜度及扭曲度Tab.3 Inclination and twist of each layer of pagoda

塔層傾斜度相鄰塔層扭曲度129'40.76″121'5.36″1316'5.16″237'59.68″146'44.22″347'28.77″154'18.4″455'41.06″165'11.6″567'23.70″173'5.17″675'1.34″

4 結(jié)語

為了精確檢測(cè)古塔變形，本文利用全站儀對(duì)古塔的每一層每一面進(jìn)行點(diǎn)位測(cè)量，保證每一面都有足夠數(shù)量均勻分布的測(cè)量點(diǎn)以擬合平面方程，然后擬合計(jì)算各層中軸點(diǎn)坐標(biāo)，進(jìn)而計(jì)算古塔的傾斜量、偏移量、各層傾斜度和扭曲度等變形參數(shù)，并通過振風(fēng)塔變形檢測(cè)實(shí)驗(yàn)驗(yàn)證了本文方法的有效性.相比傳統(tǒng)方法僅用塔底幾何中心和塔頂?shù)倪B線來表征古塔傾斜變形，通過多層中軸點(diǎn)計(jì)算及古塔空間軸線擬合，獲取了更加詳盡的古塔變形數(shù)據(jù)及規(guī)律，為其進(jìn)一步的修繕和維護(hù)提供了更詳細(xì)的依據(jù).因此，本文方法對(duì)于古塔等建筑物變形檢測(cè)工程實(shí)踐具有一定的指導(dǎo)和參考借鑒作用.

參考文獻(xiàn):

[1] 梁?？?古塔變形測(cè)量方法探討[J].城市勘測(cè)，2011(3)：113.

LIANG Haikui. The discussion of monitoring ancient buddhist pagoda deformation [J].Urban Geotechnical Investigation & Surveying, 2011(3): 113.

[2] 黃強(qiáng).古塔變形監(jiān)測(cè)探討[J]. 測(cè)繪與空間地理信息，2013(6)：218.

HUANG Qiang. Discussion on deformation monitoring of old towers [J].Geomatics & Spatial Information Technology, 2013(6):218.

[3] 劉琦, 惠小健, 李永新. 古塔變形問題研究[J]. 黑河學(xué)院學(xué)報(bào), 2015, 6(3):120.

LIU Qi, HUI Xiaojian, LI Yongxin. Study on the deformation of pagoda [J]. Journal of Heihe University, 2015, 6(3):120.

[4] 周千. 基于數(shù)學(xué)模型的古塔變形問題研究[J]. 計(jì)算機(jī)與數(shù)字工程, 2014, 42(8):1346.

ZHOU Qian. Study on pagoda deformation based on mathematical model [J]. Computer and Digital Engineering, 2014, 42(8):1346.

[5] TANG Yongjing, SHAO Zhendong, KUJAWSKI E. Zhenfeng tower maintenance and its durability[C]∥Proceedings of Third International Workshop on Civil Structural Health Monitoring. Vancouver: The University of British Columbia, 2010: 421-428.

[6] 王穗輝.誤差理論與測(cè)量平差[M].上海:同濟(jì)大學(xué)出版社,2010.

WANG Suihui, Error theory and surveying adjusting[M]. Shanghai: Tongji University Press, 2010.

[7] 王安怡,陶本藻. 顧及自變量誤差的回歸分析理論和方法[J].勘測(cè)科學(xué),2005(3)：29.

WANG Anyi, TAO Benzao. Theory and method of regression analysis for error of independent variable [J]. Site Investigation Science And Technology, 2005(3)：29.

[8] 魯鐵定,陶本藻,周世健.基于整體最小二乘法的線性回歸建模和解[J].武漢大學(xué)學(xué)報(bào):信息科學(xué)版,2008,33(5):504.

LU Tieding, TAO Benzao, ZHOU Shijian. Modeling and algorithm of linear regression based on total least squares [J]. Geomatics and Information Science of Wuhan University, 2008, 33(5):504.

[9] XU P L, LIU J N, SHI C. Total least squares adjustment in partial errors-in-variables model: algorithm and statistical analysis [J].Journal of Geodesy, 2012, 86(8):661.

[10] RYAN A G, MONTGOMERY D C, PECK E A,etal. Introduction to linear regression analysis, solutions manual to accompany [M].5th ed. Hoboken: Wiley, 2013.