これがGPUの力！three.jsによる“リアルタイム”なレイトレーシング 〜宝石編〜

これは[WebGL Advent Calendar 2015](http://qiita.com/advent-calendar/2015/webgl) 23日目の記事（の代打[^1]）です。

また、WebGL Advent Calendar 2015 15日目の記事である「[これがGPUの力！Three.jsによる“リアルタイム”なレイトレーシング](http://qiita.com/gam0022/items/03699a07e4a4b5f2d41f)」の続編です。

# はじめに
今夜は雪だそうですね。こんな寒い日にはGPUをぶん回して暖をとりましょう！

光の屈折をシミュレートすることで、輝く宝石をWebGLでレンダリングする「gem」という作品（技術デモ）をつくりました。レイトレーシングをGLSLのフラグメントシェーダで実装することで、GPUの並列計算を利用したリアルタイムな描画を実現しています。

WebGLで動くので、次のリンクを開くとブラウザ上でそのまま動作します。PCだけでなくiPhone6でも動作を確認しています。

- [webgl Raytracing example - gem](http://gam0022.net/webgl/#raytracing_gem)

宝石の屈折率、色（拡散反射・鏡面反射）、宝石の形状、自動カメラ、ガンマ補整など、各種パラメータを変更できます。負荷が重い場合は一番下の「Resolution（解像度）」を128にしてください。また一番上からパラメータのプリセットを選択できます。

![160130-0001.png](https://qiita-image-store.s3.amazonaws.com/0/17718/37598279-fad0-4f29-9db3-cb6e325842c1.png)

# 宝石と屈折
宝石のように透明な物体を描画するためには、屈折の計算が必要になります。1回までの屈折であれば、ポリゴンのラスタライズ方式でも[環境マッピング](https://wgld.org/d/webgl/w046.html)を用いて高速にレンダリングできます。

次の図は宝石内部での光の反射と屈折の様子です。宝石の輝きは宝石の内部での光の反射と屈折により生じます。そのため、複数回の屈折をシミュレートしなければ、宝石の輝きは再現できません。

![dia cut.png](https://qiita-image-store.s3.amazonaws.com/0/17718/358670fe-87ef-c756-43fd-a5dad4e36d1f.png)
この図は[jewelry-tanigawa.com](http://jewelry-tanigawa.com/y-diamond4c.html)から引用させていただきました。

複数回の屈折を行うためには、ラスタライズ方式ではなく、レイトレーシングが必要になります[^2]。レイトレーシングをするためには膨大な計算が必要ですが、GLSLのフラグメントシェーダで実装することで、GPUの並列計算を利用したリアルタイムに近い速度のレンダリングを実現しました。

# 実行結果
あらためて最初の節で紹介した「gem」の実行結果を見てみましょう。

![round_brilliant.png](https://qiita-image-store.s3.amazonaws.com/0/17718/91fc3718-a1f4-6fe9-d9e2-20b8f230461b.png)
これは「Round Brilliant」という形状です。宝石の内部で光が反射と屈折を繰り返すことで、このような結果になります。この例だと、MacBookAirでも256x256の解像度で60fps出ました。手元のiPhone6でも20fps前後でした。

![barion.png](https://qiita-image-store.s3.amazonaws.com/0/17718/6a188dbb-7568-5ab6-4014-0ffabbfb765d.png)
「Barion」と呼ばれる形状です。エメラルドグリーンを意識した色（#084d04/#acff9b）に調節しました。

![oval_brilliant.png](https://qiita-image-store.s3.amazonaws.com/0/17718/ec5f1cb0-17c6-7f43-0c57-303853f7086f.png)
「Oval Brilliant」という形状です。トパーズを意識した色（#4f5906/#ff6800）に調整しました。

![round_brilliant2.png](https://qiita-image-store.s3.amazonaws.com/0/17718/3689ce3d-b7ab-b743-c3e5-769930751c85.png)
これも「Round Brilliant」ですが、1つ目よりも上部のクラウンと呼ばれる箇所が薄いバージョンになります。

![ramiel.png](https://qiita-image-store.s3.amazonaws.com/0/17718/dfa094a0-9559-499f-73be-1a14c548cfaf.png)
「Regular Octahedron（正八面体）」です。色を#0b024c/#ffffff にして、Refractive Index（屈折率）は1.1と低めにしました。エヴァンゲリオンのラミエルっぽいですね。面数が8と少ないため、MacBookAirでもフルスクリーンで60fpsでうごきます。

# 実装の概要
詳細な実装について解説する前に、概要を述べておきます。

- WebGLのライブラリであるthree.jsで実装
- three.jsのシーンとしてはPlaneが1枚あるだけ
    - Planeのシェーディング（色の塗り方）処理として、レイトレーシングを実装
- objファイルを用意すれば様々な形状をレンダリング可能
    - CPU（JavaScript)からGPU(GLSL)に動的に形状情報を橋渡し
    - CPUからGPUに渡せるデータのサイズには厳しい制限がある
        - 頂点情報を浮動小数データテクスチャとして渡すことで制限を回避
- three.jsの資産を活かして、objファイルのパースやカメラのコントロールはthree.jsに丸投げ

# フラグメントシェーダでレイトレを実装するとは
GLSLのフラグメントシェーダとしてレイトレーシングを実装するとはどういうことでしょうか？

実装の概要で述べたとおり、three.jsのシーンとしてはPlaneが1枚があるだけです。Planeの塗り方を計算する処理として、フラグメントシェーダ（GPUで動作するポリゴンのピクセルの塗り方を決めるプログラム）でレイトレーシングを実装しました。

Planeの塗り方を計算する処理として、3D空間を描画するというアプローチは、前回の記事と全く一緒です。`THREE.RawShaderMaterial`を適用してフラグメントシェーダでゴリゴリとレイトレーシングを実装しています。前回の記事の繰り返しになってしまうので、今回の記事ではこのアプローチについては深く解説しません。

レイトレーシングって何だろう？という方は、先にdoxasさんの「[WebGL と GLSL で気軽にレイトレーシングに挑戦してみよう！](http://qiita.com/doxas/items/477fda867da467116f8d)」という記事を読むことをオススメします。レイトレーシングのGLSL実装について、基礎から解説されています。

前回の記事はレイトレーシングの中でもレイマーチング（スフィアトレーシング）という方法を利用しました。レイマーチングでは、シーンを距離関数で定義して、レイの先端が物体表面にたどり着くまでループを繰り返して衝突判定します。

今回は一般的に知られるレイトレーシングのように、解析的な衝突判定にしました。宝石のobjファイルを読み込んで表示するという都合上、距離関数でシーンを定義するのは難しいからです。

レイトレーシングでもレイマーチングでもシーンとレイ（半直線）の衝突判定が肝になります。

今回はTomas Mollerの手法を用いて、ポリゴン（三角形）とレイ（半直線）の衝突判定を実装しました。

Tomas Mollerの手法は、ポリゴンとレイの衝突判定としては高速で一般的に使われるアルゴリズムです。自分がGLSLで実装する際に参考にした資料を紹介しておきます。

- [[3D] 三角形と線分の交差判定（Raycast）](http://qiita.com/edo_m18/items/2bd885b13bd74803a368)
- [レイトレーシングによるコンピュータグラフィクス入門- レイと三角形の交差判定 -](http://kanamori.cs.tsukuba.ac.jp/jikken/inner/triangle_intersection.pdf)

次のコードはTomas Mollerの手法をGLSL実装したものです。

```glsl:GLSLによるポリゴン（三角形）とレイ（半直線）の衝突判定
// 3次正方行列の行列式をクラメルの公式で計算する
float det( vec3 a, vec3 b, vec3 c ) {
	return (a.x * b.y * c.z)
			+ (a.y * b.z * c.x)
			+ (a.z * b.x * c.y)
			- (a.x * b.z * c.y)
			- (a.y * b.x * c.z)
			- (a.z * b.y * c.x);
}

// ポリゴン（三角形）とレイ（半直線）の衝突判定する
void rayIntersectsTriangle( vec3 origin, vec3 ray, vec3 v0, vec3 v1, vec3 v2, inout Intersect nearest ) {

	vec3 invRay = -ray;
	vec3 edge1 = v1 - v0;
	vec3 edge2 = v2 - v0;

	float denominator =  det( edge1, edge2, invRay );
	if ( denominator == 0.0 ) return;

	float invDenominator = 1.0 / denominator;
	vec3 d = origin - v0;

	float u = det( d, edge2, invRay ) * invDenominator;
	if ( u < 0.0 || u > 1.0 ) return;

	float v = det( edge1, d, invRay ) * invDenominator;
	if ( v < 0.0 || u + v > 1.0 ) return;

	float t = det( edge1, edge2, d ) * invDenominator;
	if ( t < 0.0 || t > nearest.distance ) return;

	nearest.isHit    = true;
	nearest.position = origin + ray * t;
	nearest.distance = t;
	nearest.normal   = normalize( cross( edge1, edge2 ) ) * sign( invDenominator );
	nearest.material = DIAMOND_MATERIAL;
	nearest.isFront  = invDenominator > 0.0;
	
}

```

1つ目の資料のJavaScriptによる実装をGLSLに移植したものになりますが、3つ工夫をしました。

1つ目の工夫について。除算コストは乗算コストより大きいことから、`denominator`の除算を、事前に`denominator`の逆数を計算することで、乗算に置き換えて，コストを削減しました。

2つ目の工夫について。今回は物体を透過させて描画するので、背面からの衝突も検知する必要があります。そのため、`denominator`が負の場合も衝突したとみなすようにしました。表裏の判定は最後に`denominator`の符号をチェックするようにしました。

3つ目の工夫について。背面の場合は法線を反転させる処理が必要です。一般的にGPUで計算する場合にはifなどの分岐のコストが大きいので、条件分岐を避けるために、組み込みの`sign`という組み込みメソッドを利用するようしました。

# レイトレーシングにおける屈折の実装
今回のテーマは「屈折」です。本節ではレイトレーシングにおける屈折の実装について述べます。
## 屈折の方向の計算には組み込み関数の`refract`をつかう
GLSLでは`refract`という屈折のベクトルを計算する組み込み関数があるので、これを利用しました。

- https://www.opengl.org/sdk/docs/man/html/refract.xhtml

入射ベクトル、法線ベクトル、屈折率を引数にして、屈折ベクトルを計算します。全反射の判定も可能で、全反射の場合は0ベクトルが返ります。

## 衝突後のレイの原点を修正する方向の注意点
レイトレーシングで反射をする場合、レイと物体の交点を反射用のレイの原点として、再度レイを飛ばします。

しかし、交点をそのまま原点にしてしまうと、交点の位置で衝突していると誤判定されて、レイが進まないことがあります。これを回避するために、法線方向にむかってレイの原点をずらす必要があります。

ここで注意点なのですが、屈折の場合には、法線の逆方向、つまり物体の内側にずらさなくてはいけません。屈折の場合は物体の内部を通過するレイとなるので、当然といえば当然ですね。

## 屈折率の反転
レイが空気中から宝石内部に突入する場合と、宝石内部から空気中に脱出する場合では、屈折率を反転する必要があります。

「レイが脱出するとき」の屈折率は「レイが突入するとき」の屈折率の逆数にしなくてはなりません。

レイが突入した OR 脱出した の判定は、`rayIntersectsTriangle`関数で計算した`isFront`を利用します。

今回は使用した宝石の形状はソリッドモデルなので、表からポリゴンに衝突したときにはレイの突入、裏からポリゴンに衝突したときにはレイの脱出であるという法則が成り立ちます。

## GLSLで再帰が使えない制約について
GLSLでは再帰が使えないという制約があります。反射しか考慮しない場合にはループだけで処理できるのですが、反射と屈折の両方を実現するためには、再帰的に反射方向のレイと、屈折方向のレイをそれぞれトレーシングする必要があるので、ループだけでは処理できません。

苦肉の策として、宝石に関しては反射の計算は省いてしまい、全反射のときだけ反射をするようにしました。反射を考慮していませんが、ある程度はそれらしい結果になったと思っています。

# 浮動小数データテクスチャでCPUからGPUに大量の頂点情報を送る
CPU（JavaScript）からGPU（GLSL）へデータを渡すには`uniform`を使います[^3]。

モデルの頂点は3次元のベクトル（`vec3`）です。宝石のモデルの形状情報をシェーダに渡すためには、モデルのデータをトライアングルリスト[^4]に変換して`vec3`の1次元配列である`v3v`を`uniform`で渡すというのが最もシンプルな実装です。

しかし、`uniform`で渡せるデータサイズには厳しい制限があります。今回もまずは`v3v`で頂点を渡すような実装にしたのですが、環境によっては面数100くらいで`too many uniforms`というエラーになってしました。

そこで`v3v`ではなく、テクスチャの画素値として頂点情報を詰めることにしました。
つまり、画素値のRGB値をそれぞれ頂点のXYZとして、CPUからGPUに頂点情報を渡す作戦です。

通常のテクスチャでは精度が足りなかったので、浮動小数データテクスチャを利用することにしました。

通常のテクスチャは1ピクセルあたりRGBAそれぞれUint8（符号無し整数）の8bit、合計32bitが割り当てられます。浮動小数データテクスチャでは1ピクセルあたりRGBAそれぞれFloat32（浮動小数）の32bit、合計128bitが割り当てられます。

浮動小数データテクスチャはWebGLの拡張機能という位置づけで、実装されていることの保証がないAPIとなりますが、実際には多くのデバイスで利用可能のようです。

three.jsから浮動小数データテクスチャを利用する方法については、[JavaScriptoon2](https://techbooster.github.io/c89/#scriptoon2)というC89本の「SIMD.jsを超える！？WebGLで100万パーティクルへの挑戦」という記事に詳しく書かれています。

# three.jsの資産の活用
three.jsでは3D関連の役に立つコードが既にたくさん存在します。ですので、今回も利用できるものは積極的に利用しました。

例えば、objファイルの読み込みに[`THREE.OBJLoader`](http://threejs.org/docs/#Reference/Loaders/OBJLoader)を利用しました。

また、自由モード（マウスのドラッグなで視点の移動ができる）でのカメラの制御には、[`THREE.OrbitControls`](http://threejs.org/examples/#misc_controls_orbit)を利用しました。カメラの原点と向きを`uniform`を使ってフラグメントシェーダに渡す設計にしました。

<!--
なお、デフォルトでは自由カメラモードではなく、オートカメラモードになっています。自由カメラモードにするためには、`Free Camera`のチェックボックスをONにしてください。
-->

<!--
# NG集
成功した結果だけだと面白く無いので、開発途中のスクリーンショットを紹介します。

![polygon-1.png](https://qiita-image-store.s3.amazonaws.com/0/17718/9fce32cf-908c-f89f-1b27-c3ff2b85e939.png)
上記のTomas Mollerの手法で三角形の衝突判定がとりあえずできた段階です。

![polygon-2.png](https://qiita-image-store.s3.amazonaws.com/0/17718/f5067c7f-43fa-90c3-c865-7e6697a5619f.png)
三角形を2つ並べてみた段階です。

![vertexs.png](https://qiita-image-store.s3.amazonaws.com/0/17718/b178ba2d-12f2-c17e-e463-f4c718251477.png)
objの読み込みまでができたものの、愚直に頂点を渡していた頃のものです。
すべての頂点を渡すとシェーダがコンパイルできないので、面を適当に間引いた結果です。

![bug-1.png](https://qiita-image-store.s3.amazonaws.com/0/17718/0cea6601-6bd6-1e84-f1e6-ad195f79c8f9.png)
屈折の計算などがバグっていた頃の結果です。さらに背景にゴミも混じるバグがありました。

![bug-3.png](https://qiita-image-store.s3.amazonaws.com/0/17718/0c755aeb-302d-6b36-ba67-a756180e5cba.png)
屈折バグその2。

![index.png](https://qiita-image-store.s3.amazonaws.com/0/17718/49f68d2c-1865-3daa-2aa2-95106519224d.png)
市松模様の床を出した直後で、まだパラメータの調整ができていない頃の結果です。
屈折率が小さすぎて、ただの透過のようになってしまっています。

![ruby.png](https://qiita-image-store.s3.amazonaws.com/0/17718/de141efb-1633-d857-cec8-a47804a781ab.png)
今とほぼ同じ結果。Specular Color が無しの状態です。

![ruby-gamma.png](https://qiita-image-store.s3.amazonaws.com/0/17718/0446eed3-9806-bca5-9973-9ad05a3a8ef0.png)
試しにガンマ補正をかけた結果です。現在も`Gamma`のチェックボックスをONにすると利用できます。

diamondのモデルデータについて
https://3dwarehouse.sketchup.com/user.html?id=0093174706763663929804792
-->

# おわりに
ここまでお付合いいただき、本当にありがとうございます。

本当は宝石をたくさん並べてみたかったのですが、負荷的に厳しすぎて断念しました。

GLSLのフラグメントシェーダでレイトレーシングを実装しようとすると、それなりに制約があるので大変ですが、試行錯誤の過程がとても楽しいと思います。

みなさんもGLSLでレイトレーシングしてみませんか？

[^1]: もう2016/1/29ですが、WebGL Advent Calendar 2015の23日担当の人が投稿しなかったので代理で投稿しました。本来であれば、WebGLがWebテクノロジー部門ランキングで余裕のNo.1になるはずでしたが、23日に欠番があったため、[25日間全ての参加登録がされているカレンダーを対象とした](http://qiita.com/advent-calendar/2015/ranking/subscriptions)を満たせずランキング外([ストック数](https://qiita-image-store.s3.amazonaws.com/0/17718/7506381d-7bf2-e029-1fb2-8f88d420fda0.png)、[購読数](https://qiita-image-store.s3.amazonaws.com/0/17718/c24769ee-05d8-ead6-6be0-043061ecd70b.png))となっていました。個人的にとても無念だったので、代わりに投稿するに至りました。【追記】見事WebGLがランキング1位（[ストック数](http://qiita.com/advent-calendar/2015/ranking/stocks/categories/web_technologies)、[購読数](http://qiita.com/advent-calendar/2015/ranking/subscriptions/categories/web_technologies)）になりました :clap:

<!--
（[ストック数](https://qiita-image-store.s3.amazonaws.com/0/17718/9e0356f1-3c34-b83c-c02e-05e9f1c805e5.png)、[購読数](https://qiita-image-store.s3.amazonaws.com/0/17718/fe716052-a467-1e11-ab8e-9e223b57d208.png)）
-->

[^2]: レイトレーシングを使わない[2回屈折を考慮した描画の研究](http://www.npal.cs.tsukuba.ac.jp/projects/2014/liquid_rendering/index.html)も存在します。しかし、背面と全面に分けて2回屈折を実現する手法では、今回の宝石のように内部で反射と屈折をするケースや3回以上の屈折にはおそらく対応できないと考えています。

[^3]: 他には`attribute`などもありますが、それは本筋ではないので割愛します。気になる方は[シェーダの記述と基礎](https://wgld.org/d/webgl/w008.html)を読んで下さい。

[^4]: トライアングルリストとは、ポリゴンを構成する3頂点を並べたものです。面数が100の形状であれば、100*3 = 300 要素のトライアングルリストになります。重複する頂点分はメモリ的には無駄になりますが、頂点と面を別々に持つよりも、フラグメントシェーダ側から頂点を取り出す実装が簡単になります。

<!--また、検証はしていませんが、おそらくパフォーマンス的にもトライアングルリストの方が利点があると考えています。頂点と面を別々のテクスチャとして渡す場合、相互にテクスチャを参照するフェッチコストが大きくなってしまいます。私も自信が無いので、詳しい方は教えて下さい。-->