徐立東　包雪怡

摘? 要： 采用多分辨率FFT（快速傅里葉變換）波浪方程來(lái)計(jì)算法線貼圖，在低分辨率網(wǎng)格上構(gòu)建海面高度場(chǎng)，然后使用Cuda架構(gòu)對(duì)波浪方程并行求導(dǎo)算出更高分辨率的法線貼圖，用Open GL著色器來(lái)渲染海面波動(dòng)效果，可有效模擬海面。隨著法線貼圖分辨率的提高，能更加細(xì)膩、真實(shí)地展現(xiàn)海面波浪的繪制效果。實(shí)驗(yàn)結(jié)果表明，采用基于FFT法線貼圖的海面波浪模擬可視化框架，能夠不增加海面網(wǎng)格的幾何復(fù)雜度、節(jié)省計(jì)算時(shí)間，在一定程度上提升海面細(xì)節(jié)、增強(qiáng)視覺(jué)效果。

關(guān)鍵詞： 快速傅里葉變換; 法線貼圖; 海面; 波浪; 模擬

中圖分類號(hào)：TP399? ? ? ? ? 文獻(xiàn)標(biāo)志碼：A? ? ?文章編號(hào)：1006-8228（2019）12-23-05

Research on sea wave simulation method based on FFT normal map

Xu Lidong， Bao Xueyi

（Ocean University of China，College of Information Science and Engineering，Qingdao， Shandong 266100， China）

Abstract： The multi-resolution FFT （Fast Fourier Transform） wave equation is used to calculate the normal map， the sea surface height field is constructed on the low resolution grid， and then the Cuda architecture is used to parallelize the wave equation to calculate the higher resolution normal map. Textures， using Open GL shaders to render sea wave effects， can effectively simulate the sea surface. As the resolution of the normal map is improved， the drawing effect of the sea wave can be more delicately and realistically displayed. The experimental results show that the FFT normal map based sea surface visualization framework can increase no more the geometric complexity of the sea surface mesh， save the computation time， and improve the sea surface detail and enhance the visual effect to some extent.

Key words： Fast Fourier Transform; normal map; sea surface; wave; simulation

0 引言

近年來(lái)，“智慧海洋”成為海洋數(shù)據(jù)信息化發(fā)展前沿的研究熱點(diǎn)之一。海面可視化是海洋環(huán)境虛擬現(xiàn)實(shí)可視化中的一個(gè)重要的課題，由于海洋環(huán)境復(fù)雜多變性，海面會(huì)受到風(fēng)力、潮汐、航行中船舶等諸多因素影響，產(chǎn)生不同幅度的海面波動(dòng)，為了實(shí)現(xiàn)海面波浪可視化更加逼真的渲染繪制效果，提高程序的并行性和計(jì)算效率，以及提高渲染繪制的實(shí)時(shí)性成為了國(guó)內(nèi)外學(xué)者研究的難點(diǎn)。

在海面建模和渲染方面，早期的方法是采用一張預(yù)先計(jì)算好的法線貼圖映射到海面上，并通過(guò)紋理的擾動(dòng)來(lái)得到近似波動(dòng)的效果，此法是基于二維紋理波技術(shù)，不屬于基于物理的海面建模。Gerstner wave最初在物理中用于水波的模擬，其形狀比較真實(shí)，而且計(jì)算量不大，被廣泛用于計(jì)算機(jī)圖形學(xué)中水波的模擬，被作為模擬海面的物理模型。1986年，F(xiàn)ournier和Reeves[1]首次使用此模型模擬海岸線海面波動(dòng)，生成較為真實(shí)的海場(chǎng)景。2001年，Jerry Tessendorf[2] 提出快速傅里葉變換（FFT）能夠在適中規(guī)模的網(wǎng)格上生成逼真的海面并可達(dá)到實(shí)時(shí)處理的幀率。童若鋒等[3]對(duì)頂點(diǎn)波、漣漪、紊亂的短峰波等進(jìn)行分類造型，模擬水波的大致形狀。Yann[4]給出了基于三維體元的納維爾斯托克斯方程（N-S方程）求解方法來(lái)模擬海面。徐迎慶[5]等從描述不穩(wěn)定的Saint-Venant方程組來(lái)模擬海浪。2018年，Jeschke S，Sk?ivan T，Müller-Fischer M等[6]使用小波變換對(duì)流體變化進(jìn)行離散化處理，此法拓展傳統(tǒng)的FFT波浪模型，能準(zhǔn)確地模擬流體和物體之間的相互作用。本文將FFT波浪方程計(jì)算法線貼圖和納維爾斯托克斯方程（N-S方程）相結(jié)合，構(gòu)建了一個(gè)海面波浪模擬可視化框架。

1 海面波浪模擬

真實(shí)的海面模擬效果不僅需要?jiǎng)澐趾Ｃ婢W(wǎng)格，也需要逼真的海面波浪模型及真實(shí)感繪制，本文所用的渲染技術(shù)是可編程著色技術(shù)（shader），對(duì)N-S流體動(dòng)力學(xué)方程的求解和FFT變換均在CUDA上進(jìn)行，渲染工作均由Open GL著色器編寫。

1.1 輕微波浪模型

對(duì)于海面波浪模擬的簡(jiǎn)單方法是將其看作一組正弦曲線的疊加：

[H（x，z，t）=（Ai×sin（Di.x，z）×ωi+t×?i）]? ⑴

其中，Ai是每個(gè)波浪的振幅，Di是一個(gè)二維向量指向波峰傳播方向，[ωi]= 2π/L（這里，L為每個(gè)波浪的波長(zhǎng)）。

此模型模擬的波浪浪尖圓滑，適合于平靜海面。為了模擬波峰尖銳的波浪，可在波的傳播方向增加位移，得到Gerstner wave波動(dòng)函數(shù)，Gerstner wave表示為：

[X=X+（Qi Ai×X×cos（ωi（x?，z）+ψi t））Y=（Ai×sin（ωi Di?（x，z）+ψi t））]? ?⑵

在三角函數(shù)中添加一個(gè)π/2的相位偏移，Gerstner wave可表示成如下形式：

[X=X-（Qi Ai×X×cos（ωi?（x，z）+ψi t））Y=（Ai×sin（ωi Di?（x，z）+ψi t））]? ?⑶

可得，Gerstner wave的每個(gè)波浪的振幅是常數(shù) ，Gerstner wave是隨時(shí)間t變化的三角函數(shù)波形模擬波浪運(yùn)動(dòng)。海洋統(tǒng)計(jì)學(xué)中，波浪的振幅?？醋魇且粋€(gè)隨時(shí)間變化的高度函數(shù)h（x，t），X是一個(gè)二維坐標(biāo)（x，z）。將式⑶中與時(shí)間t的部分去掉，使用h（x，t）替換常數(shù) 可得：

[Y=h（x，t）×cos（ωi Di?（x，z））]? ? ? ? ⑷

在式⑷中增加一個(gè)虛部函數(shù)，設(shè)波浪波長(zhǎng)為L(zhǎng)，可得：

[ ⑸]

根據(jù)歐拉公式[eiα=cosα+isinα]，⑸可以改寫為：

[Y=h（x，t）eiωiDiX]? ? ? ? ? ? ? ? ⑹

由此得到傅里葉逆變換的海面模型的初步形式，使用快速傅里葉逆變換可快速對(duì)大量cos和sin波形求和，最終基于FFT的波浪模型為：

[h（X，t）=h（K，t）eiKX]? ? ? ? ? ? ⑺

[h（X，t）]是一個(gè)復(fù)數(shù)，實(shí)部表示波浪高度，虛部表示相位。X是水平面二維坐標(biāo)，該公式結(jié)果為X=（x，z）坐標(biāo)下的高度。該表達(dá)式在[s×s]的范圍內(nèi)的二維網(wǎng)格生成一個(gè)周期的FFT波浪，s尺度范圍10米～2000米，參與傅里葉變換的網(wǎng)格點(diǎn)數(shù)可達(dá)2048×2048，在一定視野范圍使用較大的s生成的波浪沒(méi)有較明顯的周期性。為創(chuàng)建初始時(shí)刻海面高度場(chǎng)的傅里葉頻譜，Tessendorf[6]使用了一個(gè)如下的菲利普頻譜函數(shù)：

[Ph（k）=? Aexp（-1/（kl）2）k4|K?w|2]? ? ⑻

這里，[l= v2g]是在風(fēng)速v下的最大波浪，k是波浪向量的模，向量[K]是波浪的傳播方向，向量[w]是波浪的傳播風(fēng)向。時(shí)刻的高度由式⑼得出：

[h0（k）=12（εr+iεi）Ph（K）]? ? ? ? ⑼

[εr]和[εi]是均值為0方差為1的高斯隨機(jī)數(shù)。根據(jù)t時(shí)刻的頻譜，本文可得t時(shí)刻的頻譜函數(shù)：

[h（-k，t）=h0（k）eiω（k）t+h0（-k）e-iω（k）t]? ? ? ⑽

從該式可得向相反方向傳播的波浪的高度頻譜同高，海面高度場(chǎng)有對(duì)稱性。

1.2 洶涌波浪模型

由式⑹FFT波浪模型可模擬出平滑的浪尖和波紋，此模型適用于較為平靜的海面模擬。但是在現(xiàn)實(shí)中由于風(fēng)力等強(qiáng)外力影響，海面波浪在其頂端會(huì)形成浪尖，波浪的輪廓（波長(zhǎng)和波高）受環(huán)境影響。足夠高度的波浪會(huì)在頂端發(fā)生碰撞，生成復(fù)雜的物理現(xiàn)象，比如泡沫和水花。對(duì)于波浪的模擬需研究海面運(yùn)動(dòng)的流體動(dòng)力學(xué)方程。此類方法在海面高程和動(dòng)力勢(shì)能上描述水體的運(yùn)動(dòng)規(guī)律，都源自納維爾斯托克斯方程（N-S方程）。Creamer[7]使用“米氏變換”的數(shù)學(xué)方法產(chǎn)生一系列海面高程和速度勢(shì)能的規(guī)范化變換。使用這種變換方法可將海面高程和速度勢(shì)能轉(zhuǎn)換成簡(jiǎn)單的表示形式。米氏變換可在二維空間中得出精確解，在三維空間中難以應(yīng)用。根據(jù)二維空間中精確解，利用插值形式可求出三維空間中的解：一個(gè)波浪的水平位移項(xiàng)，該位移隨著波浪和時(shí)間變化。

本文在FFT的海浪模型的基礎(chǔ)上使用了兩種波浪模型分別用來(lái)模擬風(fēng)速較小時(shí)和風(fēng)速較大時(shí)的波浪。

原始的FFT波浪方程如下：

[hX，t=hk，teikx]? ? ? ? ? ⑾

為得出風(fēng)速較大時(shí)的波浪，不修改參與傅里葉變換的高度函數(shù)，在原來(lái)的FFT表現(xiàn)形式下增加一個(gè)位移項(xiàng)，將原始波浪在水平方向上位移一定距離，使用下式：

[X=X+λDX，t]? ?⑿

λ是一個(gè)位移系數(shù)，用此法模擬海浪需要注意：波浪的頂端可能會(huì)發(fā)生重疊，即海面會(huì)經(jīng)過(guò)自身并發(fā)生卷曲，振幅較大的波浪會(huì)產(chǎn)生相互重疊的位移，解決此問(wèn)題需將λ設(shè)定為較小的值。

[DX，t]可再通過(guò)快速傅里葉變換得到如下公式：

[D（X，t）=-iKkh（K，t）eiKx]? ? ⒀

[DX，t]是一個(gè)二維偏移向量，使用式⑿可得當(dāng)前網(wǎng)格點(diǎn)經(jīng)過(guò)平移后的對(duì)應(yīng)的網(wǎng)格點(diǎn)的位置，為獲得連續(xù)的波浪，使用一次雙線性插值計(jì)算偏移位置對(duì)應(yīng)的海面高度。

1.3 基于FFT的法線貼圖

渲染海面，計(jì)算法線是關(guān)鍵。相鄰網(wǎng)格點(diǎn)的高度差距計(jì)算量小，但海面網(wǎng)格的點(diǎn)數(shù)少且波長(zhǎng)小時(shí)，此法所得近似值不準(zhǔn)確，波長(zhǎng)小于單位網(wǎng)格間距時(shí)難以看到海面光照波動(dòng)效果。

使用式⑹計(jì)算表面次法線B與正切向量T，B和T可對(duì)式⑹求x和z方向的偏導(dǎo)數(shù)求出，可再用一次傅里葉變換，（6）求導(dǎo)可得：

[?h（x，t）=ikh（k，t）eikx]? ? ? ? ? ?⒁

海面參數(shù)方程表示為：

[P（x，z，t）=（x，H（x，z，t），z）]? ? ? ? ? ?⒂

次法線B可表示為：

[B（x，z）=（?x?x，?H（x，z，t）?x，?z?x）=（1，?H（x，z，t）?x，0）]? ⒃

切向量T可表示成：

[T（x，z）=（?x?z，?H（x，z，t）?z，?z?z）=（0，?H（x，z，t）?x，1）] ⒄

頂點(diǎn)法線：

[N=T×B=（-?H（x，z，t）αx，1，-?H（x，z，t）αz）]? ?⒅

把式⑾帶入，得法線表達(dá)式：

[N=（-ik.xh（k，t）eikx，1，-ik.yh（k，t）eikx）]? ⒆

其中，k.x表示k的x分量，k.y表示k的y分量。

網(wǎng)格分辨率的影響海面的渲染效果，借助GPU的強(qiáng)大性能，可較快實(shí)現(xiàn)大規(guī)模網(wǎng)格的FFT變換，通過(guò)在式⒁中二維波浪向量k的x方向以及k的z方向添加小的偏移[Δx]、[Δz]，計(jì)算任意二維波浪向量k[+Δk]的法線值，此法可實(shí)時(shí)計(jì)算出任意分辨率的法線貼圖。將法線貼圖傳遞到片段著色器中進(jìn)行光照計(jì)算，可增加波浪的高頻部分和細(xì)節(jié)。

2 海面光照

海面真實(shí)感是光與海面復(fù)雜交互的結(jié)果。光照模型是用于計(jì)算光源經(jīng)過(guò)物體反射、散射后返回視點(diǎn)的數(shù)學(xué)表示。隨著GPU硬件性能的提升，光照模型的復(fù)雜性和硬件的限制逐漸減弱。利用Cuda架構(gòu)可對(duì)光照模型進(jìn)行有效快速的計(jì)算，構(gòu)建海面光照可有效提升海面波浪細(xì)節(jié)，增強(qiáng)視覺(jué)效果。

海面的光照模型可看作一個(gè)理想鏡面，但此模型不具真實(shí)感，需對(duì)其進(jìn)行改進(jìn)。假設(shè)到達(dá)海面的光線主要分為兩部分：穿透海面經(jīng)過(guò)折射的光線和被反射的光線，這兩部分的光強(qiáng)受到反射和折射系數(shù)的影響。

2.1 反射光線

在理想的鏡面反射中，入射角和反射角相等。本文構(gòu)建一個(gè)反射光線的簡(jiǎn)單描述，海面上的三維點(diǎn)可表示成水平方向x與波浪高度函數(shù)[h（x，t）]的線性組合：

[r（x，t）=x+yh（x，t）]? ? ? ? ? ? ? ⒇

其中，[y]是垂直向上的方向向量。點(diǎn)r處的法線由波浪方程求導(dǎo)得出，[ε（x，t）=?h（x，t）]，法線由下式得出：

[ns（x，t）=y-ε（x，t）1+ε2（x，t）]? ? ? ? ? ? 21

對(duì)于方向[ni]與入射表面交于r的入射光線，其反射光線的夾角只和入射光線和反射表面在r處的法線有關(guān)。由入射角和反射角相等，可得反射向量：

[nr（x，t）=ni-2ns（x，t）（nr（x，t）?ni）]? ? ? 22

2.2 折射光線

折射光線的方向與入射光線和折射表面的法線有關(guān)，遵循Snell定律（入射角和折射角的比值等于入射介質(zhì)和折射介質(zhì)的折射率的比）。設(shè)入射光線從折射率為ni的介質(zhì)入射，折射光線從折射率為ns的介質(zhì)出射，對(duì)于入射角有下面公式：

[sinθi=1-（ni?ns）2=|ni×ns|]? ? 23

根據(jù)Snell定律入射角和折射角有以下關(guān)系：

[ni sinθi=nt sinθt]? ? ? ? ? ? 24

折射向量由下式給出：

[nt（x，t）=±{1-（nint）2|ni×ns（x，t）|2}12]? ? 25

2.3 菲涅耳反射和折射

反射系數(shù)和透射系數(shù)決定光線在介質(zhì)反射和投射的效率。在沒(méi)有光強(qiáng)損失的限定條件下，反射系數(shù)R和T存在以下關(guān)系：

[R+T= 1]? ? ? ? ? ? ? ? ?26

R和T的推導(dǎo)結(jié)果

[R（ni，nr）=12{sin2（θt-θi）sin2（θt+θi）+tan2（θt-θi）tan2（θt+θi）}]? 27

海面上不同的反射率構(gòu)成海面的紋理。反射率是入射角與法線的函數(shù)，由此可知，海面高程梯度是影響海面紋理的關(guān)鍵因素，而浪高度起到波浪間相互遮擋作用，海面高程梯度圖和海面高程圖相比與渲染后的海面具有更加相似的紋理。

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

3.1 實(shí)驗(yàn)環(huán)境與實(shí)驗(yàn)方法

本文使用的CUDA版本為8.0，實(shí)驗(yàn)采用的機(jī)器配置為：Intel Core i5-4430，雙核四線程，主頻3.0GHz，內(nèi)存8GB，NVIDIA GTX960顯卡。

海面模擬和可視化框架由OpenGL著色語(yǔ)言實(shí)現(xiàn)，本文框架對(duì)FFT法線貼圖實(shí)現(xiàn)參數(shù)化控制，理論上其可以支持無(wú)限分辨率海面渲染，但是由于內(nèi)存與顯存容量限制，實(shí)際計(jì)算中并不能得出無(wú)限分辨率的結(jié)果。本文的渲染框架中計(jì)算開銷最大的部分是FFT法線貼圖的計(jì)算，CUDA完成全部計(jì)算過(guò)程。

框架的FFT海面高度場(chǎng)均使用128×128的二維網(wǎng)格生成，為達(dá)到真實(shí)的視覺(jué)效果并節(jié)省計(jì)算開銷，低分辨率的海面網(wǎng)格上使用高分辨率法線貼圖來(lái)增強(qiáng)視覺(jué)效果。當(dāng)法線貼圖的分辨率增加時(shí)，呈現(xiàn)更多高頻波浪及更多海面細(xì)節(jié)。

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

3.2.1 效果圖

所有實(shí)驗(yàn)結(jié)果均在1200×900像素的窗口展示，效果圖使用分辨率為512×512像素的法線貼圖。

3.2.2 測(cè)試圖表

隨著法線貼圖分辨率的提高，得出計(jì)算時(shí)間與幀率趨勢(shì)圖。從表1、圖2可以看出，計(jì)算時(shí)間隨著法線貼圖分辨率的提高而逐漸升高，計(jì)算時(shí)間增長(zhǎng)速度也隨著法線貼圖分辨率的提高而加快，幀率隨著法線貼圖分辨率的提高而逐漸下降，幀率下降速度無(wú)明顯加快的趨勢(shì)。幀率雖然下降，但在所測(cè)試的分辨率的情況之下均可保持在60FPS以上，能夠保證交互感和逼真感。將本文方法與投影網(wǎng)格方法進(jìn)行對(duì)比，本文方法的幀率有所下降，但不影響流暢程度，并且計(jì)算時(shí)間與投影網(wǎng)格方法相比降低10%左右。

4 結(jié)束語(yǔ)

本文構(gòu)建一個(gè)海面模擬可視化框架，采用了FFT波浪方程計(jì)算法線貼圖和納維爾斯托克斯方程（N-S方程）相結(jié)合來(lái)構(gòu)建海面，采用CUDA對(duì)N-S方程的求解和FFT變換，并用Open GL著色器進(jìn)行海面渲染，實(shí)現(xiàn)海面波浪模擬?？刹辉黾雍Ｃ婢W(wǎng)格復(fù)雜度而得到逼真的視覺(jué)效果，節(jié)省計(jì)算時(shí)間，滿足實(shí)時(shí)交互的幀率要求，為海面波浪模擬提供了一定的參考。

參考文獻(xiàn)（References）：

[1] Fournier A，Reeves W T.A simple model of ocean waves[J].ACM SIGGRAPH Computer Graphics，1986.20（4）：75-84

[2] Tessendorf J.Simulating ocean water[J].Simulating nature：realistic and interactive techniques.SIGGRAPH， 2001.1（2）：5.

[3] 童若鋒，汪國(guó)昭.用于動(dòng)畫的水波造型[J].計(jì)算機(jī)學(xué)報(bào)， 1996.8：594-599

[4] Yann L.Realistic water rendering[J].2003）.http：//www. andyc，org/lecture/vewlog. php，2003.

[5] 徐迎慶，蘇成，李華，et al.基于物理模型的流水及波浪模擬[J].計(jì)算機(jī)輔助設(shè)計(jì)與圖形學(xué)學(xué)報(bào)，1997.2：191-192

[6] Creamer D B，Henyey F，Schult R，et al.Improved linear representation of ocean surface waves[J].Journal of Fluid Mechanics，1989.205（-1）：135.

[7] Tessendorf J.Simulating ocean water[J].Simulating nature：realistic and interactive techniques.SIGGRAPH， 2001.1（2）：5.

[8] Creamer D B，Henyey F，Schult R，et al.Improved linear representation of ocean surface waves[J].Journal of Fluid Mechanics，1989.205：135-161