摘  要： 提出了一種基于離散曲率估計(jì)和kd-tree簡(jiǎn)化人臉點(diǎn)云的并行EM-ICP配準(zhǔn)算法。首先建立人臉點(diǎn)云的三維空間kd-tree，并結(jié)合離散高斯曲率對(duì)點(diǎn)云進(jìn)行了保留幾何特征的簡(jiǎn)化；然后基于CUDA對(duì)EM-ICP算法進(jìn)行并行加速，對(duì)簡(jiǎn)化的人臉點(diǎn)云進(jìn)行配準(zhǔn)。該算法能夠避免局部配準(zhǔn)等缺陷，同時(shí)EM-ICP算法并行保證了配準(zhǔn)工作的高效。實(shí)驗(yàn)證實(shí)了本文算法的健壯性和穩(wěn)定性。
關(guān)鍵詞： 點(diǎn)云配準(zhǔn)； EM-ICP； kd-tree； CUDA

    人臉識(shí)別是身份認(rèn)證等領(lǐng)域的重要技術(shù)，一直受到眾多研究者的關(guān)注。隨著三維掃描技術(shù)的發(fā)展及三維掃描的普及應(yīng)用，三維人臉的識(shí)別和匹配成為生物特征識(shí)別領(lǐng)域的研究熱點(diǎn)之一。在三維人臉識(shí)別與匹配過(guò)程中，有兩個(gè)步驟可應(yīng)用點(diǎn)云配準(zhǔn)過(guò)程：(1)使用三維掃描儀獲取人臉深度圖像，但要得到完整三維人臉點(diǎn)云數(shù)據(jù)往往需要多次掃描完成，因?yàn)槊看螔呙璧玫降狞c(diǎn)云數(shù)據(jù)往往只是部分人臉表面數(shù)據(jù)，所以需要對(duì)掃描得到的人臉深度圖像進(jìn)行配準(zhǔn)； (2)可將待匹配人臉模型與參考庫(kù)中人臉模型或兩者對(duì)應(yīng)特征進(jìn)行配準(zhǔn)，通過(guò)比較配準(zhǔn)誤差來(lái)判斷匹配程度。
    隨著點(diǎn)云數(shù)據(jù)規(guī)模的擴(kuò)大和配準(zhǔn)精度要求的提高，傳統(tǒng)串行點(diǎn)云配準(zhǔn)過(guò)程效率較低；而圖形處理單元GPU(Graphic Processing Unit)可應(yīng)用于并行計(jì)算，適合大規(guī)模數(shù)據(jù)處理，尤其是Nvidia公司推出的統(tǒng)一計(jì)算架構(gòu)CUDA(Compute Unified Device Architecture)能為算法并行提供更直觀(guān)的編程模型和優(yōu)化原則，提高點(diǎn)云配準(zhǔn)的運(yùn)算效率。
    本文基于GPU的高性能并行計(jì)算功能以及EM-ICP算法實(shí)現(xiàn)了對(duì)較大規(guī)模的人臉點(diǎn)云數(shù)據(jù)的并行配準(zhǔn)，實(shí)驗(yàn)證明了本文算法的配準(zhǔn)精度和速度上的顯著提高。
1 相關(guān)工作
    點(diǎn)云配準(zhǔn)通常采用迭代最近點(diǎn)ICP（Iterative Closest Point)算法[1]，是通過(guò)迭代計(jì)算使兩點(diǎn)云上對(duì)應(yīng)點(diǎn)對(duì)或點(diǎn)面的均方誤差最小，以實(shí)現(xiàn)點(diǎn)云的精確配準(zhǔn)。ICP算法的不足之處是對(duì)初始對(duì)齊敏感，故通常采用主元分析PCA(Principal Component Analysis)等方法進(jìn)行初始配準(zhǔn)，即粗對(duì)齊，再以初始配準(zhǔn)結(jié)果為條件進(jìn)行ICP配準(zhǔn)，即精確對(duì)齊。參考文獻(xiàn)[2]將最大期望算法EM(Expectation Maximization algorithm)算法[3]應(yīng)用到ICP算法中，提出了EM-ICP算法，從而避免了初始配準(zhǔn)的步驟。
    由于ICP和EM-ICP算法均含有大規(guī)模矩陣運(yùn)算，串行的配準(zhǔn)工作效率較低。參考文獻(xiàn)[4]基于CUDA對(duì)ICP算法進(jìn)行了并行加速，實(shí)現(xiàn)了深度圖像的實(shí)時(shí)配準(zhǔn)。參考文獻(xiàn)[5]基于GPU實(shí)現(xiàn)了EM-ICP的并行計(jì)算，配準(zhǔn)效率較高，但由于顯存空間的限制，該工作對(duì)規(guī)模較大的點(diǎn)云模型進(jìn)行簡(jiǎn)單的隨機(jī)采樣簡(jiǎn)化，其配準(zhǔn)精度存在明確的損失，本文實(shí)驗(yàn)也證明了該工作在應(yīng)用于人臉點(diǎn)云配準(zhǔn)時(shí)存在局部配準(zhǔn)的缺陷。
    解決因隨機(jī)采樣導(dǎo)致局部配準(zhǔn)缺陷的關(guān)鍵在于簡(jiǎn)化人臉點(diǎn)云的同時(shí)盡可能地保留點(diǎn)云的幾何特征。曲率是表示形狀的基本特征，能夠反映人臉三維表面的凹凸變化程度，且對(duì)旋轉(zhuǎn)、平移和縮放等變換具有幾何不變性，依據(jù)離散的高斯曲率估計(jì)[6]和kd-tree[7]進(jìn)行人臉點(diǎn)云簡(jiǎn)化可以保留足夠的人臉幾何特征。
    本文首先基于離散高斯曲率估計(jì)和kd-tree對(duì)三維人臉點(diǎn)云進(jìn)行簡(jiǎn)化，然后基于CUDA和EM-ICP算法對(duì)簡(jiǎn)化點(diǎn)云進(jìn)行并行加速配準(zhǔn)。實(shí)驗(yàn)證明，本文算法避免了局部配準(zhǔn)的問(wèn)題，提高了人臉點(diǎn)云配準(zhǔn)的效率和精確度。
2 利用kd-tree進(jìn)行點(diǎn)云簡(jiǎn)化
    本文在簡(jiǎn)化點(diǎn)云過(guò)程中，將點(diǎn)云分為關(guān)鍵點(diǎn)與非關(guān)鍵點(diǎn)。三維人臉原始點(diǎn)云數(shù)據(jù)為X，首先求解點(diǎn)云中每點(diǎn)的高斯曲率，并與指定閾值進(jìn)行比較，若高斯曲率大于指定閾值則判定此點(diǎn)為關(guān)鍵點(diǎn),直接復(fù)制到簡(jiǎn)化點(diǎn)云中；若小于指定閾值，判定此點(diǎn)是非關(guān)鍵點(diǎn)，再通過(guò)kd-tree進(jìn)行k鄰域搜索建立考查球，考查球內(nèi)的點(diǎn)密度，若點(diǎn)密度大于某閾值，則標(biāo)記考查球內(nèi)所有點(diǎn)的平均點(diǎn)為關(guān)鍵點(diǎn)，否則將球內(nèi)所有點(diǎn)標(biāo)記為關(guān)鍵點(diǎn)，據(jù)此實(shí)現(xiàn)對(duì)非關(guān)鍵點(diǎn)的簡(jiǎn)化。
    算法1 利用高斯曲率和kd-tree進(jìn)行點(diǎn)云簡(jiǎn)化[7]
　    輸入：原始點(diǎn)云X
　    輸出：簡(jiǎn)化點(diǎn)云XR
    具體步驟如下。
　　(1) 輸入原始點(diǎn)云X；
　　(2) 對(duì)X建立kd-tree，設(shè)定鄰域半徑，在鄰域內(nèi)計(jì)算每個(gè)點(diǎn)的高斯曲率；
　    (3) 設(shè)定閾值,高斯曲率大于閾值的點(diǎn)標(biāo)記為關(guān)鍵點(diǎn)，并復(fù)制到XR；
　    (4) 對(duì)于非關(guān)鍵點(diǎn)，通過(guò)隨機(jī)采樣選取n個(gè)隨機(jī)點(diǎn)；
　　(5) for i=0; i<n; i++ do
　　    以第i點(diǎn)為中心，利用kd-tree搜索離中心點(diǎn)最
        近的k個(gè)點(diǎn)，(a)如果找到的k個(gè)點(diǎn)均為非關(guān)鍵
        點(diǎn)，以中心點(diǎn)為球心,以第k個(gè)點(diǎn)到中心點(diǎn)的歐
        氏距離為半徑，建立考查球； (b)如果找到的離
        中心點(diǎn)最近的第l(l<k)個(gè)點(diǎn)是已標(biāo)記的關(guān)鍵點(diǎn)，
        則以中心點(diǎn)為球心， 以第l-1個(gè)點(diǎn)到中心點(diǎn)的
        歐氏距離為半徑， 建立考查球； 計(jì)算考查球內(nèi)
        點(diǎn)密度,若大于指定閾值，標(biāo)記考查球內(nèi)所有點(diǎn)
        的平均點(diǎn)為關(guān)鍵點(diǎn)，復(fù)制到XR；若小于指定閾
        值，將球內(nèi)所有點(diǎn)標(biāo)記為關(guān)鍵點(diǎn)，復(fù)制到XR；
　　    i←i+1    
　    end

5 實(shí)驗(yàn)結(jié)果及分析
    實(shí)驗(yàn)所用三維掃描儀的分辨率為640×480,幀頻為24 f/s。實(shí)驗(yàn)程序運(yùn)行硬件配置為：Intel Celeron 2.66 GHz處理器，1 GB內(nèi)存，GeForce GTS 250顯卡，128個(gè)CUDA處理器核心，1.1 GHz顯存頻率，1 GB顯存容量。系統(tǒng)環(huán)境：Gentoo Linux，CUDA 4.1，GCC4.5.3。
5.1 三維人臉點(diǎn)云簡(jiǎn)化
    實(shí)驗(yàn)結(jié)果如圖1所示?？梢钥闯?本文基于kd-tree的點(diǎn)云簡(jiǎn)化算法對(duì)人臉點(diǎn)云中高斯曲率值較大的區(qū)域有很好的保留，并且可以看到，兩個(gè)簡(jiǎn)化點(diǎn)云中曲率較大的區(qū)域基本一致，主要集中在眼、鼻、口等區(qū)域，而參考文獻(xiàn)[5]隨機(jī)采樣由于選點(diǎn)的隨機(jī)性難以得到一致性對(duì)應(yīng)。

5.2  EM-ICP算法應(yīng)用于人臉點(diǎn)云配準(zhǔn)效果
    由于參考文獻(xiàn)[5]隨機(jī)采樣的人臉點(diǎn)云中臉部點(diǎn)所占比重較大，應(yīng)用EM-ICP算法進(jìn)行點(diǎn)云配準(zhǔn)可能會(huì)出現(xiàn)局部配準(zhǔn)的缺陷，即：對(duì)臉部實(shí)現(xiàn)較好配準(zhǔn)，而對(duì)眼、鼻、口等部位未完全對(duì)齊。圖2(a)是隨機(jī)采樣點(diǎn)云對(duì)齊后的截面圖，可見(jiàn)在對(duì)兩采樣點(diǎn)云的鼻、口等部位配準(zhǔn)出現(xiàn)一定偏差；圖2(b)和圖2(c)是采用本文點(diǎn)云簡(jiǎn)化算法配準(zhǔn)結(jié)果的截面圖，效果均比圖2(a)的配準(zhǔn)好。

     本文點(diǎn)云簡(jiǎn)化算法能夠保留更多眼、鼻、口等部位的點(diǎn)， 增加這些部位的比重， 因此可以避免圖2(a)所示的局部配準(zhǔn)情形。整體配準(zhǔn)效果如圖3 所示，可以看出，最終配準(zhǔn)結(jié)果均較為理想，兩人臉點(diǎn)云實(shí)現(xiàn)很好對(duì)齊，下巴等點(diǎn)云連接處過(guò)渡平滑，圖2(b)和圖2(c)也證明了本文算法實(shí)現(xiàn)在口、鼻處的精確配準(zhǔn)?；诤?jiǎn)化點(diǎn)云的EM-ICP算法在人臉點(diǎn)云配準(zhǔn)中能夠避免局部配準(zhǔn)，提高配準(zhǔn)算法魯棒性。

5.3 EM-ICP算法并行加速
    利用參考文獻(xiàn)[5]的方法實(shí)現(xiàn)了基于CUDA的EM-ICP算法并行，加速效果明顯，針對(duì)不同點(diǎn)云規(guī)模的EM-ICP并行與串行效率對(duì)比如圖4所示，最大加速比可達(dá)近450倍。

    點(diǎn)云配準(zhǔn)精度直接影響著三維人臉識(shí)別和匹配的準(zhǔn)確度。本文提出的基于高斯曲率簡(jiǎn)化點(diǎn)云的EM-ICP并行配準(zhǔn)算法實(shí)現(xiàn)了三維人臉點(diǎn)云的有效配準(zhǔn)，改進(jìn)了局部配準(zhǔn)等不足，提高了算法的健壯性，具有較高的實(shí)際應(yīng)用價(jià)值。后續(xù)工作將考慮如何基于高斯曲率等幾何信息提取人臉點(diǎn)云的顯著性特征點(diǎn)云，以期進(jìn)一步提高配準(zhǔn)算法的計(jì)算效率和配準(zhǔn)精度。
參考文獻(xiàn)
[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] 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.
[5] 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.
[6] WOLFGANG K. Differential geometry: curves-surfaces-manifolds[M]. 2nd Edition, Kuhnel, Wolfgang: American  Mathematical Society, 2006:158-165.
[7] De Berg M, CHEONG O. Computational geometry: algo-rithms and applications[M]. 3rd Edition, New York: Springer, 2008:99-105.
[8] HORN B P. Closed-form solution of absolute orientation using unit quaternions[J]. Journal of the Optical Society of America, 1987:629-642.
[9] Nvidia. CUDA CUBLAS Library[Z]. http://cudazone.nvidia.cn/cublas/.