基于光流場動(dòng)態(tài)參與介質(zhì)繪制算法

王元龍

摘要：為了解決動(dòng)態(tài)參與介質(zhì)場景連續(xù)幀的實(shí)時(shí)繪制問題，提出了一種基于光流場動(dòng)態(tài)參與介質(zhì)場景的繪制算法。首先，采用區(qū)域匹配的方法計(jì)算關(guān)鍵幀之間的光流場；然后，通過插值的方法計(jì)算中間幀之間的光流場，采用幀間光流連貫性函數(shù)表示幀與幀之間光流一致性的度量，保證幀與幀之間介質(zhì)運(yùn)動(dòng)不會(huì)發(fā)生突變；最后，按照所產(chǎn)生的光流場，繪制連續(xù)動(dòng)態(tài)變化的參與介質(zhì)場景。在連續(xù)5幀動(dòng)態(tài)參與介質(zhì)場景繪制中，所提算法比基于徑向基函數(shù)模型（RBF）參與介質(zhì)的光子映射算法效率提高了近3倍，能夠達(dá)到連續(xù)幀的實(shí)時(shí)繪制，且繪制質(zhì)量比較高。

關(guān)鍵詞：多散射；參與介質(zhì)；光子映射；并行計(jì)算；光流場

中圖分類號(hào)：TP391 文獻(xiàn)標(biāo)志碼：A

Abstract： In order to achieve the realtime rendering of continuous frame for participating media scene， a rendering algorithm based on optical flow was proposed. First， the regional matching method was used to calculate the field of optical flow between key frames. Then the field of optical flow between intermediate frames was calculated by the interpolation method， and the optical coherence function between frames was used to denote the consistency of optical flow to guarantee that the media motion wont be suddenly changed. Finally， the dynamic scene of continuous frames was rendered according to the field of optical flow. In the participating media scene rendering for 5 continuous frames， the efficiency of the proposed algorithm increase nearly five times than that of based on Radial Basis Function （RBF） model， realtime rendering of consecutive frames is implemented and rendering quality is relatively high.

Key words：multiple scattering； participating media； photon mapping； parallel computing； optical flow field

0 引言

繪制介質(zhì)參與的場景時(shí)，通常需要考慮光照在物體表面或介質(zhì)中的散射現(xiàn)象。因此，在真實(shí)感圖像渲染領(lǐng)域，散射效果的模擬（像煙霧、云和大部分非金屬表面等等的繪制）成為越來越被研究人員關(guān)注的課題?；诠饩€跟蹤的蒙特卡羅方法能夠渲染高質(zhì)量的介質(zhì)參與場景，模擬非均勻介質(zhì)的不同屬性[1-3]。但是，由于在介質(zhì)中散射到處發(fā)生，需要跟蹤多次散射路徑直到光能被吸收或路徑離開參與介質(zhì)，因此該方法渲染效率低，渲染一幅簡單場景圖像需要幾分鐘。

研究者提出了光子映射的方法[4-7]來提高介質(zhì)參與場景的繪制效率。光子映射方法最初被應(yīng)用于全局光照下沒有參與介質(zhì)的場景渲染中[4]。Jensen等[5]把全局光照的光子映射方法擴(kuò)展到有參與介質(zhì)的場景中，與文獻(xiàn)[4]中的方法一樣，用一個(gè)球形鄰域來搜索影響該點(diǎn)的光子。該方法實(shí)現(xiàn)了各種全局光照效果，但要達(dá)到好的渲染結(jié)果需要十幾萬的光子，不但計(jì)算量比較大，而且容易出現(xiàn)搜索球疊加或存在搜索不到區(qū)域的情況。Jarosz等[6]對文獻(xiàn)[5]中的輻射度估計(jì)作了改進(jìn)，提出一種用光子束（Photon Beam， PB）來估計(jì)內(nèi)散射輻射度的方法，他們根據(jù)光子輻射度估計(jì)方法重新推出了在介質(zhì)中的光線傳輸?shù)仁剑w估計(jì)每一條從視點(diǎn)發(fā)出光線的輻射度，用可變核（Variable Kernel， VK）的方法來搜索最近的光子和用層次包圍盒（Bounding Box Hierarchy， BBH）構(gòu)建光子貼圖。Jarosz等[7]把PB的概念不僅應(yīng)用到估計(jì)輻射度的查詢操作上，還應(yīng)用到光子的數(shù)據(jù)表示中，提高了渲染圖像的質(zhì)量。但是，PB方法需要考慮光束中的所有光子，所以計(jì)算量仍然很大，且對光束寬度采用確定的值，造成了大的偏差，雖然可以通過減小寬度來減小偏差，但增加了噪聲。最近，Zhang 等[8]提出一種預(yù)處理的輻射度傳輸方法來實(shí)現(xiàn)介質(zhì)參與場景全局光照效果的實(shí)時(shí)渲染。該方法能夠快速地估計(jì)視點(diǎn)光線的輻射度，但該方法需要的預(yù)計(jì)算時(shí)間較長（73～395s）。文獻(xiàn)[9]也提出了基于徑向基函數(shù)（Radial Basis Function， RBF）模型參與介質(zhì)的光子映射算法。該算法用存儲(chǔ)RBF模型來代替存儲(chǔ)光子從而減少存儲(chǔ)空間，實(shí)現(xiàn)了實(shí)時(shí)渲染。但在整個(gè)光子映射過程中需要的預(yù)計(jì)算時(shí)間較長（0.81～5.61s）。目前，光子映射算法存在預(yù)計(jì)算時(shí)間較長的問題，難以實(shí)現(xiàn)動(dòng)態(tài)介質(zhì)參與場景連續(xù)幀的實(shí)時(shí)繪制。

本文在文獻(xiàn)[9]的基礎(chǔ)上，將光流場概念擴(kuò)展并應(yīng)用到動(dòng)態(tài)介質(zhì)參與場景連續(xù)幀的繪制中，提出一種基于光流場動(dòng)態(tài)介質(zhì)參與場景的繪制算法。首先，采用區(qū)域匹配的方法計(jì)算關(guān)鍵幀之間的光流場，通過最小化殘差找到幀與幀之間特征對的光流矢量；然后，通過插值的方法來計(jì)算中間幀之間的光流場，采用幀間光流連貫性函數(shù)來表示幀與幀之間的光流一致性的度量，保證幀與幀之間介質(zhì)運(yùn)動(dòng)不會(huì)發(fā)生突變；最后，按照所產(chǎn)生的光流場，繪制連續(xù)動(dòng)態(tài)變化的介質(zhì)參與場景。另外，在CPU+GPU的框架下實(shí)現(xiàn)了本文算法，進(jìn)一步加速動(dòng)態(tài)參與介質(zhì)的繪制。

1 相關(guān)工作

文獻(xiàn)[9]中提出了基于RBF模型參與介質(zhì)的光子映射算法。在構(gòu)建RBF模型之前，首先需要計(jì)算基函數(shù)的中心，由于聚類個(gè)數(shù)不確定，因此采用自組織映射無監(jiān)督聚類方法[10]來計(jì)算基函數(shù)的中心。在輻射度估計(jì)階段，使用構(gòu)建好的RBF模型來計(jì)算穿過參與介質(zhì)的每一條視點(diǎn)光線輻射度。用式（1）計(jì)算視點(diǎn)光線的內(nèi)散射輻射度：

2 基于光流場的動(dòng)態(tài)介質(zhì)參與場景繪制算法

大部分光流算法的研究是針對動(dòng)態(tài)場景，從圖像序列中得到運(yùn)動(dòng)目標(biāo)的運(yùn)動(dòng)矢量[11-12]。光流場是將三維動(dòng)態(tài)場景投影到二維圖像平面上，利用圖像序列的灰度變化、空域變化和相關(guān)性等信息計(jì)算得到的運(yùn)動(dòng)矢量場。本文把光流場的概念應(yīng)用到繪制介質(zhì)參與場景的光子映射算法中，來模擬介質(zhì)的動(dòng)態(tài)變化效果。首先介紹算法流程，然后分別闡述算法中采用的關(guān)鍵技術(shù)。在介紹算法之前，定義兩個(gè)關(guān)于幀的定義。

關(guān)鍵幀 在動(dòng)態(tài)顯示場景時(shí)，場景中對象變化較大的幀。

中間幀 在關(guān)鍵幀之間，對象變化比較平滑的幀。

2.1 算法流程

本文算法經(jīng)過并行設(shè)計(jì)，使用CPU+GPU的框架完成，充分利用了GPU并行計(jì)算能力，其主要步驟如下：

步驟1 分別計(jì)算關(guān)鍵幀場景的RBF模型。本文采用文獻(xiàn)[17]的方法計(jì)算參與介質(zhì)的RBF模型。

步驟2 計(jì)算關(guān)鍵幀之間光流場。將關(guān)鍵幀的RBF模型作為基本元素，通過區(qū)域匹配的方法計(jì)算關(guān)鍵幀之間光流場。本文把RBF的中心和大小作為特征進(jìn)行匹配，找到幀與幀之間特征對，根據(jù)特征對計(jì)算得到光流場。

步驟3 采用插值的方法計(jì)算中間幀之間的光流場。本文采用幀間光流連貫性函數(shù)來表示幀與幀之間的光流一致性度量。

步驟4 渲染最終的動(dòng)態(tài)介質(zhì)參與場景連續(xù)幀。

注 本文在GPU顯存中存儲(chǔ)第一幀場景的RBF模型及所有幀間光流場，第一幀以外其他幀場景的RBF模型通過幀間光流場計(jì)算得到。

2.2 關(guān)鍵幀之間光流場的計(jì)算

按照數(shù)學(xué)方法劃分，計(jì)算光流場的方法分為4種：微分法、區(qū)域匹配法、基于能量的方法和基于相位的方法，由于區(qū)域匹配法計(jì)算簡單、實(shí)時(shí)性強(qiáng)等特點(diǎn)，得到了較為廣泛的應(yīng)用[13]。本文采用區(qū)域匹配法計(jì)算關(guān)鍵幀間的光流場。區(qū)域匹配法認(rèn)為上一幀中某個(gè)區(qū)域，在下一幀中會(huì)形成對應(yīng)的模式，光流場定義為兩幀中相應(yīng)區(qū)域塊之間產(chǎn)生最佳擬合的位移矢量，實(shí)質(zhì)上是在前幀與后幀圖像之間實(shí)施的一種對應(yīng)任務(wù)，這類方法是將光流定義為不同時(shí)刻圖像區(qū)域的位移量d=（dx， dy），使得不同時(shí)刻的圖像區(qū)域匹配最佳。匹配圖像區(qū)域常常等價(jià)為最小化殘差的過程，給定兩幀圖像I和J，需要找到位移量d= （dx， dy）使殘差最?。?/p>

如果用聚類算法將像素點(diǎn)匹配轉(zhuǎn)化為聚類特征點(diǎn)會(huì)大大減少計(jì)算量[14]。本文在光子映射階段采用自組織映射無監(jiān)督聚類方法把光子的分布用RBF來模擬（詳見文獻(xiàn)[9]），因此將RBF看作區(qū)域塊來搜索相應(yīng)的位移矢量場。在連續(xù)關(guān)鍵幀I和J匹配過程中，第I、J幀分別用RBF的半徑作為匹配值，記為RI（i， j，k）、RJ（i， j，k），其中（i， j，k）為RBF的中心位置。以此中心建立一個(gè)搜索窗口W（wi，wj，wk），搜索窗口不能太小也不能太大，基于統(tǒng)計(jì)的實(shí)驗(yàn)結(jié)果，本文采用大小為6×6×6的搜索窗口。與計(jì)算圖像光流場不同（光流場通常表示物體點(diǎn)的三維運(yùn)動(dòng)在二維圖像上的表現(xiàn)，而本文計(jì)算三維空間中介質(zhì)運(yùn)動(dòng)的光流場），本文在J幀中（i， j，k）位置附近搜索與I幀（i， j，k）對應(yīng)的RBF，三維空間中計(jì)算光流場相應(yīng)的殘差定義為ε（d）=（RI（i， j，k）-RJ（i+di， j+dj，k+dk））2，通過最小化殘差的方法計(jì)算得到相應(yīng)的光流矢量，d（di，dj，dk）就是需要計(jì)算得到的光流矢量。計(jì)算光流場如算法1所示。

算法1 從關(guān)鍵幀動(dòng)態(tài)介質(zhì)中進(jìn)行光流場的計(jì)算。

2.4 GPU并行實(shí)現(xiàn)

本文算法使用計(jì)算統(tǒng)一設(shè)備架構(gòu)（Compute Unified Device Architecture， CUDA）來完成，算法流程中的步驟都經(jīng)過了并行設(shè)計(jì)，可以并行地在GPU上運(yùn)行。設(shè)計(jì)CUDA程序首先需要盡可能地發(fā)掘GPU體系結(jié)構(gòu)的特點(diǎn)：減少分支操作、減少讀寫操作、充分利用線程資源等。

在計(jì)算關(guān)鍵幀之間光流場階段，本文先將表示光子分布的RBF數(shù)據(jù)綁定到GPU上的紋理存儲(chǔ)器。因?yàn)镚PU上的紋理硬件有一個(gè)很有用的特性，它提供了一個(gè)紋理緩存，因此紋理存儲(chǔ)器的平均訪問速度比全局存儲(chǔ)器更快。存儲(chǔ)RBF采用層次包圍盒（BBH）結(jié)構(gòu)，加速光線與RBF求交的過程（詳見文獻(xiàn)[9]）。對于光流場計(jì)算，式（4）中要求在搜索窗口W（wi，wj，wk）范圍查找匹配的RBF，本文設(shè)計(jì)了在BBH結(jié)構(gòu)中搜索對應(yīng)RBF的方法，如圖3所示為在BBH中搜索具體過程的示例。首先在后一幀中找到前一幀RBF對應(yīng)的位置，該位置可能存在RBF，也可能不存在，因此在圖3中用虛線框表示，如果存在RBF，用式（4）計(jì)算前一幀與后一幀中該RBF的殘差；然后在搜索窗口 W（wi，wj，wk）范圍中搜索使殘差最小的RBF，圖3中w-1、w+1等表示向前和向后搜索，其中“-”表示向前搜索，“+”表示向后搜索。搜索過程為：1）如果位置（i， j，k）為其父節(jié)點(diǎn)的右（左）節(jié)點(diǎn)，那么向前（后）搜索其左（右）兄弟節(jié)點(diǎn)；2）如果沒有達(dá)到搜索范圍繼續(xù)搜索其父節(jié)點(diǎn)的兄弟節(jié)點(diǎn)，直到達(dá)到搜索范圍停止搜索。

在輻射度估計(jì)階段（即繪制階段），估計(jì)每條視點(diǎn)光線的輻射度，需要求得與該視點(diǎn)光線相交的所有RBF，這個(gè)過程稱為光線遍歷。由于每條視點(diǎn)光線是彼此獨(dú)立的，因此遍歷方法實(shí)現(xiàn)的一個(gè)重要特點(diǎn)是高度并行性（在GPU上并行的跟蹤每一條光線，設(shè)置的視平面大小為800×600，而光線的設(shè)置為從視點(diǎn)分別到視平面的每個(gè)像素的方向，所以跟蹤的光線也為800×600條，從而充分利用線程資源）。每一幀的輻射度估計(jì)，先從視點(diǎn)并行地向成像平面上每個(gè)像素發(fā)射視點(diǎn)光線，然后用當(dāng)前幀的RBF模型來估計(jì)輻射度，在最后顯示時(shí)，把輻射度值轉(zhuǎn)化為顏色值，再將其綁定到紋理對象上，通過OpenGL緩存技術(shù)將其最終的圖像顯示出來。但每幀得到RBF模型的方式不同，在繪制第一幀的過程中直接使用光子跟蹤階段得到的RBF模型來進(jìn)行輻射度估計(jì)，而在后續(xù)幀的繪制過程中需要先通過光流場計(jì)算得到當(dāng)前幀的RBF模型。

3 實(shí)驗(yàn)結(jié)果及討論

本文實(shí)驗(yàn)使用：PIV 2.83GHz 四核CPU，4GB主存，NVIDIA GeForce GTX470顯卡的PC機(jī)，編譯環(huán)境為VS2008和CUDA2.3，所有實(shí)驗(yàn)結(jié)果圖像的分辨率均為800×600。介質(zhì)數(shù)據(jù)與文獻(xiàn)[9]相同，介質(zhì)中每個(gè)位置上的密度值隨介質(zhì)位置的不同而不同，能很好地模擬不均勻介質(zhì)的散射。光流場在GPU紋理存儲(chǔ)器中采用一個(gè)四維數(shù)據(jù)類型float4表示（在GPU中， float4的讀取操作比4個(gè)float分開讀取獲得更高的效率），其中前三維用來表示光流場的矢量的方向，第四維用來表示RBF半徑的變化率（即矢量的大?。?jì)算光流場的平均時(shí)間大約為1s。

為了驗(yàn)證搜索窗口大小對光流場的影響，本文計(jì)算不同大小搜索窗口得到光流場，以及通過光流場繪制了后一幀介質(zhì)的散射效果，如圖4所示。圖4（a）、4（b）分別為文獻(xiàn)[9]中基于RBF光子映射算法繪制的前后兩幀介質(zhì)參與場景，圖4（c）、4（d）、4（e）分別是不同大小搜索窗口計(jì)算得到的光流場以及采用該光流場繪制后一幀的介質(zhì)散射結(jié)果。從圖4中可以看出，對于基于光流場算法繪制后續(xù)幀的介質(zhì)散射效果算法，采用大小為6×6×6的搜索窗口能夠達(dá)到與采用預(yù)處理算法[9]效果類似的繪制結(jié)果，但效率提高了大約20倍。表1為基于RBF 光子映射算法和基于光流場算法繪制后一幀介質(zhì)參與場景的運(yùn)行統(tǒng)計(jì)數(shù)據(jù)， 其中，基于RBF 光子映射算法需要在每一幀進(jìn)行預(yù)處理操作，而基于光流場算法在繪制過程中使用光流場計(jì)算后一幀相于對前一幀的RBF， 因此可以提高連續(xù)幀的繪制效率。正如圖4和表1所示，搜索窗口越大繪制效果越好，同時(shí)繪制時(shí)間越長。采用大小為6×6×6的搜索窗口就能夠達(dá)到與文獻(xiàn)[9]類似的較好效果，且能夠?qū)崿F(xiàn)連續(xù)動(dòng)畫的效果，因此本文采用大小為6×6×6的搜索窗口實(shí)現(xiàn)連續(xù)幀的實(shí)時(shí)繪制。

圖5為在關(guān)鍵幀之間插值生成中間幀光流場的示意圖，其中，第1幀和第5幀視為關(guān)鍵幀，第2、3、4幀為插值計(jì)算得到的中間幀，其中xy表示從第x幀到第y幀之間的光流場。圖6所示為通過圖5光流場產(chǎn)生的連續(xù)5幀動(dòng)態(tài)介質(zhì)參與場景的繪制結(jié)果，第1幀采用文獻(xiàn)[9]中基于RBF光子映射算法繪制，后面的連續(xù)幀采用基于光流場算法繪制。從圖中可以看出，基于光流場算法繪制的連續(xù)多幀介質(zhì)散射效果質(zhì)量比較高，而且可以實(shí)現(xiàn)幀間實(shí)時(shí)的繪制。其中基于BRF光子映射算法[9]，在每一幀都需要進(jìn)行預(yù)處理操作，每一幀平均時(shí)間為5.31s，5幀共26.55s。而基于光流場算法僅僅需要在第1幀進(jìn)行預(yù)處理操作（5.31s），后續(xù)4幀通過光流場計(jì)算各幀RBF模型（0.26×4=1.04s），總計(jì)6.35s。因此在連續(xù)5幀動(dòng)態(tài)介質(zhì)場景繪制中比基于RBF光子映射算法效率提高了近3倍。

4 結(jié)語

本文針對動(dòng)態(tài)介質(zhì)場景連續(xù)幀的繪制，提出一種基于光流場動(dòng)態(tài)參與介質(zhì)繪制算法，其借用了運(yùn)動(dòng)分析中光流場的方法進(jìn)行幀與幀之間相關(guān)性的分析，并且采用CPU+GUP的并行框架，實(shí)現(xiàn)了連續(xù)幀的實(shí)時(shí)繪制。從實(shí)驗(yàn)結(jié)果可以看出，用本文算法模擬的效果與利用基于RBF光子映射算法繪制的效果相當(dāng)，但本文算法提高了繪制連續(xù)幀的效率。

本文在插值計(jì)算中間幀光流場的過程中，需要嘗試不同的迭代值找出最佳的迭代次數(shù)，不能應(yīng)用于所有介質(zhì)場景的自動(dòng)繪制。將來要重點(diǎn)解決的問題是考慮設(shè)置自動(dòng)判別迭代次數(shù)的機(jī)制。

參考文獻(xiàn)：

[1]KAJIYA J， HERZEN B. Ray tracing volume densities[C]// Proceedings of the 11th Annual Conference on Computer Graphics and Interactive Techniques. New York： ACM， 1984： 165-174.

[2]MAX N. Efficient light propagation for multiple anisotropic volume scattering[C]// Proceedings of Eurographics Workshop on Rendering. Heidelberg： Springer， 1994： 87-104.

[3]JENSEN H， MARSCHNER S， LEVOY M， et al. A practical model for subsurface light transport[C]// Proceedings of the 28th Annual Conference on Computer Graphics and Interactive Techniques. New York： ACM， 2001： 511-518.

[4]JENSEN H W. Global illumination using photon maps[C]// Proceedings of the 7th Eurographics Workshop on Rendering Techniques. Berlin： Springer， 1996：21-30.

[5]JENSEN H W， CHRISTENSEN P H. Efficient simulation of light transport in scenes with participating media using photon maps[C]// Proceedings of the 25th Annual Conference on Computer Graphics and Interactive Techniques. New York： ACM， 1998：311-320.

[6]JAROSZ W， ZWICKER M， JENSEN H W. The beam radiance estimate for volumetric photon mapping[J]. Computer Graphics Forum， 2008， 27（2）：557-566.

[7]JAROSZ W， NOWROUZEZAHRAI D， SADEGHI I， et al. A comprehensive theory of volumetric radiance estimation using photon points and beams[J]. ACM Transactions on Graphics， 2011， 30（1）：1-19.

[8]ZHANG Y B， ZHAO D， MA LK. Realtime volume rendering in dynamic lighting environments using precomputed photon mapping[J]. IEEE Transactions on Visualization and Computer Graphics， 2013， 19（8）：1317-1330.

[9]WANG Y L， GUO P. Rendering scenes with participating media based on RBFs for photon mapping using graphics hardware[C]// Proceedings of the 12th ACM SIGGRAPH International Conference on VirtualReality Continuum and Its Applications in Industry. New York： ACM， 2013： 167-176.

[10]VESANTO J， ALHONIEMI E. Clustering of the selforganizing map[J]. IEEE Transactions on Neural Networks， 2000， 11（3）： 586-600.

[11]IKETANI A， NAGAI A， KUNO Y， et al. Realtime surveillance system detecting persons in complex scenes[J]. RealTime Imaging， 2001， 7（5）： 433-446.

[12]潘金山，蘇志勛，王偉.運(yùn)動(dòng)細(xì)節(jié)估計(jì)的光流場方法[J].計(jì)算機(jī)輔助設(shè)計(jì)與圖形學(xué)學(xué)報(bào)，2011， 23（8）：1433-1441.（PAN J S， SU Z X， WANG W. Optical flow method for motion details estimation[J]. Journal of ComputerAided Design & Computer Graphics，2011，23（8）：1433-1441.）

[13]艾海舟，蘇延超. 圖像處理、分析與機(jī)器視覺[M]. 北京：清華大學(xué)出版社， 2003： 542-548.（AI H Z， SU Y C. Image Processing Analysis and Machine Vision[M]. Beijing： Tsinghua University Press， 2003： 542-548.）

[14]林強(qiáng)，田桂蘭. 聚類匹配法光流檢測中共同運(yùn)動(dòng)匹配策略的研究[J].電子學(xué)報(bào)，1997， 25（1）：58-61.（LIN Q， TIAN G L. An investigation on the matching tactics of common motion in detecting flow by cluster algorithm[J].Acta Electronica Sinica，1997， 25（1）：58-61.）