基于信息熵的流場特征提取及可視化研究

牛 嬋，梁 猛,高瑞波,聶俊嵐,*

(1.燕山大學(xué) 信息科學(xué)與工程學(xué)院，河北 秦皇島066004；2.燕山大學(xué) 河北省計算機虛擬技術(shù)與系統(tǒng)集成重點實驗室，河北 秦皇島066004)

0 引言

流場可視化中的特征可視化是矢量場數(shù)據(jù)研究領(lǐng)域中的熱門方向，通過提取數(shù)據(jù)特征，忽略冗余數(shù)據(jù)，在保證信息準(zhǔn)確性的同時大大減少了不必要的可視化數(shù)據(jù)量，能夠很好地應(yīng)用于復(fù)雜而龐大的海量數(shù)據(jù)研究。

特征提取的流場可視化方法是指提取并繪制用戶感興趣的特征信息，當(dāng)前主要方法是分析流場拓撲結(jié)構(gòu)，1991年由Helman等首先提出，通過計算矢量場中的一階奇點來定位簡單特征點并對其進行分類，可提取出流場中的基本結(jié)構(gòu)信息[1]。2000年，Sadarjeo等在此基礎(chǔ)上實現(xiàn)了矢量場中精細特征信息的自動提取，為實現(xiàn)基于交互的特征可選可視化方式奠定了基礎(chǔ)[2]。1998年，Tang等為提高標(biāo)量場拓撲結(jié)構(gòu)信息提取的精度，在標(biāo)量場中引入矢量化操作[3]，解決了標(biāo)量場特征提取技術(shù)難以實現(xiàn)特征選擇的問題。同年，Slilver等首次提出時變矢量場可視化的理論基礎(chǔ)及處理方法，提出了一種清晰的時變矢量數(shù)據(jù)可視化的框架[4]。2014年Chaudhuri等通過分析空間流線實現(xiàn)了對矢量場特征信息的探索[5]。2015年Li等通過監(jiān)測流線的線段分割對流場進行特征提取[6]。2016年，巴振宇等基于特征信息進行種子點選取實現(xiàn)了流線的多層次可視化[7]。2017年，魯大營等提出了流場可視化中最優(yōu)視點的選擇方法[8]。此外，為了實現(xiàn)流場不同層次上的細節(jié)特征信息的可視化，近年來提出的多種流場拓撲結(jié)構(gòu)簡化繪制的理論與方法，成為了領(lǐng)域研究的熱點。

信息熵量化了矢量場中方向信息變化的程度，可以幫助研究學(xué)者有效地確定變化劇烈的流場區(qū)域。2010年，Xu等人基于信息熵概念，首次提出了一種基于矢量信息內(nèi)容引導(dǎo)的流場可視化框架[9]，先計算出矢量場各子區(qū)域的信息熵值，迭代測量信息內(nèi)容并迭代布局種子點，直至達到滿意的效果。2011年，Chen等人在此基礎(chǔ)上對流線進行了相似性分組，在保持精度的同時有效解決了遮擋問題[10]。2013年，Tao等人在此基礎(chǔ)上引入了視點概念，根據(jù)全局矢量數(shù)據(jù)創(chuàng)建不同的視點作為視點集，實現(xiàn)了不同視點下的相應(yīng)流場信息可視化[11]。2014年，Jun等人為解決3D流場的特征信息可視化問題，設(shè)計了一種內(nèi)部特征信息的自動計算方法[12]，有效解決了特征定位、流線布局、選取視點等問題。2015年，Liu等基于信息熵概念，提出了一種基于Clifford代數(shù)的流場特征信息自動檢測方法[13]，根據(jù)預(yù)設(shè)的特征模板，通過Clifford卷積方法進行自動特征匹配操作，便于用戶尋找感謝興趣的特征區(qū)域。同年，張沙等基于條件熵的概念提出了一種種子點迭代布局更新的算法[14]，反復(fù)迭代直到條件熵收斂，該方法有效刻畫了原始矢量場但計算繁瑣。2016年，劉曉帆提出了基于信息熵與Clifford代數(shù)的流場特征檢測[15]，流場信息熵的引入避免了大量無效的模版匹配，提高了特征檢測效率。2018年黃冬梅等提出了基于貪婪策略的種子點選取方法[16]，得益于流場信息熵的計算，該方法對流場中的關(guān)鍵特征具有高度敏感性。

特征信息提取和交互式可視化探索是大規(guī)模矢量場數(shù)據(jù)可視化的兩大研究重點。本文從信息熵應(yīng)用于矢量場的角度出發(fā)，提出局部最大熵值和變化規(guī)律的概念，對矢量場的特征信息進行精確提取并基于此設(shè)計合理的種子點布局策略，同時對流場基于多分辨率的交互式可視化進行了探索。

1 特征提取

特征點定位是分析流場拓撲結(jié)構(gòu)的基礎(chǔ)?，F(xiàn)有方法能夠確定特征點的大概存在區(qū)域，但缺乏足夠的精確度，只能選取熵值較大的區(qū)域進行特征提取，因此很容易遺漏掉一些特征信息。為彌補現(xiàn)有方法的不足，本文對信息熵及特征點之間的關(guān)系進行了深入探索。流場中存在很多重要的特征信息，例如旋渦、源點、匯點、鞍點、湍流等。圖1展示了不同矢量場中的矢量信息分布情況。二維矢量場中的矢量角度量化范圍為θ∈[0,360),將其劃分成60個區(qū)間段，每個區(qū)間段的長度為6。圖1(a)和圖1(b)分別為馬鞍型矢量場和中心型矢量場及其相對應(yīng)的矢量信息分布情況，圖1(c)和圖1(d)為無特征矢量場及其相對應(yīng)的矢量信息分布情況，很顯然特征區(qū)域的矢量信息比非特征區(qū)域的信息分布要均勻，即信息的變化量大。圖1(a)、1(b)、1(c)、1(d)中矢量場所對應(yīng)的信息熵值依次為5.79、5.82、4.36、2.42。信息變化越大的矢量場，其信息熵值越大，就越有可能存在特征點，這為矢量場的特征信息提取提供了可靠有效的依據(jù)。

圖1 不同矢量場中的矢量信息分布
Fig.1 The vector information distribution in different vector fields

由特征點的定義可知特征點周圍區(qū)域的矢量方向變化最為劇烈，因此特征區(qū)域范圍內(nèi)的矢量數(shù)據(jù)必定是特征點所對應(yīng)的熵值最大，即特征點處熵值必為局部最大值。但局部最大熵值不一定為特征點，因此需要進一步進行分析判斷，為此提出一種基于局部最大熵值并結(jié)合線性插值原理的特征點定位算法。

1.1 特征點定位

假設(shè)矢量場中的某點P在各個分量上的速度均為零，即滿足(up,vp,wp)=0則稱P點為特征點。設(shè)一個特征點候選網(wǎng)格區(qū)域，(xi,yi)表示頂點Vi的坐標(biāo)位置，(ui,vi)表示頂點Vi的矢量值。對特征點候選單元網(wǎng)格的4個頂點數(shù)據(jù)的某個方向分量進行正負標(biāo)記，首先分析x方向，正值標(biāo)記為1，負值標(biāo)記為0。每個候選特征單元有4個頂點，每個頂點均可能有0、1兩種標(biāo)記情況，因此理論上總共可以有16種標(biāo)記情況。

以情況0001為例進行分析，若4個頂點的u方向(u0,u1,u2,u3)標(biāo)記為狀態(tài)(0,0,0,1)的話，由于實際矢量場的變化可看成是連續(xù)線性變的，因此候選特征網(wǎng)格中的邊V2V3必定存在分量u為0的P點，邊V0V3之間必定存在分量u為0的Q點。根據(jù)線性插值原理便可求出P、Q兩點的坐標(biāo)。P點坐標(biāo)的計算公式為

(1)

采用相同的方法可求Q出點的坐標(biāo)，則線段PQ表示該候選特征網(wǎng)格中u=0的集合，即u=0的等值線段。同理可求出v=0的等值線段MN。計算線段PQ和線段MN的交點記為O，判斷O點是否在當(dāng)前的候選特征網(wǎng)格內(nèi)，若在網(wǎng)格單元中則判定該交點為特征點。

拓撲結(jié)構(gòu)信息的基礎(chǔ)構(gòu)成包括特征點和描述特征點的積分曲線，其中積分曲線一般使用幾何可視化中的流線。根據(jù)特征點周圍區(qū)域方向變化規(guī)律，即特征點周圍區(qū)域的方向變化最劇烈，距離特征點越遠，區(qū)域的方向變化劇烈程度則隨之降低。將該規(guī)律對應(yīng)于熵場，對于簡單矢量場其表現(xiàn)形式為特征點處的熵值最大，向外延伸熵值依次降低，直到熵值變化趨于零。對于復(fù)雜矢量場，即特征區(qū)域相互作用，其表現(xiàn)形式為特征點處熵值最大，向外延伸依次降低，直到熵值逐漸增大，其臨界值的位置即為特征區(qū)域的邊界。

1.2 特征點類型檢測

特征點按類型分類主要分為中心點(center)、聚點(focus)、鞍點(saddle)、交點(node)等四大類型，其中focus類型和node類型又分別可被分為排斥(repelling)類型和吸引(attracting)類型這兩種情況，因此一共細分為6類。特征點類型的判斷，可以通過列出特征點位置的雅克比矩陣J，求其特征值，然后對特征值實部和虛部進行正負判斷，特征值實部和虛部的正負狀態(tài)代表著不同類型的特征點，特征點處的Jacobian矩陣形式為

(2)

根據(jù)矩陣J的特征值λ1和λ2，用Re(λi)、Im(λi)分別表示λi的實部和虛部，特征點類型的判斷依據(jù)如下：

1) 源點：0

2) 排斥鞍點：Re(λ1)<0

3) 吸引鞍點：Re(λ1)≤Re(λ2)<0

4) 匯點：Re(λ1)≤Re(λ2)≤Re(λ3)<0。

1.3 特征區(qū)域界定

結(jié)合特征區(qū)域的信息熵值變化規(guī)律對特征區(qū)域進行界定，逐層對特征區(qū)域內(nèi)的熵值進行比較，采用內(nèi)層與外層熵值的差值進行判定，計算方式為

Variation=Entropy1-Entropy2，

(3)

其中，Entropy1表示內(nèi)層熵值，Entropy2表示外層熵值，Variation表示內(nèi)層熵值與外層熵值的差值，若差值大于零，則表示仍在該特征區(qū)域的影響范圍內(nèi)，若該差值趨于零或等于零則表示到達該特征值的影響區(qū)域邊界。在進行熵值比較時，對于大規(guī)模流場數(shù)據(jù)來說，若對特征區(qū)域各層數(shù)據(jù)內(nèi)的熵值都進行差值比較運算，則計算規(guī)模必然龐大。不同特征所對應(yīng)的影響區(qū)域均呈現(xiàn)一個圓形，因此只需對特征區(qū)域4個方向上的值進行熵值數(shù)據(jù)比較運算即可。運用差值比較方法對特征區(qū)域進行界定計算，針對特征區(qū)域的4個方向上的信息熵值進行差值比較，即減少計算量，又保持精確度。圖2給出了熵值變化規(guī)律計算的模板。

圖2 熵值變化規(guī)律計算模板
Fig.2 A template for calculating the law of entropy change

2 基于特征信息的流場交互可視化

基于流線的矢量場可視化方法具有較好的視覺表達效果，但流線刻畫流場信息的質(zhì)量嚴(yán)重依賴于種子點分布策略的選取。本文采用基于特征點的種子點布局方案，并設(shè)計一種基于多分辨率的可視化方法，完成矢量場信息的高精度繪制。在判定特征點的位置和類型后，為了更加精確地描述特征點周圍的矢量流動模式，可以在特征點影響范圍區(qū)域內(nèi)，根據(jù)其特征類型來選取相應(yīng)的流線種子點播撒模式。

種子點布局方案的主要思想為根據(jù)特征點的位置確定種子點模板的播撒位置；根據(jù)特征點的特征類型選取合適的種子點模板類型；根據(jù)特征點的影響區(qū)域的界定來確定種子點模板的放置規(guī)模。實現(xiàn)了嚴(yán)格依據(jù)矢量場特征點的三大特征信息對流場進行種子點布局。

對于矢量場的拓撲結(jié)構(gòu)信息可視化展示，結(jié)合矢量場特征點的位置、類型、影響區(qū)域范圍進行針對性的種子點播撒并進行流線積分繪制。對于流場的全局信息可視化展示，在基于拓撲結(jié)構(gòu)種子點播撒方案的基礎(chǔ)上對非特征區(qū)域采取均勻的種子點播撒模式。如圖3所示，其中圖3(a)展示的是原始矢量場數(shù)據(jù)，圖3(b)是流場特征點的位置判斷及特征影響區(qū)域的界定結(jié)果，圖3(c)是流場拓撲結(jié)構(gòu)信息可視化的種子點布局方案，圖3(d)是流場全局信息可視化的種子點布局方案。

圖3 種子點布局策略
Fig.3 The change regulation of entropy in entropy field

本文采用GPU生成積分流線。GPU的優(yōu)勢主要是其中含多個幾何著色器(Geometry Shader，GS)，能夠?qū)崿F(xiàn)多個圖元同時進行計算，很大程度上降低了計算時間，GS是存在于GPU中的頂點著色器和片段著色器之間的一種著色器，它的計算單位為圖元。相較于頂點著色器其優(yōu)勢是可以在同一時間段處理某一圖元中的所有頂點信息，并且能對這些頂點進行修改、刪除等操作，能夠?qū)崿F(xiàn)對輸入的幾何形態(tài)進行修改并輸出另外一種全新的幾何形態(tài)。利用這一原理，將積分生成流線所用到的種子點信息轉(zhuǎn)換成圖元表示形式，將矢量場數(shù)據(jù)轉(zhuǎn)換成紋理表示形式，并將其傳輸進入圖形管線中，流線的積分生成過程全部在GS中執(zhí)行，如圖4所示。以這種方式采用PK4積分器計算生成流線，很好地彌補了其復(fù)雜的計算過程要求實現(xiàn)了實時性。

圖4 基于GPU的RK4積分生成流線方法
Fig.4 Based on GPU RK4 integral to generate streamline method

3 實驗及分析

實驗在64位Win7操作系統(tǒng)下完成，中央處理器為Intel i7-4770k 3.5GHz，內(nèi)存為32G，顯卡為NVIDIA GTX770。本文實驗所使用的數(shù)據(jù)是由NOAA公開的二維全球流場數(shù)據(jù)，數(shù)據(jù)規(guī)模為1 440×720，數(shù)據(jù)格式為規(guī)則網(wǎng)格結(jié)構(gòu)。

本文利用信息熵提取流場中方向變化劇烈的區(qū)域，首先基于局部最大熵值提取出特征點，然后嚴(yán)格根據(jù)特征信息來可視化矢量場。

如圖5所示，圖5(a)為比例信息熵的計算方法，主要思想是選取信息熵值較大的區(qū)域，在其中心位置布局菱形種子點模板并積分生成流線，但沒有針對性地計算特征點位置及類型。而圖5(b)為本文算法。使用本文算法所可視化出的矢量場特征信息明顯要多，因為算法是嚴(yán)格根據(jù)特征信息布局種子點并生成流線，而不是只針對采樣點進行流線積分生成操作。譬如用紅色方框畫出的是矢量場中重要的漩渦特征，而該特征信息在圖5(a)中并沒有可視化出來，這說明圖5(a)的方法不能繪制出矢量場中的所有特征信息，如果此基礎(chǔ)上利用條件熵對流線進行刪減、添加操作，但是會使得計算量變大，尤其不適于大規(guī)模矢量場。

圖5(c)為基于區(qū)域信息熵的計算效果，圖5(d)為本文算法繪制效果，其中方框畫出的特征信息，同樣在圖5(c)中被丟失。圖5(c)采用的方法主要思想是對矢量場數(shù)據(jù)進行粒度的網(wǎng)格劃分并分別計算出網(wǎng)格區(qū)域的信息熵值，對信息熵值較大的網(wǎng)格區(qū)域進行特征點檢測，并計算特征點類型然后根據(jù)矢量場的特征點位置和類型布局種子點。該方法雖然實現(xiàn)了基于特征點的位置和類型的計算，但其最大的不足在于只選取信息熵值較大的部分流場網(wǎng)格區(qū)域進行特征信息檢測，顯然很可能會遺漏掉某些存在特征信息的網(wǎng)格區(qū)域。而本文方法則是根據(jù)特征點區(qū)域熵值,首先找出可能是特征點的所有矢量點集合，再通過零點定理對其進行篩選，以確定特征檢測的完整性和準(zhǔn)確性。繪制效果圖驗證了本文算法的有效性。

圖5 本文算法與現(xiàn)有算法效果對比圖
Fig.5 The comparison between existing algorithm and this paper′s algorithm

在特征分布較為稀疏的矢量場區(qū)域，必然會將非特征區(qū)域劃入其中，導(dǎo)致特征種子點模板的播撒規(guī)模變大，因而不能精確地提取出流場的拓撲特征。在性能方面，常用的維諾圖劃分方法主要有分治算法、增量算法、掃描線算法和基于Delaunay三角剖分的算法，其中基于Delaunay三角剖分的算法速度相對較快，但其時間復(fù)雜度也達到了o(n2)，而本文提出的基于熵值變化規(guī)律的特征區(qū)域界定算法只需對特征點4個對角方向的熵值進行變化率計算，其時間復(fù)雜度為o(n)，因此本文算法在提高了繪制效果的同時也提高了計算性能。

在繪制性能方面，選取了三種單一分辨率與本文所提出的多分辨率方法在不同的視點與矢量場距離下進行了繪制幀率的性能比較，如表1所示，其中幀率(5)、幀率(10)、幀率(15)分別表示種子點密度為每經(jīng)度5個、10個、15個的種子點集合所構(gòu)成的流場種子點布局，其中視點距離單位為百萬米。

表1 種子密度與幀率的對比Tab.1 Comparison of seed density and frame rate

分析表中數(shù)據(jù)可得出隨著距離改變積分繪制性能呈下降趨勢。表中最后一列多分辨率為采用本文算法的實驗結(jié)果，在視點距離較近時，繪制幀率較每個經(jīng)度15個種子點集合的效果有所提高，而在視點距離較遠時，運算速度大幅度提升，說明多分辨率方法在保障數(shù)據(jù)特征提取的前提下，繪制性能有很大提高。

4 結(jié)論

特征信息提取和交互式可視化探索是當(dāng)前大規(guī)模矢量場數(shù)據(jù)可視化的兩大研究重點。針對傳統(tǒng)方法對特征區(qū)域界定不準(zhǔn)確問題，提出了一種基于局部最大熵值的特征點檢測方法，將信息熵應(yīng)用于矢量場的特征提取，實現(xiàn)了特征精確提取及種子點合理布局策略；為更好地展示數(shù)據(jù)的規(guī)律，設(shè)計實現(xiàn)了基于多分辨率的交互式可視化方法，通過實驗驗證了該方法在大規(guī)模矢量場中的有效性。在未來，本文的算法可以應(yīng)用于三維矢量場數(shù)據(jù)，在更加復(fù)雜的矢量數(shù)據(jù)中進行探索，實現(xiàn)用戶對特征信息的高精度刻畫。