國立高雄大學資訊工程學系研究所
碩士論文
基於 K-D 樹以局部移動性最小平方重建法進行三維模
型修復
3D Model Surface Reconstruction Using Local Moving
Least Squares and KD Trees
研究生：洪郁媚
撰
指導教授：陳佳妍
博士
殷堂凱
博士
i
致謝
三年的研究生生活，雖然一路上跌跌撞撞，遇到許多挫折，但也因此學習了 解到許多新的事物，而且因為受到許多人的鼓勵與幫助，才有辦法完成學業論文， 很開心可以擁有這段美好的記憶，希望可以藉由這篇文章，表達對你們的感謝。 首先是指導教授陳佳妍老師，感謝老師於課業上的指導，在論文遇到問題時， 不厭其煩給予研究方向，並且帶領我們參與研討會，讓我們有不同的視野。感謝 指導教授殷堂凱老師於論文上的幫助，對於問題給於不同方向建議和啟發。感謝 陳銘志老師與彭昭暐老師特地抽空前來參加口試，並且給予許多珍貴的意見。感 謝張家泓學長於論文上給予許多的幫助，不管是研究方法或是程式撰寫。感謝實 驗室的學長、同學、學弟妹們與系助淑真，在這段時間裡，很多事情因為有你們 的幫助與鼓勵，才能解決許多大小事情。最後要感謝家人，謝謝你們支持讓我可 以完成我的學業。 最後感謝於求學過程中所有幫助過我的人，因為你們的存在豐富了我的研究 生生活，不管事課業或是為人處事都給予極大的幫忙，我會帶著你們給予的鼓勵， 繼續向未來邁進。 洪郁媚 謹致於 國立高雄大學 資訊工程學系 中華民國一零六年十月ii
基於 K-D 樹以局部移動性最小平方重建法進行三維
模型修復
指導教授：陳佳妍 博士 殷堂凱 博士 國立高雄大學資訊工程學系 學生：洪郁媚 國立高雄大學資訊工程學系研究所 摘要 近幾年許多資訊技術快速發展。三維重建也是熱門議題之一，卻在重建模型 時，可能因外在因素而導致資料缺失，而造成三維模型重建結果部分點雲缺失， 進而造成模型表面有空洞。如果重新擷取資料與三維模型重建，將會需花費過多 的時間與人力成本。 本研究提出資料修補流程，利用 K-D 樹找尋各點限制範圍內的 k 個最近鄰 居（K-Nearest Neighbors, KNN)，並且利用穩定的局部雜訊範圍偵測（Robust Local Noise Scale Estimation)計算各點所存在之最適平面，並對模型進行剃除離群值與 雙邊濾波器（Bilateral Filter)降低雜訊。接著藉由各點所存在的最適平面，偵測 空洞邊界上的點，並劃分修補範圍以進行局部移動性最小平方重建法（Local Moving Least Squares, LMLS）修補空洞。最後執行參數化攤平（Parameterization Mapping)，重建德勞內三角化（Delaunay Triangulation, DT)網格模型。iii
3D Model Surface Reconstruction Using Local Moving
Least Squares and KD Trees
Advisor(s): Dr. Chia-Yen Chen Dr. Tang-Kai Yin
Department of Computer Science and Information Engineering National University of Kaohsiung
Student: Yu-Mei Hong
Department of Computer Science and Information Engineering National University of Kaohsiung
ABSTRACT
In recent years, technologies to analyze and process the various types of acquired data are being investigated and developed. 3D reconstruction and modelling is among one of the developing technologies that received a lot of attention lately. During the process of 3D reconstruction, a multitude of factors may cause loss of data or failure to acquire correct data, resulting in missing portions of the acquired 3D point cloud. The missing or incorrect data further causes holes or artifacts within the rendered 3D model. Therefore, it is desirable to develop methods to process the 3D point cloud such that the 3D models can be rendered in a more realistic manner.
In this work, we propose a process to perform 3D surface reconstruction based on local moving least squares and KD trees. In the proposed method, K-D tree is used to locate the kth nearest neighbors within the neighbourhood of a point, robust local noise scale estimation is used to calculate the best fitting surface for each point. Outlier removal and bilateral filter are used to reduce noise on the 3D model. The best fitting surface for each point is then used to determine the points on the boundaries of holes, and the area is repaired using local moving least squares. Finally, parameterization mapping is used and Delaunay triangulation of the 3D model is performed to produce the final mesh model.
Keywords: Surface Reconstruction, Local Moving Least Square, KD-Trees, Delaunay Triangles.
iv
目錄
致謝 ... i 摘要 ... ii ABSTRACT ... iii 目錄 ... iv 圖目錄 ... v 表目錄 ... vi 符號定義 ... vii 第一章 緒論 ... 1 1.1 前言 ... 1 1.2 研究目的 ... 1 1.3 研究方法 ... 2 1.4 論文架構 ... 2 第二章 文獻探討 ... 4 第三章局部移動性最小平方重建法 ... 7 3.1 資料前處理 ... 8 3.1.1 利用 KD-Tree 找尋最近鄰居 ... 8 3.1.2 平面估測與剃除離群值 ... 9 3.1.3 雜訊平滑化 ... 11 3.2 資料修補 ... 14 3.2.1 空洞位置定義 ... 14 3.2.2 局部移動性最小平方重建法（LMLS） ... 16 3.3 德勞內三角化網格模型重建 ... 18 3.3.1 德勞內三角化 ... 19 3.3.2 參數化攤平 ... 21 第四章實驗結果與分析 ... 26 4.1 離群值剃除與雜訊平滑化 ... 27 4.2 空洞位置定義 ... 28 4.3 局部移動性最小平方重建法（LMLS） ... 31 4.4 參數化攤平網格重建結果 ... 39 第五章 結論與未來研究方向 ... 42 5.1 結論 ... 42 5.2 未來研究方向 ... 43 參考文獻 ... 44v
圖目錄
圖 2-1 研究實作流程圖 ... 5 圖 3-1 研究流程架構圖 ... 7 圖 3-2 KD-Tree 空間切割圖[10] ... 8 圖 3-3 KD-Tree 樹狀圖[10] ... 9 圖 3-4 平面偵測圖 ... 11 圖 3-5 雙邊濾波器示意圖[15] ... 12 圖 3-6 k個最近點搜尋圖[17] [18] 。 ... 14 圖 3-7 空洞邊界點集合分類圖 ... 15 圖 3-8 空洞過濾示意圖 ... 16 圖 3-9 LMLS 適應性修補圖 ... 17 圖 3-10 網格模型重建步驟圖 ... 19 圖 3-11 最大空圓原則圖 ... 20 圖 3-12 最大化最小角原則 ... 20 圖 3-13 逐點插入法步驟圖 ... 21 圖 3-14 參數化攤平初始化示意圖 ... 22 圖 3-15 參數攤平化示意圖 ... 23 圖 3-16 產生錯誤網格情況 ... 25 圖 4-1 點雲模型離群值剃除與雜訊平滑結果比較圖 ... 28 圖 4-2 曲面模型空洞定義結果比較圖 ... 29 圖 4-3 兔子點雲模型空洞定義結果比較圖 ... 30 圖 4-4 行政大樓點雲空洞定義結果比較圖 ... 31 圖 4-5 曲面模型 LMLS 與 MLS 結果比較圖 ... 33 圖 4-6 兔子模型 LMLS 與 MLS 結果比較圖 ... 35 圖 4-7 行政大樓點雲模型 LMLS 與 MLS 結果比較圖 ... 37 圖 4-8 曲面模型修補圖 ... 38 圖 4-9 三角錐模型修補圖 ... 38 圖 4-10 球體模型修補圖 ... 38vi
表目錄
表 4-1 實驗結果誤差值表 ... 39 表 4-2 實驗數據設定表 ... 39
vii
符號定義
k 最近鄰居數目 𝑛𝑛𝑢𝑢 ����⃑ 平面之單位法向量。 𝑃𝑃(𝑃𝑃𝑋𝑋, 𝑃𝑃𝑌𝑌, 𝑃𝑃𝑍𝑍) 存在於平面上的點群。 𝑑𝑑 點𝑃𝑃到平面的距離。 𝑛𝑛�⃑ 具有與單位向量𝑛𝑛����⃑同方向之法向量。 𝑢𝑢 𝛿𝛿(𝑛𝑛�⃑, 𝑃𝑃) 平面到任一點距離函數。 𝜙𝜙 空間中𝑘𝑘個點到平面距離平方的總和。 𝑃𝑃𝑖𝑖 點雲中的第 i 點。 𝑁𝑁 𝑃𝑃𝑖𝑖所存在的平面所包含的鄰居數目。 𝑀𝑀 計算點到平面距離限制所需迭代的次數。 𝜀𝜀𝑁𝑁 𝑃𝑃𝑖𝑖所存在平面中點雲可能包含雜訊的比值。 𝑃𝑃̇ 欲刪除雜訊的比例。 𝑝𝑝𝑖𝑖(𝑢𝑢, 𝑣𝑣) 點雲投影於二維平面的點。 𝛪𝛪 影像強度值。 𝛪𝛪̂(𝑝𝑝𝑖𝑖) 𝛪𝛪進行平滑化之後，獲得𝑝𝑝𝑖𝑖新的強度值𝛪𝛪̂。 𝑝𝑝𝑗𝑗 𝑝𝑝𝑖𝑖的鄰居。 𝑊𝑊𝑐𝑐(�𝑝𝑝𝑗𝑗− 𝑝𝑝𝑖𝑖�) 空間域差距權重。 𝑊𝑊𝑠𝑠(�𝛪𝛪(𝑝𝑝𝑖𝑖) − 𝛪𝛪�𝑝𝑝𝑗𝑗��) 強度域差距權重。 𝑥𝑥𝑑𝑑 𝑝𝑝𝑖𝑖與𝑝𝑝𝑗𝑗的距離。 𝑥𝑥𝐼𝐼 𝑝𝑝𝑖𝑖與𝑝𝑝𝑗𝑗的強度差距。 𝜎𝜎𝑐𝑐 空間域差距權重高斯參數。 𝜎𝜎𝑠𝑠 強度域差距權重高斯參數。 𝜌𝜌 𝑝𝑝𝑖𝑖與鄰居的距離最大範圍值。viii 𝑝𝑝� 𝚤𝚤 𝑝𝑝𝑖𝑖完成平滑化的點。 𝑛𝑛 平面法向量。 𝑑𝑑𝑗𝑗 𝑝𝑝𝑗𝑗與平面的距離。 𝑠𝑠𝑗𝑗 𝑝𝑝𝑗𝑗的權重值總和。 𝑛𝑛𝑛𝑛𝑛𝑛𝑛𝑛𝑗𝑗 𝑝𝑝𝑗𝑗的權重正規化值。 𝑠𝑠 所有𝑝𝑝𝑗𝑗的𝑠𝑠𝑗𝑗總和。 𝑛𝑛𝑛𝑛𝑛𝑛𝑛𝑛 所有𝑝𝑝𝑗𝑗的𝑛𝑛𝑛𝑛𝑛𝑛𝑛𝑛𝑗𝑗總和。 𝑅𝑅 搜尋半徑 𝑆𝑆(𝑢𝑢, 𝑣𝑣) 經過修補後的表面近似函數。 𝑓𝑓𝑖𝑖 尚未修補藉由𝑝𝑝𝑖𝑖近似的表面函數。 𝑝𝑝𝑖𝑖′(𝑢𝑢, 𝑣𝑣) 用於修補二維資料的第 i 個新取樣點。 𝐸𝐸(𝑆𝑆) 誤差函數。 𝑁𝑁𝑠𝑠 𝑝𝑝𝑖𝑖′新取樣點的數目。 𝑤𝑤𝑖𝑖 第 i 個新取樣點的權重函數。 𝑎𝑎 𝑆𝑆(𝑢𝑢, 𝑣𝑣)多項式的係數。 𝐵𝐵 利用𝑝𝑝𝑖𝑖所建立之矩陣。 𝐷𝐷𝑖𝑖 𝑝𝑝𝑖𝑖′到各𝑝𝑝𝑖𝑖的距離權重。 𝜙𝜙𝐷𝐷 𝑝𝑝𝑖𝑖′到𝑝𝑝𝑖𝑖的距離函數。 𝑛𝑛 𝑝𝑝𝑖𝑖′到𝑝𝑝𝑖𝑖的距離。 ℎ 𝑝𝑝𝑖𝑖′到𝑝𝑝𝑖𝑖的距離平均。 𝑝𝑝𝚤𝚤̇ 𝑝𝑝𝑖𝑖經由參數化攤平所獲得的新座標。 𝑤𝑤𝑖𝑖𝑗𝑗 所有鄰居權重值總和。 α 控制參數化攤平的速度。 𝜆𝜆𝑖𝑖𝑗𝑗 𝑝𝑝𝑖𝑖與𝑝𝑝𝑗𝑗距離權重值。 𝑃𝑃𝑗𝑗 𝑃𝑃𝑖𝑖的鄰居。
ix 𝑢𝑢𝑖𝑖 𝑝𝑝𝑖𝑖經由參數化攤平所得到的𝑥𝑥軸值。 𝑛𝑛𝑖𝑖 參數化攤平內部點數目。 𝑣𝑣𝑖𝑖 𝑝𝑝𝑖𝑖經由參數化攤平所得到的𝑦𝑦軸值。 𝐴𝐴 第 i 個內部點減掉其餘內部點距離權重矩陣。 𝑥𝑥 內部點座標矩陣。 𝑏𝑏 𝑝𝑝𝑖𝑖距離權重總合矩陣。
1
第一章 緒論
1.1 前言
基於對生活品質的需求和美好未來的憧憬，人類不斷地進行各種研究。電腦 是二十一世紀科技發展重要的推手之一。相較於過去因體積及執行速度造成運算 限制，現在大量資料的運算和傳輸彈在指之間就能完成，使許多資訊相關領域都 有突破性的發展。電腦視覺便是近年來不論在技術或應用上，不停有所突破的一 門主題，從無人駕駛車上的測距計算，到大範圍的環場三維重建，日新月異的創 新和拓展令人目不暇給。而在電腦視覺的三維重建中，要重建出與真實目標一樣 的模型，不僅需要高效能和精準的演算法，亦需要適度齊全的環境資訊作為輸入。 要同時滿足上述條件並不容易，故大多數三維重建的模型都具有些許瑕疵。因此， 讓這些具有瑕疵的三維模型，盡可能回到它完美的狀態，即是：三維模型修復。1.2 研究目的
如上所述，模型的完善取決於演算法的設計和輸入資料的適度齊全。三維重 建的演算法常難以完美，這是由於重建時所依靠的重要資訊：本體移動估測 （Ego-Motion）的誤差會隨尺度放大而不斷累積，當感測器失去自身在空間中的 確實位置，就會造成點雲（Point Cloud）的錯誤對齊。另一方面，對於環境資料 的蒐集，可能有一些難以拍攝到的角度平面，甚至環境中的移動物件和背景都會 影響本體移動估測的結果，造成重建瑕疵。提升重建演算法技術或對輸入資料的 要求固然可以有效改善結果，但多次重建與重錄的時間與人事成本，可能會造成 意想不到的額外負擔。因此本研究之目標在於，不經過再次重建和資料的重新蒐2 集，即達到對模型進行修飾和修補的效果。故三維模型修復通常有以下功能：去 除模型中的雜訊和離群部分、修剪重複冗餘的資訊以及修補模型中因資訊缺失造 成的破洞。本研究將提供一有效方法，對點雲模型進行妥善修復，並產生網格 （Mesh Grid）使模型具有連貫的表面後輸出。
1.3 研究方法
本研究分為三個部分。一為對資料離群值與雜訊的前處理，二為對資料遺失 部分進行資料修補，三為將修補結果再次進行網格重建。初始的點雲資料通常包 含許多離群值與雜訊，將影響到修補結果的正確性，故一開始會先進行平面偵測， 尋找各點的最適平面（Best Fitting Plane)，利用平面偵測結果對資料進行離群值 的濾除與點的平滑化。為了進行資料修補，需找出資料遺失的位置。首先以各自 點為中心，將範圍內之鄰近點依角度順序列為觀察集合。集合中依序的兩點與中 心會形成夾角，當夾角大於設定的閥值時便表示此處有空洞。最後利用倒投影將 平面資訊回復至三維空間，即完成修補。 於平面偵測步驟中，須找尋範圍內的鄰居，進行平面偵測。但遇到模型存在 銳利特徵（Sharp Features）時，可能會因點雲分部造成有錯誤。因此本研究加入 了穩定局部雜訊範圍偵測方法（Robust Local Noise Scale Estimation, RLNSE）的 距離限制。除此之外，鄰近點的選擇與資料修補邊界的搜尋息息相關。故本研究 建立 K-Ｄ樹（K-D Tree），並以局部移動性最小平方（Local Moving Least Squares, LMLS）重建法於資料修補與網格重建過程中，將三維資料修補問題簡化為二維 空間進行處理進行資料修補。1.4 論文架構
3 第一章為緒論，在前言之後介紹本研究的目的與研究方法，並簡介本論文架 構。第二章為文獻探討，將對表面重建技術發展進行討論，並著重於本論文研究 方法相關技術演化。第三章介紹本研究之三維模型修復方法。第四章針對第三章 所提出方法進行實驗結果與分析，第五章為本研究做下結論，並講述未來之研究 方向和展望。
4
第二章 文獻探討
本研究目的在於提供資料修補方法，藉此修補於三維模型重建時，因為資料 的缺失而影響到之後三維重建的品質，資料修補方法源自於三維重建中的表面重 建方法之一。
早期於 1998 年 Amenta 等人[1]提出利用沃羅諾伊圖（Voronoi Diagram）與 德勞內三角化（Delaunay Triangulation, DT）方法為基礎的表面重建方法，於存 在雜訊的情況下，表面重建的結果並不會失去原有的表面特徵，並且有平滑化的 表面的效果。2001 年 Carr 等人[2]藉由徑向基底函數（Radial Basis Functions, RBF） 與隱函數的方式表示一個高階多項式的平滑表面。2003 年 Alexa 等人[3]利用移
動最小平方方法（Moving Least Squares, MLS），將三維空間上的點投影至二維平
面進行局部性的表面近似運算。而 2006 年 Kazhdan 等人[4]利用泊松方程式 （Poisson Equation）獲得指示函數（Indicator Function）以表示表面函數。
近幾年表面重建的方法可以依據平滑程度限制，而將方法分為局部表面重建 （Local Surface Smoothness），如局部移動性最小平方（Local Moving Least Squares, LMLS），與全域表面重建（Global Surface Smoothness），如 RBF 與 Poisson Equation。上段敘述中，早期如 DT、RBF 與 Poisson Equation 等表面重建技術， 均著重於平滑表面的重建，到近幾年演變成如 2007 年 Wang 與 Oliveira[5]的 LMLS，與 2012 年 Doria 等人[6]的 Inpainting Depth Gradients 等針對資料缺失部 分進行重建。相較之下，局部表面重建方法有利於模型特徵的重建，而全域表面 重建會有過度平滑化與執行時間較久的問題。除了依平滑程度限制分類，表面重
建方法也可依據重建表面的演算法，分為：代數方法（Algebraic Methods）如 LMLS；
幾何計算方法（Computational Geometry Methods）如 Voronoi Diagram 與 DT； 以及隱函數方法（Implicit Functions Methods）像是 RBF 與 Poisson Equation。
5 雲做離群值的剃除與雜訊平滑化。在資料修補階段，先找出空洞位置與邊界，再 對合理範圍內的模型進行重建。於修補完成後再對點雲表面進行網格重建，獲得 三維模型。 圖 2-1 研究實作流程圖 而為了提升資料修補的正確性，本研究將著重於局部的資料修補方法，利用 局部移動性最小平方重建法（Local Moving Least Squares, LMLS）進行資料修補， LMLS 不僅執行效率高，並且重建結果與真實表面相似度較高。有關 LMLS 的相 關技術，最早於 1981 年 Lancaster 與 Salkauskas [7]提出三維的移動性最小平方 重建法，對點進行內插進行表面重建。2005 年 Fleishman 等人[8]先提出穩定的移 動性最小平方重建法，加入離群值的偵測，避免模型銳利特徵部分過度被平滑化。 2007 年 Wang 與 Oliveira[5] 提出使用局部移動性最小平方重建法修補表面空洞， 先利用網格資料找尋空洞邊界點，之後使用移動性最小平方重建法修補空洞。而
6 2011 年 Obermaier 等人[9]提出重建線條紋表面方法，藉由適應性的移動性最小 平方重建法，利用表面曲率（Curvature）與向量場（Vector Field）的近似，加上 DT 高效率網格重建的特性，顯示表面重建結果。 本研究以 Wang 與 Olivenian 的研究為基礎。然而，利用已經取樣過的網格 資料進行資料修補時輸入點的數量會太少，影響表面修補結果的正確性。所以本 研究改為先對於點雲資料進行修補，再產生網格。而 Wang 與 Olivenian 的研究 中，並沒有妥善剃除離群值，可能造成資料修補的誤差。故本研究將 LMLS 結 合平面偵測，並加入以 Obermaier 等人提出的適應性概念為基礎，在空洞邊界上 由外而內依序修補，藉此提升正確度。
7
第三章局部移動性最小平方重建法
本研究實作方法流程如下圖 3-1，本章節將敘述以 KD-Tree 為基礎，對資料 進 行 先 處理 ， 利 用 穩 定 局 部 雜 訊 範 圍偵 測 方 法（ Robust Local Noise Scale Estimation, RLNSE）偵測平面與剃除離群值，再使用雙邊濾波器（Bilateral Filter） 平滑雜訊。接著執行資料修補，先空洞偵測，利用偵測結果進行局部移動性最小 平方重建法（Local Moving Least Squares, LMLS）。最後利用經過參數化攤平 (Parameterization Mapping)的德勞內三角化（Delaunay Triangulation, DT）結果重 建網格表面。本章各節可分為，3.1 節說明執行資料修補時的前處理，包含利用 KD-Tree 找尋最近鄰居、平面偵測並剃除離群值與雜訊平滑化；3.2 節介紹本研 究資料修補方法；3.3 節說明利用修補點雲結果進行表面重建，藉此獲得網格模 型。
8
3.1 資料前處理
3.1.1 利用 KD-Tree 找尋最近鄰居
KD-Tree 可提供高維度空間中索引結構的搜尋，並能快速找尋範圍內 k 個最 近鄰居(K-Nearest Neighbors, KNN)[10][11]。它屬於一種二元搜尋樹（Binary
Search Tree, BST），當左子樹不為空值時，其左子樹所包含的點均會小於其根節 點之值，右子樹則為相反均會大於根節點之值，而且左右子樹也均為二元樹。 KD-Tree 建立是依據每次所挑選的維度，對該維度點進行大小比較，並依照 大小排列的結果進行切割。然而 KD-Tree 容易面臨兩個問題，一為維度順序的選 擇，二為左右子樹的建立是否達到平衡。解決維度選擇順序可以藉由計算每個維 度中點雲彼此之間的不變性（Invariance），找尋存在最大不變性的維度進行切割。 另一方面，為了使建立的 KD-Tree 接近平衡，可找尋各維度的中位數作為切割標 準。如圖 3-2[10]與圖 3-3[11]所示，將維度空間不斷區分為左右子樹，當空間無 法再次進行切割時，便將點值儲存於葉節點中。此外，為了提升搜尋速度，也會 對每個維度的切割點進行儲存。 圖 3-2 KD-Tree 空間切割圖[10]
9 圖 3-3 KD-Tree 樹狀圖[10] 當建立完 KD-Tree 之後，會以根節點為起始，依據其儲存的座標中位數與各 點位置進行比較，並向下走訪其子節點，直至訪問到葉節點時便停止搜尋。走訪 同時，也記錄途中節點彼此間的距離，依此距離判斷分支中是否存在更近的點， 若發現更近的節點，則改以走訪該子節點。 3.1.2 平面估測與剃除離群值
根據 Hesse 所提出描述線、平面空間或是超平面（Hyper Plane）的方程式[12] （Hesse Normal Form）中，利用一組單位法向量𝑛𝑛����⃑、點集合𝑃𝑃(𝑃𝑃𝑢𝑢 𝑋𝑋, 𝑃𝑃𝑌𝑌, 𝑃𝑃𝑍𝑍)與長度𝑑𝑑， 即可獲得三維空間平面：
𝑛𝑛𝑢𝑢
����⃑ ∙ 𝑃𝑃 = 𝑑𝑑 (3.1)
令𝑛𝑛�⃑ = (𝑛𝑛𝑥𝑥, 𝑛𝑛𝑦𝑦, 𝑛𝑛𝑧𝑧)為與𝑛𝑛����⃑同方向之法向量，而𝑑𝑑假設為平面與點距離，及法向量𝑛𝑛�⃑𝑢𝑢 長度，可以利用此兩項敘述獲得。
10 𝑑𝑑 = ‖𝑛𝑛�⃑‖ (3.2) 結合上式進而推導出： 𝑛𝑛�⃑ ‖𝑛𝑛�⃑‖ ∙ 𝑃𝑃 =‖𝑛𝑛�⃑‖ (3.3) 上式中描述具有單位向量之平面上，包含著所有投影長度為‖𝑛𝑛�⃑‖的投影平面。而 當找尋平面包含的點大於三個時，則需要找尋最適平面，即誤差函數最小化之平 面。此空間中平面到任一點距離定義如下： 𝛿𝛿(𝑛𝑛�⃑, 𝑃𝑃) =‖𝑛𝑛�⃑‖ ∙ 𝑃𝑃 −𝑛𝑛�⃑ ‖𝑛𝑛�⃑‖ (3.4) 利用點到平面距離平方總和可以得到下式，使其結果最小化可得到最適平面。 𝜙𝜙 = � 𝛿𝛿(𝑛𝑛�⃑, 𝑃𝑃)2 𝑘𝑘 𝑖𝑖=1 (3.5) 其中，𝜙𝜙為空間中𝑘𝑘個點到平面距離平方的總和。 如何使其式子(3.5)結果最小化，並可以得到平面法向量方法，可以利用主成 分分析（Principal Components Analysis, PCA）方式，對式子(3.5)進行特徵分解， 藉此獲得最小的距離平方總和。PCA 主要用於對於大量高維度數據的處理，由 Pearson 於 1901 年提出，利用降低數據維度方式，來提升處理效率。PCA 也可 視為一種降低維度的方法，藉由共變異數矩陣來進行特徵分解，藉此獲得特徵值 與特徵向量，而特徵向量即為該平面法向量。 然而，如下圖 3-4 所示，由於模型上存在如兔耳一般銳利的特徵，當利用 KNN 執行平面偵測時，可能會包含到不屬於該平面卻又位於搜尋範圍內的點， 進而影響到平面偵測結果。因此，可再藉由平面距離閥值限制以濾除離群值，提 升平面偵測的準確性。
11 圖 3-4 平面偵測圖 對於垂直距離限制，本研究基於 Li 等人[13]的 RLNSE 為基礎，進行可保留 銳利特徵的平面偵測。Li 等人的方法假設點點雲中任意點𝑃𝑃𝑖𝑖有 N 個鄰居，在鄰 居中任意挑選三個非共線的點產生候選平面，計算各點到候選平面距離並記錄距 離。上述方法執行 M 次之後將產生 M 個候選平面，挑選含有最小距離的平面， 即為𝑃𝑃𝑖𝑖的垂直距離限制。 方法執行次數 M 可以藉由如下式(3.6)計算，先假設 N 個鄰居中，可能蘊含 雜訊的比值𝜀𝜀𝑁𝑁，與最後剃除雜訊的比例𝑃𝑃̇去計算執行次數 M。 𝑃𝑃̇ = 1 − (1 − (1 − 𝜀𝜀̇ 𝑁𝑁)3)𝑀𝑀 (3.6) 𝑀𝑀 =log[1 − (1 − 𝜀𝜀log(1 − 𝑃𝑃̇) 𝑁𝑁)3] (3.7) 3.1.3 雜訊平滑化 除了離群值，雜訊也可會造成之後資料修補的結果有所偏差，故本研究使三 維的雙邊濾波器，對點雲資料進行平滑化，三維的雙邊濾波器是藉由修改二維的
12
雙邊濾波器方法獲得。二維的雙邊濾波器是由高斯濾波器演化而來的線性濾波器， 於 1998 年 Tomasic 與 Manduchi[14]提出，具有平滑化雜訊又可以保留邊緣特徵
的優點。圖 3-5 (a)中原始影像藉由使用高斯平滑（Gaussian Smoothing），於空間
域（Spatial Domain）與強度域（Intensity Domain）計算權重值圖 3-5 (b)，調整 鄰居對平滑結果的影響權重大小圖 3-5 (c)。 (a) (b) (c) 圖 3-5 雙邊濾波器示意圖[15] (a)為原始輸入影像、(b)為空間域與強度域濾波核而(c)為輸出影像 如下式(3.8)中𝛪𝛪̂(𝑢𝑢)為對於影像上的𝑝𝑝𝑖𝑖(𝑢𝑢, 𝑣𝑣)強度𝛪𝛪中點進行平滑化之結果： 𝛪𝛪̂(𝑝𝑝𝑖𝑖) = ∑𝑝𝑝𝑗𝑗∈𝑁𝑁𝑊𝑊𝑐𝑐(�𝑝𝑝𝑗𝑗− 𝑝𝑝𝑖𝑖�)𝑊𝑊𝑠𝑠(�𝛪𝛪(𝑝𝑝𝑖𝑖) − 𝛪𝛪�𝑝𝑝𝑗𝑗��)𝛪𝛪�𝑝𝑝𝑗𝑗� ∑𝑝𝑝𝑗𝑗∈𝑁𝑁𝑊𝑊𝑐𝑐(�𝑝𝑝𝑗𝑗− 𝑝𝑝𝑖𝑖�)𝑊𝑊𝑠𝑠(�𝛪𝛪(𝑝𝑝𝑖𝑖) − 𝛪𝛪�𝑝𝑝𝑗𝑗��) (3.8) 𝑁𝑁為鄰居數目，並包含點集合�𝑝𝑝𝑗𝑗�。𝑊𝑊𝑐𝑐(�𝑝𝑝𝑗𝑗− 𝑝𝑝𝑖𝑖�)為空間域差距權重計算， �𝑝𝑝𝑗𝑗− 𝑝𝑝𝑖𝑖�為𝑝𝑝𝑖𝑖到鄰居𝑝𝑝𝑗𝑗的距離，而𝑊𝑊𝑠𝑠(�𝛪𝛪(𝑝𝑝𝑖𝑖) − 𝛪𝛪�𝑝𝑝𝑗𝑗��)則為強度域差距權重計算， �𝛪𝛪(𝑝𝑝𝑖𝑖) − 𝛪𝛪�𝑝𝑝𝑗𝑗��為𝑝𝑝𝑖𝑖與鄰居𝑝𝑝𝑗𝑗影像強度差距，𝑊𝑊𝑐𝑐與𝑊𝑊𝑠𝑠計算如下式(3.11)與(3.12)： 𝑥𝑥𝑑𝑑 = �𝑝𝑝𝑗𝑗− 𝑝𝑝𝑖𝑖� (3.9) 𝑥𝑥𝐼𝐼 = �𝛪𝛪(𝑝𝑝𝑖𝑖) − 𝛪𝛪�𝑝𝑝𝑗𝑗�� (3.10) 𝑊𝑊𝑐𝑐(𝑥𝑥𝑑𝑑) = 𝑒𝑒−𝑥𝑥𝑑𝑑2⁄2𝜎𝜎𝑐𝑐2 (3.11) 𝑊𝑊𝑠𝑠(𝑥𝑥𝐼𝐼) = 𝑒𝑒−𝑥𝑥𝐼𝐼2⁄2𝜎𝜎𝑠𝑠2 (3.12)
13 �𝑝𝑝𝑗𝑗− 𝑝𝑝𝑖𝑖� < 𝜌𝜌 = ⌈2𝜎𝜎𝑐𝑐⌉ (3.13) 𝑊𝑊𝑐𝑐為標準的高斯濾波器（Gaussian Filter），藉由參數𝜎𝜎𝑐𝑐調整，當𝜎𝜎𝑐𝑐越大，距離越 遠的鄰居所獲得權重越大，而𝜎𝜎𝑐𝑐得設定須符合距離限制，及兩倍的𝜎𝜎𝑐𝑐大小需大於 距離差距如式(3.13)所述。而𝑊𝑊𝑠𝑠為特徵保留權重函數(Feature-preserving Weight Function)，藉由參數𝜎𝜎𝑠𝑠調整，當𝜎𝜎𝑠𝑠越大，當與鄰居的強度差距越大，所獲得的權 重越大。𝜎𝜎𝑐𝑐與𝜎𝜎𝑠𝑠兩者參數值越大時，所獲得的權重曲線結果會越平滑，影像會變 模糊，反之權重曲線結果則會變尖銳。 本研究所使用的三維的雙邊濾波器方法為 Fleishman 等人所提出的方法[16]， 由於局部的點雲資料中雜訊的占有比例可能過高，若直接利用 KNN 決定鄰居， 執行時會因為雜訊影響到平滑化結果。故本研究結合 3.1.2 所提的平面偵測結果， 將所獲得的各點鄰居帶入執行雜訊平滑化。 Fleishman 等人所提出的方法如下： 𝑝𝑝� = 𝑝𝑝𝚤𝚤 𝑖𝑖+ 𝑑𝑑 ∙ 𝑛𝑛 (3.14) 𝑝𝑝�為𝑝𝑝𝚤𝚤 𝑖𝑖進行完平滑化的點，Fleishman 等人的方法藉由將二維的雙邊濾波器中計算 影像強度方法，改成計算位移值𝑑𝑑，即平面與點距離。即新的表面上𝑝𝑝�存在之新𝚤𝚤 的平面與原有表面距離差值。利用計算鄰居與平面距離，並藉由上式(3.11)與下 式(3.15)計算𝑊𝑊𝑐𝑐與𝑊𝑊𝑠𝑠。 𝑊𝑊𝑠𝑠(𝑑𝑑𝑗𝑗) = 𝑒𝑒−𝑑𝑑𝑗𝑗2�2𝜎𝜎𝑠𝑠2 (3.15) 於三維座標中，除了𝑊𝑊𝑐𝑐與影像平滑化式子中，同樣計算與鄰居的距離差距，𝑊𝑊𝑠𝑠改 為計算點到平面距離如式(3.16)。利用𝑊𝑊𝑐𝑐(𝑡𝑡)與𝑊𝑊𝑠𝑠(ℎ)計算每個鄰居𝑝𝑝𝑗𝑗的權重值𝑠𝑠𝑗𝑗式 (3.17)與正規化值𝑛𝑛𝑛𝑛𝑛𝑛𝑛𝑛𝑗𝑗式(3.18)。
14 𝑑𝑑𝑗𝑗 = 〈𝑛𝑛, 𝑝𝑝𝑖𝑖− 𝑝𝑝𝑗𝑗〉 (3.16) 𝑠𝑠𝑗𝑗 = (𝑊𝑊𝑐𝑐× 𝑊𝑊𝑠𝑠) × 𝑑𝑑𝑗𝑗 (3.17) 𝑛𝑛𝑛𝑛𝑛𝑛𝑛𝑛𝑗𝑗 = 𝑊𝑊𝑐𝑐× 𝑊𝑊𝑠𝑠 (3.18) 𝑑𝑑 = (𝑠𝑠/𝑛𝑛𝑛𝑛𝑛𝑛𝑛𝑛) (3.19) 計算完權重值與正規化值後各自加總並相除，就可以獲得位移值𝑑𝑑，再將𝑑𝑑帶入 式子(3.14)即可獲得新的點𝑝𝑝�。 𝚤𝚤
3.2 資料修補
3.2.1 空洞位置定義 資料修補需先定義空洞位置。本研究以 Bendels 等人[17]與 Truong-Hong 等 人[18]提出方法為基礎，利用 KNN 與範圍限制獲得鄰居群，尋找可能在空洞邊 界上的點。因為單使用 KNN 在大量雜訊狀況時，會造成空洞定義錯誤。 圖 3-6(a)因為未限制搜尋範圍而使得點群包含到不同分類的點群，圖 3-6 (b)因為 加入了範圍限制而有正確的分類結果。 (a) (b) 圖 3-6 𝑘𝑘個最近點搜尋圖[17] [18] 。15
(a)為沒有限制範圍，𝑘𝑘個最近點搜尋結果、(b)為有限制範圍，限𝑘𝑘個最近點搜尋
結果。
完成鄰居搜尋之接著以角度偵測空洞(Hole detection by angle)將點雲中之一 點 P 與其鄰居投影至二維平面之後，以 P 之平面投影位置為中心並展開一限制 範圍，在範圍內尋找平面上所有大於閥值之角度視為空洞。 各點在範圍內執行完角度檢查後，需判別不同一個空洞邊界上的點。圖 3-7 中，藍點為可能在空洞邊界上，此時將紅點作為起始中心向外展開半徑為 R 搜 尋範圍，紅色圓形內的點，均判斷為同一個空洞的邊界點。接著再依據新搜尋到 的點，如圖 3-7 中的橘點，同樣再次以半徑為 R 的紅色圓形為搜尋範圍，直到沒 有再搜尋到新的點，依此區分多個空洞點群。接著計算凸包(Convex Hull)獲得空 洞邊界順序，劃分修補範圍，移除可能為物體表面邊界的邊界點，如圖 3-8(a)中 空洞位於凸包內即為空洞邊界，圖 3-8(b)當位於凸包上時則判斷為空洞。 圖 3-7 空洞邊界點集合分類圖
16 (a) (b) 圖 3-8 空洞過濾示意圖 (a)為物體邊界、(b)為空洞邊界 3.2.2 局部移動性最小平方重建法（LMLS） LMLS 主要被用於計算點雲的表面近似函數𝑆𝑆(𝑢𝑢, 𝑣𝑣)，將點群與平面之垂直距 離視為誤差，透過最小化此誤差來獲得表面函數。與傳統的最小平方法比較， LMLS 加入權重限制計算不同遠近距離的鄰居影響程度，改善複雜表面造成的缺 失。此外，在計算 LMLS 過程中宜加入[9]中提到的適應性概念，將前幾次重建 時新增的點作為最新一次重建的輸入，讓重建的表面更符合原有曲面結構。如圖 3-9 中，以凸包座標的中心計算與各點的距離為順序，由外往內展開一範圍 R（紅 色圓形）進行修補，納入新增的修補點，即圖中的綠點，能使計算表面時有更豐 富的參考資訊。
17 圖 3-9 LMLS 適應性修補圖 移動最小平方計算方式如下。為了使原有的表面函數𝑓𝑓𝑖𝑖與表面重建之後的 𝑆𝑆(𝑢𝑢, 𝑣𝑣)表面函數誤差結果最小，即 最小化將𝑆𝑆(𝑢𝑢, 𝑣𝑣)代入𝑝𝑝′(𝑢𝑢, 𝑣𝑣)時得到的距離d。 𝑆𝑆(𝑢𝑢, 𝑣𝑣)可由誤差函式獲得，此誤差函式𝐸𝐸(𝑆𝑆)可表示為： 𝐸𝐸(𝑆𝑆) = � 𝑤𝑤𝑖𝑖(𝑝𝑝𝑖𝑖′)(𝑆𝑆(𝑝𝑝𝑖𝑖) − 𝑓𝑓𝑖𝑖)2 𝑁𝑁𝑠𝑠 𝑖𝑖=1 (3.20) 為了得到最小誤差之結果，對表面係數進行偏微分解方程式，而𝑆𝑆(𝑢𝑢, 𝑣𝑣)可以藉由 多項式近似獲得下式。 𝑆𝑆(𝑢𝑢, 𝑣𝑣) = 𝑎𝑎0+ 𝑎𝑎1𝑢𝑢 + 𝑎𝑎2𝑣𝑣 + 𝑎𝑎3𝑢𝑢2+ 𝑎𝑎4𝑣𝑣2+ 𝑎𝑎5𝑢𝑢𝑣𝑣 (3.21) 其中係數 a 可以分解如下式： 𝑎𝑎(𝑝𝑝𝑖𝑖′) = (𝐵𝐵𝑤𝑤𝑖𝑖(𝑝𝑝𝑖𝑖′)𝐵𝐵𝑇𝑇)−1𝐵𝐵𝑤𝑤𝑖𝑖(𝑝𝑝𝑖𝑖′)𝑓𝑓𝑖𝑖 (3.22) 𝐵𝐵為利用輸入點所建立之矩陣，可以表示為下式：
18 𝐵𝐵 = ⎣ ⎢ ⎢ ⎢ ⎢ ⎡ 𝑢𝑢1 ⋯ 1 1 ⋯ 𝑢𝑢𝑛𝑛 𝑣𝑣1 ⋯ 𝑣𝑣𝑛𝑛 𝑢𝑢12 ⋯ 𝑢𝑢𝑛𝑛2 𝑣𝑣12 ⋯ 𝑣𝑣𝑛𝑛2 𝑢𝑢1𝑣𝑣1 ⋯ 𝑢𝑢𝑛𝑛𝑣𝑣𝑛𝑛⎦ ⎥ ⎥ ⎥ ⎥ ⎤ (3.23) 𝑤𝑤𝑖𝑖(𝑝𝑝𝑖𝑖′)為第 i 個新取樣點，於移動最小平方中的權重函數，𝑤𝑤𝑖𝑖(𝑝𝑝𝑖𝑖′)可以表示為對 角矩陣形式，如下式其式子(3.24)中的值𝐷𝐷𝑖𝑖，即為第𝑖𝑖個新取樣點到𝑝𝑝𝑖𝑖距離函數𝜙𝜙。 𝑤𝑤𝑖𝑖(𝑝𝑝𝑖𝑖′) = 𝑑𝑑𝑖𝑖𝑎𝑎𝑑𝑑(𝐷𝐷1, 𝐷𝐷2, 𝐷𝐷3, ⋯ 𝐷𝐷𝑛𝑛) (3.24) 𝐷𝐷𝑖𝑖 = 𝜙𝜙𝐷𝐷(‖𝑝𝑝𝑖𝑖′− 𝑝𝑝𝑖𝑖‖) (3.25) 本研究將使用[19]中的距離函數𝜙𝜙𝐷𝐷，來表示權重函數來計算輸入點與新取樣點彼 此影響程度，影響程度會因點群之間的距離越遠，對於新取樣點 d 的計算結果影 響越小。距離函數𝜙𝜙𝐷𝐷可以表示為下式，𝑛𝑛為𝑝𝑝𝑖𝑖′到輸入點的距離，ℎ為𝑝𝑝𝑖𝑖′到所有輸入 點距離的平均距離。 𝜙𝜙𝐷𝐷(𝑛𝑛) = 𝑒𝑒𝑥𝑥𝑝𝑝 (−𝑛𝑛 2 ℎ2 ) (3.26)
3.3 德勞內三角化網格模型重建
在網格化的階段，本研究先對點雲進行二次的平面偵測，將網格重建問題由 三維轉換為二維重建，再對各個點群進行參數化攤平，參數化攤平時，會利用 DT 找尋用於計算權重值的鄰居，最後對參數化攤平結果使用 DT 執行網格重建。 DT 是由 Delaunay Boris 於 1934 年提出網格三角化的方法，DT 根據[20]所描述具 有下列原則與特性，下圖 3-10 為網格模型重建流程。19 圖 3-10 網格模型重建步驟圖 3.3.1 德勞內三角化 DT 的產生主要有兩大原則，即最大空圓原則與最大化最小角原則。最大空 圓原則限制三角形所形成之外接圓內，除了三角形點之外不可以包含其他點。如 圖 3-11(a)中，由紅色點的三角形所形成的外接圓包含綠色的第四點，並沒有符 合最大空圓原則，而 圖 3-11 (b)符合了最大空圓原則。為了符合最大空圓原則，
可以利用局部優化方法（Local Optimization Procedure, LOP），LOP 為 Lawson 提 出，確保形成之網格符合 DT 最大空圓原則如圖 3-11。LOP 先找尋具有共同邊之 三角形，將兩個三角形相連形成多邊形，並找尋多邊形中可存在之最大空圓，最 後檢查圓內是否存在第四個點，以調整多邊形對角線位置。另外，最大化最小角 原則如圖 3-12 中則是確認四邊形中，是否會因為對角線互換，使六個內角最小 值不會再增加。因為產生的三角形不會有內角過小問題，而使網格模型結構更加 穩定。除上述 DT 的原則外，DT 本身也應具有唯一性與區域性。唯一性使每次 起始點的挑選即使不一樣，仍導致相同的網格重建的結果；而區域性則限制當三 角形新增、刪除與移動時，只有相鄰三角行會受影響。
20 (a) (b) 圖 3-11 最大空圓原則圖 (a)為未符合最大空圓原則而(b) 有符合最大空圓原則。 圖 3-12 最大化最小角原則 DT 可以利用 Lawson 於 1972 提出的遞增方法完成[21]。一開始先對未進行 DT 點群建構可包含所有點之多邊形。如 圖 3-13(a)，點群依序納入多邊形後， 依此切割出多個三角形。任選一點為控制點，並找尋控制點位於哪一個三角形中。 連接控制點與三角形後，便形成新的三角形如圖 3-13(b)。之後再檢查三角形外 接圓是否符合規則，未符合規則者進行局部優化，並對對角線進行交換動作，直 到全部外接圓均符合規則為止如圖 3-13(c)。重複執行上述動作，直到所有點均 納入多邊形，可產生如圖 3-13(d)圖中的三角網格。最後移除掉多邊形邊界，即 可以獲得 DT 網格結構。
21 (a) (b) (c) (d) 圖 3-13 逐點插入法步驟圖 (a)為多邊形中割出多個三角形情況、(b)為控制點與附近三角形相連情況、(c) 為執行局部優化結果而(d)為逐點插入法執行完成結果。 3.3.2 參數化攤平 本研究利用 Chen[22]所採用的網格重建方法，先進行參數化攤平方法，再對 參數化攤平所得到的二維平面使用 DT 重建網格。其中參數化攤平方法，意在先 確保平面點群倒投影至三維空間時能與真實特徵相符，才正式執行 DT 重建。參 數化攤平方法主要利用 Floater 與 Reimers[23]和 Hormann 與 Reimers[24]的做法，
將點群依據邊界點所形成的單位圓如，並將所有內部點𝑝𝑝7、𝑝𝑝8與𝑝𝑝9沿著黑色箭頭
22 序調整內部點位置，直到內部點位置不再改變，如圖 3-14 所示。 圖 3-14 參數化攤平初始化示意圖 調整內部點方式如下式： 𝑝𝑝𝚤𝚤̇ = 𝑝𝑝𝑖𝑖 + 𝛼𝛼 � 𝑤𝑤𝑖𝑖𝑗𝑗�𝑝𝑝𝑖𝑖 − 𝑝𝑝𝑗𝑗� 𝑁𝑁 𝑗𝑗∈{𝑝𝑝𝑖𝑖}=𝑁𝑁 (3.27) 𝑝𝑝𝚤𝚤̇ 為𝑝𝑝𝑖𝑖經由迭代所獲得的新的內部點，假設各內部點𝑝𝑝𝑖𝑖有𝑁𝑁個鄰居𝑝𝑝𝑗𝑗，藉由計算 三維座標中𝑝𝑝𝑖𝑖與鄰居𝑝𝑝𝑗𝑗距離並給予權重值𝑤𝑤𝑖𝑖𝑗𝑗進行位移，其中α為控制位移快慢的 參數。𝑤𝑤𝑖𝑖𝑗𝑗權重為各鄰居權重值𝜆𝜆𝑖𝑖𝑗𝑗加總，並正規化後獲得的數值。權重計算如下 式： 𝑤𝑤𝑖𝑖𝑗𝑗 = ∑ 𝜆𝜆𝑖𝑖𝑗𝑗 𝜆𝜆 𝑖𝑖𝑗𝑗 𝑁𝑁 𝑗𝑗∈{𝑝𝑝𝑖𝑖} (3.28) 𝜆𝜆𝑖𝑖𝑗𝑗 = 1 �𝑃𝑃 𝑖𝑖− 𝑃𝑃𝑗𝑗� � (3.29) 參數化攤平中利用 DT 決定鄰居，只需計算各點於 DT 中有連接邊的點，其 餘點的權重均為零。中，𝑝𝑝7所需要計算權重的鄰居，並不包含𝑝𝑝3與𝑝𝑝4，如此可避
23 免過遠的鄰居影響到參數化攤平結果。而為了加速參數化的計算，[23]中提出利 用解聯立方程式的方法計算內部點座標如下式： 𝑢𝑢𝑖𝑖 = ∑𝑗𝑗∈𝑁𝑁𝑖𝑖𝜆𝜆𝑖𝑖𝑗𝑗𝑢𝑢𝑗𝑗, 𝑖𝑖 = 1, … . , 𝑛𝑛𝑖𝑖 (3.30) 𝑣𝑣𝑖𝑖 = ∑𝑗𝑗∈𝑁𝑁𝑖𝑖𝜆𝜆𝑖𝑖𝑗𝑗𝑣𝑣𝑗𝑗, 𝑖𝑖 = 1, … . , 𝑛𝑛𝑖𝑖 (3.31) 聯立方程式各對𝑝𝑝𝚤𝚤̇ (𝑢𝑢, 𝑣𝑣)x 軸與 y 軸進行計算，各軸計算方法如上式(3.30)與(3.31)。 假設共有𝑛𝑛𝑖𝑖個內部點，與 N 個邊界點，以 x 軸為例，𝑢𝑢𝑖𝑖為內部點與邊界點權重總 合，所以可以藉由𝑢𝑢𝑖𝑖減去由內部點計算的權重，來獲得邊界點計算的權重。根據 式(3.33)可將𝑢𝑢𝑖𝑖減去由內部點計算的權重結果，重設為𝑛𝑛𝑖𝑖× 𝑛𝑛𝑖𝑖矩陣矩陣𝐴𝐴與 𝑛𝑛𝑖𝑖 × 1 矩陣𝑥𝑥，利用已知數𝑏𝑏邊界點權重總合，解聯立方程式得到𝑢𝑢。 𝑢𝑢𝑖𝑖 − ∑𝑛𝑛𝑗𝑗=1𝜆𝜆𝑖𝑖𝑗𝑗𝑢𝑢𝑗𝑗 = ∑𝑁𝑁𝑗𝑗=𝑛𝑛+1𝜆𝜆𝑖𝑖𝑗𝑗𝑢𝑢𝑗𝑗, 𝑖𝑖 = 1, … . , 𝑛𝑛 (3.32) 𝐴𝐴𝑥𝑥 = 𝑏𝑏 (3.33) 圖 3-15 參數攤平化示意圖 以圖 3-15 為例，紅點為邊界點，黑點則是未知的內部點，根據上述計算方法與 (3.30)、(3.31)與(3.32)式子，可將式子區分為式子(3.36)、(3.37)與(3.38)，其中𝑝𝑝7、
24 𝑝𝑝8與 𝑝𝑝9未知數藉由解聯立方程式獲得。 �𝑝𝑝𝑝𝑝78 = 𝑝𝑝= 𝑝𝑝11𝑤𝑤𝑤𝑤11+ 𝑝𝑝+ 𝑝𝑝22𝑤𝑤𝑤𝑤22+ 𝑝𝑝+ 𝑝𝑝33𝑤𝑤𝑤𝑤33+ 𝑝𝑝+ 𝑝𝑝44𝑤𝑤𝑤𝑤44+ 𝑝𝑝+ 𝑝𝑝55𝑤𝑤𝑤𝑤55+ 𝑝𝑝+ 𝑝𝑝66𝑤𝑤𝑤𝑤66+ 𝑝𝑝+ 𝑝𝑝78𝑤𝑤𝑤𝑤78+ 𝑝𝑝+ 𝑝𝑝99𝑤𝑤𝑤𝑤99 𝑝𝑝9 = 𝑝𝑝1𝑤𝑤1 + 𝑝𝑝2𝑤𝑤2+ 𝑝𝑝3𝑤𝑤3+ 𝑝𝑝4𝑤𝑤4 + 𝑝𝑝5𝑤𝑤5+ 𝑝𝑝6𝑤𝑤6+ 𝑝𝑝7𝑤𝑤7+ 𝑝𝑝8𝑤𝑤8 , (3.34) �𝑝𝑝𝑝𝑝78− 𝑝𝑝− 𝑝𝑝78𝑤𝑤𝑤𝑤78− 𝑝𝑝− 𝑝𝑝99𝑤𝑤𝑤𝑤99= 𝑝𝑝= 𝑝𝑝11𝑤𝑤𝑤𝑤11+ 𝑝𝑝+ 𝑝𝑝22𝑤𝑤𝑤𝑤22+ 𝑝𝑝+ 𝑝𝑝33𝑤𝑤𝑤𝑤33+ 𝑝𝑝+ 𝑝𝑝44𝑤𝑤𝑤𝑤44+ 𝑝𝑝+ 𝑝𝑝55𝑤𝑤𝑤𝑤55+ 𝑝𝑝+ 𝑝𝑝66𝑤𝑤𝑤𝑤66 𝑝𝑝9− 𝑝𝑝7𝑤𝑤7− 𝑝𝑝8𝑤𝑤8 = 𝑝𝑝1𝑤𝑤1 + 𝑝𝑝2𝑤𝑤2+ 𝑝𝑝3𝑤𝑤3+ 𝑝𝑝4𝑤𝑤4 + 𝑝𝑝5𝑤𝑤5+ 𝑝𝑝6𝑤𝑤6 (3.35) A = �−𝑤𝑤17 −𝑤𝑤18 −𝑤𝑤−𝑤𝑤99 −𝑤𝑤7 −𝑤𝑤8 1 � (3.36) x = [𝑝𝑝7 𝑝𝑝8 𝑝𝑝9]𝑇𝑇 (3.37) b = �𝑝𝑝𝑝𝑝11𝑤𝑤𝑤𝑤11 + 𝑝𝑝+ 𝑝𝑝22𝑤𝑤𝑤𝑤22+ 𝑝𝑝+ 𝑝𝑝33𝑤𝑤𝑤𝑤33+ 𝑝𝑝+ 𝑝𝑝44𝑤𝑤𝑤𝑤44+ 𝑝𝑝+ 𝑝𝑝55𝑤𝑤𝑤𝑤55+ 𝑝𝑝+ 𝑝𝑝66𝑤𝑤𝑤𝑤66 𝑝𝑝1𝑤𝑤1 + 𝑝𝑝2𝑤𝑤2+ 𝑝𝑝3𝑤𝑤3+ 𝑝𝑝4𝑤𝑤4+ 𝑝𝑝5𝑤𝑤5+ 𝑝𝑝6𝑤𝑤6 � (3.38) 本研究將依據上述所執行各點群參數攤平化結果，產生二維的 DT 網格，再 將所獲得的各網格索引點，再直接帶入三維座標來獲得網格模型。但當執行 DT 網格重建時，靠近邊界部分的重建結果可能會有錯誤，所以將以各點群產生網格 改為以各點產生，只保留影連接該點的網格。如圖 3-16(b)，剔除圖中包含紅線 的網格，即利用 3.3.1 所提到的區域性規範，當對 DT 的三角形進行刪除動作時， 只會影響鄰居，而不是整體結構的特性。藉由移除不合適表面的網格，提升網格 模型平滑程度，使網格模型更接近原本物體表面如圖 3-16(a)。
25 (a) (b) 圖 3- 16 產生錯誤網格情況 (a)為當平面具有大幅度彎曲情況而(b)為平面於二維空間重建情況。 本章節之資料修補方法，如前段所述，將點雲經前處理優化後，進行空洞定 義與修補，最後以參數攤平化再輸出為網格即完成。下一章節將呈現參數與資料 修補的實驗結果。
26
第四章實驗結果與分析
本章節將利用第三章所提出的方法，進行實驗與分析。使用的資料來自 The Stanford 3D Scanning Repository 實驗室點雲資料，與利用 Agisoft PhotoScan 重建 實驗室空拍機所拍攝高雄大學影像點雲模型。於本章節實驗結果中，會比對基礎 資料修補流程與加入第三章方法結果做比較。基礎資料修補方法藉由找尋各點鄰 居，並對所存在之平面進行表面近似，接著對表面重新取樣，最後將取樣倒投影 回三維空間，藉此完成資料修補動作。而實驗結果實作平台為使用 Intel i7-2600 CPU，記憶體為 16GB，並利用 Matlab2015b 軟體執行實驗結果。 而本研究除了對於資料修補方法進行流程整合之外，並於各步驟加入方法提 升重建結果的完整性。於 3.1.3 節加入 3.1.2 的平面偵測結果，避免只利用 k 個最 近鄰居(K-Nearest Neighbors, KNN)決定鄰居時，當局部的點雲雜訊的占有比例= 過高，而造成平面偵測錯誤。而在 3.2.1 節中加入空洞過濾方法，解決只單一使 用角度辨別空洞邊界，可能造成的物體邊界的誤判。最後於 3.3.2 節中加入網格 修剪方法，使網格重建結果更貼近原有模型表面。 實驗結果依據加入方法的先後順序呈現各節實驗結果。於 4.1 節中加入剃除 離群值與雜訊平滑化，4.2 節加入空洞偵測與判斷是否存在多個空洞，4.3 節為利 用局部移動性最小平方重建法（Local Moving Least Squares, LMLS）進行空洞修 補，最後利用參數攤平化來執行德勞內三角化（Delaunay Triangulation, DT）網 格重建。因為本研究主要針對資料修補方法進行整合，並提出資料修補流程，藉 由不同方法之間互相彌補各自可能存在的缺失，所以實驗結果呈現方式，才會藉 由比較各方法加入基礎資料修補後所帶來的效果。
27
4.1 離群值剃除與雜訊平滑化
如果點雲本身就存有雜訊，在資料修補過程中，便會因為參考到雜訊位置而 造成最後修補結果錯誤。下圖 4-1 中為高雄大學行政大樓牆面點雲模型，有些窗 戶邊框因為雜訊關係，而凹凸不平。經過離群值剃除與雜訊平滑化後，任意挑選 圖 4-1(a)中的一點，如圖 4-1(a)中的紅點，計算離群剃除與雜訊平滑化執行前後 結果。計算結果為圖 4-1(b)，其中左圖藍色點為原始點群，右圖中的藍色點則為 經過剃除離群值與雜訊平滑化點群，紅色點均為點群投影結果，當未執行前點到 平面距離較執行後遠。 (a)28 (b) 圖 4-1 點雲模型離群值剃除與雜訊平滑結果比較圖 (a)為行政大樓牆面點雲模型、(b)為未執行離群值剃除與雜訊平滑化點群結 果，與有執行離群值剃除與雜訊平滑化點群結果。
4.2 空洞位置定義
當找尋空洞邊界上的點時，如果誤把建築物邊界視為空洞邊界時，可以藉由 定義空洞邊界上的點，區分不同空洞範圍，藉此避只使用 KNN 近似平面進行資 料修補，造成資料修補犯錯誤。下圖 4-2 中，為了計算重建結果誤差值而產生的 曲面模型，會於 4.3 節表示於 LMLS 修補結果的誤差值，當進行空洞邊界上的點 偵測時，如圖 4-2 邊界的部分被誤判為空洞邊界上的點，修補完會出現如圖 4-2(b) 中多餘的點雲資料，這是因為修補範圍依據鄰居所在位置劃分，所以造成修補錯 誤。29 (a) (b) 圖 4-2 曲面模型空洞定義結果比較圖 (a) 執行空洞過濾結果、(b) 為只針對各點所存在 KNN 之平面進行資料修補。 於下圖 4-3 中兔子點雲模型，當進行空洞邊界上的點偵測時，如圖 4-3 靠近 耳朵的地方被判斷為空洞邊界上的點，造成修補範圍錯誤，所以耳朵的部分會出 現多餘的點雲資料圖 4-3(b)。
30 (a) (b) 圖 4-3 兔子點雲模型空洞定義結果比較圖 (a) 執行空洞過濾結果、(b) 為只針對各點所存在 KNN 之平面進行資料修補。 下圖 4-4 為行政大樓牆面點雲模型空洞偵測結果，其中會判定圖 4-4 (a)中為 經過空洞過濾的修補結果。如果只單用 KNN 近似平面進行資料修補，會產生如 圖 4-4 (b)中，建築物邊界因為錯誤的修補於呈現鋸齒狀的結果，但藉由空洞過濾 區分不同空洞並劃分修補範圍，可以剃除避免於建築物邊界的地方誤判為需要修
31 補的範圍。 (a) (b) 圖 4-4 行政大樓點雲空洞定義結果比較圖 (a)執行空洞過濾結果、(b)為只針對各點所存在 KNN 之平面進行資料修補。
4.3 局部移動性最小平方重建法（LMLS）
圖 4-5 為曲面模型進行資料修補結果，圖 4-5 (a)為曲面初始未修補資料，經 過 LMLS 後，可以獲得圖 4-5(b)的結果，而圖 4-5(c)則是利用 MLS 進行修補的 結果。相較於圖 4-5(b)，圖 4-5(c)中因為附近曲線弧度不一致，所以當 MLS 進行32
修補時一次參考所有邊界點，就可能帶入錯誤的點群位置資料，造成於 MLS 於 近似曲面時而有誤差。
(a)
33 (c) 圖 4-5 曲面模型 LMLS 與 MLS 結果比較圖 (a)為原始未進行修補影像、(b) 為利用 LMLS 進行資料修補結果而(c)為利用 MLS 近似平面藉此修補空洞。 圖 4-6 為兔子模型底部區域進行資料修補結果，圖 4-6 (a)為初始未修補資料， 經過 LMLS 後，可以獲得圖 4-6(b)的結果，而圖 4-6(c)則是利用 MLS 進行修補 的結果。相較於圖 4-6(b)，圖 4-6(c)底部的修補結果並未保有原本曲線，因為 MLS
34
進行修補時一次參考所有空洞上的邊界點，所以可能參考到附近點，造成於 MLS 於近似曲面時而有誤差。
(a)
35 (c) 圖 4-6 兔子模型 LMLS 與 MLS 結果比較圖 (a)為原始未進行修補影像、(b)為利用 LMLS 進行資料修補結果而(c)為利用 MLS 近似平面藉此修補空洞。 下圖 4-7 為行政大樓點雲模型進行資料修補結果，圖 4-7 (a)為初始未修補資 料，經過 LMLS 後，可以獲得圖 4-7(b)的結果，而圖 4-7(c)則是利用 MLS 進行 修補的結果。相較於圖 4-7(b)，圖 4-7(c)錯誤修補資料，並且超過原有模型邊界。
36
(a)
37 圖 4-7 行政大樓點雲模型 LMLS 與 MLS 結果比較圖 (a)為原始未進行修補影像、(b)為利用 LMLS 進行資料修補結果而(c)為利用 MLS 近似平面藉此修補空洞 本研究利用產生三個標準模型，並於模型上產生空洞，再執行 LMLS 修補 空洞，藉此計算修補結果的誤差值。因為兔子與建築物模型本身，可能因為於模 型重建時的誤差累積，而造成計算修補誤差的錯誤，所以本研究改以產生標準模 型進行計算，如下圖 4-8、圖 4-9 與圖 4-10 中的曲面、三角錐與球體模型，利用 刪除已知的位置，並進行資料修，最後對入修補結果計算誤差，誤差計算結果為 下表 4-1 中。
38 (a) (b) (c) 圖 4-8 曲面模型修補圖 (a)為原始模型、(b)為產生空洞之後的模型而(c)為利用 LMLS 修補結果 (a) (b) (c) 圖 4-9 三角錐模型修補圖 (a)為原始模型、(b)為產生空洞之後的模型而(c)為利用 LMLS 修補結果 (a) (b) (c) 圖 4-10 球體模型修補圖 (a)為原始模型、(b)為產生空洞之後的模型而(c)為利用 LMLS 修補結果
39 表 4-1 實驗結果誤差值表 曲面模型 三角錐模型 球體模型 最小誤差 0.00028 0.000034 0.0025 最大誤差 0.0309 0.0538 0.0243 平均誤差 0.0139 0.0217 0.01 中位數誤差 0.0127 0.0192 0.0085 下表 4-2 設定為本研究於各模型搜尋鄰居，與找尋空洞邊界點所需要設定的 閥值，平面偵測中於 KNN 步驟的所需的最大鄰居數目，並限制最大搜尋半徑與 角度，藉此判斷平面上所存在最大角度是否符合空洞邊界點條件。 表 4-2 實驗數據設定表 點雲數目 鄰居數目 搜尋半徑 搜尋角度 球體 2436 32 32 60 三角錐 18001 32 16 90 曲面 9009 32 0.8 90 兔子 360325 128 256 90 高雄大學行政大樓 176712 128 64 60
4.4 參數化攤平網格重建結果
如下圖 4-11(a)擷取行政大樓紅色圓圈區域進行網格重建結果，對比參數化 攤平對於網格重建結果的影響，圖 4-11(b)中因為未執行參數化攤平，使得點群 投影至二維時部分位置錯誤，而於產生網格時連接錯誤的點，導致倒投影回三維 空間網格重建產生大量錯誤。相對於圖 4-11(c)的部分較為整齊，與原有模型表40
面較為相符合。
(a)
41
(c)
圖 4-11 參數化攤平選擇結果比較圖
(a)為行政大樓點雲模型、(b)為直接利用最適平面進行網格重建而(c)為利用參數 化攤平進行網格重建
42
第五章 結論與未來研究方向
本研究目的在於利用資料修補方法重建模型，減少因遺失於需重新蒐集資料 可能需要花費的時間與人力成本。本方法不只侷限於物體，也能利用於建築模型 的修復。故本研究提出保留物體特徵之資料修補流程，將各方法做結合，使模型 資料修補結果更加完整。5.1 結論
本研究以點雲模型資料作為輸入，並先對資料進行前處理:以 k 個最近鄰居 (K-Nearest Neighbors, KNN)並限制搜尋範圍，對各點進行平面偵測。為了避免銳 利特徵而造成的偵測錯誤，找尋最適平面時，加入點到平面限制，此距離限制藉 由穩定局部雜訊範圍偵測方法（Robust Local Noise Scale Estimation, RLNSE）計 算獲得並利用於離群值剃除。除了平面偵測與離群值剃除，也利用雙邊濾波器 （Bilateral Filter）平滑降低表面雜訊。 在修補流程中，要確認每個點於平面上與鄰居相連之角度，是否有超過設定 的罰值，依此判斷該點是否為空洞邊界點。接著利用空洞過濾方法，先設定搜尋 範圍區分不同空洞點群，並結合凸包(Convex Hull)，區分是空洞還是物體邊界點。 尋找空洞位置之後，使用 LMLS 資料修補方法，增加新的點群來填補空洞。為 了提升修補資料的完整，提出適應性的方式，將每次新增的修補點納入參考選擇， 讓修補表面更符合原有曲面。 完成資料修補之後，最後以先處理所偵測平面進行網格重建，以各點的為中 心產生網格，先於各平面計算凸包判別外部點與內部點，再加入參數化攤平方式， 使內部點群可以重新調整二維資料點位置。使點群再次產生 DT 網格時，可以更 符合於三維空間中的相對位置，讓產生的網格模型與原有表面更加契合。但於實43 驗結果中對於模型空洞位置定義部分，但當模型存在類似空洞表面結構時，可能 會因為本研究所使用空洞判斷方法，而造成修補錯誤。另外於空洞修補步驟時， 因為是依據空洞邊界點作為空洞修補依據，所以當空洞範圍涵蓋曲面過於複雜時， 雖然可以填補原有表面使其與模型契合，但與原有表面圖案仍存在誤差。
5.2 未來研究方向
之後可嘗試藉由二維影像資料，處理空洞偵測可能造成判斷錯誤的問題，利 用也於二維影像偵測空洞，藉此再次確認於三維空間中的空洞偵測結果是否有誤 。也可以利用偵測二維影像中物體曲線，輔助於空洞修補時，遇到表面過愈複雜 時可能造成修補不夠精確。除了上述的問題，為了執行適應性空洞修補方法，而 造成執行時間過久，也是之後需要改善的重點之一。 二維影像資料的使用，也可以用於提升各點相對顏色的正確性，藉而影像得 知正確的模型色彩紋路，對網格模型進行材質拼貼，提升模型完成度。而資料修 補方法也可以跟相關重建系統做結合，並銜接 3D 模型巢狀結構技術，延伸應用 於 3D 列印輸出，未來除了對於現有問題深入討論，也會延伸該方法的應用，使 本研究可以有廣泛的應用，不只侷限於三維模型資料修補方法研究。44
參考文獻
[1] N. Amenta, M. Bern, and M. Kamvysselis, “A New Voronoi-Based Surface Reconstruction Algorithm,” in Proceedings of the 25th annual conference on
Computer graphics and interactive techniques, 1998, pp.415-421.
[2] C. J. Carr, K. R. Beatson, B. J. Cherrie, J. T. Mitchell, R. W. Fright, C. B. McCallum, and R. T. Evans, “Reconstruction and Representation of 3D Objects with Radial Basis Functions,” in Proceedings of the 28th annual conference on
Computer graphics and interactive techniques, 2001, pp. 67-76.
[3] M. Alexa, J. Behr, D. Cohen-Or, S. Fleishman, D. Levin, and T. C. Silva, “Computing and rendering point set surfaces,” IEEE Transactions on
visualization and computer graphics, vol. 9, no. 1, pp. 3-15, 2003.
[4] M. Kazhdan, M. Bolitho and H. Hoppe, “Poisson Surface Reconstruction,” in
Proceedings of the fourth Eurographics symposium on Geometry processing, vol.
7, 2006, pp. 61-70.
[5] J. Wang and M. M. Oliveira, “Filling holes on locally smooth surfaces reconstructed from point clouds,” Image and Vision Computing, vol. 25, no. 1, pp. 103-113, 2007.
[6] D. Doria and RJ. Radke, “Filling Large Holes in Lidar Data by Inpainting Depth Gradients," in Proceeding of IEEE Computer Vision and Pattern Recognition
Workshops (CVPRW), 2012, pp. 65-72.
[7] P. Lancaster and K. Salkauskas, “Surfaces Generated by Moving Least Squares Methods,” Mathematics of computation, vol. 37, no. 155, pp.141-158, 1981. [8] S. Fleishman, D. Cohen-Or and T. C. Silva, “Robust moving least-squares fitting
45
544-552, 2005.
[9] H. Obermaier, M. Hering-Bertram, J. Kuhnert, and H. Hagen, “Generation of adaptive streak surfaces using moving least squares,” Proceedings of Dagstuhl
Follow-Ups Scientific Visualization Seminar, vol. 2, 2011.
[10] L. J. Bentley, “Multidimensional binary search trees used for associative searching,” Communications of the ACM, vol. 18, no. 9, pp. 509-517, 1975. [11] I. Wald and V. Havran, “On Building Fast KD-Trees for Ray Tracing, and on
Doing that in O (N log N),” in Proceeding of IEEE Symposium on Interactive
Ray Tracing, 2006, pp. 61-69.
[12] S. W. Huang, “Integration of LIDAR and Vision Based Approaches for Textured 3D Scene Reconstruction.” M. S. Thesis, Department of Computer Science and Information Engineering, National University of Kaohsiung, 2012.
[13] B. Li, R. Schnabel, R. Klein, Z. Cheng, G. Dang, and S. Jin, “Robust normal estimation for point clouds with sharp features,” Computers & Graphics, vol. 34, no. 2, pp. 94-106, 2010.
[14] C. Tomasi and R. Manduchi , “Bilateral Filtering for Gray and Color Images,” in
Proceedings of the International Conference on Computer Vision (ICCV) , 1998,
pp. 839-846.
[15] C. Pal, A. Chakrabarti, and R. Ghosh, “A brief survey of recent edge-preserving
smoothing algorithms on digital images,” Procedia Computer Science, pp. 1-40, 2015.
[16] S. Fleishman, I. Drori, and D. Cohen-Or, ”Bilateral mesh denoising.” ACM
transactions on graphics (TOG), vol. 22, no. 3, pp. 950-953, 2003.
[17] H. G. Bendels, R. Schnabel and R. Klein, “Detecting Holes in Point Set Surfaces,”
46
[18] L. Truong‐Hong, F. D. Laefer, T. Hinks, and H. Carr, “Combining an Angle Criterion with Voxelization and the Flying Voxel Method in Reconstructing Building Models from LiDAR Data,” Computer‐Aided Civil and Infrastructure
Engineering, vol. 28, no. 2, pp. 112-129, 2013.
[19] D. Levin, “The approximation power of moving least-squares,” Mathematic of
Computation, vol. 67, no. 224, pp. 1517-1531, 1998.
[20] C. Y. Chen, “A Hierarchical Spatial Clustering Algorithm Based on Delaunay Triangulation.” M. S. Thesis, Department of Computer Science and Information Engineering, Feng Chia University, 2005.
[21] C. L. Lawson, “Generation of a triangular grid with applications to contour plotting,” California Institute of Technology Jet Propulsion Laboratory Technical Memorandum, no. 299, 1972.
[22] B. C. Chen, “Unorganized Point Cloud Reconstruction Using Mean-Shift,” Ph. D. Thesis, Department of Geomatics, National Cheng Kung University, 2007.
[23] M. S. Floater and M. Reimers, “Meshless parameterization and surface
reconstruction,” Computer Aided Geometric Design, vol. 18, no. 2, pp. 77–92, 2001.
[24] K. Hormann and M. Reimers, “Triangulating point clouds with spherical topology,” Curve and Surface Design, pp. 215–224, 2002.
