陳濱，虞鴻，吳哲夫

(1.浙江工業(yè)大學(xué) 藝術(shù)學(xué)院，浙江 杭州 310023;2.浙江工業(yè)大學(xué) 信息工程學(xué)院，浙江 杭州 310023)

在光線追蹤繪制真實場景過程中，加速結(jié)構(gòu)是用來減少光線-圖元相交檢測所需要的時間［1］.早期提出的加速結(jié)構(gòu)包括Bentley提出的kd樹［2］，Rubin和Whited提出的 BVH［3-4］，以及 Fujimoto提出的均勻柵格法［5-6］，Sachidhar等人已經(jīng)把這種算法在GPU上實現(xiàn)［7］.對于均勻柵格法的加速結(jié)構(gòu)，它是將空間細(xì)分為具有相同大小的長方體形式的三維柵格.這種空間組織方式創(chuàng)建簡單，主要就是把空間中的圖元分配到關(guān)聯(lián)的網(wǎng)格當(dāng)中，因此它的建立速度要快于kd樹和BVH，但是它的光線遍歷并求交的速度要慢于kd樹和BVH.對于動態(tài)真實場景的繪制，每種加速結(jié)構(gòu)的效率是由兩部分組成:加速結(jié)構(gòu)建立時間和光線-圖元求交時間.在真實復(fù)雜場景繪制的時候，由于場景中的圖元分布未知，Wald等人認(rèn)識到很難選擇一種較為平衡的方法，最終提出了采用建立速度較快的相對比較平衡的grids方法來替代建立速度較慢的kd樹和BVH等方法［8］.Shirley提出只對每個擁有圖元的柵格構(gòu)造對象［9］，擁有圖元的柵格只存一個圖元指針，該優(yōu)化方法大大降低了該方法所需要的內(nèi)存.然而，均勻柵格法在光線遍歷的過程中還沒有很好的提升，主要因為在類似由一個較小的物體和較大的房間組成的場景中，有許多沒有涉及到實體的空間也會進(jìn)行同樣數(shù)量的柵格劃分［10］.因此，這種方法就增加了許多冗余的光線與加速結(jié)構(gòu)求交計算.

本文描述了一種改進(jìn)后的均勻柵格方法，它根據(jù)每個柵格所包含的圖元數(shù)量動態(tài)分割柵格，在那些比較空曠的空間就不進(jìn)行柵格劃分了，除去了因為沒有必要的柵格劃分而帶來的多余計算量.最后，在CPU上實現(xiàn)動態(tài)柵格加速結(jié)構(gòu)的光線追蹤方法.

1 柵格數(shù)據(jù)結(jié)構(gòu)

建立柵格數(shù)據(jù)結(jié)構(gòu)過程包括把整個空間中所有物體的包圍盒分割成大小相同的柵格(本文把空間中的一個長方體體元當(dāng)作一個柵格)，并把組成物體的圖元(一般為三角形)放入這些柵格中.在均勻柵格法當(dāng)中，柵格的數(shù)量是根據(jù)整個空間中圖元的個數(shù)決定的，每個維度的柵格數(shù)量nx、ny、nz可進(jìn)行如下計算［11］:

式中:wx、wy、wz為空間中所有物體包圍盒 x、y、z方向的跨度，m用來改變柵格數(shù)量的系數(shù)，通常取2，s為每個柵格的邊長.由于它均勻地分割成相同大小的柵格，可將這些柵格根據(jù)對應(yīng)的坐標(biāo)序號放入一個一維數(shù)組當(dāng)中，如圖1.柵格數(shù)組中的每個柵格對象保存著指向與該柵格關(guān)聯(lián)的三角形圖元指針，若柵格沒有關(guān)聯(lián)的三角形圖元，那么可以把這個柵格對象設(shè)置為空.

圖1 均勻柵格法的空間劃分和數(shù)據(jù)結(jié)構(gòu)表示Fig.1 Uniform grids space partition and data struct presentation

動態(tài)柵格法的空間劃分也是基于均勻柵格劃分策略，但是動態(tài)柵格劃分并不是一下子均勻地劃分成相等大小的柵格，如圖2.柵格樹的根部代表整個空間所有物體的包圍盒，每個子節(jié)點(diǎn)表示父節(jié)點(diǎn)進(jìn)行一次均勻劃分4個部分的子?xùn)鸥?和均勻柵格的數(shù)據(jù)表示一樣，每個柵格對象保存著指向與該柵格關(guān)聯(lián)的三角形圖元指針，若沒有關(guān)聯(lián)圖元對象可以設(shè)置為空.柵格樹有2個重要的性質(zhì):1)涉及多個柵格圖元，它的指針可能會保存在多個柵格節(jié)點(diǎn);2)樹左側(cè)節(jié)點(diǎn)中圖元離視點(diǎn)距離一定小于樹右側(cè)節(jié)點(diǎn)中圖元離視點(diǎn)的距離，并且隨著樹的深度變化，節(jié)點(diǎn)中的圖元離視點(diǎn)的距離也越遠(yuǎn).這樣在光線遍歷的過程中，如發(fā)現(xiàn)與某個圖元相交，那么就結(jié)束該光線的遍歷.

圖2 動態(tài)柵格法的空間劃分和數(shù)據(jù)結(jié)構(gòu)表示Fig.2 Dynamic grids space partition and data struct presentation

2 基于視點(diǎn)距離的柵格建立

為了提高光線追蹤中光線與柵格遍歷速度，許多方法在建立加速數(shù)據(jù)結(jié)構(gòu)的時候都進(jìn)行了排序［12-13］.根據(jù)從頂層到底層的順序建立柵格樹，這個從頂層到底層的順序?qū)?yīng)著柵格節(jié)點(diǎn)中的圖元距離從近到遠(yuǎn).先把根節(jié)點(diǎn)初始化，對空間中所有物體的包圍盒做一次均勻柵格劃分，此次劃分在x、y、z方向一分為二，即在三維空間中把物體包圍盒分成8個柵格.然后根據(jù)離視點(diǎn)的距離，從左到右依次作為根節(jié)點(diǎn)的子節(jié)點(diǎn).在本文的實現(xiàn)當(dāng)中，視點(diǎn)放在z軸的負(fù)半軸上，因此柵格的z軸序列號較小的就是離視點(diǎn)較近的.最后判斷這些子節(jié)點(diǎn)中的圖元是否到達(dá)柵格的飽和狀態(tài)，若沒有到達(dá)柵格的飽和狀態(tài)，那么用相同的方法繼續(xù)劃分子?xùn)鸥?，反之則停止.圖3描述了一個柵格的飽和狀態(tài)，左圖中的柵格中有2個三角形圖元a和b，如果對這個柵格再進(jìn)行每個維度二分的柵格劃分，那么它的每個子?xùn)鸥襁€是關(guān)聯(lián)著2個圖元三角形，這樣的重復(fù)操作將會消耗很大的遞歸?？臻g，因此在實現(xiàn)當(dāng)中，針對遇到這種情況就停止柵格劃分.

圖3 柵格的飽和狀態(tài)Fig.3 Saturation state of the grid

如上所述的柵格數(shù)據(jù)結(jié)構(gòu)構(gòu)造方法，顯然比較適合一些非規(guī)則分布的模型，例如當(dāng)視點(diǎn)對著廣闊的天空，天空中排布這零星的物體，這些物體相對天空就比較小.因此在天空這部分做柵格劃分的時候，很快就到達(dá)了柵格飽和狀態(tài)，而劃分物體的時候會劃分地比較密，算法1給出了本文的動態(tài)柵格劃分方法.

算法1 動態(tài)柵格建立

輸入:空間中三角形圖元數(shù)組objects，這些圖元的包圍盒bbox，柵格根節(jié)點(diǎn)recursivecells

1)if recursivecells equal null then

2)初始化根節(jié)點(diǎn)的子?xùn)鸥駭?shù)，numcells=2×2×2，并把根節(jié)點(diǎn)recursivecells的8個子?xùn)鸥穸荚O(shè)置為null;

3)for each objectsi，do

4)objbbox=objectsi-＞bbox;

5)根據(jù)單個圖元包圍盒objbbox的位置，把它分配到父節(jié)點(diǎn)包圍盒bbox其中8個子?xùn)鸥癞?dāng)中，得到這個圖元在父節(jié)點(diǎn)當(dāng)中的3個維度最大最小序列號 cellminx，cellminy，cellminz，cellmaxx，cellmaxy，cellmaxz

6)for iz from cellminz to cellmaxz

for iy from cellminy to cellmaxy

for ix from cellminx to cellmaxx

7)獲得子?xùn)鸥裥蛄刑朿ellindex=ix+iy×2+iz×2;

8)if recursivecells［cellindex］equal null then

9)構(gòu)造一個類型和根節(jié)點(diǎn)一樣的柵格對象recursivecell，把 objectsi分配到 recursivecell

10)else then

直接把objectsi分配到recursivecell

11)for each recursivecelli

if recursivecelli has reached saturate state then break;

else then繼續(xù)調(diào)用該函數(shù);

3 動態(tài)柵格建立算法分析

算法1在實現(xiàn)動態(tài)柵格建立的過程當(dāng)中采用函數(shù)遞歸的方法，必須考慮到算法效率和調(diào)用?？臻g開銷中尋求平衡.首先有以下假設(shè):整個空間中圖形元組成的包圍盒是一個正方體，這個包圍盒中有N(8的冪次數(shù))個三角形圖元;每個柵格的飽和狀態(tài)為容納一個三角形圖形元.

算法1消耗的時間主要由均勻柵格劃分時間、分配圖元時間組成.無論柵格大小、均勻柵格劃分時間Tgrid都是一樣的.分配圖元的時間與圖元的數(shù)量成正比，Tprimitive為分配單位圖元的時間.算法1中，最壞的情況就是每次柵格都分割為8個子?xùn)鸥瘢脠D4中的遞歸樹分析可得，算法1在最壞情況下的消耗時間T為

式中:i表示在遞歸樹中的深度，根節(jié)點(diǎn)的深度為0.深度為i的子節(jié)點(diǎn)個數(shù)為8i，每個子節(jié)點(diǎn)中的圖元數(shù)目為n/8i-1.進(jìn)一步可得到化簡式:

從式(3)可得算法1的時間復(fù)雜度為O(n)，與均勻柵格劃分方法的時間復(fù)雜度是一個數(shù)量級別的［14］.

圖4 算法1的遞歸樹Fig.4 Recursive tree of algorithm 1

為了驗證算法1的時間復(fù)雜度和均勻柵格劃分的時間復(fù)雜度一致，在Intel i5 CPU，2 GB內(nèi)存的PC機(jī)上進(jìn)行實驗，取4個不同大小的Stanford大學(xué)網(wǎng)站上的模型Bunny，其圖元數(shù)量大約為4k、10k、16k、69k，分別測量柵格結(jié)構(gòu)的建立時間.測量結(jié)果如表1所示.從表1中可以看出，本文方法的運(yùn)行時間略慢于均勻柵格法，由前文分析可得，該方法需要對每個柵格判斷是否處于飽和狀態(tài)，然而均勻柵格法是不需要這些判斷的.本文方法和均勻柵格法的運(yùn)行時間增長趨勢是一致的，說明該方法在建立柵格結(jié)構(gòu)的運(yùn)行時間是合理的.

表1 不同模型下均勻柵格法和本文方法的建立時間Table 1 Building time of uniform grid and our method by different models

4 基于動態(tài)柵格法的光線遍歷

和均勻柵格光線遍歷方法一樣，每根光線只與它穿過的柵格內(nèi)部圖元做相交運(yùn)算［15］.這比沒有加速結(jié)構(gòu)的時候能節(jié)省許多時間，特別是當(dāng)場景中有大量物體的時候，光線只與其中一部分物體做相交運(yùn)算.本文的光線遍歷方法也繼承了以上的優(yōu)點(diǎn)，并對場景中存在較大空曠區(qū)域的情況進(jìn)行了優(yōu)化.

圖5 一條光線的動態(tài)柵格結(jié)構(gòu)遍歷Fig.5 Dynamic grids struct traversal of a ray

圖5所表示的就是一條光線穿越本文提出的動態(tài)柵格結(jié)構(gòu)的情景.斜杠標(biāo)記的柵格即光線穿越的柵格，光線將根據(jù)柵格與視點(diǎn)距離進(jìn)行求交運(yùn)算，順序依次為(c，0)- ＞ (d，3)- ＞ (g，0)- ＞(h，2)- ＞(f，0).括號中的第1項為柵格編號，第2項為該柵格中的圖元編號.第2章節(jié)中已經(jīng)提到，本文的柵格劃分方法就是基于柵格與視點(diǎn)的距離建立的.把這個順序?qū)?yīng)到圖2中去，這個操作順序是非常方便的，柵格樹同一深度的子樹，操作順序一定是從左往右，不同深度的子樹，操作順序一定是從上往下，這在一定程度上增強(qiáng)了算法性能.

本文動態(tài)柵格結(jié)構(gòu)光線遍歷方法關(guān)鍵在于如何沿著光線方向依次遍歷柵格，如圖5，即如何沿著(c，0)- ＞ (d，3)- ＞ (g，0)- ＞ (h，2)- ＞ (f，0)該路徑遍歷柵格.觀察到圖5的動態(tài)柵格結(jié)構(gòu)是把根節(jié)點(diǎn)均勻劃分為a、b、c、d 4個柵格，然后再把b節(jié)點(diǎn)劃分為e、f、g、h 4個柵格.因此可以對根節(jié)點(diǎn)采用均勻柵格的光線遍歷方法，當(dāng)檢測到節(jié)點(diǎn)b時，再一次采用均勻柵格的光線遍歷方法，即動態(tài)柵格結(jié)構(gòu)光線遍歷方法就是一種均勻柵格結(jié)構(gòu)光線遍歷方法的遞歸模式.

均勻柵格結(jié)構(gòu)中為了沿著光線方向計算遍歷步進(jìn)，采用以下比較巧妙的方法.即使光線和柵格邊的交點(diǎn)不是均勻的，如圖6(a)，但是在x方向和y方向上的交點(diǎn)是均勻的，如圖6(b)和圖6(c).這樣允許分別從x、y方向計算光線參數(shù)t，每個方向的光線參數(shù)t的步進(jìn)為

圖6 每條光線的動態(tài)柵格結(jié)構(gòu)遍歷Fig.6 Dynamic grids struct traversal of every ray

針對每個柵格，必須算出下一個光線經(jīng)過的柵格，在二維柵格的情況下，通過x方向步進(jìn)和通過y方向的步進(jìn)得到的下一個柵格是不一樣的，圖7說明了這2種不同情況，黑色的正方形就是應(yīng)該選擇的下一個柵格.在圖7中，假設(shè)光線剛進(jìn)入初始柵格時，在x方向最小的光線參數(shù)為tminx，在y方向最小的光線參數(shù)為tminy，橫斜杠標(biāo)記的柵格為當(dāng)前進(jìn)入的柵格，豎斜杠標(biāo)記的柵格為下一個進(jìn)入的柵格，計算txnext和tynext，它們分別為x方向和y方向的下一個光線經(jīng)過柵格的光線參數(shù).如果txnext＜tynext，就把x方向的下一個柵格作為光線下一個經(jīng)過的柵格，反之把y方向的下一個柵格作為光線下一個經(jīng)過的柵格.

圖7 決定下一個光線經(jīng)過的柵格過程Fig.7 Process of next grid a ray pass through

本文的動態(tài)柵格光線遍歷方法也是基于以上方法，它決定下一個光線經(jīng)過的柵格之后，必須還要判斷這柵格是否還有子?xùn)鸥?，如果有子?xùn)鸥竦脑掃€需要重復(fù)上述過程，反之則檢測這個柵格和圖元的求交情況.算法2給出了本文光線遍歷算法，為了簡化算法描述，假設(shè)光線的起點(diǎn)位于空間所有物體的外面.

算法2 動態(tài)柵格光線遍歷

輸入:當(dāng)前光線對象ray，光線參數(shù)tmin，當(dāng)前父柵格的包圍盒bbox，光線圖元求交之后的信息hitrecord，柵格樹 recursivecells.

1)初始化光線進(jìn)入空間柵格區(qū)域每個方向的柵格序號ix，iy，iz;每個方向的柵格光線參數(shù)步進(jìn)dtx=(txmax-txmin)/nx，dty=(tymax-tymin)/ny，dtz=(tzmax-tzmin)/nz;

2)計算出每個方向下一個柵格的光線參數(shù)，每次步進(jìn)的柵格個數(shù)，已經(jīng)遍歷的最后的柵格序號.

算法2遞歸次數(shù)由某個空間區(qū)域的圖元數(shù)量決定，數(shù)量少的區(qū)域很快完成遞歸，數(shù)量多的區(qū)域會進(jìn)行多次遞歸.

5 動態(tài)柵格法光線遍歷算法分析

與第3節(jié)分析動態(tài)柵格建立算法一樣，采取相同的假設(shè):整個空間中圖形元組成的包圍盒是一個正方體，這個包圍盒中有N(8的冪次數(shù))個三角形圖元;每個柵格的飽和狀態(tài)為容納一個三角形圖形元.

由算法2可知，光線遍歷算法消耗時間由2部分組成:確定初始化柵格具體序號的時間、光線與它經(jīng)過的柵格檢測相交的時間.每條光線只需進(jìn)行一次初始化柵格具體序號的確定，把這段時間定為Tinit，第2部分消耗時間則與空間內(nèi)部具體圖元分布有關(guān)，每個柵格求交時間為Tintersect.考慮最壞情況，即空間的圖元是分布均勻的，那么空間動態(tài)柵格建立的子?xùn)鸥褚彩蔷鶆蚍植嫉模饩€在這種情況下沿著對角線經(jīng)過子?xùn)鸥?，并且光線到達(dá)最后一個子?xùn)鸥癫艡z測到與圖元相交.對角線上的柵格數(shù)量也是與N有關(guān)，如圖4空間中物體包圍盒每個方向可以分成8log8N個柵格，那么物體包圍盒最長對角線的數(shù)量Nray可以由式(5)得到.那么在最壞情況下，每條光線的遍歷時間Ttraversal如式(6)所示.因此，本文算法計算復(fù)雜度為O((log8N)2/3).這個數(shù)量級和均勻網(wǎng)格遍歷方法是一樣的，但是在實際場景中，特別是針對空曠的場景性能是有明顯提升的.式(6)只是給出了本文算法的估計的一個運(yùn)行數(shù)量級時間，與實際往往有比較大的偏差，主要原因在于Tintersect這個變量在每個柵格中并不是都是一樣的，可能在某些柵格它的飽和圖元數(shù)量為c1，其他柵格的飽和圖元數(shù)量為c2.這會影響光線與當(dāng)前柵格的求交時間.

6 實驗及結(jié)果分析

在Intel i5 CPU，2GB內(nèi)存的PC機(jī)上進(jìn)行了實驗.實驗以如下方式進(jìn)行:選擇若干典型模型，它們都處于由三角形圖元組成的地面.分別對這些場景用均勻柵格法和本文提出的動態(tài)柵格法進(jìn)行光線追蹤，最后把這2種方法進(jìn)行比較.

實驗實用的模型包括多個不同大小的Bunny、Horse、Fish、Hand，為了體現(xiàn)空間的空曠性，在它們的下方放著由三角形圖元組成的地面，圖8為使用本文方法光線追蹤繪制的模型圖.實驗一共分為2組，第1組使用均勻柵格方法并記錄運(yùn)行時間和劃分的柵格數(shù)目.第2組使用本文提出的動態(tài)柵格光線追蹤方法并記錄相同的數(shù)據(jù)項目.表2給出了這些模型具體圖元數(shù)目，采用2種方法的柵格劃分?jǐn)?shù)目以及這2種方法的運(yùn)行時間.為了突出本文算法在空曠場景的作用，本次實驗的模型包括三角形圖元的地面，這樣更加凸顯出場景的空曠性.

圖8 用本文方法光線追蹤繪制的場景Fig.8 Scenes drawn by ray tracing through our method

實驗結(jié)果如表2所示，采用本文提出的動態(tài)柵格光線追蹤方法繪制場景確實提高了繪制速度.注意到Horse、Fish、Bunny4k的圖元個數(shù)雖然遠(yuǎn)遠(yuǎn)小于其余的模型，但是繪制時間還是比較長，這是由于放大了渲染模型.可以看到，本文方法所需要的柵格遠(yuǎn)遠(yuǎn)小于均勻柵格所需要的柵格，這從算法的初衷就可以猜測到:采用均勻柵格結(jié)構(gòu)光線追蹤算法沒有對模型進(jìn)行任何判斷，直接對場景按照式(1)進(jìn)行均勻柵格劃分，然后本文算法是根據(jù)場景的疏密程度動態(tài)地進(jìn)行柵格劃分，因此所需要的柵格會大大減少.表2中還可以看到不同圖元的Bunny隨著圖元數(shù)量的增加本文方法的運(yùn)行時間并沒有減少，這和柵格劃分?jǐn)?shù)有關(guān)，在圖元集中的區(qū)域柵格數(shù)量的增加是可以加速場景繪制速度的，從另一個角度說明本文算法隨著場景圖元疏密程度增加性能就越高.

表2 不同模型下均勻柵格法和本文方法繪制場景運(yùn)行時間Table 2 Scene drawing time of uniform grid and our method by different models

7 結(jié)束語

對于光線追蹤繪制場景的加速均勻柵格結(jié)構(gòu)，本文提出了一種新的動態(tài)柵格加速結(jié)構(gòu).與已有的方法相比，本文的方法對需要繪制的場景根據(jù)不同分布特征進(jìn)行按需柵格劃分.在建立動態(tài)柵格加速結(jié)構(gòu)的時候，不同深度從根節(jié)點(diǎn)到葉子節(jié)點(diǎn)，同一深度從左邊節(jié)點(diǎn)到右邊節(jié)點(diǎn)，都是隨著離視點(diǎn)的距離逐漸增大，這有利于光線的遍歷.實驗結(jié)果表明，本文方法在處理空曠場景的時候，優(yōu)于同類的已有工作，本文的方法還能推廣到基于GPU光線追蹤加速結(jié)構(gòu)上，這也是未來將擴(kuò)展的工作.

［1］ANDREW S.An introduction to ray tracing［M］.London:Academic Press，1989:201-204.

［2］BENTLEY J L.Multidimensional binary search trees used for associative searching［J］.Communications of the ACM，1975，18(9):509-517.

［3］RUBIN S M，WHITTED T.A 3-dimensional representation for fast rendering of complex scenes［C］//Proceedings of the 7th Annual Conference on Computer Graphics and Interactive Techniques.New York，USA，1980:509-517.

［4］KIRILL G，SIMON P，ALEXANDER B，et al.Grid-based SAH BVH construction on a GPU［J］.The Visual Computer，2011，27(6/7/8):697-706.

［5］TANAKA T，IWATA K.Arts:accelerated ray-tracing system［J］.IEEE Computer Graphics and Applications，1986，6(4):16-26.

［6］童星，袁道華.基于GPU和均勻柵格法的光線追蹤算法研究［J］.計算機(jī)工程與設(shè)計，2011，32(10):3499-3502.

TONG Xing，YUAN Daohua.Research of ray-tracing algorithm based on GPU and uniform grid method［J］.Computer Engineering and Design，2011，32(10):3499-3502.

［7］GUNTURY S，NARAYANAN P J.Raytracing dynamic scenes on the GPU using grids［J］.IEEE Transactions on Visualization and Computer Graphics，2012，18(1):5-16.

［8］WALD I，IZE T，ANDREW K，et al.Ray tracing animated scenes using coherent grid traversal［C］//Proceedings of ACM SIGGRAPH.New York，USA，2006，485-493.

［9］SHIRLEY P，ASIKHMIN M，GLEICHER M，et al.Fundamentals of computer graphics［M］.2nd ed.Wellesley:AK Peters，2005:224-226.

［10］KALOJANOV J，BILLETER M，SLUSALLEK P.Two-level grids for ray tracing on GPUs［J］.Computer Graphics Forum，2011，30(2):307-314.

［11］KALOJANOV J，SLUSALLEK P.A parallel algorithm for construction of uniform grids［C］//Proceedings of the Conference on High Performance Graphics.New York，USA，2009:23-28.

［12］GARANZHA K，LOOP K.Fast ray sorting and breadthfirst packet traversal for GPU ray tracing［J］.Computer Graphics Forum，2010，29(2):289-298.

［13］PANTALEONI J，LUEBKE D.HLBVH:hierarchical LBVH construction for real-time ray tracing of dynamic geometry［C］//Proceedings of the Conference on High Performance Graphics.Switzerland，2010:87-95.

［14］VACLAV S.An efficient space partitioning method using binary maps［C］//Proceedings of the 11th International Conference on Telecommunications and Informatics.Wisconsin，USA，2012:121-124.

［15］PERROTTE L，SAUPIN G.Fast GPU perspective grid construction and triangle tracing for exhaustive ray tracing of highly coherent rays［J］.International Journal of High Performance Computing Applications，2012，26(3):192-202.