基于改進三次B樣條曲線的奶牛點云缺失區域修復方法

何東健 牛金玉 張子儒 郭陽陽 譚 益

(1.西北農林科技大學機械與電子工程學院， 陜西楊凌 712100； 2.農業部農業物聯網重點實驗室， 陜西楊凌 712100； 3.陜西省農業信息感知與智能服務重點實驗室， 陜西楊凌 712100)

0 引言

奶牛養殖業在我國國民經濟中具有重要的戰略地位[1-2]。奶牛的規模化標準養殖是提升我國奶業競爭力和品質的必經之路，也是我國奶牛產業發展的必然要求[3]。規模化標準養殖需要優良的奶牛品種，而奶牛的體型評定對選育適應機械化擠奶的優良奶牛品種具有重要意義。

奶牛的三維點云數據包含奶牛身體尺寸和外形等幾何信息，因此，獲取奶牛的三維點云數據對奶牛的體型評價、體尺測量以及健康狀況評定和行為模型研究均具有重要的應用價值[4]。目前，獲取三維點云數據的途徑主要有雙目立體視覺系統、深度相機和激光掃描儀等[5]，由于環境和設備精度影響，特別是奶牛養殖場圍欄等復雜環境因素，導致獲取的點云存在孔洞現象，點云孔洞不僅使模型實現正確的可視化受限，而且嚴重影響點云的后期處理工作[6-7]。

國內外學者針對點云孔洞修復問題進行了許多研究。WANG等[8]使用移動最小二乘法并依據模型的幾何形狀和陰影信息對孔洞進行擬合重建；CHALMOVIANSKY等[9]首先計算模型中每個點的鄰近點，以此確定孔洞邊界，然后使用邊界的鄰近點構造代數曲面片填充孔洞；MARCHANDISE等[10]通過正交梯度法求解徑向基函數在任意曲面上的偏微分方程，較好地修補了網格模型中存在的孔洞。NGO等[11]將復雜的孔洞分成較小的簡單子孔洞，再用三角網格算法修補子孔洞，最終完成復雜孔洞修補。QUINSAT等[12]以數值模型的先驗知識作為名義網格，識別出數字化孔洞后，計算名義網格和點云之間的差異，然后確定最小化的變形能量進行網格變形，對點云模型上的孔洞進行修補。蔣剛[13]先對空間采樣點曲線在3個坐標平面上投影，然后用支持向量機分別對3組二維投影數據進行回歸運算，計算得到孔洞的三維坐標，完成孔洞修補；陳相等[6]利用三角形網格特征尋找到孔洞邊界，然后依據孔洞邊界夾角進行三角形網格填充，得到修補效果較好的點云模型。晏海平等[14]提出一種基于徑向基函數的孔洞修復算法，該算法提取孔洞多邊形邊界特征并用最小二乘法填充孔洞，再建立孔洞邊界點和鄰域信息的徑向基隱式曲面，結合梯度下降法調整孔洞填充點的位置，實現點云孔洞的修復。上述算法均針對區域相對較小的封閉孔洞進行修復，而去除遮擋欄桿后的奶牛點云，不僅點云缺失區域大，且為非封閉孔洞，上述算法難以有效對非封閉較大孔洞進行修復。

為解決奶牛點云非封閉較大孔洞的修復問題，本文提出一種基于改進三次B樣條曲線的修復方法，以提高奶牛點云模型精度。

1 奶牛點云數據獲取及預處理

1.1 奶牛點云數據獲取

在陜西省楊凌科元克隆股份有限公司奶牛養殖場，用美國微軟公司Kinect v2體感設備(簡稱Kinect)及配套軟件，采集健康荷斯坦奶牛深度數據。

安裝在支架上的Kinect設置于活動場和擠奶間欄桿通道(長約4 m，寬約1.4 m)一側，Kinect傳感器有效檢測范圍為0.5～4.5 m，為有效拍攝到奶牛多幀全身點云且降低環境點云影響，故視野范圍應覆蓋此段欄桿通道，經測試，將Kinect放置于其中心距最近欄桿(直徑0.05 m)軸線2.15 m處，距離地面高度為0.9 m。當奶牛整個身體完全進入Kinect的視野范圍時，開始采集深度碼流數據，采集幀率30 f/s，深度圖像分辨率為512像素(水平)×424像素(垂直)，數據保存在便攜式計算機固態硬盤(三星850 EVO SATA3，120 G)中。

于2017年3月24—30日20:00—21:30，采集45頭奶牛全身側視深度數據，每頭奶牛采集5幀，共計225幀。

1.2 點云數據預處理

在臺式計算機(CPU為Intel Core i5-4590 3.3 GHz、8 GB內存)上，用Microsoft Visual Studio 2013、PCL點云庫、Matlab R2017a進行數據處理。

由于拍攝環境復雜，獲取的點云中除奶牛外，還包括地面、水平欄桿、豎直欄桿等干擾，故需要去除非奶牛點云數據，如圖1所示，點云數據坐標系以紅外攝像頭感光元件中心為坐標原點o，x軸水平向左，y軸豎直向上，z軸沿拍攝方向垂直xoy面。

圖1 Kinect獲取的奶牛點云Fig.1 Dairy cow’s point cloud obtained by Kinect

奶牛點云預處理流程如圖2所示，具體步驟為：①用直通濾波方法去除部分背景點云，如圖3a所示。②用查詢點與鄰域點集之間平均距離估計方法[15]過濾環境中的離群點。③用帶法線的平面模板匹配法去除地面點云數據，如圖3b所示。④用與給定軸線平行的直線模板匹配法去除水平欄桿點云。⑤用查詢點與鄰域點集之間平均距離估計的方法去除余留下的少量欄桿邊緣噪聲點，如圖3c所示。⑥對于豎直欄桿和與之相連的少量地面點云數據，用歐氏聚類方法[16]予以去除，最終得到奶牛個體的點云數據，如圖3d所示。

圖2 預處理流程圖Fig.2 Flow chart of preprocessing

圖3 奶牛點云預處理結果Fig.3 Results of preprocessing dairy cow’s point cloud

2 奶牛缺失區域修復方法

經過預處理后的奶牛點云存在被水平欄桿遮擋的缺失較大的區域，將對后續體尺參數提取產生重要影響，故本文提出了基于改進的三次B樣條曲線的點云修復方法。先去除環境背景后的奶牛點云進行切片投影，再對每個切片點列中間距較大的相鄰點進行補點，最后對補點后得到的新點列進行三次B樣條曲線擬合，即得到修復后奶牛點云。

2.1 奶牛點云切片投影

奶牛點云缺失區域基本沿x軸方向且上下側均有點云，故可用上下側點云信息修復缺失區域。先將點云沿x軸進行切片投影處理，步驟如下：

(1)確定奶牛點云中點與點之間最小空間距離dmin，沿x軸方向以dmin等間隔將點云切片采樣為N份。由于部分切片采樣點云中存在斷連間隔，故對切片寬度前后分別擴展dmin。

(2)考慮到擴展切片采樣區間后點數增加較多，對曲線擬合有較大影響，故對每個切片點云沿y軸以dmin等間隔將切片點云分為Ni份，對每份取其所有點的中心作為該份的代表點。

(3)將沿x軸方向切片的點云向yoz面投影，再對投影的點進行三次B樣條曲線擬合，將獲得的新點列還原到三維空間時，將x坐標設置為切片采樣區間中點。對N份點云切片做同樣處理。

奶牛前后肢處切片投影相對復雜且具有代表性，故選取前肢處一切片投影作為擬合的典型代表，其點云切片投影如圖4所示。

圖4 奶牛前肢點云切片投影Fig.4 Slice projection of dairy cow’s forelimb point cloud

2.2 奶牛點云切片投影的擬合方法

不同部位的奶牛點云切片投影形狀和點數各不相同，且相鄰點之間密集程度差異較大，故無法用最小二乘法[17-19]確定通用擬合函數以擬合不同切片投影，使用樣條曲線[20-21]易出現過擬合現象，而三次B樣條曲線[22-25]具有局部性、凸包性等性質，可由4個點擬合一段曲線，且不會出現過擬合現象，故本文用三次B樣條曲線對奶牛點云切片投影進行擬合。

2.2.1三次B樣條曲線擬合

對于任意相鄰的4個頂點Pi、Pi+1、Pi+2、Pi+3，三次B樣條曲線段表達式為

(1)

對三次B樣條方法擬合后的曲線進行擬合取點時，取決于參數t的取值，為了使擬合點列均勻，可將t均勻取值。但因奶牛點云中存在水平欄桿遮擋，切片投影中存在間隔較大的相鄰點，因此，用三次B樣條曲線對切片投影進行擬合時，若參數t均勻取值，則擬合后在間隔較大的相鄰點處會出現點列稀疏問題。如圖5所示，對奶牛肢蹄處切片投影點云分別進行三次B樣條曲線擬合，由局部放大圖可以看出，在原點列間隔較大處擬合得到的新點列均勻性較差。

圖5 三次B樣條曲線擬合Fig.5 Cubic B-spline curve fitting

2.2.2改進的三次B樣條擬合方法

針對三次B樣條曲線擬合出現的點列稀疏問題，一種方法是再次用三次B樣條曲線擬合，但這會使擬合的曲線更加偏離原曲線位置，而且擬合的形狀與原曲線差異更大；另一種方法是先對其進行補點，再用三次B樣條曲線擬合，該方法可以使擬合的曲線更加逼近原曲線，如文獻[25]中增加型值點的方法。但本文中奶牛切片投影存在較大的間隔區域，文獻[25]的補點方法會使增加的點超出補點區域，導致擬合曲線局部偏差較大。為此，本文對文獻[25]方法進行改進，提出一種對較大間隔區域補點的方法，然后再用三次B樣條曲線擬合，以使擬合的曲線更加逼近原曲線。

改進方法原理如圖6所示，具體步驟如下：

圖6 補點示意圖Fig.6 Diagram of adding points

(1)已知Pi、Pi+1、Pi+2、Pi+3為平面上的點列，過Pi+1點作PiPi+2的平行線li+1，過Pi+2點作Pi+1Pi+3的平行線li+2。

(2)設補點間距閾值為L，若|Pi+1Pi|≥L、|Pi+1Pi+2|≥L，則分別做li+1與Pi+1Pi和Pi+1Pi+2夾角的角平分線Pi+1Pi+1,0和Pi+1Pi+1,1，使

(2)

式中h——比例系數

若相鄰點間距小于L，則不進行補點。

同理，在Pi+2點處可得到Pi+2,0和Pi+2,1點。

(3)將Pi、Pi+1,0、Pi+1、Pi+1,1、Pi+2,0、Pi+2、Pi+2,1、Pi+3按順序合并成新點列，再返回步驟(1)、(2)，直到所有相鄰點間距均小于L。

(4)將補點后得到的新點列進行B樣條曲線擬合，即可得到平滑的點列。

3 試驗結果與分析

為了驗證本文方法對奶牛點云欄桿遮擋較大缺失區域的修復性能，驗證以下方面：

(1)奶牛典型肢體部位切片投影擬合后點列的均勻性以及擬合點列對原點列的逼近性能。點列的均勻性用相鄰點間距均值和標準差評價；擬合點列對原點列的逼近性能用切片逼近誤差評價，即計算原點列的每個點與擬合點列中點的最小間距，將所有最小間距的均值作為切片逼近誤差。

(2)擬合后的每幀奶牛點云對原點云的逼近性能。用幀逼近誤差，即每幀奶牛點云所有切片逼近誤差的均值，評價對原點云的逼近性能。

3.1 參數h和L取值分析

本文算法中直接影響點云切片投影均勻性和逼近性能的參數是比例系數h和間距閾值L，因此需對二者的較優取值進行試驗。

對圖4所示的前肢處切片投影，分別取h為1/6、1/5、1/4和1/3，L/dmin在1～20間取整數進行擬合，統計相鄰點間距均值、標準差以及擬合后點數和切片逼近誤差，如圖7所示。

分析圖7可知，當L/dmin=1時，相鄰點間距均值、標準差和逼近誤差達到最小，但均值小于擬合間距dmin(該幀點云為7.114 mm)的0.5倍，同時，擬合后點數也大于三次B樣條方法獲得點數的2倍，點數的成倍增加會對后期奶牛點云體尺提取效率等產生重要影響。當h取值為1/4，L/dmin取值為5時，可使標準差和逼近誤差達到較小值。

取h=1/4，L/dmin=5，對奶牛前肢處點云切片投影進行擬合，結果如圖8所示，擬合后相鄰點間距均值為7.552 mm，標準差為1.305 mm，逼近誤差為0.720 mm，點數由擬合前150個增加到擬合后的256個，與圖5僅用三次B樣條曲線擬合相比，本文方法擬合后點列較為均勻，且更逼近原點列。

圖7 奶牛前肢處點云切片投影擬合結果Fig.7 Fitting results of slice projection constructed by point cloud of dairy cow’s forelimb

圖8 本文方法擬合結果Fig.8 Fitting result of proposed method

3.2 奶牛擬合點云逼近性能分析

取h=1/4，L=5dmin，用本文方法對45頭奶牛的225幀點云進行欄桿遮擋缺失區域修補，計算出幀逼近誤差和平均幀逼近誤差，如圖9所示。從圖9可知，僅用三次B樣條曲線擬合得到的平均幀逼近誤差為2.214 mm，最大和最小幀逼近誤差為3.918 mm和0.819 mm，而本文方法得到的平均幀逼近誤差為1.623 mm，最大和最小幀逼近誤差為2.547 mm和0.737 0 mm，平均幀逼近誤差降低了26.7%，更加逼近原點云。

圖9 奶牛點云幀逼近誤差Fig.9 Fitting error of dairy cow’s point cloud by frame

圖10為三次B樣條方法和本文方法修復效果圖，圖10a～10d為未修復的整體和局部細節圖，圖10e～10h為三次B樣條曲線修復得到的整體和局部細節圖，圖10i～10l為使用本文算法修復得到的整體和局部細節圖。可以看出，圖10f與10g在奶牛點云欄桿遮擋區域邊緣有明顯的空缺，而本文方法得到的圖10j與10k在欄桿遮擋區域邊緣的修復效果較好。從奶牛點云的背部俯視圖可以看出，三次B樣條方法和本文方法對較小間隔均能取得較好的修復效果，但在背部有較大間隔處，本文方法修復效果更好。

圖10 奶牛點云修復效果對比圖Fig.10 Repairing effect comparison of dairy cow’s point cloud

從圖10中還可以看出，在奶牛頭部點云缺失區域，仍存在未得到修復的區域，這是因為該處點云缺失區域僅一側有點云信息，而另一側無點云信息。

3.3 噪聲敏感性試驗

由于奶牛身體存在黑白塊區域，Kinect發出的紅外光會被黑色區域吸收一部分，導致獲取的點云存在噪聲。為考察噪聲對本文方法擬合效果的影響，對含有噪聲的奶牛點云用本文方法進行擬合試驗，結果如圖11所示。

圖11 噪聲對本文方法的影響Fig.11 Effect of noise on proposed method

圖11b為對圖11a包含噪聲的肢體部分擬合的結果，可見在肢體邊緣會產生不屬于奶牛的多余點。圖11e為對圖11d所示存在噪聲的腹部進行擬合的結果，擬合后的點云噪聲點數增大。因此，在點云擬合前應先進行光順和降噪處理。對圖11a和圖11d點云進行高斯濾波[16]，用本文方法進行擬合的結果如圖11c和圖11f所示，可見肢體邊緣和腹部擬合效果得到改善，但由于采用統一的高斯濾波參數，仍存在肢蹄變細、腹部有少量噪聲問題，故可利用局部平滑方法進一步改善。

4 結論

(1)提出了一種基于三次B樣條曲線的奶牛點云較大缺失區域修復方法，參數h=1/4，L/dmin=5時，修復后的點列均勻性和逼近性能均優于三次B樣條方法。

(2)對45頭奶牛共225幀點云修復試驗，結果表明，與三次B樣條方法相比，本文方法平均幀逼近誤差降低了26.7%，更加逼近原點云。

(3)在奶牛肢體和腹部存在噪聲時，修復效果會受到較大影響，用高斯濾波進行光順降噪后可以取得較好的平滑效果，但肢蹄變細、腹部仍存在少量噪聲問題，需對局部光順降噪算法做進一步研究。

1 李順.中國畜牧業發展歷程分析及趨勢預測[J].中國畜牧雜志,2010,46(12):25-28.

LI Shun. Analysis on China livestock development process and trend forecast[J]. Chinese Journal of Animal Science,2010,46(12):25-28. (in Chinese)

2 劉長全, 楊洋.中國奶業產業政策的發展及成效[J].中國奶牛,2017(10):58-64.

3 莫佳蓓, 劉慶友.中國和日本奶牛養殖業規模化發展的比較分析[J].江蘇農業科學,2017,45(9):303-307.

4 何東健, 劉冬, 趙凱旋.精準畜牧業中動物信息智能感知與行為檢測研究進展[J/OL].農業機械學報,2016,47(5):231-244.http:∥www.j-csam.org/jcsam/ch/reader/view_abstract.aspx?file_no=20160532&flag=1.DOI:10.6041/j.issn.1000-1298.2016.05.032.

HE Dongjian, LIU Dong, ZHAO Kaixuan. Review of perceiving animal information and behavior in precision livestock farming[J/OL]. Transactions of the Chinese Society for Agricultural Machinery,2016,47(5):231-244. (in Chinese)

5 邵小寧.基于Kinect的植物三維點云獲取與重建方法研究[D].楊凌:西北農林科技大學,2016.

SHAO Xiaoning. Research on 3D point cloud acquisition and reconstruction method of plant based on Kinect[D].Yangling:Northwest A&F University,2016. (in Chinese)

6 陳相, 童小華.基于三角格網的點云空洞修補算法及精度研究[J].測繪通報,2013(4):1-3.

CHEN Xiang, TONG Xiaohua. Research on points cloud hole filling algorithm and accuracy in triangle mesh[J]. Bulletin of Surveying and Mapping,2013(4):1-3. (in Chinese)

7 王果, 崔希民, 俞世煒, 等.應用掃描光線重建法的建筑立面遮擋點云恢復[J].測繪通報,2014(5):63-65.

WANG Guo, CUI Ximin, YU Shiwei, et al. Recovery of occluded building facade point cloud based on the reconstructed scan light and planar segmentation[J]. Bulletin of Surveying and Mapping,2014(5):63-65. (in Chinese)

8 WANG J, OLIVEIRA M M. Filling holes on locally smooth surfaces reconstructed from point clouds[J]. Image and Vision Computing, 2007, 25(1): 103-113.

9 CHALMOVIANSKY P, JUTTLER B. Filling holes in point clouds[J]. Lecture Notes in Computer Science, 2003, 2768: 196-212.

10 MARCHANDISE E, PIRET C, REMACLE J F. CAD and mesh repair with radial basis functions[J]. Journal of Computational Physics,2012, 231(5): 2376-2387.

11 NGO T M, LEE W S. Feature-first hole filling strategy for 3D meshes[M]. Springer Berlin Heidelberg, 2013: 53-68.

12 QUINSAT Y, LARTIGUE C. Filling holes in digitized point cloud using a morphing-based approach to preserve volume characteristics[J]. International Journal of Advanced Manufacturing Technology,2015, 81(1-4): 411-421.

13 蔣剛.基于SVM和空間投影的點云空洞修補方法[J].計算機工程,2009,35(22):269-271.

JIANG Gang. Point cloud hole filling method based on SVM and space projection[J]. Computer Engineering, 2009,35(22):269-271. (in Chinese)

14 晏海平,吳祿慎,陳華偉.基于徑向基函數的散亂點云孔洞修復算法[J].計算機工程與設計,2014,35(4):1253-1257.

YAN Haiping, WU Lushen, CHEN Huawei. Hole repairing algorithm in scattered point cloud based on radial basis function[J]. Computer Engineering and Design,2014,35(4):1253-1257. (in Chinese)

15 RUSU R B, MARTON Z C, BLODOW N, et al. Towards 3D point cloud based object maps for household environments[J]. Robotics and Autonomous Systems,2008, 56(11): 927-941.

16 朱德海. 點云庫PCL學習教程[M]. 北京: 北京航空航天大學出版社, 2012.

17 劉俊焱, 云挺, 周宇, 等.基于最小二乘的點云葉面擬合算法研究[J].西北林學院學報,2014,30(5):70-77.

LIU Junyan, YUN Ting, ZHOU Yu, et al. Leaf surface fitting of point cloud based on least square method[J]. Journal of Northwest Forestry University,2014,30(5):70-77. (in Chinese)

18 閆蓓, 王斌, 李媛.基于最小二乘法的橢圓擬合改進算法[J].北京航空航天大學學報,2008,34(3):295-298.

YAN Bei, WANG Bin, LI Yuan. Optimal ellipse fitting method based on least-square principle[J]. Journal of Beijing University of Aeronautics and Astronautics,2008,34(3):295-298. (in Chinese)

19 田垅, 劉宗田.最小二乘法分段直線擬合[J].計算機科學,2012,39(6A):482-484.

TIAN Long, LIU Zongtian. Least-squares method piecewise linear fitting[J]. Computer Science,2012,39(6A):482-484. (in Chinese)

20 苗莎, 鄭曉薇.三次插值樣條曲線擬合多核并行算法[J].計算機應用,2010,30(12):3194-3196.

MIAO Sha, ZHENG Xiaowei. Multi-core parallel algorithm for cubic spline curve fitting[J]. Journal of Computer Applications,2010,30(12):3194-3196. (in Chinese)

21 海嘯, 朱志杰.最小二乘法和三次樣條曲線擬合的比色測溫誤差修正對比分析[J].激光雜志,2015,36(6):72-76.

HAI Xiao, ZHU Zhijie. The comparison between the least square method and cubic spline interpolation in errors of colorimetry temperature measurement[J]. Laser Journal,2015,36(6):72-76. (in Chinese)

22 喻德生, 程程.基于離散曲率的三次均勻B樣條的局部光順算法[J].浙江大學學報:理學版,2011,38(5):511-517.

YU Desheng, CHENG Cheng. On local fairing algorithm for cubic B-spline with discrete curvature[J]. Journal of Zhejiang University:Science Edition,2011,38(5):511-517. (in Chinese)

23 韓江, 江本赤, 夏鏈, 等.基于輪廓關鍵點的B樣條曲線擬合算法[J].應用數學和力學,2015,36(4):423-431.

HAN Jiang, JIANG Benchi, XIA Lian, et al. A B-spline curve fitting algorithm based on contour key points[J]. Applied Mathematics and Mechanics,2015,36(4):423-431. (in Chinese)

24 王文珂, 李思昆.無序B樣條曲線的曲面擬合算法[J].計算機輔助設計與圖形學學報,2013,25(5):679-685.

WANG Wenke, LI Sikun. B-spline surface fitting to unorganized curves[J]. Journal of Computer-Aided Design & Computer Graphics,2013,25(5):679-685. (in Chinese)

25 段振云, 王寧, 楊旭, 等.一種改進B樣條曲線擬合算法研究[J].機械設計與制造,2016(5):17-19.

DUAN Zhenyun, WANG Ning, YANG Xu, et al. An algorithm of improved B-spline curve fitting[J]. Machinery Design & Manufacture,2016(5):17-19. (in Chinese)