趙新，李杰，岳丹丹，韓重陽，唐婷，吳偉斌

1.華南農(nóng)業(yè)大學(xué)工程學(xué)院，廣州 510642；2.華南農(nóng)業(yè)大學(xué)南方農(nóng)業(yè)機械與裝備關(guān)鍵技術(shù)教育部重點研究室/國家柑橘產(chǎn)業(yè)體系機械化研究室/廣東省山地果園機械創(chuàng)新工程技術(shù)研究中心，廣州 510642

路面不平度作為車輛行駛時振動的主要激勵來源，對車輛行駛平順性、乘坐舒適性、操縱穩(wěn)定性、零部件疲勞壽命、運輸效率、油耗等均有影響，對路面結(jié)構(gòu)本身也有較大的影響［1］。農(nóng)用機械在果園非結(jié)構(gòu)化路面上工作時，受到的沖擊載荷會大于在普通路面上受到的沖擊。因此，對于農(nóng)用機械而言，果園路面不平度的研究對農(nóng)用機械的設(shè)計有重要意義。采集真實果園路面不平度信息，可以提高農(nóng)業(yè)機械農(nóng)田地面功率譜室內(nèi)模擬試驗準(zhǔn)確度［2-3］，但是傳統(tǒng)的路面不平度采集裝置存在使用范圍局限性的問題，不能滿足果園道路路面不平度的測量。隨著激光雷達(dá)和三維點云處理技術(shù)的應(yīng)用［4-5］，為路面不平度的測量提供了新思路［6］。

結(jié)構(gòu)化路面的不平度測量裝置目前相對成熟，但3 m 直尺、水平儀-標(biāo)桿、多輪儀等直接測量路面不平度的裝置由于效率低而逐漸被淘汰［7-9］。BPR 拖車式顛簸累積儀以及車載式路面平整度測試儀［10-11］由于應(yīng)用場景的限制，不適用于路面狀況較差的非結(jié)構(gòu)化路面不平度測量。三維激光雷達(dá)及點云處理技術(shù)應(yīng)用日益廣泛，國內(nèi)外研究者對激光雷達(dá)技術(shù)的運用做了大量的研究工作。如李永強等［12］利用激光雷達(dá)點云的處理方法對道路旁路燈的數(shù)據(jù)采集提取進(jìn)行了研究；武鵬［13］依據(jù)三維激光雷達(dá)點云平面特征，完成了點云的配準(zhǔn)。還有應(yīng)用在路面高程信息的研究報道，如鄔佳琪［14］采用激光雷達(dá)與相機兩種傳感器對路面數(shù)據(jù)進(jìn)行了采集，研究了全地面起重機路面高程信息；梁冠群等［15］提出利用長短期記憶網(wǎng)絡(luò)和時序輪心加速度的方法進(jìn)行實時路面不平度等級辨識。Liu 等［16］提出利用指定功率譜密度和相干函數(shù)的左右輪路徑的路面不平度的構(gòu)造方法。Yu 等［17］提出一種有效的紋理映射方法來重構(gòu)路面模型，利用激光掃描儀獲取被測物體的數(shù)據(jù)從而得到路面的細(xì)節(jié)信息，通過對重構(gòu)后的路面研究得到路面不平度；Kumar 等［18］采用固定激光掃描儀對道路進(jìn)行掃描以獲取路表信息，采用point cloud library 對采集的點云進(jìn)行運算處理，構(gòu)建了路面不平度。

隨著精準(zhǔn)農(nóng)業(yè)布局的展開，路面不平度對智能農(nóng)業(yè)裝備的精準(zhǔn)作業(yè)［19］、無人農(nóng)場的平穩(wěn)構(gòu)建有著重要影響，基于激光雷達(dá)的點云數(shù)據(jù)采集的方式精度比傳統(tǒng)方法更高，但是在農(nóng)業(yè)機械的研發(fā)過程中較少的研究者采取激光雷達(dá)的方式對山地果茶園的路面不平度進(jìn)行研究。本研究采用三維激光雷達(dá)測量與點云處理算法相結(jié)合的方法進(jìn)行路面不平度采集，旨在為山地林果茶園農(nóng)用機械的優(yōu)化設(shè)計提供參考。

1 材料與方法

1.1 系統(tǒng)組成

路面不平度采集系統(tǒng)主要由履帶運輸車、三維激光雷達(dá)、計算機、移動電源和安裝架組成（圖1）。三維激光雷達(dá)由安裝架固定在履帶運輸車上，三維激光雷達(dá)可以調(diào)節(jié)相應(yīng)角度進(jìn)行采集點云數(shù)據(jù)，采集到的點云數(shù)據(jù)保存至計算機中的VeloView 軟件。對比不同測量原理、不同維數(shù)、不同工作介質(zhì)的激光雷達(dá)，選擇了精度測量更高、適合于三維空間測量的VLP-16型激光雷達(dá)。

圖1 不平度采集系統(tǒng)示意圖Fig.1 Schematic diagram of unevenness acquisition system

1.2 路面選擇

測試路面需滿足以下要求：①覆蓋果園中常見的農(nóng)機行駛路面，如水泥路面、泥土路面、植被覆蓋路面以及砂石路面；②路面盡量不出現(xiàn)積水，路況良好。依據(jù)農(nóng)用車輛經(jīng)常行走的果園路面，選擇果園路面中常見的水泥路面、泥土路面以及砂石路面進(jìn)行試驗（圖2）。在實地考察后，分別在廣東省現(xiàn)代農(nóng)業(yè)裝備研究所黃埔基地以及華南農(nóng)業(yè)大學(xué)樹木園選擇了滿足本次試驗要求的路面。

開始測試前首先安裝激光雷達(dá)，如圖2B 所示，使激光雷達(dá)與地面成45°。然后打開上位機，選擇通訊接口配置菜單，正確設(shè)置各參數(shù)，觀察激光雷達(dá)是否正常，記錄數(shù)據(jù)。一切正常之后，遙控啟動履帶運輸車，在不同路面運行采集數(shù)據(jù)，再停止運輸車，利用Veloview 軟件記錄保存數(shù)據(jù)，即單次測量結(jié)束。路面分別選擇20 m的水泥、泥土和砂石路段，重復(fù)10次對該20 m路段進(jìn)行測試。計算平均速度，保留3組最接近2 km/h的結(jié)果。

圖2 3種果園模擬路況Fig.2 Three kinds of orchards to simulate road cjonditions

1.3 坐標(biāo)標(biāo)定與轉(zhuǎn)換

坐標(biāo)系的建立主要分為激光雷達(dá)坐標(biāo)系CL、車輛坐標(biāo)系CV和地面坐標(biāo)系CG的建立。由于CL與CV均為空間直角坐標(biāo)系，而且激光雷達(dá)與車身是固連的，其位置關(guān)系并不會隨著車輛的行駛而產(chǎn)生變化。標(biāo)定時，標(biāo)定物采用一種帶有直桿的紙板，選用直桿與紙板的交點作為被測量點。首先測量得出車體坐標(biāo)系下標(biāo)定點的測量值。然后通過三維激光雷達(dá)測得標(biāo)定點在三維激光雷達(dá)坐標(biāo)系下的坐標(biāo)值，再根據(jù)轉(zhuǎn)換公式計算得出標(biāo)定點在車體坐標(biāo)系下的計算值。激光雷達(dá)坐標(biāo)CL與車輛坐標(biāo)CV轉(zhuǎn)換關(guān)系如公式1所示。

車體坐標(biāo)系下標(biāo)定點的測量值與計算值以及誤差如表1 所示。通過表1 可以看出，最大誤差小于0.03 m，這樣的標(biāo)定誤差在激光雷達(dá)系統(tǒng)誤差0.03 m內(nèi)，即可認(rèn)定標(biāo)定結(jié)果產(chǎn)生的誤差不影響路面不平度裝置的測量。

表1 標(biāo)定誤差Table 1 Calibration errorm

2 結(jié)果與分析

2.1 點云數(shù)據(jù)去噪

采用VLP-16 激光雷達(dá)進(jìn)行數(shù)據(jù)采集時，雷達(dá)每秒可獲得數(shù)十萬原始點云數(shù)據(jù)。但是得到的原始點云數(shù)據(jù)會含有噪聲和異常數(shù)據(jù)，需要從原始點云數(shù)據(jù)中將噪聲和異常數(shù)據(jù)剔除，只保留關(guān)鍵特征。

果園路面不平度采集裝置在有植被覆蓋的路面進(jìn)行數(shù)據(jù)采集時，為獲取路面高程信息，需先將路表面覆蓋的植被（非地面點）剔除。為去除這部分多余的植被噪聲點，選用CloudCompare 中的布模擬濾波（cloth simulation filter，CSF）對每幀采集的點云上植被數(shù)據(jù)進(jìn)行濾波（圖3），首先對激光雷達(dá)點云進(jìn)行反向處理，然后設(shè)置合適的布的柔軟程度和網(wǎng)格分辨率生成布網(wǎng)。通過分析布節(jié)點與相應(yīng)激光雷達(dá)點之間的相互作用，可以確定布節(jié)點的位置，從而生成地面的近似。最后，通過對比原始激光雷達(dá)點和生成的表面，從激光雷達(dá)點云中提取出地面點，迭代完成所有的布粒子，如果雷達(dá)點云與布模擬粒子的距離小于閾值，則判斷為地面點，否則即是非地面點。經(jīng)布模擬濾波處理后地表覆蓋的植被即被去除，處理后地面點云圖如圖4C 所示，被去除的植被等非地面點則如圖4D所示。

圖3 布模擬濾波算法示意圖Fig.3 Schematic diagram of the cloth simulation algorithn

圖4 植被覆蓋路面布模擬濾波結(jié)果Fig.4 Simulation filtering result map of vegetation covered pavement

果園路面中部分被雜草樹木覆蓋的路面經(jīng)過布模擬濾波處理得到了無植被覆蓋的點云數(shù)據(jù)。當(dāng)路面不平度采集裝置在采集無植被覆蓋的路面點云信息時，激光雷達(dá)同樣采集了大量離群地面點云信息，為了對此部分無植被覆蓋的果園路面進(jìn)行精簡化處理，以及對布模擬濾波處理后的數(shù)據(jù)進(jìn)行離群點去除，選用統(tǒng)計濾波對每幀點云進(jìn)行濾波處理。

結(jié)合CloudCompare中的統(tǒng)計分析濾波器（statistical outlier removal filter），依據(jù)點到鄰近點的距離分布對點云進(jìn)行濾波處理。計算每個點與所有鄰近點的平均距離，通過分析，該結(jié)果滿足高斯分布，其均值和標(biāo)準(zhǔn)差決定高斯分布的形狀，若計算的平均距離在給定標(biāo)準(zhǔn)范圍之外，該點即為離群點，計算公式為：

圖5 布模擬濾波結(jié)合統(tǒng)計分析濾波處理結(jié)果Fig.5 Distributed simulation filtering combined with statistical analysis of filtering processing results

2.2 點云數(shù)據(jù)配準(zhǔn)

通過激光雷達(dá)獲取大量的點云數(shù)據(jù)，對點云數(shù)據(jù)進(jìn)行配準(zhǔn)處理，使之拼接為三維激光雷達(dá)點云數(shù)據(jù)是獲取路面高程信息的前提。利用激光雷達(dá)掃描過程中多幀點云數(shù)據(jù)重疊的信息進(jìn)行特征融合，從而提高激光雷達(dá)的數(shù)據(jù)量，進(jìn)而提高最終的識別效果。

通過CloudCompare 中ICP 模塊對2 幀點云進(jìn)行配準(zhǔn)，配準(zhǔn)后的路面點云信息量得到了倍增，提高了信息密度，有利于下一步進(jìn)行路面的高程信息提取和計算路面功率譜密度。

圖6 ICP配準(zhǔn)后的點云Fig.6 Point cloud after ICP registration

2.3 點云數(shù)據(jù)分割

搭載激光雷達(dá)的路面不平度采集裝置工作時，雷達(dá)測量到的點云數(shù)據(jù)包含360°范圍內(nèi)的點云信息。針對路面點云進(jìn)行了研究，圖7為多幀點云數(shù)據(jù)經(jīng)過濾波與配準(zhǔn)得到的道路配準(zhǔn)示意圖。采用Cloud-Compare中分割模塊對試驗路面進(jìn)行分割，即可分割出滿足要求的點云數(shù)據(jù)組。在得到的多幀路面配準(zhǔn)點云圖中，保留車輛行駛前方寬度為4 m、長20 m 的矩形范圍內(nèi)的部分點云用于路面高程研究。此時分割算法的設(shè)計，考察每個激光測量點的坐標(biāo)值，當(dāng)其坐標(biāo)值在公式（3）所示的范圍外時，將其濾除。留下的區(qū)域即為后續(xù)路面識別的區(qū)域，后續(xù)的路面高程信息的提取以及路面不平度的計算均在此區(qū)域內(nèi)進(jìn)行。

圖7 多幀點云配準(zhǔn)圖Fig.7 Multi-frame point cloud registration map

經(jīng)過了直通濾波分割后的點云如圖8 所示。用一種較為理想的直通濾波的方法，直接留下車輛最有可能駛過4 m×20 m的路面，進(jìn)行后續(xù)的路面高程信息提取。

圖8 點云分割后路面Fig.8 Pavement after point cloud segmentation

2.4 基于點云數(shù)據(jù)的路面高程提取

通過點云去噪、點云配準(zhǔn)與點云分割提取了4 m×20 m 的車輛可能行駛路面如圖9 所示。采用數(shù)字高程模型（DEM）對路面功率譜進(jìn)一步進(jìn)行計算，該模型按照一定的結(jié)構(gòu)，采用數(shù)字形式來表示地面高程。克里金（Kriging）算法能夠較好反映出各地形的變化，故為保證路面識別的準(zhǔn)確性，采用克里金模型提取路面的高程序列。

圖9 4 m ×20 m路面分割圖Fig.9 4 m ×20 m road segmentation map

對路面不平度進(jìn)行計算，從而提高對路面不平度判斷的準(zhǔn)度。如圖10所示，履帶車輛前進(jìn)時，提取左輪前進(jìn)方向與地面點云相交序列q（l）、右輪前進(jìn)方向與地面點云相交序列q（r）以及車輛左右對稱中心線處與地面點云相交序列q（z），得到的水泥路面二維高程序列如圖11所示。

圖10 點云數(shù)據(jù)提取示意圖Fig.10 Schematic diagram of point cloud data extraction

圖11 二維路面高程序列Fig.11 Two-dimensional road height program sequence

2.5 路面不平度計算

路面不平度是典型的平穩(wěn)隨機過程，可通過建立AR 模型來表征其統(tǒng)計特性。其中AR 模型中AR（P）功率譜密度估計方法如式（4）所示。

通過求解Yule-Walker 方程可得AR 模型中AR（P）的參數(shù)φk和σα2，從而求解路面不平度q(n)序列的功率譜密度估計。采用Matlab 中的基于AR 模型的功率譜密度估計函數(shù)Pyulear求解試驗路面的功率譜，果園水泥路面3條縱向序列數(shù)據(jù)的路面不平度功率譜密度曲線如圖12所示。

由圖12 可看出，水泥路面在高頻區(qū)出現(xiàn)了豐富的頻率分量，這是由于采用定帶寬分析來計算功率譜密度，突出了所謂真實功率譜分布或由噪聲引起的功率譜密度的起伏變化。根據(jù)國家標(biāo)準(zhǔn)GB/T 7031―2005 的要求，為了得到可觀測和可分類的數(shù)據(jù)，應(yīng)將譜值分解到各倍頻程中進(jìn)行平滑處理。以水泥路面為例，在Matlab 中將原始功率譜密度分解到各倍頻程中進(jìn)行平滑處理，處理后的平滑功率譜密度和空間頻率在對數(shù)坐標(biāo)下的關(guān)系如圖13A 所示。根據(jù)上述方法對泥土路面（圖13B）以及砂石路面（圖13C）高程信號進(jìn)行處理得到位移功率譜密度。

圖12 水泥路面功率譜密度Fig.12 Power spectral density of cement pavement

根據(jù)比例分析法對果園水泥路面、泥土路面以及砂石路面試驗結(jié)果進(jìn)行評價，計算q(l)、q(r)、q(z)3 個方向路面不平度各自占比的平均值（表2）。由平滑后的功率譜密度值（圖13）可知，同一路況下不同方向路面不平度曲線走勢大致相似。由表2 可知，水泥路面主要集中在B 級，B 級占比82.33%；泥土路面主要集中在D、E 兩個等級，D 級路面占比48.67%，E 級占比31.00%；砂石路面主要集中再C級，C級占比84.00%。

表2 不同路面不平度計算值Table 2 Calculated value of different road roughness

圖13 平滑后位移功率譜密度圖Fig.13 Smoothed power spectral density map

2.6 路面功率譜密度測試與驗證

履帶車輛行駛道路路面譜測試原理是將加速度傳感器放置于履帶車輛的負(fù)重輪處，測量相應(yīng)輪軸的振動加速度信號，采用傅里葉變換以及積分運算得路面功率譜密度，其原理流程如圖14所示。

圖14 測試原理流程圖Fig.14 Test principle flow chart

本研究采用Dytran4400B3 加速度振動記錄儀對路面針對車輛的振動進(jìn)行采集，其橫向靈敏度為5%，加速度振動記錄儀安裝圖如圖15 所示，框線標(biāo)識處為傳感器。將Dytran4400B3 加速度振動記錄儀依靠4 個磁體吸附安裝在負(fù)重輪軸支撐梁上。

圖15 加速度傳感器安裝位置Fig.15 Acceleration sensor installation location

為了驗證激光雷達(dá)果園路面不平度采集裝置的可行性，選用加速度振動儀采集相同試驗路面不平度信息進(jìn)行對比驗證。各路面加速度測量結(jié)果以時間為橫坐標(biāo)進(jìn)行繪圖，即得到加速度時域曲線，可以直觀的比較振幅大小，而疏密程度也從一定程度上反映了頻率高低（圖16）。試驗對加速度傳感器采集的振動數(shù)據(jù)進(jìn)行功率譜密度計算，采用DytranVibracoder 以及Matlab 對采集的數(shù)據(jù)進(jìn)行分析處理，可得3 個路面功率譜密度，確定路面不平度等級。依據(jù)AR 模型參數(shù)估計法，應(yīng)用Matlab 軟件對插值后的幅值譜序列進(jìn)行處理并以對數(shù)坐標(biāo)進(jìn)行繪圖，參考國標(biāo)GB/T7031―2005 中給出的不同倍頻程帶寬的平滑處理方法對功率譜密度進(jìn)行平滑處理，得到處理后的平滑曲線（圖17）。按照路面等級比例分析法，對路面等級進(jìn)行評價，可得3種路面等級（表3）。

圖16 加速度測量結(jié)果Fig.16 Acceleration measurement result

表3 不同路面不平度等級占比Table 3 Proportion of different road roughness levels%

圖17 平滑后路面功率譜密度圖Fig.17 Smoothed road power spectrum density map

比較激光雷達(dá)與加速度振動法測試的路面不平度等級，結(jié)果顯示，激光雷達(dá)的檢測結(jié)果為水泥路面主要集中在B 級，B 級路面占比82.33%；泥土路面主要集中在D、E 兩個等級，D 級路面占比48.67%，E 級占比31.00%；砂石路面主要集中在C 級，C 級占比84.00%。加速度振動法的檢測結(jié)果為水泥路面主要集中在B 級，B 級占比70.00%；泥土路面主要集中在D、E 兩個等級，D 級路面占比51.00%，E 級路面占比35.00%；砂石路面主要集中在C 級，C 級路面占比79.00%。結(jié)果表明，2種方法對同一試驗路面的測試的路面不平度等級一致。

2種方法的測量結(jié)果對比如圖18所示，水泥路面位于B 級，兩者誤差為15%，水泥砂石路面集中在C級，兩者誤差為6%，泥土路面集中在D、E 級，D 級兩者誤差為5%，E 級兩者誤差為13%，由此可見，2 種方法的誤差控制在15%以內(nèi)。其誤差來源大致如下：①在多次試驗的過程中，雖然盡量保持的相同路面，但是在車輛行駛過程中難以保證路面的完全一致。②基于3D 激光雷達(dá)和加速度振動傳感器的試驗精度的差別也會導(dǎo)致出現(xiàn)一定的誤差。雖然兩者存在一定誤差，但是誤差較小且路面評級一致可以驗證基于激光雷達(dá)的對果園路面不平度采集具有可行性。因此，基于激光雷達(dá)的果園路面不平度采集系統(tǒng)可以用于果園路面不平度的采集。

圖18 2種方法路面不平度分級對比圖Fig.18 Comparison chart of road roughness classification between two methods

3 討 論

本研究針對果園典型路面，采取基于三維激光雷達(dá)結(jié)合點云處理技術(shù)，搭建了基于三維激光雷達(dá)果園路面不平度信息采集系統(tǒng)平臺。通過提取果園路面的高程信息，計算果園路面的功率譜密度，采用比例分析法確定了路面不平度等級評定，通過加速度振動記錄儀開展了采集系統(tǒng)驗證試驗，試驗結(jié)果顯示所測水泥路面中A 級路面占比2.67%，B 級路面占比82.33%，C 級路面占比15.00%；泥土路面中C級路面占比4.33%，D 級路面占比48.67%，E 級路面占比31.00%，F(xiàn) 級路面占比16.00%；砂石路面中B級路面占比4.33%，C 級路面占比84.00%，D 級路面占比8.67%，E 級路面占比3.00%。試驗結(jié)果表明，基于激光雷達(dá)果園路面不平度的采集系統(tǒng)的測試結(jié)果與加速度振動記錄儀測試結(jié)果一致。同時，本研究采用單個VLP-16 激光雷達(dá)采集路面點云信息，對采集到的點云信息只能一幀一幀手動進(jìn)行迭代完成點云構(gòu)圖，當(dāng)數(shù)據(jù)量過大時，會導(dǎo)致三維點云配準(zhǔn)失真，測量誤差增大的情形。今后可以采用多傳感器融合技術(shù)來進(jìn)行測量，以提高路面不平度信息采集的精度。