劉存泰,鄭德華,王 浩,程宇翔,胡 創(chuàng)
(河海大學(xué) 地球科學(xué)與工程學(xué)院,江蘇 南京 210098)
三維激光掃描技術(shù)[1]被廣泛地應(yīng)用于測繪及相關(guān)領(lǐng)域,獲取的點(diǎn)云數(shù)據(jù)在拼接處理時(shí)可采用公共特征目標(biāo)實(shí)現(xiàn)多測站點(diǎn)云的空間直角變換[2]。目前公共特征目標(biāo)主要有球形目標(biāo)和平面目標(biāo)兩種類型[3],其中球形目標(biāo)因球心固定且唯一、不受激光束的位置束縛、識別距離遠(yuǎn)等優(yōu)點(diǎn),廣泛地應(yīng)用于激光點(diǎn)云配準(zhǔn)中[4]。
球形目標(biāo)球心坐標(biāo)的精度對點(diǎn)云配準(zhǔn)精度和建模與應(yīng)用的質(zhì)量有重要的影響[5]。目前,國內(nèi)外學(xué)者主要從擬合方法、模型改進(jìn)和定權(quán)方法三個(gè)方面對球形目標(biāo)球心坐標(biāo)擬合開展研究:①擬合方法的適應(yīng)性。主要有最小二乘法[1]、整體最小二乘球形目標(biāo)定位方法[8]、最小中值方差一致性估計(jì)算法[9]等,最小二乘法難以抵御粗差及異常點(diǎn)的干擾,導(dǎo)致定位精度較低;整體最小二乘方法引入系數(shù)矩陣誤差提高了擬合精度;方差一致性估計(jì)方法采用隨機(jī)采樣一致性改進(jìn)算法提高估計(jì)模型參數(shù)精度,具有結(jié)果穩(wěn)定性較低的不足。②改進(jìn)擬合模型方面主要有多項(xiàng)式擬合和幾何球面模型的方法[10]、外部綜合影響參數(shù)因子[11,12]等擬合模型,提高原有模型的抗差能力。③定權(quán)方法優(yōu)化方面主要采用鄰域定權(quán)[13]、穩(wěn)健加權(quán)[14]以及構(gòu)建大型稀疏矩陣[17]等方法,得到較整體最小二乘方法更高的精度和效率。
本文針對球形目標(biāo)掃描點(diǎn)云中存在的偶然誤差、系統(tǒng)誤差和粗差三類基本誤差影響精確配準(zhǔn)[18]的問題,設(shè)計(jì)球形目標(biāo)激光掃描實(shí)驗(yàn);在充分了解濾波算法的基礎(chǔ)上[19],結(jié)合統(tǒng)計(jì)濾波較半徑濾波及直通濾波[20]方法在保留原始點(diǎn)云細(xì)節(jié)幾何特征信息和避免過度去噪[21]方面的優(yōu)勢,設(shè)計(jì)基于殘差分布定權(quán)的M估計(jì)球形目標(biāo)點(diǎn)云擬合方法,以解決不同掃描距離條件下采集球形目標(biāo)點(diǎn)云三類誤差變化的球形目標(biāo)擬合的抗差能力和適用性問題。
球形目標(biāo)點(diǎn)云具有各方向分布均勻、表面特征易識別等特點(diǎn),由于各站掃描均只能得到半球點(diǎn)云,且邊緣部分點(diǎn)云誤差較大,因此,噪聲剔除環(huán)節(jié)尤為重要。針對偏離球形目標(biāo)主體點(diǎn)云較遠(yuǎn)噪點(diǎn)稀疏且不均勻的空間分布特征,采用統(tǒng)計(jì)濾波方法對其進(jìn)行剔除。統(tǒng)計(jì)濾波是對輸入的點(diǎn)云數(shù)據(jù)的查詢點(diǎn)與鄰域點(diǎn)集之間的距離進(jìn)行統(tǒng)計(jì)分析,將目標(biāo)點(diǎn)云中粗差和較大的系統(tǒng)性誤差作為噪聲濾除。首先,遍歷點(diǎn)云,計(jì)算每個(gè)點(diǎn)與其最近的K個(gè)鄰域點(diǎn)之間的平均距離;其次,計(jì)算所有平均距離的均值μ0與標(biāo)準(zhǔn)差σ0,設(shè)置距離閾值R進(jìn)行判別;最后,再次遍歷點(diǎn)云,剔除與K個(gè)鄰域點(diǎn)的平均距離大于R的點(diǎn)。算法具體原理如下:
首先,對點(diǎn)云中所有數(shù)據(jù)點(diǎn)si(Xi,Yi,Zi),其中,i=1, 2, 3, …,n(n為點(diǎn)云數(shù)量)進(jìn)行距離統(tǒng)計(jì)分析,則si到其任一鄰近點(diǎn)hm(Xm,Ym,Zm)的距離rm為:
(1)
按式(2)計(jì)算si至其所有鄰近點(diǎn)的平均距離r0:
(2)
該平均距離的集合滿足正態(tài)分布,由數(shù)學(xué)期望μ0和標(biāo)準(zhǔn)差σ0決定,該分布的概率密度函數(shù)為:
(3)
則可按式(4)計(jì)算平均距離閾值R:
R=μ0+σ0·std
(4)
式(4)中,std為標(biāo)準(zhǔn)差倍數(shù)閾值。若點(diǎn)si到其鄰域內(nèi)的平均距離小于R,記其為有效點(diǎn),反之為噪聲點(diǎn),予以剔除。針對不同掃描距離下std的取值效果進(jìn)行實(shí)驗(yàn)設(shè)計(jì),對比分析不同取值條件下的濾波效果,得到最適合該距離的std值作為最優(yōu)取值。
統(tǒng)計(jì)濾波流程如圖1所示。
圖1 統(tǒng)計(jì)濾波流程Fig.1 Statistical filtering flow chart
在統(tǒng)計(jì)濾波的基礎(chǔ)上,球形目標(biāo)點(diǎn)云中較為明顯的噪聲得到去除,但由于系統(tǒng)誤差的存在,一些小尺度的噪聲點(diǎn)依然不能完全去除干凈,通過改善數(shù)學(xué)模型的方法減小系統(tǒng)誤差,由此引入M估計(jì)擬合模型,進(jìn)一步提高球形目標(biāo)點(diǎn)云質(zhì)量。
M估計(jì)是Huber提出的一種穩(wěn)健估計(jì)方法,該方法將殘差平方和用另一種關(guān)于殘差的函數(shù)替代,其目標(biāo)函數(shù)為:
(5)
其中,ri為第i點(diǎn)的殘差值,對于一組數(shù)量為n的點(diǎn)云(xi,yi,zi),其中,i=1,2,…,n。誤差方程式為
(6)
其中,vi為掃描觀測時(shí)的偶然誤差,(xi,yi,zi)為球形目標(biāo)的表面觀測數(shù)據(jù),(x0,y0,z0)為球心三維坐標(biāo),r為球形目標(biāo)半徑。將上式經(jīng)過線性化處理后可得誤差方程模型為:
(7)
其中,
根據(jù)球形目標(biāo)點(diǎn)云擬合殘差繪制頻數(shù)直方圖,如圖2所示。無粗差的點(diǎn)云,其擬合殘差為正態(tài)分布特征,可得到殘差分布的數(shù)學(xué)期望μ和標(biāo)準(zhǔn)差σ。按式(8)確定權(quán)函數(shù)ωi,權(quán)函數(shù)ωi形式如圖3所示。
(8)
圖2 殘差頻數(shù)分布直方圖Fig.2 Residual frequency distribution histogram
圖3 權(quán)函數(shù)圖像Fig.3 Weight function image
圖4 M估計(jì)算法流程Fig.4 M-estimation algorithm flow chart
式(8)中,σ為殘差標(biāo)準(zhǔn)差。該權(quán)函數(shù)將殘差小于2.5σ的點(diǎn)權(quán)重賦為1;若某點(diǎn)殘差在2.5σ~3.75σ之間,則根據(jù)其值大小進(jìn)行定權(quán);若某點(diǎn)殘差大于3.75σ,該點(diǎn)不參與參數(shù)估計(jì)。
則球形目標(biāo)擬合參數(shù)的解為:
X=(BTωiB)-1BTωiL
(9)
由于權(quán)函數(shù)確定要根據(jù)估計(jì)參數(shù)的殘差迭代計(jì)算得到,因此,設(shè)計(jì)基于殘差分布定權(quán)的M估計(jì)擬合算法流程,如圖4所示。
為了分析不同掃描距離的球形目標(biāo)點(diǎn)云中3類誤差對擬合精度的影響,驗(yàn)證本文方法擬合精度和抗差能力,實(shí)驗(yàn)設(shè)計(jì)了一個(gè)移動框架結(jié)構(gòu)裝置,3個(gè)原裝進(jìn)口球形目標(biāo)設(shè)置在支撐框架位置,如圖5所示。實(shí)驗(yàn)數(shù)據(jù)采集采用德國的Z+F IMAGER 5016三維激光掃描儀,激光波長為1 550 nm。3個(gè)球形目標(biāo)半徑值為72.5 mm,編號從左到右依次為1號、2號、3號。實(shí)驗(yàn)選取60 m范圍實(shí)驗(yàn)場地,掃描距離從10 m到60 m,每隔10 m采集一次球形目標(biāo)點(diǎn)云,共獲取了3個(gè)球形目標(biāo)在6種等間隔掃描距離采集的點(diǎn)云數(shù)據(jù),如圖6所示。
圖5 實(shí)驗(yàn)裝置Fig.5 Experimental device
圖6 實(shí)驗(yàn)整體點(diǎn)云Fig.6 The overall point cloud of the experiment
實(shí)驗(yàn)所得1號、2號、3號球形目標(biāo)點(diǎn)云,如圖7所示。
圖7 各球形目標(biāo)點(diǎn)云Fig.7 Each spherical target point cloud
針對掃描獲取的3個(gè)球形目標(biāo)點(diǎn)云,分別采用本文方法、最小中值方差一致性估計(jì)算法(Least Median of Squares,LMedS)以及最小二乘法擬合處理,3個(gè)球形目標(biāo)半徑擬合誤差與掃描距離關(guān)系如圖8及表1所示。
圖8 3種方法半徑擬合誤差比較Fig.8 Comparison of radius fitting errors of three methods
由圖8及表1可知,以1號球、掃描距離10 m為例,本文方法半徑擬合誤差較LMedS算法減小0.08 mm,較最小二乘法減小0.14 mm;以1號球、掃描距離40 m為例,本文方法半徑擬合誤差較LMedS算法增加0.15 mm,較最小二乘法減小1.43 mm。3個(gè)球形目標(biāo)在10~60 m掃描距離條件下,本文方法半徑擬合誤差均小于最小二乘法半徑擬合誤差;在10~30 m掃描距離條件下,本文方法半徑擬合誤差均小于LMedS算法半徑擬合誤差,在40~60 m掃描距離條件下,1號球和3號球采用LMedS算法所得半徑擬合誤差小于本文方法。
針對掃描距離40 m之后半徑擬合誤差較大的情況,本文對掃描點(diǎn)云進(jìn)行統(tǒng)計(jì)濾波,鄰域點(diǎn)數(shù)目(K)設(shè)為50,標(biāo)準(zhǔn)差閾值(std)的選取根據(jù)掃描距離確定,各距離最優(yōu)取值見表2;以掃描距離10 m時(shí)獲取的點(diǎn)云為例,得到濾波前后點(diǎn)云殘差分布情況,如圖9、圖10所示。
由圖9、圖10可知,在掃描距離為10 m的條件下,統(tǒng)計(jì)濾波對于殘差在±1.5 mm以上的噪聲點(diǎn)有明顯的過濾效果。對去噪后的點(diǎn)云采用本文方法進(jìn)行球形目標(biāo)擬合,3個(gè)球形目標(biāo)去噪前后的半徑擬合誤差與掃描距離關(guān)系如圖11所示。
圖11 濾波前后半徑擬合誤差Fig.11 Radius fitting error before and after filtering
由圖11可知,以1號球、掃描距離10 m為例,濾波后半徑擬合誤差較濾波前減小0.03 mm;以1號球、掃描距離40 m為例,濾波后半徑擬合誤差較濾波前減小1.21 mm。3個(gè)球形目標(biāo)在10~60 m掃描距離條件下,濾波后半徑擬合誤差均小于濾波前半徑擬合誤差,且隨掃描距離增加,減小幅度逐漸變大。
根據(jù)點(diǎn)云擬合殘差計(jì)算中誤差,建立實(shí)驗(yàn)中3個(gè)球形目標(biāo)點(diǎn)云擬合中誤差與掃描距離關(guān)系,如圖12及表3所示。
圖12 濾波前后點(diǎn)云擬合中誤差Fig.12 Error in point cloud fitting before and after filtering
由圖12及表3可知,以1號球、掃描距離10 m為例,濾波后點(diǎn)云擬合中誤差較濾波前減小0.06 mm,精度提高15.4 %;以1號球、掃描距離40 m為例,濾波后點(diǎn)云擬合中誤差較濾波前減小0.27 mm,精度提高24.8 %。3個(gè)球形目標(biāo)在10~60 m掃描距離條件下,濾波后點(diǎn)云擬合中誤差均小于濾波前點(diǎn)云擬合中誤差。
1)針對現(xiàn)有球形目標(biāo)球心擬合方法精度不足的問題,提出基于殘差分布定權(quán)的M估計(jì)球形目標(biāo)點(diǎn)云擬合方法,半徑擬合誤差較LMedS算法減小21.1 %,較最小二乘法減小48.1 %。本文方法既克服了最小二乘法抗粗差能力不足的缺點(diǎn),又解決了LMedS算法穩(wěn)健性不足的問題。
2)引入統(tǒng)計(jì)濾波剔除粗差,濾波后球形目標(biāo)半徑擬合誤差較濾波前減小49.6 %,點(diǎn)云擬合精度提高24.8 %。統(tǒng)計(jì)濾波的引入,提高了球形目標(biāo)點(diǎn)云擬合精度,特別是改善了40~60 m掃描距離球形目標(biāo)點(diǎn)云擬合精度稍低的不足。
3)設(shè)計(jì)的等距離間隔條件下掃描球形目標(biāo)實(shí)驗(yàn),獲取了3個(gè)球形目標(biāo)在不同掃描距離下的點(diǎn)云數(shù)據(jù),通過不同掃描距離條件下球形目標(biāo)點(diǎn)云的半徑擬合精度對比,驗(yàn)證了本文方法的穩(wěn)健性和適用性。