應用激光雷達與相機信息融合的障礙物識別

黃 興，應群偉

(上海宇航系統工程研究所，上海 201109)

0 引言

近年來，NASA宣布重啟探月任務(EM)，計劃首先開始實施無人繞月飛行(EM-1階段)，技術成熟后實施載人探月(EM-2階段)，同時其火星、木星、太陽等地外天體的探測活動持續開展。俄羅斯在2007提出了自己的“探月三階段”發展戰略，預計2025年前實現載人登月，2032年前在月球建立常駐考察基地，同時計劃與歐空局合作開展火星著陸巡視探測及采樣返回[1]。日本則在小行星探測領域一枝獨秀，“隼鳥一號”、“隼鳥二號”相繼成功在小行星表面采樣返回，令人矚目[2]。歐空局、印度等機構或國家也都規劃了其地外天體探測計劃[3]，深空探測將逐漸成為未來航天領域發展的一個重要方向。

以中國探月工程與美國火星探測為代表，地外天體探測將向著陸巡視、采樣返回及載人探測返回的方向發展，軟著陸過程中對星體表面的障礙物精確識別尤為關鍵[4-5]，現有工程應用上，尚未采用信息融合進行障礙物識別。

為此，本文提出了一種融合三維點云與灰度圖像進行精確障礙物識別的方法。

1 方法概述

地外天體探測任務中，探測器軟著陸的過程一般分為主減速段、接近段和垂直著陸段[6]，如圖1所示。

圖1 星體表面軟著陸過程

本文所述方法應用在著陸器垂直著陸段，實現探測器在垂直著陸過程中提供精確障礙物識別功能，同時為后續的著陸點選取及相對導航提供基礎。

垂直著陸段一般距離星體表面在百米量級，經過姿態控制，探測器此時已經垂直于星體表面。距離與姿態都滿足導航相機與激光雷達的工作要求，導航相機開機工作，獲取探測器下方星體表面的光學圖像信息；激光雷達開機工作，獲取探測器下方，雷達視場內環境的三維激光點云數據。

圖2 方法流程圖

如圖2所示，探測器視覺導航計算機首先進行兩類傳感器數據的融合，該環節實質是獲取相機三維坐標系與激光雷達三維坐標系之間的轉換關系，從而將激光雷達三維坐標系下的空間點坐標轉換到相機三維坐標系下，最終獲得視場環境內各點的灰度信息與空間位置信息。其次，分別對灰度圖像與三維點云數據進行障礙物識別獲取光學障礙圖與點云障礙圖。最后通過雷達到相機坐標系下的轉換矩陣，將點云障礙圖與光學障礙圖進行歸一化，形成最終的障礙圖，實現精確障礙物識別。

2 三維點云數據與圖像數據的融合

本方案采用坐標歸一化的方法來完三維點云數據與灰度圖像數據的融合，通過求解激光雷達坐標系與相機坐標系的轉換關系進行數據的歸一化處理。

如圖3所示，導航相機與激光雷達的安裝位置是固定的，分別定義兩者的坐標系，則相機坐標系與激光雷達坐標系存在固定的轉換關系。

圖3 相機與激光雷達安裝示意圖

2.1 坐標系定義及相機內參

相機坐標系Oc-XYZ：相機光心為坐標原點，相機光軸指向月面為Z軸，OXY平面同相機像平面，形成右手坐標系。

激光雷達坐標系Ol-XYZ：以掃描中心為坐標原點，Z軸指向月面，X、Y軸與相機坐標系保持一致，形成右手系。

相機內參主要描述相機三維坐標系下空間點到像平面的投影矩陣，一般使用小孔成像的模型對相機的投影關系進行模擬，圖4為相機模型[7]。

圖4 相機模型

相機內參矩陣的形式一般采用齊次坐標形式，如式1所示。

(1)

其中：f為相機的焦距；dx*dy為相機感光器件每個像素對應的物理尺寸；u0、v0為對應像平面像素原點與像平面中心的偏移量，單位為像素。

相機內參矩陣為相機的固有屬性，一般出廠后能夠通過相關參數的計算獲取。工程應用上會再通過棋盤圖像進行標定，修正內參矩陣。

2.2 激光雷達坐標系轉換矩陣

本方案選用線陣掃描式激光雷達，它一般有16線、32線和64線等規格，通過安裝在平行于激光雷達Z軸的轉動軸上的16/32/64個激光收發單元進行距離計算(基于激光Time of Flight原理)，收發單元在轉動軸縱向按固定角度分布，稱為激光雷達的垂直角分辨率，轉軸的轉動角范圍為激光雷達的水平視場角[8]。

激光雷達轉動軸水平方向以一定的頻率持續旋轉，掃描計算得到的數據以數據流的形式通過網絡傳輸至視覺計算機。在現有的技術水平下，工業級的激光雷達可達到120萬點/秒的測量速度，相當于百萬級像素的相機，轉動頻率在5～30 Hz的范圍內，能夠滿足距離維度測量的實時性[9]。

圖5 激光雷達掃描示意圖

采集數據的處理上，首先根據激光點對應的激光線ID可以查詢激光線傾角ω，根據測距值r計算激光點在掃描面平面坐標系下的坐標X0。

(2)

當激光雷達進行掃描時，輸出的掃描角α為掃描面與激光雷達坐標系XOZ平面的夾角。將掃描面平面坐標轉換至激光雷達三維坐標，旋轉矩陣Rx為:

(3)

最后可得到空間測量點在激光雷達坐標系下的坐標：

(4)

通過導航相機坐標系與激光雷達坐標系的建立，可以得到相機和激光雷達兩種傳感器在各自視場內空間點的坐標。而相機與激光雷達安裝在同一水平支架的兩處，相機坐標系與激光雷達坐標系存在固定的轉換關系：

(5)

3 基于灰度圖像的障礙物識別

綜合技術特點和工程應用兩方面考慮，采用基于改進K均值聚類的光學障礙物識別算法。該算法首先通過對灰度圖片進行中值濾波預處理，濾除成像過程中產生的噪聲；然后采用K均值聚類算法獲取區分背景與非背景的最佳分離像素值，并進行二值化；最后采用邊緣檢測與輪廓提取方法，分離出障礙物，生成光學障礙物圖。

3.1 圖像預處理

通過分析發現星體表面灰度圖像中的噪聲主要包含加性噪聲和隨機干擾噪聲，為此首先采用中值濾波器對圖像預處理，濾除圖像的噪聲。同時為更好地適應光照對圖像采集的影響，采用直方圖均衡化的方式，對灰度圖像進行增強[10]。

3.2 基于改進K均值聚類的灰度圖像分割

傳統K均值聚類是基于誤差平方和最小準則對樣本數據進行聚類分析，算法對樣本的區分依據是同類樣本數據的相似性原則。本方案中，它根據灰度信息相似性度量原則對圖像進行分割，獲取能最優區分背景與非背景的分割像素值，依據背景像素值占據圖像絕大部分，而非背景物體則占據少數，同時星表灰度圖像中，背景暗，像素值偏低，非背景物體亮，像素值偏高的規律，對非背景最有可能出現的圖像灰度中值至最大值的像素范圍進行均值分類，選取相似性距離測度最大的像素值作為最優區分值，實現圖像分割。算法的收斂性好，實際效果理想[11]。

定義誤差度量函數Jk為：

(6)

mj=E(Cj)

(7)

定義相似性距離測度Lj為:

(8)

(9)

其中:Cj(j=1,2,3…k)表示k個聚類，mj為k個聚類的中心，這里采用聚類像素值的平均值作為聚類中心。

算法輸入為預處理后圖像，輸出為分割后的二值化圖。算法步驟如下。

步驟1：首先初始化k個聚類中心，本方法采用圖像平均像素值到最大像素值進行k階等分形成k個聚類中心；

步驟2：利用直方圖對平均像素到最大像素值之間的所有像素進行誤差度量計算，根據誤差度量函數Jk最小的原則，將每個像素劃分至相應的聚類中；

步驟3：根據聚類重新劃分后的結果，重新計算每個聚類的像素均值作為聚類中心，更新k個聚類中心；

步驟4：循環運行步驟3，直至更新后的k個聚類中心與更新前的k個聚類中心無變化，退出循環；

步驟5：根據式(8)計算步驟4得到的k個聚類中心的距離測度Lj，選取距離測度最大的Lx，然后根據式(9)計算最優分割像素值T；

步驟6：以最優分割像素值T為閾值對灰度圖進行二值化，得到分割二值化圖，算法結束。

3.3 邊緣檢測輪廓提取

利用3.2節得到的二值化圖像可以使背景與非背景物體進行最優分割。然后對二值化圖像進行Canny算子邊緣檢測，并將二值化圖與Canny邊緣檢測圖相加，最后對相加的圖像進行輪廓提取，對每個輪廓進行分析，根據輪廓大小及面積決定是否判定為障礙物。下圖為本算法在月貌仿真實驗室進行檢測的實際效果圖，算法效果較為理想。

圖6 基于K均值聚類的灰度圖像分割效果

4 基于激光三維點云的精確障礙物識別

激光測距能夠有效避免光照影響，實時性也能滿足垂直著陸段精確避障的要求，在星體軟著陸任務中有良好的應用前景。

激光雷達測量數據傳輸至視覺計算機后，經過坐標轉換，形成激光雷達坐標系下的三維點云。本方案中借助開源點云數據處理庫PCL(point cloud library, PCL) 對點云數據進行處理[12]，算法流程如圖7。

圖7 點云數據處理流程

4.1 點云數據預處理

對激光雷達數據進行Kalman濾波的主要作用是修正雷達的系統誤差以及探測器振動帶來的測量誤差[13]。

激光雷達由于其垂直角分辨率的原因，對應的空間分辨率會隨高度的升高的降低，在200 m的距離下，點云的空間分辨率在幾十個公分的量級，為提高其空間分辨率，依據預選著陸區地理環境一般都是平坦連續的考慮，采用插值的方法提高點云密度。

4.2 地形水平面提取

地形水平面提取是指從三維點云中提取星體表面，并從點云中移除，這樣剩下的點云數據即為高于星表的突起物或者低于星表的凹坑，以此來實現對星體表面突起物、凹坑及斜坡等障礙物的提取。

采用隨機采樣一致性算法(Random Sample Consensus，RANSAC)提取環境三維點云的地形水平面[14]。隨機采樣一致性算法是針對特定模型的最優估計算法，如平面、圓柱等模型，算法假設觀測數據集是一組包含模型局外點和局內點的集合，它通過重復迭代達到獲取包含最多局內點的模型參數之目的，算法與最小二乘法的區別在于隨機采樣一致性算法區分局內點與局外點，并只是期望找到能夠包含盡可能多局內點的具體模型；而最小二乘法則是試圖找到一種模型能夠對所有觀測數據點進行最優匹配[15]。針對環境三維點云來說，地形水平面與突起物、凹坑分別屬于不同的類型，使用最小二乘法來找尋一個平面包含所有物體是不合理的，而隨機采樣一致性算法則最有可能從包含少量突起物和凹坑的環境三維點云中獲取地形平面，并把突起物和凹坑與水平面區分開來。

算法的關鍵參數主要包括算法迭代次數kc、距離閾值Tn。其中kc用于約束算法的最大迭代次數，RANSAC算法從理論上講只有一定的概念保證找到正確的模型參數，迭代次數越高得到正確估計的概率就越大。算法采用歐式距離作為局內點與局外點的度量，對每個點計算其與每輪迭代得到的平面的歐式距離，小于距離閾值Tn屬于局內點，否則屬于局外點。

算法的具體計算流程如下。

步驟1：隨機選取點云中3個不共線樣本點(xa,xb,xc)，可以唯一地確定平面：

n·x+d=0

(10)

其中：平面法向量為n=(xb-xa)×(xc-xa)，平面截距為d=-n×xa。

步驟2：利用步驟1 中得到的平面，計算點云中每個樣本點xi到該平面的距離：

ri=(n×xi+d)2

(11)

步驟3：對步驟2計算得到的ri，如果ri≤Tn，則該點屬于局內點，否則屬于局外點，記錄本輪迭代局內點的個數Nj及模型參數(n、d)；

Nj=NUM(xi),ri≤Tn

(12)

步驟4：算法退出條件：

j=kc,迭代次數達到設定值；or

Nj=Num(xi),Num(xi)表示樣本點個數

(13)

步驟5：重復步驟1、2、3，直到滿足步驟4的算法退出條件，取Nj最大值所對應的平面參數(nbest,dbest)作為當地地形的最佳擬合水平面。

借助PCL點云庫可實現上述算法，下圖為實際效果圖(kc=1000、Tn=0.15m)。

圖8 三維點云平面提取示意圖

4.3 點云分割

提取點云中的地形水平面后，剩余點云中就只包含突起物、凹坑等障礙物。采用歐式聚類作為度量，基于距離最小原理對點云進行聚類，分割出一個個獨立的子點云。分割后的子點云根據其在地形水平面法向量方向上的投影判斷該聚類是高于地形平面還是低于地形水平面，相應地區分子點云為突起物或者為凹坑。

最后對分割后的子點云進行輪廓提取，提取各子點云的輪廓，對輪廓進行尺寸分析，結合障礙物的定義，確定該輪廓是否確定為障礙物，從而去除一些較小的、不影響探測器著陸的物體。同時，通過激光雷達坐標系與相機坐標的轉換關系可以將點云障礙物轉換至相機坐標系中，并投影至像平面上，從而實現障礙物的精確識別，生成最終障礙圖。

基于三維點云的障礙物識別可以有效避免圖像紋理信息對障礙物識別的影響，且光照影響對激光雷達數據的生成影響較小。在垂直著陸段的精確障礙物識別應用環境下，效果良好。

5 實驗結果與分析

通過以下實驗對本文算法進行了驗證。

5.1 實驗設備

固定相機與激光雷達，對室內模擬場景進行成像。

實驗結果如圖10～12。

圖10 實驗現場圖

圖11 點云及分割圖

圖12 激光雷達點云處理后生成的障礙圖

5.2 結果分析

實驗中，使用矩形紙盒、圓形鐵盒及一個易拉罐模擬障礙物。同時，實驗場景的光學圖像中包含地磚縫隙、日光燈倒影等難以使用光學障礙物算法剔除的干擾物，但點云數據可以有效避免這些干擾的影響。通過點云分割、輪廓提取及對像平面的投影，成功生成基于激光雷達點云處理后生成的障礙圖，實驗效果比較理想。

不足之處主要有幾點：

1) 16線激光雷達的角分辨率(1～2°)較低，當工作距離比較高時，存在漏檢的可能；

2) 物體尺寸較小時，點云數據稀疏，歸一化、投影過程中可能漏檢；

3) 算法的運算量大，實時性不足。

6 結束語

本文提出了一種融合激光三維點云數據與光學圖像的精確障礙物識別方法。與嫦娥三號探測器垂直著陸段采用導航避障相機+激光三維成像敏感器接力避障的方式相比具有以下優勢[16-18]：通過融合光學灰度信息與三維點云信息，可以同時獲取視場環境內各點的灰度與空間坐標位置信息，為后續的最佳落點選取和相對導航解算提供數據基礎，拓展視覺系統的導航功能，而在嫦娥三號任務中，相機和激光三維成像儀僅用作避障，不提供導航信息；此外，融合光學灰度信息與三維點云信息進行障礙物識別，可以有效解決采用單一灰度信息進行粗避障處理中圖像紋理對障礙物識別的干擾，提高障礙物識別的準確度，避免單一信息源存在的缺陷，豐富避障信息來源，提高了算法的魯棒性和可靠性。經過仿真運算，本文方法具備可行性，效果較為理想，但算法的運算量較大，對視覺計算機的運算能力要求較高。

隨著航天技術的發展，地外星體的巡視探測將越來越普遍，人類航天器走向深空的步伐也會越來越遠，對探測器的自主性要求更高，本文提出的方法可以很好地適應深空探測的應用需求。