基于持久同調的三維模型檢索方法

況立群，李 麗，幸嘉誠，諶鐘毓，韓 燮

(1.中北大學 大數據學院，山西 太原 030051；2.南昌大學 軟件學院，江西 南昌 330047；3.華東交通大學 軟件學院，江西 南昌 330013)

0 引 言

面對數量眾多、種類繁雜的三維模型，如何及時有效地獲取到用戶所需模型，成為三維建模領域的研究熱點[1]。根據特征描述子的不同提取方式，三維模型檢索方法可分成統計特征獲取、投影視圖計算和拓撲結構描述3類。基于統計特征的檢索方法研究每個個體特征組成的總體分布情況，在描述模型內部結構方面有所欠缺，使得檢索效率不高[2]。基于投影視圖的檢索方法將三維模型投影為多個二維視圖[3]，在降維時會丟失三維模型的部分信息[4]。以上兩類檢索方法主要是對模型的幾何結構進行分類，忽略了模型的本質拓撲特征，而基于拓撲結構的檢索方法則通過描述三維模型的骨架結構來獲取模型的拓撲特征[5]，但難以應用于無明顯骨架結構的三維模型檢索。

近年來，國際上興起了基于持久同調理論的拓撲數據分析方法，適用于研究任意形狀物體的拓撲特征，已廣泛應用于大數據[6]、腦影像特征分析[7]等領域，然而在三維模型檢索中的應用不多，且國內鮮有相關研究報道。因此，本文提出一種基于持久同調的三維模型檢索方法。該方法關注三維模型在多尺度下的拓撲結構，從本質上獲取三維模型內部的拓撲特征，可刻畫一些無明顯骨架結構的三維模型的拓撲特征，具有旋轉平移和尺度不變性。

1 相關工作

1.1 持久同調

持久同調(persistent homology)是拓撲數據分析中的一個重要研究方向，關注的是三維模型結構中點與點之間的拓撲不變量，用于研究多個尺度下的定性特征[8]。本文運用持久同調理論提取三維模型的特征描述子，其中的重要一步是為三維模型構建復形過濾流(filtered complex)，其主要過程為：在三維模型中，每個點被一個直徑為d的球體所包圍，該球的直徑d從0開始逐步增大，當三維模型中任意兩點的間距小于該時刻的直徑d時，便將這兩點連在一起。在該過程中使用單純復形(simplicial complex)去填充整個變化的結構，最后形成一系列嵌套的單純復形，如圖1所示。

圖1 復形過濾流構建

在構建復形過濾流過程中獲取到三維模型在多尺度范圍上的拓撲結構，一些拓撲結構短暫出現之后被填充，被認為是噪音和擾動；另外一些拓撲結構長時間持續存在，被認為是穩定的拓撲特征，這些穩定的拓撲特征描述了三維模型內部本質的特征，將其作為特征描述子來表征三維模型，這便是持久同調的全過程[9]。

常構建的復形過濾流有Cech復形、Alpha復形、Witness復形、Vietoris-Rips復形。本文中構建的是Vietoris-Rips復形，即在高維空間中有一點集P，在過濾尺度λ下生成變化的P的子集σ，且子集中任意兩點之間距離都小于λ[10]，Vietoris-Rips復形的形式化如式(1)所示，其構建過程如圖2所示[11]

VλP=du,v≤λ,?u≠v∈σ

(1)

圖2 Vietoris-Rips構建過程

1.2 特征描述子

為了通過比較三維模型間的差異來實現檢索，需要為三維模型生成一個特征描述子。在持久同調過程中，產生了一些持續時間長的穩定拓撲特征，通常用貝蒂數(Betti number)的變化來量化這些特征[12]。

令K為一個原始的單純復形，可在某一時刻添加其子單純復形，構建變化的復形過濾流，建立單純復形間的嵌套關系

K0?K1?K2…?Kn=K

(2)

因為Ki-1?Ki，在單純復形之間引入同態f:Hp(Ki-1)→Hp(Ki)，使嵌套的復形序列與同態連接的同調群序列一一對應，其中每個維度對應一個不同的p，則同調群秩的序列如下[8]

Hp(K0)→Hp(K1)→Hp(K2)…→Hp(Kn)=Hp(K)

(3)

同調群Hp(K)的秩叫作單純復形K的p維貝蒂數，它是代數拓撲中定義在同調群上的拓撲不變量，0維貝蒂數表示連通分支的個數，1維貝蒂數表示一維下孔洞的個數，可用條形碼[11]和持久性[13]來表示貝蒂數的變化，如圖3和圖4所示。在條形碼圖中直線的起始點表示孔洞的產生時間，直線的終止點表示孔洞的消失時間；持久性圖中橫坐標表示孔洞的產生時間，縱坐標表示孔洞的消失時間。條形碼圖與持久性圖可以相互轉換，本文采用持久性圖來表征三維模型特征，作為三維模型的特征描述子。

圖3 條形碼

圖4 持久性

1.3 Wasserstein距離

由于持久性圖是二維散點的集合，難以在一維空間直接度量持久性圖之間的相似性。最傳統的度量方法當屬瓶頸距離法，近年來有研究表明Wasserstein距離[14]具有更好的度量性能。Wasserstein距離最早被用來比較兩個直方圖的差異，之后被用來刻畫兩個分布之間的距離，可以通俗理解為將一堆物體移動到另一個位置的最小成本，它的定義為

(4)

其中，X和Y分別為滿足μ和ν概率分布的d維空間集合，p表示Lp范數，inf表示下界，本文方法中p值取1。

2 基于持久同調的三維模型檢索方法

本文提出一種基于持久同調的三維模型檢索方法，整個檢索框架如圖5所示。首先，采用Ripser軟件包來構建Vietoris-Rips復形，計算得到三維模型的拓撲特征(1維貝蒂數)，并將其表示在持久性圖中，形成特征描述子；然后利用分段Wasserstein距離計算兩個持久性圖之間的相似性；最后，檢索相應模型，通過獲取到的評價參數來評價檢索效率。

圖5 總體檢索框架

2.1 數據預處理

本文采用Vietoris-Rips復形構建三維點云模型的拓撲結構，該結構比傳統的三角網格化模型的結構復雜得多，嵌套復形的信息量異常龐大，難以對其直接拓撲化，需預先對三維模型數據集進行簡化。本文利用meshlab中的插件meshlabserver進行批處理降采樣，通過減少面的個數來減少點的個數，同時在降采樣前設置“保留拓撲特征”等參數。這樣不僅保留了三維模型內部的拓撲特征，又減少了模型的稠密度，縮短了計算時間。實驗中考慮到時間和設備的限制，將模型降為800個面，400多個點。

2.2 提取特征描述子

構建Vietoris-Rips復形的軟件包有Javaplex、PHAT、GUDHI、DIPHA、Ripser等，其中Ripser是用于計算Vietoris-Rips持久性條形碼的精簡C++代碼，與其它軟件包相比，Ripser的運行速度超出40倍，內存效率超過15倍，且沒有外部依賴，可以大大提高計算的時間和效率[11]。因此，本文選用Ripser軟件包來構建復形結構，編程環境為Python語言的Pycharm。

在讀取簡化后三維模型點的信息后，調用Ripser庫中的ripser函數為每個三維模型構建Vietoris-Rips復形，獲取到復形過濾流下的一維貝蒂數，它比零維貝蒂數能更好地表征三維模型內部本質的拓撲特征。圖6為對三維模型數據集中的一個杯子模型構建Vietoris-Rips復形后，得到零維和一維下的貝蒂數，將其表示在條形碼圖(左圖)和持久性圖(右圖)中，形成特征描述子。

圖6 杯子的特征描述子圖

2.3 相似性度量

Wasserstein距離在數值上只能計算一維高斯分布的情況，而本文計算的是二維持久性圖間的相似性，故此提出一種分段的Wasserstein距離來度量持久性圖之間的相似性。分段的Wasserstein距離是對Wasserstein距離的一個改進，其基本原理為：讓通過原點的多條線去劃分持久性圖的平面，將二維的點投影到每條線的方向上，從而將二維點集投影為一維點集，然后再計算兩個一維點集間的Wasserstein距離，最后將各個方向上Wasserstein距離累加求和再平均，即為最終距離。本方法將提取到的拓撲特征以點的形式記錄在持久性圖中，在第一象限內的多個方向上進行多次迭代，得到分段Wasserstein距離。計算過程分為對角線投影、Wasserstein距離計算、迭代求平均這3個步驟。

(1)投影到對角線

假設有兩個持久性圖dg1和dg2，其中dg1為m個點構成的特征集合，dg2為n個點構成的特征集合，將dg1的點投影到y=x這條對角線上，并將投影后的點加入到dg2中。同理，將dg2的點投影到y=x這條對角線上，將投影后的點加入到dg1中。投影到對角線保證了兩個持久性圖在進行相似性度量時特征向量個數一致，并且仍保留各自的特征分布，原始的持久性圖表示如下

(5)

投影到對角線后得到

(6)

(2)計算Wasserstein距離

計算Wasserstein距離的示意圖如圖7所示。從原點射出的某條直線與X正半軸成θ角，且滿足sinθ2+cosθ2=1，由式(7)依次計算dg1中每個點和這個方向的點積

(cosθ,sinθ)·(xi,yi)=xicosθ+yisinθ=Xi,i=1…(m+n)

(7)

圖7 計算Wasserstein距離

將結果計入V1數組并排序，同理計算dg2中每個點與該方向的點積，將結果計入V2數組并排序

(8)

(3)迭代求平均

將第一象限平均分成M個方向，即迭代次數為M。為了加快計算過程，可以設置較小的迭代次數，實驗結果表明迭代次數較少時也能得到較好的檢索效果。從0開始以s=π/2/M的角度遞增，在每個方向上計算出一維數據的Wasserstein距離，之后累加求和再平均得到分段Wasserstein距離(式(9))。以此作為兩個持久性圖之間的差異值，該值越小表示兩個持久性圖的分布越近似，說明兩個三維模型越相似。迭代的目的是使兩個持久性圖在多個方向上計算差異，通過迭代求平均使實驗結果更準確，更具有穩定性

(9)

3 實驗結果及分析

3.1 實驗環境及數據

本研究實驗硬件環境為Intel(R)Xeon(R)CPU E5-1603 v3 @ 2.80 GHz 2.80 GHz，內存為16 GB，軟件環境為64位的windows 10、matlab2013b、PyCharm Community Edition 2018.3.2。實驗數據集來源于2011年的三維模型檢索競賽(SHREC TRACK 2011)的標準數據庫[15]，本數據集包括600個三維模型，共分為30個類，每類20個三維模型。

3.2 評價標準

為了衡量檢索效果的好壞，本文使用檢索領域常用的評價參數NN(nearest neighbor)、F-T(first tire)、S-T(second tire)、E(e-measure)、DCG(discounted cumulative gain)來進行評估。假設目標類模型個數為C，K為需檢索出的模型個數，則NN表示當K=1時正確檢索到的目標模型的比例；F-T表示當K=C-1時正確檢索到的目標模型的比例；S-T表示當K=2*C-1時正確檢索到的目標模型的比例；DCG表示根據檢索結果的排列順序來評價檢索優劣；E為查全率和查準率的綜合評價參數。為了更好地表征三維模型檢索的效果，通常還采用P-R圖表示數據集整體的查全率和召回率。

3.3 實驗結果

本文實驗結果分為兩部分，第一部分是本文方法的檢索結果，第二部分是本方法與其它方法的實驗數據對比。

(1)采用本文方法對“螞蟻”、“鉗子”三維模型進行檢索，這兩類三維模型在本數據集中各有20個，所以本實驗檢索出與目標模型最相近的前20個三維模型，以評價檢索的準確性，檢索結果如圖8所示。由圖8可以看出，除個別三維模型外，其余模型與目標模型同屬一類，由此得出，本文方法可以有效地檢索出目標模型，準確率較高。

圖8 檢索結果

(2)本文實驗與Osada等的形狀概率分布方法(shape distribution，D2)、Sun等的熱核特征方法(heat kernel signature，HKS)[16]、Sipiran等的角點檢測方法(Harris 3D geodesic map，Harris3DGeoMap)[17]和焦世超等的多特征融合流行排序方法(multi-features fused manifold ran-king，M-F FMR)[18]進行了對比。D2方法首先計算所有點之間的距離，然后將所有距離值作為特殊分布表示在直方圖中。HKS方法對模型上每個點的熱量變化曲線離散化，之后獲取全部點的特征，進而得到整個模型的熱核特征。Harris3DGeoMap方法采用3D Harris算法選取原始模型小部分的顯著興趣點，用特征分布直方圖表示興趣點間測地線距離的差異。Harris3DGeoMap16和Harris3DGeoMap32分別表示直方圖區間數目為16和32條件下的Harris3DGeoMap實驗。M-F FMR方法采用融合特征對二維視圖提取特征，然后利用深度學習對草圖進行語義標注，再使用流行排序算法實現二維草圖到三維模型的檢索。

這幾種方法都是在SHREC TRACK 2011的三維數據集上進行實驗，實驗結果具有可比性。對NN、F-T、S-T、E、DCG參數進行實驗對比，結果見表1，這些值越接近于1表明檢索效果越好，從表中看出本文的全部參數優于所有對比方法。同時繪制P-R曲線圖進行對比，結果如圖9所示，該圖走向越趨向于右上方說明檢索效果越好，從圖中可以看出本文方法的曲線均高于其它方法。

表1 6種檢索方法的評價參數

圖9 6種檢索方法的P-R曲線

綜上，實驗結果表明，本文提出的基于持久同調的三維模型檢索方法具有較高的準確率，在檢索結果(圖8)中可以看到準確檢索到了待查詢的模型；此外，本文方法的評價參數和P-R曲線圖均優于其它方法。

4 結束語

本文將持久同調理論應用于三維模型檢索，提出了一種基于持久同調的三維模型檢索方法，能夠在不降維的情況下提取出三維模型的全局特征描述子，該描述子具有旋轉、平移和尺度不變性，在數據擾動的情況下具有很好的穩定性，可全面有效地表征任意拓撲結構的三維模型的特征。另外提出了區分性更好的分段Wasserstein距離，對持久性圖之間進行了更為準確的相似性度量，進一步提高了三維模型檢索的查準率和查全率。本文研究結果對今后三維模型檢索提供了一種新的研究思路和設計方法。鑒于時間與空間上的計算資源限制，本文未能將三維模型的全部點計算在內，而是進行了降采樣處理，客觀上降低了檢索效率。在今后的研究中將利用并行計算資源來完善實驗環節，提高檢索準確率。