[WebGL] GPU でベジェ曲面をものすごい勢いで描く

[WebGL Advent Calendar] (http://qiita.com/advent-calendar/2015/webgl) 14日目の記事です
***

# モチベーション
WebGL こそ曲面が活躍するプラットフォームだと思うわけです。なぜなら少ないデータで綺麗な面が描画できるからです。100Mバイトもするような大量の三角形をサーバからいちいち送ってる場合ではありませんし、それをアニメーションで変形させるとなれば大変な話です。

そこでまず、どうやって曲面を WebGL で描画するのか考えてみます。曲面と言ってもいろいろあるわけですが、モデリングの自由度などを考えるとやはり最終的にはサブディビジョンサーフェスなりを描画したいという目標があります。そのためにはベジェパッチを描画できることがひとつ有力な選択肢としてあると思います（詳しい説明は [CEDEC2015 サブディビジョンサーフェスのすべてがわかる](http://www.slideshare.net/takahitotejima/cedec2015-52266249) をごらんください）

完成した動くデモはこちら
http://takahito-tejima.github.io/tessTutorial
コードはこちらから
https://github.com/takahito-tejima/takahito-tejima.github.io/tree/master/tessTutorial

# ベジェパッチ
双３次ベジェパッチは、次のような16個のコントロールポイントで定義されます。

![bezier.png](https://qiita-image-store.s3.amazonaws.com/0/64544/649e24d4-fee3-c487-70ed-1cd842ca32b4.png)

ラスタライズでこれを描くためには、適当な密度で三角形にテセレーションしてやる必要があります。OpenGLES は 3.1 でテセレーションが入りましたが、2015年末現在では WebGL はまだテセレーションが自由に使える状況ではありませんので、VertexShader のみでどうにか実装します。基本的には Instanced Tessellation と呼ばれる比較的歴史ある手法にそって実装します。
（参考文献: [Instanced Tessellation, GDC2008] (http://developer.download.nvidia.com/presentations/2008/GDC/Inst_Tess_Compatible.pdf) )

まず１パッチだけ描画する方法を考えてみましょう。幸い WebGL でも VertexShader でテクスチャがフェッチできますので、１６個のコントロールポイントの座標をテクスチャに入れて、VertexShader にテセレーションする密度のトポロジのメッシュ（**リファインメントパターン**と呼びます）を与えてあげて、頂点座標を計算する、という流れが使えます。図にするとこんな感じになります。

![VSflow.png](https://qiita-image-store.s3.amazonaws.com/0/64544/638aae6f-81ee-ec71-5137-a7b46f3bcd93.png)

いやちょっと待て、それなら最初っから uv 空間でメッシュつくらずにベジェパッチを計算したら簡単でいいんじゃない、と思われるかもしれません。まあその通りなんですが、こうすることで

- 高密度のテセレーション（ベジェ計算）もGPUで効率的にできる
- uv メッシュは静的に作りおきでき、アニメーションごとに再計算はいらない
- インスタンス描画 (ANGLE_instanced_arrays) を使って、一つの uv メッシュで複数のパッチを一気に描画できる

などのメリットがあるのです。

# リファインメントパターン（UVメッシュ）生成
VertexShader ではパッチ上のパラメータ座標さえわかれば良いので、どんなメッシュを与えてもよいのですが、今回は簡単のため上図のような uv 空間で二次元のグリッドを作ることにします。これは初期化時に一度だけ行いますので、適当に作っても大丈夫です。たとえば 16x16 分割であれば、

```js:initialize_refinement_pattern_mesh
var div = 16;
var ibo = [], vbo = [];
for (var iu = 0; iu < div; iu++) {
    for (var iv = 0; iv < div; iv++) {
        vbo.push(iu/(div-1));
        vbo.push(iv/(div-1));
        if (iu != 0 && iv != 0) {
            ibo.push(vid - div - 1);
            ibo.push(vid - div);
            ibo.push(vid);
            ibo.push(vid - 1);
            ibo.push(vid - div - 1);
            ibo.push(vid);
            numTris += 2;
        }
        ++vid;
    }
}
```

こんな感じで、0.0 - 1.0 の (u, v) グリッド頂点と、それをつなぐ三角形のインデクスバッファを作ってやります。これらはいつもどおり、それぞれ ARRAY_BUFFER と ELEMENT_ARRAY_BUFFER にバインドしてやります。

```js:bind_vbo_and_ibo
var iboBuffer = gl.createBuffer();
gl.bindBuffer(gl.ELEMENT_ARRAY_BUFFER, iboBuffer);
gl.bufferData(gl.ELEMENT_ARRAY_BUFFER, new Uint16Array(ibo), gl.STATIC_DRAW);
gl.bindBuffer(gl.ELEMENT_ARRAY_BUFFER, null);

var vboBuffer = gl.createBuffer();
gl.bindBuffer(gl.ARRAY_BUFFER, vboBuffer);
gl.bufferData(gl.ARRAY_BUFFER, new Float32Array(vbo), gl.STATIC_DRAW);
gl.bindBuffer(gl.ARRAY_BUFFER, null);
```

# コントロールポイントテクスチャ作成
ベジェパッチのコントロールポイントはテクスチャとして与えてやります。ここはちょっと手抜きして、OES_texture_float 拡張を使わせてください。コンテキストの初期化時に
`gl.getExtension('OES_texture_float'))`
を呼んでおきます。

今回は１パッチしかありませんので、16個の XYZ だけ必要です。シンプルに 16x1 の float RGB な TEXTURE_2D を構築することにします。R, G, B にそれぞれ X, Y, Z の座標を入れます。

```js:populate_control_points
// create 16 control points.
var cp = new Float32Array(16*3);  // 16 * xyz
// todo: ここで cp を設定します。

// テクスチャとして送ります。
// フィルタリングはいらないので NEAREST 設定で。
gl.activeTexture(gl.TEXTURE0);
var tex = gl.createTexture();
gl.bindTexture(gl.TEXTURE_2D, tex);
gl.pixelStorei(gl.UNPACK_ALIGNMENT, 1);
gl.texParameteri(gl.TEXTURE_2D, gl.TEXTURE_MAG_FILTER, gl.NEAREST);
gl.texParameteri(gl.TEXTURE_2D, gl.TEXTURE_MIN_FILTER, gl.NEAREST);
gl.texParameteri(gl.TEXTURE_2D, gl.TEXTURE_WRAP_S, gl.CLAMP_TO_EDGE);
gl.texParameteri(gl.TEXTURE_2D, gl.TEXTURE_WRAP_T, gl.CLAMP_TO_EDGE);
gl.texImage2D(gl.TEXTURE_2D, 0, gl.RGB, 16, 1, 0, gl.RGB, gl.FLOAT, cp);
```

# シェーダー
さてリソースは全部揃ったので、vertex shader でベジェパッチの評価をしてやります。fragment shader はなんでも構いません。もしシェーディングに法線が必要であれば、法線はベジェパッチの評価時に同時に計算が可能です。

まず頂点データのテクスチャを引くテクスチャ座標を計算する関数を作ります。

```glsl:get_vertex
vec2 getVertex(float cpIndex)
{
  return vec2((cpIndex + 0.5)/16.0, 0.5);
}
```

怪しげなコードですが、cpIndex は 0 ~ 15 の整数を float で取ります（すいません）。
横幅 16 ピクセルのテクスチャなので、0.5 足して正規化してやるとテクスチャ座標にできます。

つぎに３次ベジェ基底を計算をする関数を書きます。

```glsl:eval_bezier
void evalCubicBezier(in float u, out float B[4], out float D[4]) {
    float t = u;
    float s = 1.0 - u;
    float A0 = s * s;
    float A1 = 2.0 * s * t;
    float A2 = t * t;

    B[0] = s * A0;
    B[1] = t * A0 + s * A1;
    B[2] = t * A1 + s * A2;
    B[3] = t * A2;

    D[0] =    - A0;
    D[1] = A0 - A1;
    D[2] = A1 - A2;
    D[3] = A2;
}
```
入力パラメータ u について、基底 B および一次微分 D を出力する関数になります。この辺はコンパイラがいろいろ最適化してくれることに期待して、わかりやすく書きます。

最後にこれらを使って vertex shader の main 関数を作ります。u 方向、v 方向 の順に計算していきます。パラメータ座標は頂点アトリビュートとして inUV で与えられ、コントロールポイントは texCP からフェッチできます。

```glsl:vertex_shader

attribute vec2 inUV;
uniform sampler2D texCP;
varying vec3 normal;

void main() {
    float B[4], D[4];
    vec3 cp[16];
    vec3 BUCP[4], DUCP[4];

    // 16個のコントロールポイントを取得
    for (int i = 0; i < 16; ++i) {
      vec2 st = getVertex(float(i));
      cp[i] = texture2D(texCP, st).xyz;
    }

    // U 方向の基底を計算
    evalCubicBezier(inUV.x, B, D);

    // BUCP, DUCP にそれぞれ 4 つの頂点を計算して一旦保存
    BUCP[0] = cp[0]*B[0] + cp[4]*B[1] + cp[ 8]*B[2] + cp[12]*B[3];
    BUCP[1] = cp[1]*B[0] + cp[5]*B[1] + cp[ 9]*B[2] + cp[13]*B[3];
    BUCP[2] = cp[2]*B[0] + cp[6]*B[1] + cp[10]*B[2] + cp[14]*B[3];
    BUCP[3] = cp[3]*B[0] + cp[7]*B[1] + cp[11]*B[2] + cp[15]*B[3];

    DUCP[0] = cp[0]*D[0] + cp[4]*D[1] + cp[ 8]*D[2] + cp[12]*D[3];
    DUCP[1] = cp[1]*D[0] + cp[5]*D[1] + cp[ 9]*D[2] + cp[13]*D[3];
    DUCP[2] = cp[2]*D[0] + cp[6]*D[1] + cp[10]*D[2] + cp[14]*D[3];
    DUCP[3] = cp[3]*D[0] + cp[7]*D[1] + cp[11]*D[2] + cp[15]*D[3];

    // V 方向の基底を計算
    evalCubicBezier(inUV.y, B, D);

    // BUCP を V の基底で混合して座標を得る
    vec3 Pos       = B[0]*BUCP[0] + B[1]*BUCP[1] + B[2]*BUCP[2] + B[3]*BUCP[3];
    // DUCP を使うとU、V方向それぞれの接線が計算できる
    vec3 Tangent   = B[0]*DUCP[0] + B[1]*DUCP[1] + B[2]*DUCP[2] + B[3]*DUCP[3];
    vec3 BiTangent = D[0]*BUCP[0] + D[1]*BUCP[1] + D[2]*BUCP[2] + D[3]*BUCP[3];

    // 法線は二つの接線の外積で求める。
    vec3 n = normalize(cross(BiTangent, Tangent));

    // ビュー変換（法線側は手抜きなので注意...）
    vec3 p = (modelViewMatrix * vec4(Pos.xyz, 1)).xyz;
    normal = (modelViewMatrix * vec4(n, 0)).xyz;

    // プロジェクション行列をかけて、gl_Position に出力する
    gl_Position = projectionMatrix * vec4(p, 1);
}
```


# 描画
リファインメントパターンの uv 値を頂点アトリビュートに、コントロールポイントをテクスチャをバインドして上記のシェーダを動かせば、無事ベジェパッチが描画できます。

```js:draw
gl.useProgram(program);

gl.activeTexture(gl.TEXTURE0);
gl.bindTexture(gl.TEXTURE_2D, tex);

gl.bindBuffer(gl.ELEMENT_ARRAY_BUFFER, iboBuffer);
gl.bindBuffer(gl.ARRAY_BUFFER, vboBuffer);
gl.enableVertexAttribArray(0);
// リファインメントパターンは2要素 (u, v)
gl.vertexAttribPointer(0, 2, gl.FLOAT, false, 0, 0);

gl.drawElements(gl.TRIANGLES, numTris*3, gl.UNSIGNED_SHORT, 0);
```

![image](https://qiita-image-store.s3.amazonaws.com/0/64544/c6dcbde6-3e1e-1efd-8e4a-63e2d45b9487.png)

# インスタンス描画
さて１パッチは描画できたので、次は複数パッチの描画です。これには [ANGLE_instanced_array 拡張] (https://www.khronos.org/registry/gles/extensions/ANGLE/ANGLE_instanced_arrays.txt) を使います。大まかな流れは同様ですが、以下の図の様になります。

![VSIflow.png](https://qiita-image-store.s3.amazonaws.com/0/64544/78a9ac0b-9b50-b3ae-d02c-3e16e556110b.png)


まずコントロールポイントのテクスチャを、パッチ数分まで拡張します。ここでは単純に縦に並べることにしました。128 個パッチがあれば 16*128 のテクスチャを作ることになります。

```js:populate_control_points
// create 16 control points * numPatches
var cp = new Float32Array(numPatches*16*3);  // numPatches * 16 * xyz
// todo: ここで cp を設定します。

var tex = gl.createTexture();
gl.activeTexture(gl.TEXTURE0);
gl.bindTexture(gl.TEXTURE_2D, tex);
gl.pixelStorei(gl.UNPACK_ALIGNMENT, 1);
gl.texParameteri(gl.TEXTURE_2D, gl.TEXTURE_MAG_FILTER, gl.NEAREST);
gl.texParameteri(gl.TEXTURE_2D, gl.TEXTURE_MIN_FILTER, gl.NEAREST);
gl.texParameteri(gl.TEXTURE_2D, gl.TEXTURE_WRAP_S, gl.CLAMP_TO_EDGE);
gl.texParameteri(gl.TEXTURE_2D, gl.TEXTURE_WRAP_T, gl.CLAMP_TO_EDGE);
// 幅は 16 のままですが、高さが numPatches になりました
gl.texImage2D(gl.TEXTURE_2D, 0, gl.RGB, 16, numPatches, 0, gl.RGB, gl.FLOAT, cp);
```

ANGLE_instanced_array 拡張を使うために、コンテキスト初期化時に
`var ext = gl.getExtension('ANGLE_instanced_arrays')`
とします。この戻り値は関数の呼び出しで使いますので、取っておいてください。

シェーダ内ではインスタンスごとにコントロールポイント座標を読み分ける必要があります。標準GLSLだと gl_InstanceID を使えるんですが、やはり GLES では 3.0 以降でないと使えません。ただ AttribDivisor を使うと gl_InstanceID と同じようにインスタンスごとにパラメータ化するができます。このためにまず、補助的なVBOを用意します。ここではこのVBOは描画したいインスタンスの数だけ 0, 1, 2, ... のように整数を並べたものを作ります。つまり、

```js:patch_index_buffer
var indices = [];
for (var i = 0; i < numPatches; ++i) indices.push(i);
var patchIndexBuffer = gl.createBuffer();
gl.bindBuffer(gl.ARRAY_BUFFER, patchIndexBuffer);
gl.bufferData(gl.ARRAY_BUFFER, new Float32Array(indices), gl.STATIC_DRAW);
```

こんな感じです。numPatches には描画したいパッチの数を入れておきます。

これをシェーダに渡すときに、先ほどのGL拡張を使って

```js:instanced_draw_call
gl.useProgram(program);

gl.activeTexture(gl.TEXTURE0);
gl.bindTexture(gl.TEXTURE_2D, tex);

gl.bindBuffer(gl.ELEMENT_ARRAY_BUFFER, iboBuffer);
gl.bindBuffer(gl.ARRAY_BUFFER, vboBuffer);
gl.enableVertexAttribArray(0);
gl.vertexAttribPointer(0, 2, gl.FLOAT, false, 0, 0);

// ここまでは前回と同じ

// もう一つ頂点アトリビュートを与える
gl.enableVertexAttribArray(1);
gl.bindBuffer(gl.ARRAY_BUFFER, patchIndexBuffer);
gl.vertexAttribPointer(1, 1, gl.FLOAT, false, 0, 0);

// ただし、1 インスタンスごとに 1 アトリビュート (この場合 float 1個）進める
ext.vertexAttribDivisorANGLE(1, 1);

// numPatches 個だけインスタンス描画
ext.drawElementsInstancedANGLE(gl.TRIANGLES, numTris*3, gl.UNSIGNED_SHORT, 0, numPatches);
```

シェーダ側ではパッチ番号も使ってテクスチャ座標を生成します。numPatches は uniform 変数として渡してやりましょう。

```glsl:get_vertex
uniform float numPatches;
vec2 getVertex(float patchIndex, float cpIndex)
{
  return vec2((cpIndex + 0.5)/16.0,
              (patchIndex + 0.5)/numPatches);
}
```
あとは２つめの AttribDivisor が効いた頂点アトリビュートからパッチ番号を取得します。頂点アトリビュートが複数になるので、シェーダのコンパイル時にそれぞれバインディングポイントの辻褄が合うように気をつけてください。ちゃんとgl.bindAttribLocation を呼んだほうが安全です。

```glsl:vertex_shader
attribute vec4 inUV;
attribute float patchIndex;
void main() {
    vec3 cp[16];
    for (int i = 0; i < 16; ++i) {
      vec2 st = getVertex(patchIndex, float(i));
      cp[i] = texture2D(texCP, st).xyz;
    }

    //以下略
```

これで複数パッチを一つのUVメッシュで一気に描画することができるようになりました。パッチごとに色分けしてみました。
![image](https://qiita-image-store.s3.amazonaws.com/0/64544/df4c122d-3cbc-8027-1e8e-e89e4c8a2e31.png)

このサンプルの見た目はパッとしませんが、100万頂点くらい、わりと動いてます。
ワイヤーフレームの図
![image](https://qiita-image-store.s3.amazonaws.com/0/64544/97ba9188-8fd6-9368-cf25-57e55e869179.png)

# おわりに
OES_texture_float と ANGLE_instanced_arrays はほとんどのブラウザがサポートしているようですし、いまどきのモバイルだとこれくらいの vertex shader なら問題なく動くので、実用的に使えるんじゃないかと思います。

デモを見るとわかるかと思いますが、このサンプルではパッチ間で法線が連続してません。これはベジェのコントロールポイントをそれぞれ適当に動かしているからで、これをなめらかにするにはパッチ間で辻褄が合うようにちゃんと計算してやるか、あるいはBスプラインサーフェスの出番になります。まあ PN でもいいんですけど・・・

ここまでできれば、パッチ法線にそって頂点を動かすディスプレイスメントなども出来るようになりますね。