張伯雄， 謝伙生

(福州大學(xué)數(shù)學(xué)與計(jì)算機(jī)科學(xué)學(xué)院， 福建 福州 350116)

基于特征分布的三維流線相似性研究

張伯雄， 謝伙生

(福州大學(xué)數(shù)學(xué)與計(jì)算機(jī)科學(xué)學(xué)院， 福建 福州 350116)

針對(duì)利用流線的形狀特征對(duì)流線進(jìn)行分類和選取， 可以方便用戶洞察三維流場(chǎng)的特征, 提出一種高效的基于特征分布的三維流線相似性比較方法. 該方法在Lu Kewei等流線相似性比較方法的基礎(chǔ)上， 引入曲折度和速度方向熵兩個(gè)全局屬性. 首先將流線均勻分段， 然后計(jì)算每段的曲率直方圖、 扭率直方圖、 曲折度直方圖、 速度方向熵直方圖， 構(gòu)建相應(yīng)的二維直方圖， 最后利用堆土機(jī)距離(EMD)及k-means聚類方法進(jìn)行流線相似度計(jì)算和分類. 實(shí)驗(yàn)結(jié)果表明， 該方法在引入全局幾何屬性后能夠產(chǎn)生更魯棒性的查詢和分類結(jié)果.

流線相似性； 特征分布； 二維直方圖； 流線聚類

0 引言

流場(chǎng)可視化是科學(xué)計(jì)算可視化的一個(gè)重要課題, 在許多科學(xué)和工程領(lǐng)域， 可視化矢量場(chǎng)在理解潛在的流場(chǎng)特征和模式方面起著重要的作用. 流場(chǎng)可視化研究方法包括基于圖標(biāo)的方法[1]、 基于紋理的方法[2]、 基于幾何圖形的方法[3]和基于特征結(jié)構(gòu)的可視化方法[4]. 其中因?yàn)榱骶€便于計(jì)算， 可在不同分辨率下以交互的速率進(jìn)行渲染， 所以流線可視化方法是最普遍使用的方法之一. 但當(dāng)成千上百條流線被渲染描繪一個(gè)流場(chǎng)時(shí)， 由于雜亂和遮擋問題， 用戶很難發(fā)覺潛在的流場(chǎng)特征， 而且用戶往往只想查看具有特定形狀特征的流線. 因此， 將流線基于形狀特征進(jìn)行分類， 展示和用戶目標(biāo)相關(guān)的流線是非常重要的.

流線相似性度量和分類工作中， 有許多是基于點(diǎn)和點(diǎn)之間的距離. Chen等[5]把流線間局部相似性度量的方法引入到流線生成中， 這種度量方法涉及流線形狀和方向的統(tǒng)計(jì)信息， 以及每對(duì)流線的歐式距離； Yu等[6]計(jì)算流線間點(diǎn)對(duì)的歐式距離， 使用平均最短點(diǎn)距離度量流線相似性， 最后結(jié)合單鏈接方法進(jìn)行層次聚類. 但這些相似性度量方法不具有旋轉(zhuǎn)、 平移、 縮放不變性， 兩條相似的流線可能因?yàn)榉较虿煌?距離過遠(yuǎn)、 長(zhǎng)短不一而導(dǎo)致度量距離過大不能劃分為同一類. 魯大營等[7]提出了一種基于迭代最鄰近點(diǎn)(ICP)算法與k-means聚類方法的流線生成算法， 首先利用ICP算法實(shí)現(xiàn)流線間輪廓特征上的配準(zhǔn)并根據(jù)幾何相似性進(jìn)行排序， 然后利用k-means聚類算法對(duì)流線分組， 這種方法雖然具有旋轉(zhuǎn)、 平移不變性， 但計(jì)算量過大.

近年來， 基于流線特征分布的相似性度量方法越來越受到重視. Mc Loughlin等[8]計(jì)算流線每點(diǎn)的各種特征值， 構(gòu)建1D直方圖， 并用χ2測(cè)試比較相似性， 但由于形狀不相近的流線也可能具有相似的1D直方圖， 文獻(xiàn)[8]的方法將不相似的流線歸為一類. Lu等[9]考慮特征的空間分布信息， 將流線分段， 構(gòu)建2D直方圖， 解決了文獻(xiàn)[8]方法中由于空間信息丟失導(dǎo)致的流線相似性度量錯(cuò)誤， 但他的相似度量方法只考慮流線的曲率和扭率兩種局部屬性， 當(dāng)流線選取的比例增大時(shí)， 新增流線的整體形狀可能與目標(biāo)流線有較大的偏差. Li等[10]受特征包想法的啟發(fā)， 選擇k個(gè)具有代表性的特征矢量作為詞匯， 在詞匯組合中考慮兩點(diǎn)的弧長(zhǎng)信息從而構(gòu)建空間敏感的特征包， 通過加權(quán)的曼哈頓距離進(jìn)行流線間的相似性計(jì)算， 并考慮流線的曲折度和速度方向熵兩個(gè)全局屬性， 使得查詢的流線結(jié)果與目標(biāo)流線在整體上更加貼近. 但由于文獻(xiàn)[10]中加權(quán)曼哈頓距離不一定能正確表達(dá)人對(duì)流線相似性的理解， 在流線聚類實(shí)驗(yàn)中， 可能使一條螺旋上升的流線被錯(cuò)誤地劃到與其形狀相差甚遠(yuǎn)的流線組中. 相比之下， 文獻(xiàn)[9]中基于2D直方圖的方法仍能正確地聚類， 具有更好的魯棒性.

因此， 采用2D直方圖比較流線相似性， 保證流線相似性比較和分類的魯棒性， 并引入文獻(xiàn)[10]中所用的曲折度和速度方向熵兩個(gè)全局屬性， 使得流線相似性比較在全局屬性上也能更加接近目標(biāo)流線， 產(chǎn)生更準(zhǔn)確的相似性度量結(jié)果.

1 三維流線特征的選擇

選取合適的特征描述子是構(gòu)建直方圖重要的一步. 為了選取形狀相類似的流線, 本研究選取能夠描述流線形狀的相關(guān)特征描述子. 首先選擇曲率和扭率， 因?yàn)樗鼈兎謩e表示了曲線的彎曲程度和脫離密切面的變化程度. 曲率和扭率屬于流線的局部幾何屬性， 為了能結(jié)合流線的全局幾何屬性進(jìn)行比較， 接著引入曲折度和速度方向熵. 曲折度反映了曲線整體的彎曲程度； 速度方向熵反映了速度方向的均勻程度， 亦即流線的平坦程度.

1.1 曲率

曲率是正切矢量關(guān)于弧長(zhǎng)的變化率. 曲線上某點(diǎn)p的曲率， 可以通過下式進(jìn)行計(jì)算：

(1)

1.2 扭率

扭率的計(jì)算公式如下：

(2)

其中:det表示矩陣的行列式.

1.3 曲折度

曲折度是曲線長(zhǎng)度和始末兩點(diǎn)最短距離之比. 曲線的長(zhǎng)度通過組成流線的直線段的累計(jì)長(zhǎng)度來近似， 始末兩點(diǎn)的最短距離采用歐幾里得距離進(jìn)行計(jì)算.

對(duì)于流線上點(diǎn)p的曲折度， 先算出起點(diǎn)和該點(diǎn)間的弧長(zhǎng)s(p)， 然后除以兩點(diǎn)間的歐幾里得距離， 其計(jì)算公式如下所示：

(3)

其中:p0是起始點(diǎn)位置.

1.4 速度方向熵

對(duì)于流線上每一點(diǎn)， 使用該點(diǎn)小范圍領(lǐng)域的N點(diǎn)來計(jì)算熵值. 為了計(jì)算熵值， 使用文獻(xiàn)[11]的算法把單位球體分成50份相同面積的區(qū)域， 并判斷每個(gè)矢量指向哪一塊. 對(duì)于流線上點(diǎn)p的速度方向熵可通過下面的公式進(jìn)行計(jì)算：

(4)

其中:di是領(lǐng)域中樣本點(diǎn)落入i方向塊的頻率.

2 基于特征分布的流線相似性查詢與分類

2.1 流線特征的2D直方圖

流線特征分布可選1D直方圖與2D直方圖. 對(duì)于1D直方圖， 存在兩條形狀不相同的流線并可能有著相似的特征分布的問題. 例如， 圖1(d)曲線由圖1(a)曲線經(jīng)過分割后拼接生成. 顯然， 它們有著相同的一維特征分布， 但是兩者的形狀并不相同. 兩條曲線的1D曲率直方圖分別如圖1(b)和圖1(e)所示.

本研究采用2D直方圖， 對(duì)于給定的流線， 將其分成M段， 每段具有相同數(shù)目的樣本點(diǎn)數(shù)， 接著對(duì)每一段建立一個(gè)1D直方圖， 從而構(gòu)建2D直方圖. 由于攜帶了空間信息， 2D直方圖可彌補(bǔ)由于丟失空間信息而導(dǎo)致的流線相似性度量錯(cuò)誤的不足.

圖1(c)和圖1(f)分別為曲線1(a)和曲線1(d)的2D直方圖； 在把曲線分成4等分后， 曲線各段的特征分布情況便體現(xiàn)在2D直方圖中.

2.2 基于k-means的三維流線分類

本文中k-means新的聚類中心計(jì)算公式為：

(5)

其中:C(j， t)為第j簇新中心的第t個(gè)特征直方圖；H(i， t)為流線i(或聚類中心i)的第t個(gè)特征直方圖；nj為第j簇所包含的流線下標(biāo)集合.

流線(或聚類中心)i，j的EMD相似性距離度量公式為：

(6)

k-means聚類的目標(biāo)函數(shù)為：

(7)

計(jì)算具體步驟如下：

Step1: 選擇k條流線作為初始聚類中心；

Step2: 使用公式(6)計(jì)算出流線到各簇類域間的距離， 并把流線劃分到距離最短的一類中；

Step3: 使用公式(5)更新各簇類中心；

Step4: 重復(fù)Step2和Step3直到目標(biāo)函數(shù)(7)小于一定的誤差或者達(dá)到一定的迭代次數(shù)， 終止計(jì)算.

2.3 交互式流線相似性查詢與分類方法

首先通過密集播種， 生成大量的流線； 然后根據(jù)每條流線的弧長(zhǎng)， 均勻地插值采樣相同數(shù)量的N個(gè)點(diǎn)； 接著計(jì)算流線上每點(diǎn)的特征值并建立2D直方圖； 在最后的用戶交互里， 根據(jù)用戶選擇的流線， 利用EMD計(jì)算每條流線與目標(biāo)流線的直方圖距離， 根據(jù)指定的數(shù)量比例顯示最相似的流線組； 然后使用k-means聚類方法進(jìn)行流線分類， 按用戶要求顯示感興趣的聚類結(jié)果. 在交互期間用戶可以調(diào)整流線段數(shù)和特征值等分?jǐn)?shù)； 也可以采用不同屬性組合或者賦予不同權(quán)值進(jìn)行相似性比較試驗(yàn).

3 實(shí)驗(yàn)結(jié)果與分析

本研究使用兩組數(shù)據(jù)進(jìn)行實(shí)驗(yàn). 第一組數(shù)據(jù)由Roger Crawfis的臺(tái)風(fēng)模擬算法生成， 第二組數(shù)據(jù)是來自文獻(xiàn)[12]提供的五臨界點(diǎn)數(shù)據(jù)集. 下面對(duì)這兩組數(shù)據(jù)分別進(jìn)行流線相似性查詢和流線聚類實(shí)驗(yàn).

3.1 運(yùn)行時(shí)間測(cè)試

實(shí)驗(yàn)采用CPU-GPU混合解決方案. 硬件配置為： Intel Core i5-3470 CPU， 頻率 3.2 GHz； 16 G內(nèi)存； GT630顯卡. 使用CPU進(jìn)行流線的生成， 采用CUDA技術(shù)對(duì)特征以及直方圖進(jìn)行并行計(jì)算. 表1列出了特征值估計(jì)、 直方圖構(gòu)建及相似性比較所用的時(shí)間.

表1 實(shí)現(xiàn)時(shí)間Tab.1 Implementation time

從表1中可以看出， 特征值估計(jì)和直方圖構(gòu)建所用時(shí)間都非常短， 用戶可以調(diào)整流線劃分的段數(shù)M和特征值的劃分?jǐn)?shù)N. 由于調(diào)整是對(duì)單個(gè)特征的直方圖進(jìn)行， 所以在1～2 s內(nèi)便可獲得調(diào)整效果反饋.

3.2 流線相似性查詢實(shí)驗(yàn)

試驗(yàn)中臺(tái)風(fēng)數(shù)據(jù)集和五臨界點(diǎn)數(shù)據(jù)集分別生成600條和500條流線， 每條流線采樣1 000點(diǎn). 流線均分成10段， 四種特征值均劃分為30塊. 速度方向熵中領(lǐng)域的取點(diǎn)數(shù)左右各取20點(diǎn). 在選取完目標(biāo)流線后， 加入各種屬性進(jìn)行相似性比較. 實(shí)驗(yàn)對(duì)比效果如圖2～4所示.

實(shí)驗(yàn)顯示結(jié)果可以明顯地看出：

1) 圖2(b)中, 當(dāng)臺(tái)風(fēng)流線最相似流線增加到50%時(shí)， 新增流線雖然也是漩渦狀， 但已開始向外伸展， 而圖2(c)中在加入全局屬性后新增流線依然貼近目標(biāo)流線.

2) 圖3～4中, 五臨界點(diǎn)數(shù)據(jù)流線在考慮全局屬性后整體也更向目標(biāo)流線靠攏.

3) 圖4(d)中, 不加入扭率， 選取的流線更為平坦.

由此可見, 不同屬性體現(xiàn)流線的不同形狀特征， 簡(jiǎn)單地將各特征直方圖通過EMD計(jì)算出距離后相加可能不足以突顯流線的形狀特征. 不同屬性間的組合， 賦予各特征直方圖相似性距離以不同的權(quán)值可以更好地體現(xiàn)流線形狀特征.

3.3 流線分類實(shí)驗(yàn)

實(shí)驗(yàn)中兩組數(shù)據(jù)均用k-means聚類方法， 在相同初始聚類中心的條件下聚成3類， 不同顏色的流線組表示不同的類簇.

圖5～6為聚類效果對(duì)比圖. 從結(jié)果中可以看出， 圖5(c)中紅色簇下端流線部分更接近于圖5(a)中藍(lán)色的流線簇. 考慮了全局屬性的聚類結(jié)果將紅色部分歸為藍(lán)色簇， 聚類后整體在形狀上更加貼合.

圖6為五臨界點(diǎn)數(shù)據(jù)流線聚類效果. 從結(jié)果可以看出， 加入全局屬性后黃色簇流線都較順直， 紅色簇流線中都是兩端螺旋. 從兩組數(shù)據(jù)的聚類結(jié)果可以看出， 加入全局屬性后簇內(nèi)流線更加相似， 簇間差別更加明顯.

4 結(jié)論

本研究通過特征2D直方圖距離進(jìn)行流線相似性比較， 克服了1D直方圖空間信息丟失的缺陷， 同時(shí)引入曲折度和速度方向熵兩個(gè)全局屬性， 使得流線相似性查詢和聚類結(jié)果都更具魯棒性， 整體形狀更加緊湊. 該方法采用CUDA技術(shù)對(duì)特征值和直方圖進(jìn)行并行計(jì)算， 消耗的時(shí)間都非常短， 為用戶提供了實(shí)時(shí)的交互查詢體驗(yàn). 該方法中直方圖如何量化， EMD距離如何計(jì)算， 是影響相似性比較結(jié)果的關(guān)鍵. 在流線特征值分布極不均勻時(shí)， 如何以較小的塊數(shù)進(jìn)行劃分， 如何設(shè)計(jì)EMD距離， 不僅影響相似性比較的計(jì)算量， 也影響其比較效果. 這些問題將是下一步的研究方向.

[1] KLASSEN R V， HARRINGTON S J. Shadowed hedgehogs: technique for visualizing 2D slices of 3D vector fields[C]//Proceedings of the 2nd Conference on Visualization′91. [S.l.]: IEEE Computer Society Press， 1991: 148-153.

[2] SUNDQUIST A. Dynamic line integral convolution for visualizing streamline evolution[J]. IEEE Transactions on Visualization and Computer Graphics， 2003， 9(3): 273-282.

[3] MA J， WANG C， SHENE C K. Coherent view-dependent streamline selection for importance-driven flow visualization[C]//Proceedings of SPIE - The International Society for Optical Engineering. San Francisco: SPIE， 2013， 8 654(5 755): 1-15.

[4] OZEN C A ， ROCKWELL D. Flow structure on a rotating plate[J]. Experiments in fluids， 2012， 52(1): 207-223.

[5] CHEN Y， COHEN J， KROLIK J. Similarity-guided streamline placement with error evaluation[J]. IEEE Transactions on Visualization and Computer Graphics， 2007， 13(6): 1 448-1 455 .

[6] YU H F， WANG C L， SHENE C K，etal. Hierarchical streamline bundles[J]. IEEE Transactions on Visualization and Computer Graphics， 2012， 18(8): 1 353-1 367.

[7] 魯大營， 朱登明， 王兆其. 三維流場(chǎng)的流線提取算法[J]. 計(jì)算機(jī)輔助設(shè)計(jì)與圖形學(xué)學(xué)報(bào)， 2013(5): 666-673.

[8] MCLOUGHLIN T， JONES M W， LARAMEE R S，etal. Similarity measures for enhancing interactive streamline seeding[J]. IEEE Transactions on Visualization and Computer Graphics， 2013， 19(8): 1 342-1 353.

[9] LU K， CHAUDHURI A， LEE T Y，etal. Exploring vector fields with distribution-based streamline analysis[C]//2013 IEEE Pacific Visualization Symposium. Sydney: IEEE， 2013: 257-264.

[10] LI Y F， WANG C L， SHENE C K. Streamline similarity analysis using bag-of-features[C]//Proceedings of SPIE - The International Society for Optical Engineering. Francisco: SPIE， 2013， 9017: 628-637.

[11] LEOPARDI P. A partition of the unit sphere into regions of equal area and small diameter[J]. Electronic Transactions on Numerical Analysis， 2006， 25(12): 309-327.

[12] YE X， KAO D， PANG A. Strategy for seeding 3D streamlines[C]//16th IEEE Visualization Conference VIS 2005. Minneapolis: IEEE Computer Society， 2005: 471-478.

(責(zé)任編輯： 林曉)

3D stremlines similarity analysis based on distribution of measurements

ZHANG Boxiong， XIE Huosheng

(College of Mathematics and Computer Science， Fuzhou University， Fuzhou， Fujian 350116， China)

Considering the feature of shape， classification and selection of streamlines advance the understanding of the features in 3D flow field.In this paper， we propose an efficient 3D streamlines similarity comparison method based on distribution of features. This method is based on the streamline similarity comparison method proposed by Lu Kewei et al. and introduces two global geometric properties: tortuosity and velocity direction entropy. At first， we divide each streamline into segments evenly， then we construct curvature histogram、 torsion histogram、 tortuosity histogram and velocity direction entropy histogram for each segment， and form corresponding 2D histograms.In the end， we use the Earth Mover’s Distance(EMD) and k-means clustering method to measure the similarity between streamlines and classify these streamlines. The final experiment showed that after combining the two global features， this method can produce more robust query and clustering results.

streamline similarity; feature distribution; 2D histogram; streamline clustering

10.7631/issn.1000-2243.2016.05.0633

1000-2243(2016)05-0633-06

2014-09-24

謝伙生(1964-)， 副教授， 主要從事智能圖形圖像處理、 數(shù)據(jù)挖掘、 大規(guī)模數(shù)據(jù)可視化、 機(jī)器學(xué)習(xí)等方面研究， xiehs@sina.com

福建省自然科學(xué)基金資助項(xiàng)目(2014J01229)

TP391.41

A