專利名稱：擴展的矢量場線卷積積分方法
技術(shù)領(lǐng)域：
本發(fā)明涉及一種擴展的矢量場線卷積積分LIC方法，屬于科學計算可視化領(lǐng)域。
背景技術(shù)：
在科學計算可視化研究領(lǐng)域，矢量場可視化的研究目的在于運用計算機提供的圖形、圖像、視頻等多種手段將科學實驗、仿真、測量得到的大量矢量場數(shù)據(jù)準確、有效地展現(xiàn)出來，并支持人機間的多種交互處理。矢量場可視化的方法首先必須保證矢量場信息描述的準確性和完整性；其次還應(yīng)具備較高的計算效率和較強的針對矢量場結(jié)構(gòu)及幅度的表現(xiàn)力；矢量場可視化方法中的交互處理主要指多角度及多分辨處理。為了達成上述目標，矢量場可視化的發(fā)展大致經(jīng)歷了三個階段，第一個階段是基于圖標(例如小的箭頭)的方法，此類方法的主要缺陷在于無法有效地表現(xiàn)數(shù)據(jù)密集的矢量場；第二個階段是基于流線，流面和流體的方法，此類方法的主要缺陷在于流線、流面、流體的形成預(yù)先選定的種子點關(guān)系密切，有可能忽略矢量場中重要的特征信息；第三個階段是基于紋理映射的方法，主要有點噪聲紋理合成的方法和線卷積積分方法LIC。基于紋理映射的可視化方法尤其是基于LIC的方法，由于能夠較好的克服前兩類方法的主要缺陷，因而得到了快速的發(fā)展，其應(yīng)用的領(lǐng)域十分廣泛。
目前，市面上有一些針對可視化的商業(yè)或開源軟件，其中有一部分提供了矢量場可視化功能，但是大多采用基于幾何形狀或基于顏色、光學特征的矢量場映射方法，而且針對的都是結(jié)構(gòu)化矢量場。其中，AVS/Express使用流線、矢量、圖形符來表示矢量；IDL使用箭頭、顏色等表示矢量；OpenDX也使用顏色和箭頭來表示矢量；VTK使用流線等幾何形狀表示矢量，因?qū)儆陂_源項目，目前已有人為VTK開發(fā)了基于LIC方法的矢量場可視化方法，但針對的是結(jié)構(gòu)化矢量場；ParaView本身是基于VTK的，因此其映射方法與VTK相同。
非結(jié)構(gòu)矢量場的最大特點是采樣點位置隨機，將傳統(tǒng)的結(jié)構(gòu)化的線卷積積分方法LIC直接應(yīng)用于非結(jié)構(gòu)矢量場的LIC處理存在的問題主要有兩個方面(1)結(jié)構(gòu)化的LIC方法會因分辨率等因素導致矢量場信息的失真。以圖2所示的非結(jié)構(gòu)矢量場為例，有兩種方法可以將其轉(zhuǎn)換為結(jié)構(gòu)化矢量場，一種方法是在網(wǎng)格內(nèi)對隨機采樣點上的矢量進行加權(quán)處理，獲得一個位于網(wǎng)格中心的矢量(如圖4所示)，通過對比圖4b和圖2b可以看出，這種加權(quán)處理由于降低了采樣點的分辨率，導致矢量場信息失真。另一種方法是通過加密網(wǎng)格使非結(jié)構(gòu)矢量場的隨機采樣點均落在單一網(wǎng)格的中心，從而在保證分辨率的前提下將非結(jié)構(gòu)矢量場轉(zhuǎn)化成為結(jié)構(gòu)化矢量場，但由于隨機采樣點的位置有可能非常接近，這樣的處理方法人為提高了網(wǎng)格密度，導致計算效率低下。此外，結(jié)構(gòu)化LIC的DDA曲線(矢量線，下同)形成過程是以網(wǎng)格單元為基本單位的，這就意味著結(jié)構(gòu)化LIC處理的最高分辨率為網(wǎng)格單元的大小，分辨率的限制也會導致矢量場信息的失真。
(2)結(jié)構(gòu)化LIC方法的矢量幅度表現(xiàn)力不強。在矢量場的研究領(lǐng)域，人們已經(jīng)習慣于用矢量線的疏密來表示矢量場幅度的大小(如圖5所示)，而傳統(tǒng)的結(jié)構(gòu)化LIC方法并不支持這種幅度表現(xiàn)方式。
本發(fā)明經(jīng)過國內(nèi)外專利文獻和非專利文獻的檢索尚未發(fā)現(xiàn)針對非結(jié)構(gòu)矢量場的線卷積積分方法。

發(fā)明內(nèi)容
本發(fā)明的技術(shù)解決問題是克服現(xiàn)有技術(shù)的不足，提供一種擴展的矢量場線卷積積分LIC方法，該方法既可以面向非結(jié)構(gòu)矢量場，又可面向結(jié)構(gòu)矢量，增強矢量幅度的表現(xiàn)力，提高了預(yù)處理計算速度。
本發(fā)明的技術(shù)解決方案擴展的矢量場線卷積積分LIC方法，其特點在于通過以下步驟實現(xiàn)(1)從矢量場中尋找下一個未被標記的采樣點，如果沒有找到，則對矢量場的LIC值進行歸一化處理，進入步驟(5)；如果找到，則進入下個步驟(2)；(2)以上述步驟(1)找到的采樣點為起點，采用基于數(shù)值計算的DDA曲線形成方法形成DDA曲線，該方法是首先從起點出發(fā)進行曲線的前向預(yù)測，然后還是從點出發(fā)進行曲線的后向預(yù)測，最后將前向預(yù)測得到的曲線和后向預(yù)測得到的曲線合并，形成最終的DDA曲線；(3)對此輪循環(huán)形成的DDA曲線進行有色噪聲序列直接賦值，計算DDA曲線上各點的LIC結(jié)果；(4)依據(jù)此輪循環(huán)形成的DDA曲線對整個矢量場進行區(qū)域標記，依據(jù)該曲線上的矢量場采樣點的幅度大小確定覆蓋區(qū)域，標記該曲線區(qū)域中的所有矢量場采樣數(shù)據(jù)點，并返回步驟(1)；(5)對矢量場的LIC值進行歸一化處理。
本發(fā)明與現(xiàn)有技術(shù)相比的有益效果是(1)現(xiàn)有的LIC方法不適用于非結(jié)構(gòu)矢量場的LIC處理，而本發(fā)明給出的方法和流程即能處理各類結(jié)構(gòu)化矢量場，也能處理各類非結(jié)構(gòu)矢量場，這就為各類科學計算、仿真、及實測獲得的矢量數(shù)據(jù)(如醫(yī)學檢測數(shù)據(jù)、地質(zhì)勘探測量數(shù)據(jù)、氣象數(shù)據(jù)、計算流體力學數(shù)據(jù)等)的可視化提供了一個更為方便和實用的手段。
(2)現(xiàn)有的LIC方法一方面不能很好地表現(xiàn)矢量場的幅度信息，另一方面會導致三維情況下處理結(jié)果因密集而互相遮擋，影響矢量場場結(jié)構(gòu)的表現(xiàn)效果。而本發(fā)明采用的是區(qū)域標記技術(shù)，該技術(shù)在矢量場幅度方面具有更強的表現(xiàn)力，體現(xiàn)為紋理稀疏的位置幅度小，紋理密集的位置幅度大，這一方法合理地實現(xiàn)了LIC結(jié)果的稀疏化，因而在矢量場結(jié)構(gòu)信息的方面也具有更好的表現(xiàn)力，有利于各領(lǐng)域的分析人員更好地把握矢量場的結(jié)構(gòu)，從而找出最有價值的數(shù)據(jù)信息。
(3)現(xiàn)有的LIC方法采用隨機紋理累加計算的方式進行曲線紋理映射，而本發(fā)明提出可以通過有色噪聲序列直接進行曲線紋理映射，一方面適用于結(jié)構(gòu)化的數(shù)據(jù)場，另一方面解決了非結(jié)構(gòu)矢量場中無結(jié)構(gòu)化輸入紋理情況下的紋理映射問題，同時可以通過預(yù)處理提高計算速度。


圖1為本發(fā)明流程圖；圖2為本發(fā)明DDA曲線前向預(yù)測示意圖，其中P0為當前DDA曲線點的位置，A0為P0點處的矢量，P為前向預(yù)測的下一點位置，A為P點處的矢量，Pi，i＝1，2，3，4為預(yù)測區(qū)域內(nèi)的矢量場采樣點位置，Ai，i＝1，2，3，4為Pi處的矢量；圖3為本發(fā)明DDA曲線后向預(yù)測示意圖，其中P0為當前DDA曲線點的位置，A0為P0點處的矢量，P為后向預(yù)測的下一點位置，A為P點處的矢量，Pi，i＝1，2，3，4為預(yù)測區(qū)域內(nèi)的矢量場采樣點位置，Ai，i＝1，2，3，4為Pi處的矢量；圖4為本發(fā)明區(qū)域標記示意圖，其中，c為DDA曲線上一點矢量幅度的大小，w為該點相應(yīng)的區(qū)域?qū)挾龋粓D5-7為傳統(tǒng)LIC方法與本發(fā)明方法獲得的處理結(jié)果的對比圖，其中圖5為放射形矢量場的LIC結(jié)果，左圖為傳統(tǒng)LIC的結(jié)果，右圖為本發(fā)明方法的結(jié)果，圖像分辨率均為210*210象素；圖6為圓形矢量場的LIC結(jié)果，左圖為傳統(tǒng)LIC的結(jié)果，右圖為本發(fā)明方法的結(jié)果，圖像分辨率均為210*210象素；圖7為雷達矢量場的LIC結(jié)果，左圖為傳統(tǒng)LIC的結(jié)果，右圖為本發(fā)明方法的結(jié)果，圖像分辨率均為210*210象素。
具體實施例方式
本發(fā)明的具體實現(xiàn)方法如下(1)從矢量場中尋找下一個未被標記的采樣點，如果沒有找到，則對矢量場的LIC值進行歸一化處理，結(jié)束；如果找到，則進入步驟(2)；(2)以步驟(1)找到的采樣點為起點，采用基于數(shù)值計算的DDA曲線形成方法形成DDA曲線，首先從起點出發(fā)進行曲線的前向預(yù)測，然后還是從點出發(fā)進行曲線的后向預(yù)測，最后將前向預(yù)測得到的曲線和后向預(yù)測得到的曲線合并，形成最終的DDA曲線。
如圖2所示，DDA曲線前向預(yù)測方法如下前向預(yù)測體現(xiàn)為一個迭代循環(huán)，如為第一次循環(huán)，則以步驟(1)獲得的采樣點為當前的DDA曲線點，以該采樣點處的矢量作為當前DDA曲線點處的矢量；否則以上一次前向預(yù)測得到的點作為當前DDA曲線點，以該點處的矢量作為當前DDA曲線點處的矢量。
每一次循環(huán)的處理方法如下為表述方便，設(shè)當前DDA曲線點的位置為P0＝[P0，x，P0，y，P0，z]，P0點處的矢量為A0＝[A0，x，A0，y，A0，z]，以P0為原點，以A0為基準，根據(jù)事先設(shè)定的半徑R(例如取10.0)和前向角θ(例如取180°)，可以得到一個扇形區(qū)域(三維情況下為一個圓錐體)，則矢量場在該區(qū)域中存在一組隨機采樣點Pi＝[Pi，x，Pi，y，Pi，z]，i＝1，2，…N，Pi點上的矢量為Ai＝[Ai，x，Ai，y，Ai，z]，Pi和Ai為給定的矢量場中的已知值。
沿前向計算兩鄰近點位置事先設(shè)定一個Δk(例如取0.1)，可以計算出兩點位置PΔk＝[PΔk，x，PΔk，y，PΔk，z]和P2Δk＝[P2Δk，x，P2Δk，y，P2Δk，z]PΔk,x=P0,x+ΔkA0,x|A0|,PΔk,y=P0,y+ΔkA0,y|A0|,PΔk,z=P0,z+ΔkA0,z|A0|]]>P2Δk,x=P0,x+2ΔkA0,x|A0|,P2Δk,y=P0,y+2ΔkA0,y|A0|,P2Δk,z=P0,z+2ΔkA0,z|A0|]]>根據(jù)兩鄰近點位置計算兩鄰點處的矢量獲得PΔk和P2Δk兩點位置后，可以計算兩點PΔk和P2Δk處的矢量AΔk＝[AΔk，x，AΔk，y，AΔk，z]，A2Δk＝[A2Δk，x，A2Δk，y，A2Δk，z]AΔk,x=Σi=1Nexp{-(PΔk,x-Pi,x)2+(PΔk,y-Pi,y)2+(PΔk,z-Pi,z)2Σj=1Nexp{-(PΔk,x-Pj,x)2+(PΔk,y-Pj,y)2+(PΔk,z-Pj,z)2Ai,x]]>AΔk,y=Σi=1Nexp{-(PΔk,x-Pi,x)2+(PΔk,y-Pi,y)2+(PΔk,z-Pi,z)2Σj=1Nexp{-(PΔk,x-Pj,x)2+(PΔk,y-Pj,y)2+(PΔk,z-Pj,z)2Ai,y]]>
AΔk,z=Σi=1Nexp{-(PΔk,x-Pi,x)2+(PΔk,y-Pi,y)2+(PΔk,z-Pi,z)2Σj=1Nexp{-(PΔk,x-Pj,x)2+(PΔk,y-Pj,y)2+(PΔk,z-Pj,z)2Ai,z]]>以及A2Δk,x=Σi=1Nexp{-(P2Δk,x-Pi,x)2+(P2Δk,y-Pi,y)2+(P2Δk,z-Pi,z)2Σj=1Nexp{-(P2Δk,x-Pj,x)2+(P2Δk,y-Pj,y)2+(P2Δk,z-Pj,z)2Ai,x]]>A2Δk,y=Σi=1Nexp{-(P2Δk,x-Pi,x)2+(P2Δk,y-Pi,y)2+(P2Δk,z-Pi,z)2Σj=1Nexp{-(P2Δk,x-Pj,x)2+(P2Δk,y-Pj,y)2+(P2Δk,z-Pj,z)2Ai,y]]>A2Δk,z=Σi=1Nexp{-(P2Δk,x-Pi,x)2+(P2Δk,y-Pi,y)2+(P2Δk,z-Pi,z)2Σj=1Nexp{-(P2Δk,x-Pj,x)2+(P2Δk,y-Pj,y)2+(P2Δk,z-Pj,z)2Ai,z]]>根據(jù)兩鄰點處的矢量計算誤差的一階導數(shù)獲得PΔk和P2Δk的位置及相應(yīng)的矢量后，可以計算出誤差的一階導數(shù)e1′=(AΔk,x-A0,x)2+(AΔk,y-A0,y)2+(AΔk,z-A0,z)2Δk]]>e2′=(AΔk,x-A0,x)2+(AΔk,y-A0,y)2+(AΔk,z-A0,z)2Δk]]>根據(jù)誤差的一階導數(shù)計算誤差的二階導數(shù)獲得誤差的一階導數(shù)后，可以計算誤差的二階導數(shù)e′′=e2′-e1′Δk]]>根據(jù)誤差的一階及二階導數(shù)計算預(yù)測步長獲得誤差的一階及二階導數(shù)后，可以計算出預(yù)測步長k=-e1′+(e1′)2-2e′′|A0|μe′′,]]>其中μ為事先設(shè)定的允許誤差比例，例如取0.1。將k與最大允許步長v(v為事先設(shè)定的值，例如取1.0)加以比較，如果k大于v則令k＝v，如果k小于等于v，則k不變。
根據(jù)預(yù)測步長計算此次循環(huán)的下一點預(yù)測位置獲得此輪循環(huán)的預(yù)測步長k后，可以計算出前向預(yù)測的下一點位置P＝[Px，Py，Pz]
Px=P0,x+kA0,x|A0|,Py=P0,y+kA0,y|A0|,Pz=P0,z+kA0,z|A0|]]>根據(jù)下一點預(yù)測位置計算該點處的矢量獲得P點位置后，可以計算P點處的矢量A＝[Ax，Ay，Az]Ax=Σi=1Nexp{-(Px-Pi,x)2+(Py-Pi,y)2+(Pz-Pi,z)2Σj=1Nexp{-(Px-Pj,x)2+(Py-Pj,y)2+(Pz-Pj,z)2Ai,x]]>Ay=Σi=1Nexp{-(Px-Pi,x)2+(Py-Pi,y)2+(Pz-Pi,z)2Σj=1Nexp{-(Px-Pj,x)2+(Py-Pj,y)2+(Pz-Pj,z)2Ai,y]]>Az=Σi=1Nexp{-(Px-Pi,x)2+(Py-Pi,y)2+(Pz-Pi,z)2Σj=1Nexp{-(Px-Pj,x)2+(Py-Pj,y)2+(Pz-Pj,z)2Ai,z]]>至此完成了前向預(yù)測一次循環(huán)過程。
不斷以預(yù)測得到的點P及該點矢量A作為新的當前DDA曲線點及相應(yīng)的矢量，并按上述步驟迭代循環(huán)，直到P超出矢量場的空間范圍或矢量A的模為0，則結(jié)束前向預(yù)測，此時可以得到一串點構(gòu)成的前向預(yù)測曲線及各點相應(yīng)的矢量。
如圖3所示，DDA曲線后向預(yù)測方法如下與前向預(yù)測類似，后向預(yù)測也體現(xiàn)為一個迭代循環(huán)，如為第一次循環(huán)，則以步驟(1)獲得的采樣點為當前的DDA曲線點，以該采樣點處的矢量作為當前DDA曲線點處的矢量；否則以上一次后向預(yù)測得到的點作為當前DDA曲線點，以該點處的矢量作為當前DDA曲線點處的矢量。
每一次循環(huán)的處理方法如下為表述方便，設(shè)當前DDA曲線點的位置為P0＝[P0，x，P0，y，P0，z]，P0點處的矢量為A0＝[A0，x，A0，y，A0，z]，以P0為原點，以-A0為基準，根據(jù)事先設(shè)定的半徑R(例如取10.0)和后向角θ(例如取180°)，可以得到一個扇形區(qū)域(三維情況下為一個圓錐體)，則矢量場在該區(qū)域中存在一組隨機采樣點Pi＝[Pi，x，Pi，y，Pi，z]，i＝1，2，…N，Pi點上的矢量為Ai＝[Ai，x，Ai，y，Ai，z]，Pi和Ai為給定的矢量場中的已知值。
沿后向計算兩鄰近點位置事先設(shè)定一個Δk(例如取0.1)，可以計算出兩點位置PΔk＝[PΔk，x，PΔk，y，PΔkz]和P2Δk＝[P2Δk，x，P2Δk，y，P2Δk，z]PΔk,x=P0,x-ΔkA0,x|A0|,PΔk,y=P0,y-ΔkA0,y|A0|,PΔk,z=P0,z-ΔkA0,z|A0|]]>P2Δk,x=P0,x-2ΔkA0,x|A0|,P2Δk,y=P0,y-2ΔkA0,y|A0|,P2Δk,z=P0,z-2ΔkA0,z|A0|]]>根據(jù)兩鄰近點位置計算兩鄰點處的矢量獲得PΔk和P2Δk兩點位置后，可以計算兩點PΔk和P2Δk處的矢量AΔk＝[AΔk，x，AΔk，y，AΔk，z]，A2Δk＝[A2Δk，x，A2Δk，y，A2Δk，z]AΔk,x=Σi=1Nexp{-(PΔk,x-Pi,x)2+(PΔk,y-Pi,y)2+(PΔk,z-Pi,z)2Σj=1Nexp{-(PΔk,x-Pj,x)2+(PΔk,y-Pj,y)2+(PΔk,z-Pj,z)2Ai,x]]>AΔk,y=Σi=1Nexp{-(PΔk,x-Pi,x)2+(PΔk,y-Pi,y)2+(PΔk,z-Pi,z)2Σj=1Nexp{-(PΔk,x-Pj,x)2+(PΔk,y-Pj,y)2+(PΔk,z-Pj,z)2Ai,y]]>AΔk,z=Σi=1Nexp{-(PΔk,x-Pi,x)2+(PΔk,y-Pi,y)2+(PΔk,z-Pi,z)2Σj=1Nexp{-(PΔk,x-Pj,x)2+(PΔk,y-Pj,y)2+(PΔk,z-Pj,z)2Ai,z]]>以及A2Δk,x=Σi=1Nexp{-(P2Δk,x-Pi,x)2+(P2Δk,y-Pi,y)2+(P2Δk,z-Pi,z)2Σj=1Nexp{-(P2Δk,x-Pj,x)2+(P2Δk,y-Pj,y)2+(P2Δk,z-Pj,z)2Ai,x]]>A2Δk,y=Σi=1Nexp{-(P2Δk,x-Pi,x)2+(P2Δk,y-Pi,y)2+(P2Δk,z-Pi,z)2Σj=1Nexp{-(P2Δk,x-Pj,x)2+(P2Δk,y-Pj,y)2+(P2Δk,z-Pj,z)2Ai,y]]>A2Δk,z=Σi=1Nexp{-(P2Δk,x-Pi,x)2+(P2Δk,y-Pi,y)2+(P2Δk,z-Pi,z)2Σj=1Nexp{-(P2Δk,x-Pj,x)2+(P2Δk,y-Pj,y)2+(P2Δk,z-Pj,z)2Ai,z]]>根據(jù)兩鄰點處的矢量計算誤差的一階導數(shù)獲得PΔk和P2Δk的位置及相應(yīng)的矢量后，可以計算出誤差的一階導數(shù)e1′=(AΔk,x-A0,x)2+(AΔk,y-A0,y)2+(AΔk,z-A0,z)2-Δk]]>e2′=(AΔk,x-A0,x)2+(AΔk,y-A0,y)2+(AΔk,z-A0,z)2-Δk]]>根據(jù)誤差的一階導數(shù)計算誤差的二階導數(shù)獲得誤差的一階導數(shù)后，可以計算誤差的二階導數(shù)e′′=e2′-e1′-Δk]]>根據(jù)誤差的一階及二階導數(shù)計算預(yù)測步長獲得誤差的一階及二階導數(shù)后，可以計算出預(yù)測步長k=-e1′+(e1′)2-2e′′|A0|μe′′,]]>其中μ為事先設(shè)定的允許誤差比例，例如取0.1。將k與最大允許步長v(v為事先設(shè)定的值，例如取1.0)加以比較，如果k大于v則令k＝v，如果k小于等于v，則k不變根據(jù)預(yù)測步長計算此次循環(huán)的下一點預(yù)測位置獲得此輪循環(huán)的預(yù)測步長k后，可以計算出前向預(yù)測的下一點位置P＝[Px，Py，Pz]Px=P0,x-kA0,x|A0|,Py=P0,y-kA0,y|A0|,Pz=P0,z-kA0,z|A0|]]>根據(jù)下一點預(yù)測位置計算該點處的矢量獲得P點位置后，可以計算P點處的矢量A＝[Ax，Ay，Az]Ax=Σi=1Nexp{-(Px-Pi,x)2+(Py-Pi,y)2+(Pz-Pi,z)2Σj=1Nexp{-(Px-Pj,x)2+(Py-Pj,y)2+(Pz-Pj,z)2Ai,x]]>Ay=Σi=1Nexp{-(Px-Pi,x)2+(Py-Pi,y)2+(Pz-Pi,z)2Σj=1Nexp{-(Px-Pj,x)2+(Py-Pj,y)2+(Pz-Pj,z)2Ai,y]]>Az=Σi=1Nexp{-(Px-Pi,x)2+(Py-Pi,y)2+(Pz-Pi,z)2Σj=1Nexp{-(Px-Pj,x)2+(Py-Pj,y)2+(Pz-Pj,z)2Ai,z]]>
至此，完成了后向預(yù)測一次循環(huán)過程，不斷以預(yù)測得到的點P及該點矢量A作為新的當前DDA曲線點及相應(yīng)的矢量，并按上述步驟迭代循環(huán)，直到P超出矢量場的空間范圍或矢量A的模為0，則結(jié)束后向預(yù)測，此時可以得到一串點構(gòu)成的后向預(yù)測曲線及各點相應(yīng)的矢量。
將前向預(yù)測得到的曲線和后向預(yù)測得到的曲線合并的方法如下從后向預(yù)測的最后一點開始，按逆序連接后向預(yù)測得到的曲線上各點，再按順序連接前向預(yù)測得到的曲線上各點，即可形成一條完整的DDA曲線。
(3)對此輪循環(huán)形成的DDA曲線進行有色噪聲序列直接賦值，計算DDA曲線上各點的LIC結(jié)果。
為表述方便，設(shè)此輪循環(huán)生成的DDA曲線有M個點，記為Ci＝[Ci，x，Ci，y，Ci，z]，i＝1，2，…，M，并設(shè)到此輪循環(huán)開始前已經(jīng)獲得的結(jié)果點集為Rj＝[Rj，x，Rj，y，Rj，z]，j＝1，2，…，N，各結(jié)果點相應(yīng)的LIC值為Lj，j＝1，2，…，N，各結(jié)果點相應(yīng)的相關(guān)點數(shù)為nj，j＝1，2，…，N。
具體的計算步驟如下生成有色噪聲序列并為DDA曲線點賦LIC值具體方法為，先生成一個隨機白噪聲序列Sk，k＝1，2，…，M+K，K為附加長度(例如可以取10)，再由白噪聲序列構(gòu)造有色噪聲序列Sci，i＝1，2，…，M，Sci=Σk=ii+KSk.]]>Sci即為DDA曲線點Ci對應(yīng)的LIC值。
判斷DDA曲線上各點與已獲得的結(jié)果點集中點的關(guān)系，并做相應(yīng)處理將DDA曲線各點Ci，i＝1，2，…，M分別與以往循環(huán)獲得的結(jié)果集中的點Rj＝[Rj，x，Rj，y，Rj，z]，j＝1，2，…，N加以比較，如果Ci與所有Rj，j＝1，2，…，N的距離(Ci,x-Rj,x)2+(Ci,y-Rj,y)2+(Ci,z-Rj,z)2>G]]>(G為給定門限值，例如取1.0)，表明Ci與以往循環(huán)獲得的結(jié)果集中的任何一點均不相同同，則將Ci做為新的結(jié)果加入到結(jié)果集合中，并以Ci對應(yīng)的有色噪聲序列中的值Sci作為新結(jié)果的LIC值，同時將該新結(jié)果的相關(guān)次數(shù)置1；如果Ci與以往循環(huán)獲得的結(jié)果集中某Rj的距離(Ci,x-Rj,x)2+(Ci,y-Rj,y)2+(Ci,z-Rj,z)2≤G]]>(G為給定門限值，例如取1.0)，表明Ci與Rj相同，則不生成新的結(jié)果，而將Ci對應(yīng)的有色噪聲序列中的值Sci加到Rj的LIC值上，即新的Lj等于原Lj+Sci，并將Rj的相關(guān)次數(shù)加1，即新的nj等于原nj+1。
(4)依據(jù)此輪循環(huán)形成的DDA曲線對整個矢量場進行區(qū)域標記，依據(jù)該曲線上的矢量場采樣點的幅度大小確定覆蓋區(qū)域，標記該曲線區(qū)域中的所有矢量場采樣數(shù)據(jù)點。具體的計算步驟如下如圖4所示，以此輪循環(huán)獲得的DDA曲線為中心線，根據(jù)曲線上各點的矢量幅度c的大小決定相應(yīng)的區(qū)域?qū)挾萕，c越大則W越小(可令W＝(cmax-c)/K，K為某預(yù)先設(shè)定的常數(shù)，cmax為矢量場中采樣點矢量幅度的最大值)，從而形成一個標記區(qū)域(在三維空間，這樣的區(qū)域是一個旋轉(zhuǎn)體圍成的柱空間)，對落入該標記區(qū)域的矢量場中的隨機采樣點進行標記，被標記的點不在作為以后各輪循環(huán)的DDA曲線形成的起點。
(5)對矢量場的LIC值進行歸一化處理。具體的計算步驟如下對最終形成的LIC結(jié)果集{Rj，j＝1，2，…，N}，其中每一個Rj均對應(yīng)有一個LIC值Lj和一個相應(yīng)的相關(guān)點數(shù)nj，歸一化即將Rj的LIC值Lj除以相關(guān)點數(shù)nj作為Rj的最終LIC值。
矢量場數(shù)據(jù)預(yù)處理得到的最終結(jié)果為結(jié)果點集的位置序列{Rj，j＝1，2，…，N}，其中Rj＝[Rj，x，Rj，y，Rj，z]，及相應(yīng)的LIC值序列{Lj，j＝1，2，…，N}如圖5-7為傳統(tǒng)LIC方法與本發(fā)明方法獲得的處理結(jié)果的對比，左圖為傳統(tǒng)LIC的結(jié)果，右圖為本發(fā)明方法的結(jié)果，圖像分辨率均為210*210象素、以三種典型的矢量場為例進行對比實驗，在每組實驗中，隨機選擇40000個采樣點的位置，并給出采樣點上位置上的矢量，從而構(gòu)造出相應(yīng)的非結(jié)構(gòu)矢量場。再以非結(jié)構(gòu)化的矢量場為基礎(chǔ)，采用矢量插值的方法構(gòu)造出結(jié)構(gòu)化的矢量場。分別運用傳統(tǒng)的LIC方法和本文的LIC方法進行LIC處理。從兩種方法的結(jié)果對比可以看出，傳統(tǒng)的LIC方法因分辨率的原因不能準確地給出矢量線的走向，而本發(fā)明的方法由于采用了基于數(shù)值的DDA曲線形成方法，因而能夠更準確地描繪出矢量線走向的細節(jié)信息。此外，和傳統(tǒng)的結(jié)構(gòu)化LIC方法相比，本發(fā)明方法采用的區(qū)域標記方法在矢量場幅度方面具有更強的表現(xiàn)力，體現(xiàn)為紋理稀疏的位置幅度小，紋理密集的位置幅度大，這一方法合理地實現(xiàn)了LIC結(jié)果的稀疏化，因而在矢量場結(jié)構(gòu)信息方面也具有更好的表現(xiàn)力。
權(quán)利要求
1.擴展的矢量場線卷積積分方法，其特征在于通過以下步驟實現(xiàn)(1)從矢量場中尋找下一個未被標記的采樣點，如果沒有找到，則對矢量場的LIC值進行歸一化處理，結(jié)束；如果找到，則進入步驟(2)；(2)以步驟(1)找到的采樣點為起點，采用基于數(shù)值計算的DDA曲線形成方法形成DDA曲線，即首先從起點出發(fā)進行曲線的前向預(yù)測，然后還是從起點出發(fā)進行曲線的后向預(yù)測，最后起將前向預(yù)測得到的曲線和后向預(yù)測得到的曲線合并，形成最終的DDA曲線；(3)對此輪循環(huán)形成的DDA曲線進行有色噪聲序列直接賦值，計算DDA曲線上各點的LIC結(jié)果；(4)依據(jù)此輪循環(huán)形成的DDA曲線對整個矢量場進行區(qū)域標記，依據(jù)該曲線上的矢量場采樣點的幅度大小確定覆蓋區(qū)域，標記該曲線區(qū)域中的所有矢量場采樣數(shù)據(jù)點；(5)對矢量場的LIC值進行歸一化處理。
2.根據(jù)權(quán)利要求1所述的擴展的矢量場線卷積積分方法，其特征在于所述步驟(2)中的DDA曲線前向預(yù)測方法如下前向預(yù)測體現(xiàn)為一個迭代循環(huán)，如為第一次循環(huán)，則以步驟(1)獲得的采樣點為當前的DDA曲線點，以該采樣點處的矢量作為當前DDA曲線點處的矢量；否則以上一次前向預(yù)測得到的點作為當前DDA曲線點，以該點處的矢量作為當前DDA曲線點處的矢量，每一次循環(huán)的處理方法如下(1)沿前向計算兩鄰近點位置(2)根據(jù)兩鄰近點位置計算兩鄰點處的矢量(3)根據(jù)兩鄰點處的矢量計算誤差的一階導數(shù)(4)根據(jù)誤差的一階導數(shù)計算誤差的二階導數(shù)(5)根據(jù)誤差的一階及二階導數(shù)計算預(yù)測步長(6)根據(jù)預(yù)測步長計算此次循環(huán)的下一點預(yù)測位置(7)根據(jù)下一點預(yù)測位置計算該點處的矢量(8)將此次循環(huán)獲得的預(yù)測位置及矢量作為新的循環(huán)的當前DDA曲線點及矢量迭代循環(huán)，直至預(yù)測位置超出矢量場空間范圍或預(yù)測位置處的矢量為0矢量，終止前向預(yù)測，從而獲得由一串點構(gòu)成的前向曲線及相應(yīng)點處的矢量。
3.根據(jù)權(quán)利要求1所述的擴展的矢量場線卷積積分方法，其特征在于所述步驟(2)中的DDA曲線后向預(yù)測方法如下后向預(yù)測也體現(xiàn)為一個迭代循環(huán)，如為第一次循環(huán)，則以步驟(1)獲得的采樣點為當前的DDA曲線點，以該采樣點處的矢量作為當前DDA曲線點處的矢量；否則以上一次后向預(yù)測得到的點作為當前DDA曲線點，以該點處的矢量作為當前DDA曲線點處的矢量，每一次循環(huán)的處理方法如下(1)沿后向計算兩鄰近點位置(2)根據(jù)兩鄰近點位置計算兩鄰點處的矢量(3)根據(jù)兩鄰點處的矢量計算誤差的一階導數(shù)(4)根據(jù)誤差的一階導數(shù)計算誤差的二階導數(shù)(5)根據(jù)誤差的一階及二階導數(shù)計算預(yù)測步長(6)根據(jù)預(yù)測步長計算此次循環(huán)的下一點預(yù)測位置(7)根據(jù)下一點預(yù)測位置計算該點處的矢量(8)將此次循環(huán)獲得的預(yù)測位置及矢量作為新的循環(huán)的當前DDA曲線點及矢量迭代循環(huán)，直至預(yù)測位置超出矢量場空間范圍或預(yù)測位置處的矢量為0矢量，終止后向預(yù)測，從而獲得由一串點構(gòu)成的后向曲線及相應(yīng)點處的矢量。
4.根據(jù)權(quán)利要求1所述的擴展的矢量場線卷積積分方法，其特征在于所述的步驟(2)中將前向預(yù)測得到的曲線和后向預(yù)測得到的曲線合并的方法如下從后向預(yù)測的最后一點開始，按逆序連接后向預(yù)測得到的曲線上各點，再按順序連接前向預(yù)測得到的曲線上各點，即可形成一條完整的DDA曲線。
5.根據(jù)權(quán)利要求1所述的擴展的矢量場線卷積積分方法，其特征在于所述步驟(3)的DDA曲線進行有色噪聲序列直接賦值，計算DDA曲線上各點的LIC結(jié)果的方法如下(1)生成有色噪聲序列并為DDA曲線點賦LIC值；(2)判斷DDA曲線上各點與已獲得的結(jié)果點集中點的關(guān)系，并做相應(yīng)處理，如果DDA曲線上某點不與結(jié)果點集中的任何點相同，則將該DDA曲線點作為新的結(jié)果點加入結(jié)果點集；如果DDA曲線上某點與結(jié)果點集中的一點相同，則將DDA曲線上該點的LIC值疊加到結(jié)果點集中與之相同的點的LIC值上，并對該結(jié)果點集中點的相關(guān)點數(shù)做加1處理。
6.根據(jù)權(quán)利要求1所述的擴展的矢量場線卷積積分方法，其特征在于所述步驟(4)的區(qū)域標記方法如下以DDA曲線為中心線，根據(jù)曲線上各點的矢量幅度c的大小決定相應(yīng)的區(qū)域?qū)挾葁，c越大則w越小，從而形成一個標記區(qū)域，對落入該標記區(qū)域的所有矢量場中的隨機采樣點進行標記。
全文摘要
擴展的矢量場線卷積積分LIC方法，(1)從矢量場中尋找下一個未被標記的采樣點；(2)以上述找到的采樣點為起點，采用基于數(shù)值計算的DDA曲線形成方法形成DDA曲線，首先從起點出發(fā)進行曲線的前向預(yù)測，然后還是從起點出發(fā)進行曲線的后向預(yù)測，最后將前向預(yù)測得到的曲線和后向預(yù)測得到的曲線合并，形成最終的DDA曲線；(3)對此輪循環(huán)形成的DDA曲線進行有色噪聲序列直接賦值，計算DDA曲線上各點的LIC結(jié)果；(4)依據(jù)此輪循環(huán)形成的DDA曲線對整個矢量場進行區(qū)域標記，依據(jù)該曲線上的矢量場采樣點的幅度大小確定覆蓋區(qū)域，標記該曲線區(qū)域中的所有矢量場采樣數(shù)據(jù)點；(5)對矢量場的LIC值進行歸一化處理。本發(fā)明該方法既可以面向非結(jié)構(gòu)矢量場，又可面向結(jié)構(gòu)矢量場，增強了矢量幅度及矢量場結(jié)構(gòu)的表現(xiàn)力，提高了預(yù)處理計算速度。
文檔編號G06T1/00GK1873640SQ20061008950
公開日2006年12月6日 申請日期2006年6月30日 優(yōu)先權(quán)日2006年6月30日
發(fā)明者何兵, 趙沁平, 郝愛民, 王莉莉 申請人:北京航空航天大學