, , ,

(1.武漢市國土資源和規(guī)劃信息中心，武漢 430014;2.湖北大學 資源環(huán)境學院，武漢 430062

0 引言

機載激光雷達(airborne laser radar，LiDAR)作為一種空間測量系統(tǒng)，其對地定位精度受到內(nèi)外在多種因素的綜合影響，并會有系統(tǒng)誤差產(chǎn)生[2-3]。國內(nèi)外學者對無控制條件下的航帶平差技術進行了大量研究：一種是基于數(shù)據(jù)驅動的無控制航帶平差，這類航帶平差將機載LiDAR中的各個誤差的綜合影響看做是航帶的各種形變，以單獨某條航帶作為基準或以多條航帶間共軛點的平均坐標作為基準，對點云數(shù)據(jù)坐標進行變換改正，從而減小相鄰航帶同名特征之間以及點云數(shù)據(jù)與對應的實際控制點之間的三維偏移[4-9]。這類方法受基準選取因素的影響較大，可能造成誤差的累積效應[10]。另一種是基于模型驅動的無控制航帶平差，這類航帶平差將機載LiDAR中最主要的誤差來源——安置角誤差等作為模型求解改正參數(shù)，以航線飛行軌跡作為基準，進行航帶的自檢校[11-15]。這類方法依賴航跡信息的POS精度，很難補償或改正由于POS系統(tǒng)導致的系統(tǒng)偏差。為此，本文借鑒基礎測繪領域擬穩(wěn)平差[16]思想，提出了一種無控制條件下的機載LiDAR數(shù)據(jù)航帶擬穩(wěn)平差方法，探討了該方法對提高點云數(shù)據(jù)質量精度的有效性。

1 算法流程

該方法基本思想是：在航帶平差前后，測區(qū)整體的重心穩(wěn)定不變，翻轉幅度為零，扭曲系數(shù)為零，即擬定測區(qū)整體的重心和形變程度“穩(wěn)定”。具體流程方法如圖1。

圖1 平差流程圖

1.1 航帶形變模型構建

(21)

式中，ΔX、ΔY為航帶的水平坐標改正參數(shù)，a、b、c、d為航帶高程在旁向上的線性改正及航向上的非線性改正參數(shù)。

對應航帶點云的改正方程為：

(2)

式中，x,y,z為該航帶點云的原始坐標，x改,y改,z改為點云經(jīng)過航帶形變模型糾正后的坐標，ΔH(U,V)=aU+bV2+cV+d為高程坐標改正，其中U、V為點云在該航帶坐標系中的坐標，如圖2所示。

圖2 航帶坐標系與地理坐標系間的關系

點云在航帶坐標系下的水平坐標(U,V)與在地理坐標系下的水平坐標(X,Y)之間的相互轉換關系如下：

(3)

式中，θ為航帶航向與Y軸正方向的夾角，x0、y0為航帶幾何中心的水平坐標。

1.2 航帶擬穩(wěn)約束條件構建

1.2.1 無控制航帶平差秩虧數(shù)

在僅存在連接特征作為觀測值且連接特征分布均勻、數(shù)量足夠的情況下，由于缺少必要起算條件，組成的觀測方程的系數(shù)矩陣B觀為秩虧矩陣：

V觀=B觀X-L觀

(4)

R(B觀)=t0

(5)

式中，t為系數(shù)矩陣列數(shù)。

即，該系數(shù)矩陣的秩虧數(shù)為6。為此，本文討論的擬穩(wěn)約束條件共有6組，用于擬穩(wěn)航帶平差。

1.2.2 擬穩(wěn)約束條件構建

(1)擬穩(wěn)約束條件①和②。

約束原理：平差后的測區(qū)整體在水平方向上的重心位移為0。

①和②以觀測方程的形式可表達為：

VΔX=ΔX1·P1+…+ΔXi·Pi+…+ΔXn·Pn

(6)

VΔY=ΔY1·P1+…+ΔYi·Pi+…+ΔYn·Pn

(7)

式中，VΔX,VΔY為測區(qū)整體偏移的觀測值，ΔXi，ΔYi為第i個航帶形變模型中的水平改正參數(shù)，Pi為第i個航帶的擬穩(wěn)權值(見1.2.3小節(jié))。

(2)擬穩(wěn)約束條件③。

約束原理：平差后的垂直方向上的重心位移為0。

③以觀測方程的形式可表達為。

VΔZ=ΔZ1·P1+…+ΔZi·Pi+…+ΔZn·Pn

(8)

式中，ΔZi為第i條航帶的平均垂直重心位移。

設該航帶寬lui，長lvi，則ΔZi可下式計算：

ΔZi·SD=?DaiU+biV2+ciV+didUV

(9)

即：

(10)

解得：

(11)

將上式代入公式(28)得整體測區(qū)在垂直方向上的重心位移的觀測方程。

(3)擬穩(wěn)約束條件④和⑤。

約束條件④原理：平差后的整個測區(qū)以其重心為中心饒Y軸負方向的旋轉系數(shù)為0，即平差后整個測區(qū)以x正方向為力臂方向，以重心為旋轉中心的高程偏移的一階矩Mx應為0。

Mx可表達為：

Mx=M1x·P1+…+Mix·Pi+…+Mnx·Pn

(12)

式中，Mix為第i條航帶的高程偏移一階矩。

仍設該航帶寬lui，長lvi，則Mix可下式計算：

Mix=?D(aiU+biV2+ciV+di)xdUV

(13)

式中，x為航帶上一點到測區(qū)重心在x方向上的距離。

將公式(23)代入上式，則有：

(cosθi·U-sinθi·V+x0i)dUdV

(14)

解得：

(15)

同理，約束條件⑤原理：平差后的整個測區(qū)以其重心為中心饒X軸正方向的旋轉系數(shù)為0，即平差后整個測區(qū)以y正方向為力臂方向，以重心為旋轉中心的高程偏移的一階矩My應為0。

My可表達為。

My=M1y·P1+…+Miy·Pi+…+Mny·Pn

(16)

式中,Miy為第i條航帶的高程偏移一階矩。仍設該航帶寬lui，長lvi。

解得：

(17)

再進一步消除旋轉系數(shù)觀測值與高程偏移一階矩之間的量綱差異，則公式(12)和(16)以觀測方程的形式分別表達為：

(18)

(19)

(4)擬穩(wěn)約束條件⑥。

約束原理：平差后的整個測區(qū)扭曲系數(shù)之合為0。

⑥以觀測方程的形式表達為。

VR=R1·P1+…+Ri·Pi+…+Rn·Pn

(20)

式中,Ri為航帶i的扭曲系數(shù)。

由于航帶產(chǎn)生形變的受航帶形變模型改正參數(shù)b的影響，仍設該航帶寬lui，長lvi，則Ri可按下式計算：

(21)

將上式代入公式(20)，則為擬穩(wěn)約束條件⑥的完整表達式。

根據(jù)擬穩(wěn)基準構建上述6組擬穩(wěn)約束條件后，即可用于解求基于無控制航帶擬穩(wěn)平差的最優(yōu)解。

1.2.3 航帶擬穩(wěn)權值設定

擬穩(wěn)權值反映了各個航帶對整體擬穩(wěn)基準的貢獻程度，可根據(jù)各個航帶獲取的點云質量酌情設置。

如果各個航帶由同一架次飛機獲取，那么，各個航帶對平差基準的貢獻權值Pi可按等權分配，即：

P1=P2=……=Pn=1

(22)

如果某條航帶在飛行時，出現(xiàn)GPS信號失鎖等突發(fā)狀況導致其點云獲取精度降低時，可將其從平差基準中剔除，即：

Pi=0

(23)

如果各個航帶飛行條件不同時，例如：不同架次、不同航高、不同設備等，可根據(jù)各航帶的具體航飛參數(shù)設定不同權值。

特殊的，如果將某一條航帶以外的其他航帶權值設為0，即：

Pj=1,Pi=0(i=1,2,…,n且i≠j)

(24)

此時上述情況轉變?yōu)楹綆Ｐ头ㄆ讲?，即以某一條航帶為基準改正其他航帶。

1.3 平差模型求解計算

(1)聯(lián)合觀測方程構建。

利用點云三維表面匹配技術得到同名特征對，構建同名特征觀測誤差方程組：

V觀=B觀X-L觀

(25)

利用擬穩(wěn)約束條件的構建策略，構建擬穩(wěn)約束觀測方程組：

V擬=B擬X-L擬

(26)

聯(lián)立上述公式(25)(26)構建聯(lián)合觀測方程組：

V=BX-L

(27)

(2)權值矩陣構建。

匹配的同名特征觀測值的權值Pi主要根據(jù)其匹配相似度估計，具體計算方法依匹配算法而定，但需使Pi∈[0,1]，而后構建同名特征的權值矩陣P觀：

(28)

擬穩(wěn)約束條件作為航帶平差中的基準，其觀測誤差的權值應遠大于同名特征匹配點的權值，可設定為10000，即:

P擬=10000E

(29)

聯(lián)合上述公式(28)(29)求得權值矩陣P。

(30)

(3)模型改正參數(shù)解算。

根據(jù)聯(lián)合觀測誤差方程(27)和權值矩陣(30)，利用最小二乘原理解算方程系數(shù)向量X：

X=(BTPB)-1BTPL

(31)

(4)迭代終止條件。

平差迭代的終止條件為：1)平差改正前后的航帶形變模型改正值的變化量足夠小，這里以實驗中的經(jīng)驗確定。2)最大迭代次數(shù)。本文航帶擬穩(wěn)平差模型收斂速度較快，設定的最大迭代次數(shù)為5。

2 實驗分析

2.1 對檢校場測區(qū)數(shù)據(jù)的實驗

實驗數(shù)據(jù)一來自中國敦煌市區(qū)4×0.8 km2大小的檢校場區(qū)域，該測區(qū)包含兩條東西向相向飛行的航帶數(shù)據(jù)，飛行軌跡近似重合，其相對航高約為800 m，測區(qū)整體概況如圖3所示。

圖3 敦煌市區(qū)檢校場點云數(shù)據(jù)

該實驗數(shù)據(jù)為系統(tǒng)檢校前的點云數(shù)據(jù)，具有明顯的系統(tǒng)誤差。主要表現(xiàn)為：

(1)航帶兩側存在較大高程差，航帶邊緣處高差可達0.60 m。

(2)航帶間在存在明顯水平偏移，東西偏移近2.2 m，南北偏移近0.8 m。

圖4 剖面截取位置與方向

圖5 平差前測區(qū)局部細節(jié)剖面圖

采用本文擬穩(wěn)平差模型進行平差解算，兩條航帶的擬穩(wěn)權值均設為1，得到2個航帶形變模型改正參數(shù)，如表1所示。

由表1可知，兩航帶水平坐標改正幅度相同，在旁向上的高程坐標改正系數(shù)一致，航向上的高程坐標基本不變。

利用改正參數(shù)對點云進行坐標改正，測區(qū)局部細節(jié)效果如下所示。

表1 各航帶形變模型改正參數(shù)值

圖6 平差后測區(qū)局部細節(jié)剖面圖

由圖6可知，航帶間的相對偏移被消除，同名地物吻合情況良好。

為進一步驗證上述平差結果對絕對精度的影響，又對該組數(shù)據(jù)做了嚴格的系統(tǒng)檢校，并以系統(tǒng)檢校后的點云作為參考點云與平差結果點云進行比對。兩者實驗結果的局部細節(jié)的比對如圖7所示。

圖7 擬穩(wěn)平差與系統(tǒng)檢校結果比對

由圖7可知，平差結果點云與參考點云吻合情況良好。對上述實驗結果進行精度評定，評定結果如下表。

分析表2可知，經(jīng)過擬穩(wěn)平差后的點云數(shù)據(jù)，航帶間的偏移被消除，相對精度大大提高；此外，平差結果與參考點云間相對偏移量極小，間接說明了擬穩(wěn)平差后點云的絕對精度也得到了提升。

表2 原始數(shù)據(jù)、航帶平差、系統(tǒng)檢校結果精度評定

2.2 對大范圍測區(qū)數(shù)據(jù)的實驗

為進一步驗證航帶擬穩(wěn)平差對大范圍測區(qū)數(shù)據(jù)的實用性，以及相較經(jīng)典航帶模型法平差在抑制誤差累積方面的有效性。選取了實驗數(shù)據(jù)二進行實驗，該數(shù)據(jù)來自于美國Niagara地區(qū)，包含來回飛行的五條航帶，相對航高約1 000 m，并有若干個地面高程控制檢查點，盡管該數(shù)據(jù)經(jīng)過檢校后數(shù)據(jù)質量已經(jīng)達標，但仍存在少量殘余系統(tǒng)誤差。測區(qū)整體概況如圖8所示。

圖8 Niagara測區(qū)點云數(shù)據(jù)

首先采用經(jīng)典航帶模型法平差，在無控制點條件下以航帶1為基準(即航帶1的改正參數(shù)設定為0)進行平差實驗。而后采用本文擬穩(wěn)平差模型進行平差解算，航帶擬穩(wěn)權值均設為1。最后，選用若干處于平坦地面(如馬路、停車場)上的高程控制點作為檢查點。對上述兩種實驗結果進行了精度評定，評定結果如表3所示。

表3 不同平差方法下的精度比較

由上表可知，采用經(jīng)典航帶模型法平差方法，盡管相對精度有所提高，但由于基準航帶為單挑航帶并存在翻滾角誤差，在平差過程中產(chǎn)生了累積效應，整體高程發(fā)生一定程度的漂移，絕對精度反而降低。而采用本文航帶擬穩(wěn)平差方法，受擬穩(wěn)約束條件的約束，航帶間的累積誤差得到很好的抑制，而且隨著相對精度的提高，絕對精度也有所提高。

3 結束語

文章歸納分析了現(xiàn)階段機載LiDAR數(shù)據(jù)航帶平差的主流方法的優(yōu)劣，并提出了本文研究方法的思路。首先介紹了本文航帶平差方法采用的航帶形變模型，然后提出了航帶擬穩(wěn)平差模型，推導了航帶擬穩(wěn)平差約束條件，最后利用兩組典型的測區(qū)數(shù)據(jù)進行了航帶擬穩(wěn)平差實驗，并與嚴格系統(tǒng)檢校結果、經(jīng)典航帶模型法平差結果進行了對比。實驗結果表明，該方法能夠有效減弱或消除機載LiDAR數(shù)據(jù)的系統(tǒng)誤差，同時抑制航帶平差過程中累積誤差的產(chǎn)生，從而提高點云數(shù)據(jù)的質量和精度。