摘  要： 提出了一種新型兩步式迭代最近點算法對三維人耳點云模型進行配準，該過程主要分為兩步完成：（1）采用基于CUDA并行加速的ICP" title="EM-ICP" target="_blank">EM-ICP算法進行初始配準，從而使人耳點云數據大致調整為同一姿態，并且為下一步提供良好的初始變化；（2）基于ICP算法對三維人耳點云數據進行精確配準。該方式能夠有效避免ICP算法配準過程中局部對齊等缺陷。實驗結果證明，采用兩步式迭代最近點算法配準后的三維人耳數據具有良好的配準效果與配準速度。

　　關鍵詞： EM-ICP；ICP；人耳；點云配準；CUDA

0 引言

　　在當今信息化時代，隨著科學技術的不斷發展，傳統的基于身份證、學生證、磁卡等的身份鑒別技術存在容易被偽造、被盜取以及容易遺失等問題，暴露出越來越多的缺陷。它們已經不能滿足人們對快速、便捷、有效的身份識別技術的需求。在此情況下，生物特征識別技術應時而生。人耳識別是以人耳作為識別媒介來進行身份的鑒別，是一種很有發展潛力的生物特征識別技術，受到了國內外眾多研究機構的廣泛關注。研究指出，沒有任何兩個人（即使是雙胞胎）的耳朵是完全一樣的，并且在8~70歲之間都不會有顯著的變化，可以作為個體生物識別的依據。人耳形狀特征很豐富，其表面具有大量的溝和脊，不受胡須、化妝品、年齡、表情等影響，具有更高的穩定性、唯一性和健壯性，為人耳識別技術提供了理論研究價值和實際應用前景。

　　隨著三維掃描技術的迅速發展，三維數據的獲取變得更加方便，三維模型已經成為繼數字音頻、數字圖像、數字視頻之后一種新的數字媒體形式。三維人耳模型包含的特征信息比二維圖像更為豐富，因此基于三維人耳的識別技術便逐漸發展起來。三維人耳模型不但能夠很好的反應人耳的輪廓信息，而且能夠很好地描述人耳的結構和姿態信息。采用三維人耳數據進行識別能有效解決姿態變換、陰影和光照條件改變等問題對識別率的影響，因此更適合采用三維的方式來進行采集和識別。

　　三維人耳模型識別的步驟大致如下：首先使用三維掃描儀獲取人側臉的深度圖像；其次將人耳數據從人側臉數據中準確地提取出來；最后將不同人耳數據或其特征進行配準，通過比較人耳數據之間的配準誤差，從而實現三維人耳識別。

1 相關工作

　　迭代最近點（Iterative Closest Points，ICP）算法[1]是目前最常用的三維數據配準算法，通過迭代最小化兩待配準點云上對應點間的距離誤差，獲得最佳的旋轉矩陣和平移矩陣，實現精確配準。迭代最近點算法能夠滿足大多數三維數據的配準要求，但這些算法在不知道待配準模型之間對應點的情況下都需要有一個良好的初始變換，不好的初始變換會導致三維模型局部收斂，直接影響著三維數據的配準效果。因此，為避免該缺陷，許多研究者采用了很多解決方式。2002年Granger等[2]提出EM-ICP算法，將最大期望（EM）算法[3]應用到ICP算法中，從而避免了初始配準的步驟。2005年，Hui等[4]利用兩步迭代最近點算法對人耳進行配準，首先利用ICP算法對耳輪數據進行粗配準，然后將該變換矩陣作為初始變換再次應用在ICP算法中，對第一步的匹配進行優化，提高識別效率。2007年，Yan等[5]通過主成分分析（PCA）[6]算法對待配準點云先進行初始配準，調整人耳的姿態，再對初始配準后的結果使用ICP算法進行精確配準。同年，Chen等[7]利用四元組計算初始變換進行粗對齊，再利用ICP算法進行精確匹配。

　　隨著三維掃描技術的不斷發展，三維掃描儀的掃描精度不斷提高，數據規模也隨之增大。由于ICP算法、EM-ICP算法均需要對大規模的矩陣進行運算，數據規模的增大必然導致工作效率的降低，傳統的串行配準合并算法的效率已無法滿足實時性的需求。圖形處理單元GPU進行并行計算，由多個流處理器分別進行數值運算，實現任務級和數據級的并行，能夠很好地解決上述問題。NVIDIA公司推出的統一計算架構CUDA提供了高性能的GPU并行計算環境，可用于大規模三維數據的處理。由CPU作為主機負責邏輯性強的事務處理和串行計算，GPU作為協處理器完成可并行計算的部分，高度線程化的并行處理任務則由CPU和GPU共同完成，大大提高了程序的運行效率和數據的處理速度，使由于數據規模較大、精度要求較高造成的配準及合并效率降低等問題得以解決。2008年Choi等[8]基于CUDA對ICP算法進行了并行加速，實現了對深度圖像進行實時配準。2010年Tamaki等[9]基于CUDA對EM-ICP算法進行了并行加速，配準速度明顯提高。

　　根據上述學者們的研究，本文提出了一種兩步式迭代最近點算法對三維人耳點云模型進行配準。首先采用基于CUDA加速的EM-ICP算法作為ICP算法的初始配準，使人耳點云數據大致調整為同一姿態，然后基于該算法提供的初始變化采用ICP算法對三維人耳點云數據進行精確配準，相當于進行了兩次配準，最終達到配準效果。

2 基于EM-ICP的初始配準

　　初始配準能夠有效調整三維模型的位置與姿態，為精確配準提供一個理想的初始變換。本文采用基于CUDA加速的EM-ICP算法對三維人耳模型進行初始變換，EM-ICP算法不需要建立初始對應關系，以權重表示兩點間的配準概率，迭代運算優化配準概率，最終實現點云配準。

　　已知三維人耳點云模型X={xi，i=1，…，n}與三維人耳點云模型Y={yi，i=1，…，m}，n與m分別表示人耳點云模型X與Y中點云個數，模型X上的任意一點xi與模型Y上所有點都存在一個對應關系，且用權重的大小來表示配準概率。通過求解模型的變換矩陣R與t，更改人耳點云模型Y的位置，直到點云模型間誤差函數E最小。

　　其中，ij表示xi與對應點yi的配準概率。

　　因此，點云模型間誤差函數E可重寫為：

　　EM-ICP算法[2]具體步驟如下：

　　（6）若E大于閾值τ且?滓p小于0.3，則返回到（2），否則迭代結束。

　　EM-ICP算法中，點云模型X上的每一個點都與點云模型Y上所有點存在一個對應關系，即匹配概率ij，因此計算全映射的關系矩陣A=（ij）與兩個點云模型的規模密切相關，當點云模型的規模較大時，對矩陣A運算的時間很長，對矩陣A的計算進行GPU并行加速，加快算法效率。

　　對原算法進行并行加速的關鍵問題是將運算過程分為向量與矩陣的運算和矩陣內元素間的運算兩種，利用CUBLAS（CUDA Basic Linear Algebra Subprograms）[10]對向量與矩陣間的運算進行加速，編寫CUDA kernel函數對矩陣元素間的運算進行加速[9]。具體步驟如下：

　　（1）將模型X、Y拷貝到顯存中，并且對CUDA環境及CUBLAS庫函數初始化；

　　（2）計算模型X、Y對應點之間的距離dij；

　　（6）求解旋轉矩陣R、平移矩陣t；

　　（7）更新模型Y的位置Y=RY+t；

　　（8）更新控制參數；

　　（9）若E大于閾值τ且p小于0.3，則返回到（2），否則迭代結束。

3 基于ICP的精確配準

　　本文采用ICP算法對粗配準后的三維人耳模型進行精確配準。ICP算法能夠對深度圖像進行有效的配準，是當前眾多配準算法的基礎。ICP算法不斷地更新一個點云模型的位置，直到該模型與另一個點云模型對應點之間的距離達到某閾值為止。在ICP算法中，點云模型X上的任意一點xi在點云模型Y上有且僅有一個對應點。

　　已知點云模型X={xi，i=1，…，n}與點云模型Y={yi，i=1，…，m}，n與m分別表示人耳點云模型X與Y中點云個數，尋找點云模型X上每一個點xi到點云模型Y上的最近點yi，通過求解模型的變換矩陣R與t，更改模型點云Y的位置，直到點云模型間誤差函數E最小。

　　ICP算法的主要目的是求解兩個點云模型之間的空間變換，通過這個空間變換使得兩點云模型之間的距離最小，其具體步驟如下：

　　（1）點云模型X與模型Y初始對齊；

　　（2）找到點云Y中距離點云X中xi最近的點yi；

　　（3）采用四元數方法解旋轉矩陣R，平移矩陣t，并求解；

　　（4）更新模型Y的位置X=RX+t；

　　（5）若E大于閾值τ，則返回到（2），否則迭代結束。

4 實驗結果及分析

　　實驗所用三維掃描儀的分辨率為640×480，幀頻為24 f/s。實驗程序運行硬件配置為：Intel Xeon E5-2609@2.40 GHz處理器，16 GB內存，NVIDIA Quadro 2000顯卡，192個CUDA核心，1 GB GDDR5顯存容量，計算能力2.1。系統環境：Fedora 16 Linux，CUDA6.5，GCC4.6.3。

　　4.1 數據采集

　　通過三維激光掃描儀可以得到人耳側面的掃描數據，但是得到的數據不僅包括人耳數據，還包括人耳附近的皮膚數據，需要將這些無用的數據除去，將人耳數據提取出來。

　　提取到人耳數據后，去掉其顏色信息，得到需要的三維人耳數據模型。在下面的實驗中將使用提取得到的三個人耳數據，如圖1所示，提取得到的人耳數據，方向各不相同，分別為其編號為ear_a，ear_b，ear_c。

　　4.2 配準效果

　　本文選用CUDA加速的EM-ICP算法作為ICP算法的初始配準，再使用ICP算法進行精確配準。進行兩次配準，既保證了配準速度，又保證了配準精度，最終得到了理想的配準效果。如表1所示，將ear_a作為待配準模型，對ear_a與ear_b，ear_a與ear_c分別進行EM-ICP粗配準與ICP精確配準。

　　由表1能夠清楚看出，人耳模型ear_b、ear_c經過EM-ICP粗配準后，能夠初步調整人耳模型的位置，ICP精確配準后均達到了理想的配準效果。

　　如圖2所示，待配準模型ear_a與配準后的ear_b、ear_c模型位置。可以看出，配準后的ear_b、ear_c均調整到與模型ear_a姿態一致的位置。

　　4.3 配準精度

　　將ear_a作為配準模型，對ear_a與ear_b、ear_a與ear_c進行不同方式的配準，如圖3所示為分別基于ICP[1]、EM-ICP[9]、兩步式ICP[4]以及本文提出的EM-ICP和ICP算法相結合的方式進行配準得到的配準精度。由圖可見，本文算法與其他算法相比，具有較高的配準精度，配準效果優于其他方式。

　　4.4 配準時間

　　如表2所示，將分別基于ICP[1]、EM-ICP[9]、兩步式ICP[4]以及本文提出的EM-ICP和ICP算法相結合的方式進行配準所用時間進行對比。顯然，ICP算法效率略低，EM-ICP算法具有很高的效率，采用兩種迭代方式的時間消耗比采用一種迭代方式的時間消耗高，然而在均采用兩種迭代方式前提下，本文算法的時間消耗要優于兩步式ICP算法，并且本文算法與只基于ICP算法相比，其時間消耗差距不大。

參考文獻

　　[1] BESL P J， MCKAY N D. A method for registration of 3-d shapes[J]. IEEE Transactions on Pattern Analysis and Machine Intelligence， 1992，14（2）：239-256.

　　[2] GRANGER S， PENNEC X. Multi-scale EM-ICP： a fast and robust approach for surface registration[C]. Proceedings of the 7th European Conference on Computer Vision. Copen-hagen， Denmark： Springer-Verlag， 2002：418-432.

　　[3] DEMPSTER A， LAIRD N， RUBIN D. Maximum likelihood estimation from incomplete data via EM Algorithm[J]. Journal of the Royal Statistical Society， 1977，39（1）：1-38.

　　[4] CHEN H， BHANU B. Contour matching for 3D ear recognition[C]. In：Proceedings of IEEE Workshop on Application of Computer Vision， 2005：123-128.

　　[5] YAN P， BOWYER K W， Biometric recognition using three-dimensional ear shape[J]. IEEE Trans PAMI， 2007，29（8）：1297-1308.

　　[6] ELAD M， TAL A， AR S. Content based retrieval of VRML objects-an iterative and interactive approach[C]. Proceedings of the Eurographics Workshop in Manchester， United Kingdom， 2001：107-118.

　　[7] CHEN H， BHANU B. Human ear recognition in 3D[J]. IEEE Transaction PAMI， 2007，29（4）：718-737.

　　[8] CHOI S I， PARK S Y， KIM J， et al. Multi-view range image registration using CUDA[C]. Proceedings of the 23rd International Technical Conference on Circuits/Systems， Computers and Communications， 2008：733-736.

　　[9] TAMAKI T， ABE M， RAYTCHEV B， et al. Softassign and EM-ICP on GPU[C]. Proceedings of the 2010 1st International Conference on Networking and Computing， Washington DC， USA： IEEE， 2010：179-183.

　　[10] NVIDIA. CUDA CUBLAS（CUDA Basic Linear Algebra Subprograms）Library[EB/OL].（2015-04-15）. http：//cudazone.nvidia.cn/cublas/.