李榮華，李 恒，林婷婷，王 蒙

(1. 大連交通大學(xué)機(jī)械工程學(xué)院，大連 116028；2. 上海宇航系統(tǒng)工程研究所，上海 201109)

0 引 言

在軌服務(wù)技術(shù)和空間碎片清除技術(shù)是太空探測(cè)技術(shù)的重要組成部分，如何在復(fù)雜空間環(huán)境下有效保護(hù)空間軌道資源，進(jìn)行衛(wèi)星維修、壽命延長(zhǎng)及垃圾清理是目前亟需解決的關(guān)鍵問(wèn)題[1]。要解決上述問(wèn)題，必須獲得空間非合作目標(biāo)運(yùn)動(dòng)參數(shù)及三維形貌信息，并通過(guò)質(zhì)量評(píng)估方法驗(yàn)證三維形貌精度，只有符合精度要求，服務(wù)航天器才能接近目標(biāo)航天器實(shí)施在軌服務(wù)[2]。

目前在對(duì)空間非合作目標(biāo)進(jìn)行三維形貌恢復(fù)時(shí)，重建方法主要分為兩種：第一種是采用序列圖像的重建方式[3-4]。第二種是采用激光雷達(dá)設(shè)備獲取精度較高的三維數(shù)據(jù)，再經(jīng)過(guò)目標(biāo)運(yùn)動(dòng)特性分析、點(diǎn)云預(yù)處理、點(diǎn)云配準(zhǔn)、三維重建等技術(shù)手段，獲得三維模型[5-6]。激光雷達(dá)突破了傳統(tǒng)的成像概念和模式，具有極高速度分辨率、且工作距離長(zhǎng)、受光照條件影響小，可以快速獲取目標(biāo)某個(gè)視角的三維空間信息[7]。與傳統(tǒng)的幾何建模相比，基于激光雷達(dá)的三維模型重建方法更加精確、高效，應(yīng)用領(lǐng)域也更加廣泛[8]。

三維重建完成后，重建點(diǎn)云是否可用于控制服務(wù)航天器對(duì)目標(biāo)航天器進(jìn)行在軌服務(wù)，三維模型是否可用于目標(biāo)航天器的狀態(tài)和關(guān)鍵部位識(shí)別，進(jìn)而引導(dǎo)平臺(tái)接近服務(wù)目標(biāo)等,都是未知的[9]?？梢?jiàn)，對(duì)重建結(jié)果進(jìn)行質(zhì)量評(píng)估是非常重要的環(huán)節(jié)。Javaheri等[10]通過(guò)選定的客觀質(zhì)量指標(biāo)與人類質(zhì)量評(píng)估的相關(guān)性進(jìn)行評(píng)價(jià)。陳鳳等[11]對(duì)基于序列圖像獲得的三維重建結(jié)果提出點(diǎn)云精度、重建目標(biāo)識(shí)別完整性和吻合度等評(píng)估方法。Cheng等[12]提出一種雙目標(biāo)卷積神經(jīng)網(wǎng)絡(luò)，可以估計(jì)圖像的退化類型和質(zhì)量評(píng)分。Arora等[13]提出一種模型中的圖像質(zhì)量評(píng)估方法，包括水平和垂直質(zhì)量評(píng)估以及對(duì)角線和反對(duì)角線質(zhì)量評(píng)估的混合。于康龍等[14]基于結(jié)構(gòu)相似度和尺度空間理論，提出一種針對(duì)超分辮率重建圖像的弱參考質(zhì)量評(píng)價(jià)算法。田陽(yáng)等[15]對(duì)三維重構(gòu)誤差傳播機(jī)理進(jìn)行分析，給出了三維重建誤差計(jì)算方法。文獻(xiàn)[16]從噪聲、密度、完整性、準(zhǔn)確性等方面對(duì)采集的點(diǎn)云數(shù)據(jù)進(jìn)行質(zhì)量評(píng)價(jià)。陳西江等[17]將誤差橢球引入到點(diǎn)云精度分析中，重點(diǎn)分析了誤差橢球與點(diǎn)位協(xié)方差的關(guān)系。

上述評(píng)估方法多數(shù)都是針對(duì)圖像序列三維重建或有已知模型進(jìn)行對(duì)比分析，而對(duì)于線陣激光雷達(dá)測(cè)量未知構(gòu)型的空間非合作目標(biāo)重建質(zhì)量評(píng)估研究相對(duì)較少，還沒(méi)有統(tǒng)一的評(píng)判標(biāo)準(zhǔn)及評(píng)估方法。如何基于激光雷達(dá)測(cè)量的空間非合作目標(biāo)未知構(gòu)型三維重建結(jié)果建立有效而又簡(jiǎn)便的質(zhì)量評(píng)估體系，是亟需解決的關(guān)鍵問(wèn)題。因此，本文提出空間非合作目標(biāo)未知構(gòu)型重建質(zhì)量評(píng)估方法。該方法著重解決空間非合作目標(biāo)未知構(gòu)型重建質(zhì)量評(píng)估問(wèn)題，確保服務(wù)航天器能夠精確實(shí)施在軌服務(wù)。

1 重建質(zhì)量評(píng)估流程

基于線陣激光雷達(dá)測(cè)量的空間非合作目標(biāo)三維重建結(jié)果直接影響在軌服務(wù)質(zhì)量。由于目標(biāo)是非合作的，將無(wú)法利用重建模型與真實(shí)模型進(jìn)行對(duì)比來(lái)驗(yàn)證重建算法是否達(dá)到預(yù)期效果。因此，針對(duì)非合作目標(biāo)三維重建結(jié)果，必須建立合理、完善的評(píng)估體系，確保在軌服務(wù)能夠順利進(jìn)行。因此，本文從重建點(diǎn)云密度、重建幾何性質(zhì)、重建表面完整度以及建立的多因素綜合評(píng)估對(duì)空間非合作目標(biāo)重建結(jié)果進(jìn)行分析，獲得對(duì)重建結(jié)果的滿意度，其流程如圖1所示。

圖1 質(zhì)量評(píng)估流程圖Fig.1 Quality assessment flowchart

2 3σ優(yōu)化擬合

這里主要認(rèn)定基于激光雷達(dá)測(cè)量獲得的重建點(diǎn)云模型表面多數(shù)為空間平面與直線，再結(jié)合人眼視覺(jué)特性判斷重建點(diǎn)云中的平面與直線位置，構(gòu)建數(shù)學(xué)模型，即在指定表面或直線上選取一定數(shù)量特征點(diǎn)，利用這些特征點(diǎn)進(jìn)行擬合；然后利用3σ法則不斷循環(huán)去除效果不好的數(shù)據(jù)點(diǎn)，優(yōu)化擬合結(jié)果。

2.1 空間平面擬合

根據(jù)三維重建結(jié)果，選取模型平面點(diǎn)云數(shù)據(jù)(Xi,Yi,Zi),i=1,2,3,…,N，進(jìn)行區(qū)域均勻采樣，獲得部分采樣點(diǎn)(xi,yi,zi),i=1,2,3,…,n。依據(jù)采樣點(diǎn)擬合出最優(yōu)平面，使得每一個(gè)采樣點(diǎn)到平面的距離的平方和最小，其流程如圖2所示。

圖2 平面擬合流程圖Fig.2 Flow chart of plane fitting

平面方程的一般表達(dá)式為：

Ax+By+Cz+D=0,(D≠0)

(1)

兩邊除以D得到：

則有：

ajx+bjy+cjz+1=0,j=1,2,…,n

式中：

所有采樣點(diǎn)到平面距離的平方和為：

(2)

只有當(dāng)S最小時(shí)，才能獲得最優(yōu)平面，因此應(yīng)滿足：

(3)

化簡(jiǎn)整理：

采用克拉默法則求解線性方程組得到aj,bj,cj，即求得平面方程為：

ajx+bjy+cjz+1=0,j=1,2,…,n

(4)

圖3 正態(tài)分布示意圖Fig.3 Schematic diagram of normal distribution

在圖3正態(tài)分布中，很容易發(fā)現(xiàn)：

P(μ-3σ

意味著在(μ-3σ，μ+3σ]以外的概率不到0.3%，幾乎不可能發(fā)生，稱為小概率事件。因此，要精確得到平面方程就必須將上述獲得的平面方程采用3σ法則循環(huán)剔除誤差較大的點(diǎn)，直到?jīng)]有瑕疵點(diǎn)參與擬合平面為止，獲得最優(yōu)空間平面擬合結(jié)果。

2.2 空間直線擬合

首先進(jìn)行閾值設(shè)定，選取直線點(diǎn)云數(shù)據(jù)(Xj,Yj,Zj),j=1,2,3,…,N，再?gòu)闹羞M(jìn)行區(qū)域均勻采樣獲得部分采樣點(diǎn)(xj,yj,zj),j=1,2,3,…,n。依據(jù)采樣點(diǎn)進(jìn)行空間直線擬合。

設(shè)空間直線方程為：

(5)

整理獲得XOZ與YOZ平面的直線射影方程：

(6)

式中：

因此，空間直線可以看作上述兩個(gè)平面相交的直線，如圖4所示，分別對(duì)2個(gè)方程進(jìn)行數(shù)據(jù)擬合。設(shè)式(7)表示按擬合方程x=az+b求得的近似值。一般地，它不同于實(shí)測(cè)值xi，兩者之差由式(8)所示：

(7)

(8)

同理得：

(9)

圖4 空間直線示意圖Fig.4 Schematic diagram of space line

當(dāng)滿足下列方程時(shí)Q值最小，即可獲得方程的系數(shù)a,b,c,d的值。

(10)

令：

(11)

則方程組可寫成：

FFTA=FX,FFTB=FY

式中：

根據(jù)n組數(shù)據(jù)點(diǎn)解方程組獲得a,b,c,d的值，即求得直線方程。同理，采用上述3σ法則循環(huán)剔除誤差較大的點(diǎn)，直到?jīng)]有瑕疵點(diǎn)參與擬合直線為止，獲得最優(yōu)空間直線擬合結(jié)果。

3 空間非合作目標(biāo)重建質(zhì)量評(píng)估

空間非合作目標(biāo)三維形貌結(jié)果是否滿足精度要求至關(guān)重要。因此，本文建立三個(gè)評(píng)估指標(biāo)以及多因素綜合評(píng)估對(duì)重建結(jié)果進(jìn)行評(píng)價(jià)，驗(yàn)證空間非合作目標(biāo)重建精度及滿意度。

3.1 重建點(diǎn)云密度

選取重建點(diǎn)云平面分布范圍為M×N，點(diǎn)云數(shù)量為Q，將其均勻劃分成m×n塊，用Kd表示含有點(diǎn)的分塊個(gè)數(shù)，用Kt表示總的分塊數(shù)，Md表示點(diǎn)云的密集程度。

(12)

(13)

令Nc為分割小塊的點(diǎn)云總數(shù)量，Dij為塊內(nèi)點(diǎn)云密度，如圖5所示。其中i,j為整體點(diǎn)云第i行j列的小分塊，設(shè)小分塊正方形的邊長(zhǎng)為L(zhǎng)，可由式(14)求得。因此小分塊正方形的面積Sc為：

(14)

(15)

且正方形的面積也可以由點(diǎn)云數(shù)據(jù)的X、Y坐標(biāo)值范圍求得：

Sc=(ymax-ymin)×(xmax-xmin)

(16)

因此可得塊內(nèi)點(diǎn)云密度的計(jì)算公式為：

圖5 分塊點(diǎn)云密度計(jì)算Fig.5 Block point cloud density calculation

其中:Md對(duì)整體點(diǎn)云密度也有影響?？赡艽嬖邳c(diǎn)云數(shù)量多，但Md較小，意味著有較多分塊不存在點(diǎn)云，導(dǎo)致局部密度很大，因此本文密度計(jì)算公式如下所示:

(17)

Dden受分塊閾值影響，因此，本文采用先確定閾值，再自動(dòng)確定m與n值。同時(shí)密度的標(biāo)準(zhǔn)為D標(biāo)準(zhǔn)，實(shí)際應(yīng)用中應(yīng)考慮實(shí)際情況加以選用。

3.2 重建幾何性質(zhì)

采用2.1中平面擬合結(jié)果，分析平面與平面之間的平行度與垂直度關(guān)系。由于測(cè)量存在誤差，導(dǎo)致獲得的目標(biāo)示意圖如圖6所示:

圖6 目標(biāo)示意圖Fig.6 Schematic diagram of the target

則面與面之間夾角為：

(18)

采用2.2中直線擬合結(jié)果，分析直線與直線之間的平行度與垂直度關(guān)系。則直線與直線之間夾角為：

(19)

3.3 重建表面完整度

通過(guò)下式計(jì)算出模型平面所有點(diǎn)云數(shù)據(jù)到擬合平面距離：

di=ajXi+bjYi+cjZi+1,i=1,2,…,N,j=1,2,…,n

(20)

將平面沿著Z方向投影，則殘差值的極大值為：

max|zi-(-1-a0xi-a1yi)/a2|,i=0,1,2,…,n-1

(21)

當(dāng)前觀測(cè)點(diǎn)到擬合平面的距離小于各觀測(cè)點(diǎn)與擬合平面殘差絕對(duì)值的極值,則當(dāng)前觀測(cè)點(diǎn)在選擇的表面上重構(gòu)，即為表面重構(gòu)點(diǎn)。

(22)

3.4 多因素綜合評(píng)估

假設(shè)將最優(yōu)結(jié)果100%對(duì)應(yīng)數(shù)學(xué)描述為5，即可通過(guò)計(jì)算反求得點(diǎn)云密度、幾何性質(zhì)以及表面完整度對(duì)應(yīng)的數(shù)字S密度，S幾何，S完整度。

S總=S密度×γ+S幾何×η+S完整度×μ

(23)

其中：γ+η+μ=1

采用偏大型柯西分布和對(duì)數(shù)函數(shù)作為隸屬函數(shù)，如下式所示：

(24)

其中:α,β,a,b為待定系數(shù)：

將滿意度劃分為5個(gè)等級(jí)，分別對(duì)應(yīng)數(shù)學(xué)描述為5、4、3、2、1。

當(dāng)對(duì)評(píng)估結(jié)果“很滿意”時(shí)，隸屬度設(shè)為100%，即f(5)=100%；

當(dāng)對(duì)評(píng)估結(jié)果“較滿意”時(shí)，隸屬度設(shè)為80%，即f(3)=80%；

當(dāng)對(duì)評(píng)估結(jié)果“很不滿意”時(shí)，隸屬度設(shè)為1%，即f(1)=1%。

將上述結(jié)果代入上式即可求得：

所以隸屬函數(shù)為：

將評(píng)估結(jié)果S總代入上式即可獲得多因素評(píng)估結(jié)果滿意度。

4 目標(biāo)重建實(shí)驗(yàn)評(píng)估結(jié)果

4.1 目標(biāo)三維重建結(jié)果

1)定義非合作目標(biāo)坐標(biāo)系、平臺(tái)雷達(dá)坐標(biāo)系、世界坐標(biāo)系，明確其之間的關(guān)系；建立可測(cè)點(diǎn)云提取機(jī)制，精確提取某一視場(chǎng)角數(shù)據(jù)；依據(jù)線陣激光雷達(dá)工作原理與掃描機(jī)制，采用MATLAB模擬空間環(huán)境，建立仿真系統(tǒng)，實(shí)現(xiàn)激光雷達(dá)環(huán)繞目標(biāo)做360°繞飛跟蹤掃描，如圖7所示；再利用坐標(biāo)系之間的轉(zhuǎn)換關(guān)系，獲得N幀可測(cè)點(diǎn)云數(shù)據(jù)，如圖8所示。

圖7 雷達(dá)采集數(shù)據(jù)示意圖Fig.7 Schematic diagram of radar acquisition data

圖8 可測(cè)點(diǎn)云數(shù)據(jù)Fig.8 Measurable point cloud data

2)將上述獲得的可測(cè)點(diǎn)云數(shù)據(jù)進(jìn)行粗配準(zhǔn)，再經(jīng)過(guò)點(diǎn)云聚類分割后，采用最近迭代點(diǎn)(Iterative closest point，ICP)進(jìn)行精配準(zhǔn)，即可獲得旋轉(zhuǎn)與平移矩陣,其結(jié)果如圖9所示。

圖9 點(diǎn)云配準(zhǔn)Fig.9 Point cloud registration

3)根據(jù)ICP點(diǎn)云配準(zhǔn)方法得到的差分掃描點(diǎn)云數(shù)據(jù)之間的位姿增量關(guān)系，采用逆序三維重建方式，獲得重建點(diǎn)云模型[6],如圖10所示。

圖10 三維重建結(jié)果Fig.10 Three-dimensional reconstruction result

4.2 擬合結(jié)果

根據(jù)要求，分別提取主體點(diǎn)云數(shù)據(jù)的四個(gè)面，如圖11所示，進(jìn)行評(píng)估實(shí)驗(yàn)。

圖11 主體的四個(gè)平面點(diǎn)云Fig.11 Four plane point clouds of the main body

采用2.1中平面擬合方法，獲得4個(gè)面的法向量為：

采用2.2中直線擬合方法，獲得圖12中4條直線的方向向量,見(jiàn)表2。

4.3 質(zhì)量評(píng)估結(jié)果

1)點(diǎn)云密度

本文根據(jù)需求D標(biāo)準(zhǔn)=0.02 m，設(shè)分塊大小為200 mm×200 mm，將主體各平面點(diǎn)云劃分成10×15塊小正方形，各平面點(diǎn)云數(shù)量為6864、11233、8814、9021，通過(guò)計(jì)算獲得各面點(diǎn)云密度見(jiàn)表3：

表1 各平面法向量Table 1 Normal vector of each plane

圖12 直線點(diǎn)云Fig.12 Straight line point cloud

表2 各直線的方向向量Table 2 Direction vector of each line

表3 各平面點(diǎn)云密度Table 3 Point cloud density of each plane

2)幾何性質(zhì)

根據(jù)擬合的直線公式，獲得圖13的直線擬合結(jié)果,再計(jì)算出直線之間的平行度與垂直度，見(jiàn)表4與表5。根據(jù)擬合的平面公式，獲得圖14的平面擬合結(jié)果。

圖13 直線擬合結(jié)果Fig.13 Straight line fitting result

表4 直線平行度Table 4 Straight line parallelism

表5 直線垂直度Table 5 Straight line perpendicularity

圖14 平面擬合結(jié)果Fig.14 Plane fitting result

計(jì)算出平面之間的平行度與垂直度，見(jiàn)表6與表7：

表6 平面平行度Table 6 Plane parallelism

3)表面完整度

進(jìn)行重構(gòu)點(diǎn)篩選后，將平面1與2中重構(gòu)點(diǎn)數(shù)量沿著x軸方向投影與z軸方向投影，結(jié)果如圖15(a)與圖15(b)所示；將平面3與4中重構(gòu)點(diǎn)數(shù)量沿

表7 平面垂直度Table 7 Plane perpendicularity

著y軸方向與z軸方向投影，結(jié)果如圖15(c)與圖15(d)所示。各平面的點(diǎn)云總數(shù)量如表8所示。再依據(jù)表面完整度計(jì)算公式,求得平面1與2的表面完整度如表9所示,平面3與4的表面完整度如表10所示。

圖15 重構(gòu)點(diǎn)云數(shù)量投影折線圖Fig.15 Reconstructed point cloud number projection line chart

表8 各平面重構(gòu)點(diǎn)數(shù)Table 8 Number of reconstruction points for each plane

4)多因素綜合評(píng)估結(jié)果

根據(jù)獲得的點(diǎn)云密度、幾何性質(zhì)及表面完整度，求得對(duì)應(yīng)的數(shù)字描述為：

表9 平面1與2表面完整度Table 9 Surface integrity of plane 1 and 2

表10 平面3與4表面完整度Table 10 Surface integrity of plane 3 and 4

S總=(S密度+S幾何+S完整度)/3=4.467275

將S總代入式(24)中得到滿意度為：

f(x)=95.5889%

5 結(jié) 論

綜上所述，本研究針對(duì)空間非合作目標(biāo)重建質(zhì)量評(píng)估問(wèn)題，利用3σ優(yōu)化法則進(jìn)行平面與直線擬合，從點(diǎn)云密度、幾何性質(zhì)與表面完整度三個(gè)方面進(jìn)行評(píng)估，再進(jìn)行多因素綜合評(píng)估獲得滿意度。仿真結(jié)果表明，該重建模型點(diǎn)云密度達(dá)到標(biāo)準(zhǔn)的概率為97.41%，表明重建結(jié)果較好；該重建模型幾何性質(zhì)達(dá)到標(biāo)準(zhǔn)的概率為93.72%，重建幾何性質(zhì)較好；該重建模型表面完整度達(dá)到標(biāo)準(zhǔn)的概率為76.9075%，是由于點(diǎn)云分布的比較隨機(jī)，因此重建表面完整度較差一些；多因素綜合評(píng)估滿意度為95.5889%，表明整體重建結(jié)果較好，符合要求。該研究不僅能夠驗(yàn)證基于線陣激光雷達(dá)測(cè)量空間非合作目標(biāo)重建結(jié)果的準(zhǔn)確性，且建立的多因素綜合評(píng)估方法還為未知構(gòu)型評(píng)估提供一種新的研究方向。