14. レイマーチング

前の章では、ゲームが普通やる方法で3Dを描きました。
三角形、行列、そしてラスタライザーです。
此の章では其れを全て捨て、ジオメトリを全く使わずに、1つのフラグメントシェーダーの中だけで3Dのシーンを描きます。
全てのピクセルがシーンにレイを飛ばし、何かに当たるまでレイに沿って歩きます。
其の為に第6章の距離関数を、今度は3Dで使います。
此の手法をレイマーチングと呼び、其の中でも此処で使う種類はスフィアトレーシングと呼ばれます。

Shadertoyの有名なデモの殆どが此の手法で作られていて、シリーズ全体の考え方が一緒に働く様子を見るのに最適です。
座標（第5章）、距離関数と滑らかな和集合（第6章）、変換と繰り返し（第7章）、ライティング（第9章と第10章）、そして線型代数の第11章のカメラの基底です。

第5章のシェーダープレイグラウンドに戻ります。
此の章はCの変更が全く要らないので、第13章のプロジェクトではなく、プレイグラウンドのプロジェクトから始めて下さい。
全てはshader/playground.glslの中で行い、プログラムを動かしたまま、変更する度にRキーを押せます。
最後には、柔らかい影、アンビエントオクルージョン、霧、無限の柱、周りを回るカメラ付きの此のシーンになります：

完成したレイマーチングのシーン：チェッカー模様の床の上の角の丸い箱、くっついた球、トーラス、霧に消える柱

14.1 考え方

ラスタライズでは、三角形をループして「此の三角形はどのピクセルを覆うか？」と聞きます。
レイを使う手法は其れを逆にします。
ピクセルをループして（フラグメントシェーダーが既にやってくれます）、「此のピクセルには何が見えるか？」と聞きます。
ピクセル毎にレイを作ります。
始点 $\mathbf{r}_o$ （カメラ）と単位方向ベクトル $\mathbf{r}_d$ で、レイ上の全ての点は：

\mathbf{p}(t) = \mathbf{r}_o + t\,\mathbf{r}_d, \qquad t \ge 0

$t$ はレイに沿って進んだ距離です。
問題は、 $\mathbf{p}(t)$ が表面に触れる最小の $t$ は何か、です。

球や平面なら $t$ の方程式を解けます（其れが古典的なレイトレーシングです）。
然し、複雑な、混ぜ合わせた、捻じった形には其の様な式がありません。
でも符号付き距離関数（SDF）があれば、全ての形に使える非常に簡単なループがあり、其れが此の章の秘訣です。

14.2 レイで出来たカメラ

先ず、ピクセル毎に1本のレイが必要です。
第5章の中心を原点にした座標を使います。
$y$ は−1〜1で、 $x$ は縦横比で伸ばされます（800×600なら−1.333〜1.333）：

  vec2 uv = (2.0 * gl_FragCoord.xy - uResolution) / uResolution.y;

カメラは $\mathbf{r}_o$ にあり、注視点を見ます。
其の2点から、線型代数の第11章のlookAt（と第13章のmat4LookAt）と全く同じ正規直交基底を作ります：

\mathbf{f} = \frac{\mathbf{t} - \mathbf{r}_o}{\lVert \mathbf{t} - \mathbf{r}_o \rVert}, \qquad \mathbf{r} = \frac{\mathbf{f} \times \mathbf{u}_{\text{world}}}{\lVert \mathbf{f} \times \mathbf{u}_{\text{world}} \rVert}, \qquad \mathbf{u} = \mathbf{r} \times \mathbf{f}

そしてカメラの前の距離 $F$ （焦点距離）に画像平面があると考えます。
uvのピクセルは、カメラから見て其の平面上の点 $u_x\mathbf{r} + u_y\mathbf{u} + F\mathbf{f}$ で、レイは其処を通ります：

\mathbf{r}_d = \operatorname{normalize}\left(u_x\,\mathbf{r} + u_y\,\mathbf{u} + F\,\mathbf{f}\right)

上から見た図：1つのピクセルのレイはカメラから画像平面上の其のピクセルを通る

画面の上端は $u_y = 1$ なので、縦の視野角の半分は $\arctan(1/F)$ です。
$F = 1.8$ なら $\arctan(0.5556) = 29.05^\circ$ で、縦の視野角は約58°です。
$F$ を大きくするとズームインし、小さくすると広角になります。

プレイグラウンドを此れに置き換えてRキーを押して下さい：

#version 450 core

in vec2 vUv;
out vec4 out_color;

uniform float uTime;       // 起動してからの秒数
uniform vec2 uResolution;  // 画面の大きさ（ピクセル）
uniform vec2 uMouse;       // マウスの位置（ピクセル、左下が原点）

void main() {
  // 中心が原点、縦が-1〜1の座標（第5章）
  vec2 uv = (2.0 * gl_FragCoord.xy - uResolution) / uResolution.y;

  vec3 ro = vec3(0.0, 2.2, 5.0); // レイの原点＝カメラの位置
  vec3 target = vec3(0.0, 0.6, 0.0);

  // lookAtと同じ正規直交基底（線型代数の第11章）
  vec3 forward = normalize(target - ro);
  vec3 right = normalize(cross(forward, vec3(0.0, 1.0, 0.0)));
  vec3 up = cross(right, forward);
  const float FOCAL = 1.8; // 焦点距離：大きい程視野が狭い
  vec3 rd = normalize(uv.x * right + uv.y * up + FOCAL * forward);

  out_color = vec4(rd * 0.5 + 0.5, 1.0);
}

レイの方向を色で表示：赤は右へ、緑は上へ行く程強くなる

何時もの様に、先ず色でデバッグします。
赤は右へ（レイの $x$ ）、緑は上へ行く程強くなり、全てのレイは殆ど−Zを向いているので青の部分は何処でも弱いです。
画像が左右反転していたら、外積の順番が逆です。

14.3 スフィアトレーシング

3Dの距離関数

符号付き距離関数は、任意の点 $\mathbf{p}$ について一番近い表面までの距離を返します。
外側では正、内側では負、丁度表面上では0です。
球は第6章の円に座標を1つ足しただけです：

d_{\text{sphere}}(\mathbf{p}) = \lVert \mathbf{p} - \mathbf{c} \rVert - r

float sdSphere(vec3 p, float r) {
  return length(p) - r;
}

（球を別の場所に置くには、第7章と同じく、代わりに点を動かします。
sdSphere(p - center, r)です。）

ループ

ある点のSDFの値は、何も無い事が保証された球の半径です。
シーンの中に其れより近い物はありません。
従って、レイがどの方向を向いていても、レイに沿って丁度其の距離だけ進むのは常に安全です。
距離が非常に小さくなる（何かに当たった）か、 $t$ が非常に大きくなる（外れた）まで繰り返します：

t_0 = 0, \qquad t_{i+1} = t_i + d\left(\mathbf{r}_o + t_i\,\mathbf{r}_d\right)

スフィアトレーシング：各点で空の球がどれ位進めるかを教えてくれる

全てから遠い所では歩幅が大きく、表面の近くでは小さくなり、点は表面に収束します。
箱の近くを通るレイは、当たらないのに小さな歩幅を何回か取らなければならない事に注目して下さい。

計算例

$\mathbf{r}_o = (0, 0, 0)$ から $\mathbf{r}_d = (0, 0, -1)$ 方向のレイと、半径1の2つの球を考えます。
$A$ はレイの横の $(2, 0, -3)$ 、 $B$ は真正面の $(0, 0, -8)$ です。
シーンの距離は2つの小さい方 $d = \min(d_A, d_B)$ で、 $\mathbf{p} = (0, 0, -t)$ です：

ステップ	$t$	$d_A = \sqrt{4 + (3 - t)^2} - 1$	$d_B = (8 - t) - 1$	$d$
0	0.000	2.606	7.000	2.606
1	2.606	1.039	4.394	1.039
2	3.644	1.101	3.356	1.101
3	4.745	1.654	2.255	1.654
4	6.400	2.944	0.600	0.600
5	7.000	3.472	0.000	当たり

例えばステップ1は $\sqrt{4 + (3 - 2.606)^2} - 1 = \sqrt{4.155} - 1 = 2.038 - 1 = 1.039$ です。
レイが球 $A$ の横を通る間は、 $A$ が近いので歩幅が短く（約1）、通り過ぎると歩幅は再び大きくなり、最後の1歩が丁度 $t = 7$ で $B$ の表面に着きます。
正しく $8 - 1$ です。

コード

const int MAX_STEPS = 128;   // レイ1本あたりの最大ステップ数
const float MAX_DIST = 40.0; // 此れより遠くは「何もない」
const float SURF_DIST = 0.001;

float sdSphere(vec3 p, float r) {
  return length(p) - r;
}

// シーン全体：一番近い表面までの距離
float map(vec3 p) {
  return sdSphere(p - vec3(0.0, 0.7, 0.0), 0.6);
}

// 戻り値：当たった距離t（何にも当たらなかったら-1）
float raymarch(vec3 ro, vec3 rd) {
  float t = 0.0;
  for (int i = 0; i < MAX_STEPS; i++) {
    float d = map(ro + rd * t);
    if (d < SURF_DIST * t) return t; // 十分に近い：当たり
    t += d;                          // 安全な距離だけ進む
    if (t > MAX_DIST) break;         // 遠過ぎ：外れ
  }
  return -1.0;
}

そしてmainの最後で、レイの色の代わりに：

  float t = raymarch(ro, rd);
  vec3 color = t > 0.0 ? vec3(1.0) : vec3(0.0);
  out_color = vec4(color, 1.0);

最初のレイマーチングの物体：白い円、球のシルエット

白い円だけですが、本物の3Dの球です。
roを動かしても丸いままで、どの方向から見ても球です。

mapは殆ど何処でも使われる約束事です。
シーン全体を距離場として表す1つの関数です。
ループについて3つの点：

当たりの判定はd < SURF_DISTではなくd < SURF_DIST * tです。遠くのピクセルは近くのピクセルよりシーンの広い範囲を覆うので、許容誤差を距離と共に大きくします。遠くで沢山のステップを節約出来、見た目の損はありません。
当たり以外に必ず2つの出口があります。MAX_STEPS（レイが表面に沿って這い、届かない）とMAX_DIST（空へ飛んで行った）です。此れが無いと、あるピクセルが非常に長い間ループし続ける可能性があります。
GLSLのループに特別な構文は要りませんが、隣のピクセルと繰り返し回数が違うループは、お互いを待たせます（14.9で詳しく）。

コストを見る

床を足しましょう。
平面 $y = 0$ で、其のSDFは単に高さ $p_y$ です（上にある点から平面上の一番近い点は真下です）。
2つのSDFのminは和集合です（第6章）：

// シーン全体：一番近い表面までの距離
float map(vec3 p) {
  float sphere = sdSphere(p - vec3(0.0, 0.7, 0.0), 0.6);
  return min(sphere, p.y); // 球と床（y = 0の平面）
}

便利なデバッグ表示は、表面を描く代わりにステップ数を数える事です：

// 何ステップで終わったかを数える（デバッグ用）
int countSteps(vec3 ro, vec3 rd) {
  float t = 0.0;
  for (int i = 0; i < MAX_STEPS; i++) {
    float d = map(ro + rd * t);
    if (d < SURF_DIST * t || t > MAX_DIST) return i;
    t += d;
  }
  return MAX_STEPS;
}

  float s = float(countSteps(ro, rd)) / 64.0;
  vec3 color = mix(vec3(0.0, 0.0, 0.3), vec3(1.0, 0.9, 0.2), clamp(s, 0.0, 1.0));
  out_color = vec4(color, 1.0);

ステップ数：濃い青＝少ない、黄色＝64以上。球の縁と地平線が高い。

床や球の真ん中に真っ直ぐ当たるレイは数ステップで済みます。
表面を掠めるレイは高く付きます。
球の縁が光り、地平線も光ります。
地平線ではレイが床と殆ど平行に飛び、床までの距離が長い間小さいままだからです。
レイマーチャーが遅い時、最初に見るべき画像です。

14.4 勾配から法線を求める

ライティングには法線が必要ですが、保存するメッシュがありません。
幸い、SDFが知っています。
距離は表面から真っ直ぐ離れる方向に一番速く増え、一番速く増える方向は勾配です：

\mathbf{n} = \operatorname{normalize}\left(\nabla d\right) = \operatorname{normalize}\left(\frac{\partial d}{\partial x}, \frac{\partial d}{\partial y}, \frac{\partial d}{\partial z}\right)

微分を式として持っていないので、中心差分で測ります。
各軸に沿って小さな距離 $\varepsilon$ だけ前と後ろでSDFを評価します。

\frac{\partial d}{\partial x} \approx \frac{d(\mathbf{p} + \varepsilon\mathbf{i}) - d(\mathbf{p} - \varepsilon\mathbf{i})}{2\varepsilon}

$2\varepsilon$ での割り算は3つの成分で同じで、どうせ正規化するので省略出来ます：

// 法線＝距離場の勾配（中心差分）
vec3 calcNormal(vec3 p) {
  const vec2 e = vec2(0.001, 0.0);
  return normalize(vec3(
    map(p + e.xyy) - map(p - e.xyy),
    map(p + e.yxy) - map(p - e.yxy),
    map(p + e.yyx) - map(p - e.yyx)));
}

スウィズルは第2章の便利な技です。
e = (0.001, 0)なら、e.xyyは $(0.001, 0, 0)$ 、e.yxyは $(0, 0.001, 0)$ 、e.yyxは $(0, 0, 0.001)$ です。

計算例：球の勾配は中心からの方向 $(\mathbf{p} - \mathbf{c}) / \lVert \mathbf{p} - \mathbf{c} \rVert$ になるはずです。
あたし達の球（中心 $(0, 0.7, 0)$ 、半径0.6）で、xy平面上の頂点から0.5ラジアンの点 $\mathbf{p} = (0.6\sin 0.5,\; 0.7 + 0.6\cos 0.5,\; 0) = (0.2877, 1.2265, 0)$ を取ります。
正確な法線は $(\sin 0.5, \cos 0.5, 0) = (0.4794, 0.8776, 0)$ で、中心差分は小数点以下4桁で $(0.4794, 0.8776, 0.0000)$ になります。

法線を色で表示します（床で画面が埋まらない様に、mapは再びsdSphereだけにします）：

  vec3 color = vec3(0.0);
  float t = raymarch(ro, rd);
  if (t > 0.0) {
    vec3 p = ro + rd * t;
    vec3 n = calcNormal(p);
    color = n * 0.5 + 0.5;
  }
  out_color = vec4(color, 1.0);

法線で色を付けた球

代償はピクセル毎にmapの呼び出し6回です。
複雑なシーンでは目立つので、課題にもっと安い4回版があります。

14.5 マテリアルとライティング

物体毎に違う色を付ける為に、mapは2つの数を返す様にします。
距離とマテリアル番号です。
vec2が丁度良く、minの代わりに近い方の結果を残します：

// シーン全体：x＝距離、y＝マテリアル番号
vec2 map(vec3 p) {
  // 床（y = 0の平面）
  vec2 res = vec2(p.y, 0.0);

  float sphere = sdSphere(p - vec3(0.0, 0.7, 0.0), 0.6);
  if (sphere < res.x) res = vec2(sphere, 1.0);

  return res;
}

raymarchも番号を返し（−1は「外れ」）、calcNormalは距離map(...).xだけを使います：

// 戻り値：x＝当たった距離t、y＝マテリアル番号（何にも当たらなかったら-1）
vec2 raymarch(vec3 ro, vec3 rd) {
  float t = 0.0;
  for (int i = 0; i < MAX_STEPS; i++) {
    vec2 h = map(ro + rd * t);
    if (h.x < SURF_DIST * t) return vec2(t, h.y); // 十分に近い：当たり
    t += h.x;                                     // 安全な距離だけ進む
    if (t > MAX_DIST) break;                      // 遠過ぎ：外れ
  }
  return vec2(t, -1.0);
}

// 法線＝距離場の勾配（中心差分）
vec3 calcNormal(vec3 p) {
  const vec2 e = vec2(0.001, 0.0);
  return normalize(vec3(
    map(p + e.xyy).x - map(p - e.xyy).x,
    map(p + e.yxy).x - map(p - e.yxy).x,
    map(p + e.yyx).x - map(p - e.yyx).x));
}

床にはワールド空間のチェッカー模様を付けます（第7章のfloorとmodを、当たった点の $x$ と $z$ に）：

vec3 materialColor(float id, vec3 p) {
  if (id < 0.5) {
    float checker = mod(floor(p.x) + floor(p.z), 2.0);
    return mix(vec3(0.08), vec3(0.35), checker);
  }
  return vec3(0.9, 0.25, 0.35);
}

（id == 1.0ではなくid < 0.5と比べるのは安全の為の習慣です。
番号はfloatで、floatの完全一致の比較は脆いです。）
ライティングは第9章のBlinn-Phongモデルに、第10章の空の色の半球アンビエントを加えた物です。
表面からカメラへの視線方向は単に $-\mathbf{r}_d$ なので、ハーフベクトルは $\operatorname{normalize}(\mathbf{l} - \mathbf{r}_d)$ です：

  vec3 skyColor = vec3(0.45, 0.6, 0.85) - 0.4 * rd.y;
  vec3 color = skyColor;

  vec2 hit = raymarch(ro, rd);
  if (hit.y >= 0.0) {
    vec3 p = ro + rd * hit.x;
    vec3 n = calcNormal(p);
    vec3 albedo = materialColor(hit.y, p);

    vec3 l = normalize(vec3(0.8, 0.6, 0.1)); // 光の方へ向かう方向
    vec3 h = normalize(l - rd);              // ハーフベクトル（視線は-rd）
    float diffuse = max(dot(n, l), 0.0);
    float specular = pow(max(dot(n, h), 0.0), 32.0) * diffuse;
    float sky = 0.5 + 0.5 * n.y;

    vec3 lit = albedo * vec3(1.3, 1.2, 1.0) * diffuse;
    lit += albedo * vec3(0.12, 0.15, 0.22) * sky;
    lit += vec3(0.8) * specular;
    color = lit;
  }

  out_color = vec4(pow(color, vec3(1.0 / 2.2)), 1.0); // ガンマ補正

空の色は上へ行く程暗くなり（- 0.4 * rd.y）、画像全体は第10章と同じくガンマ補正します。

チェッカー模様の床の上の照らされた球

14.6 形を増やして組み合わせる

箱

3Dの箱は第6章の2Dの箱に軸を1つ足した物です。
原点が中心で半分の大きさが $\mathbf{b}$ の箱について、 $\mathbf{q} = \lvert\mathbf{p}\rvert - \mathbf{b}$ （成分毎）とします。
正の成分は其の軸に沿って箱の外にいるという意味です。
外側では距離は正の部分の長さで、内側では一番近い面までの（負の）距離、つまり $\mathbf{q}$ の一番大きい成分です：

d_{\text{box}}(\mathbf{p}) = \lVert \max(\mathbf{q}, \mathbf{0}) \rVert + \min\left(\max(q_x, q_y, q_z), 0\right)

float sdBox(vec3 p, vec3 b) {
  vec3 q = abs(p) - b;
  return length(max(q, 0.0)) + min(max(q.x, max(q.y, q.z)), 0.0);
}

$\mathbf{b} = (1, 1, 1)$ で3つ確かめます：

$\mathbf{p} = (1.5, 0.2, 0)$ ： $\mathbf{q} = (0.5, -0.8, -1)$ なので $\lVert(0.5, 0, 0)\rVert + 0 = 0.5$ 。+Xの面の0.5単位前です。
$\mathbf{p} = (2, 2, 0)$ ： $\mathbf{q} = (1, 1, -1)$ なので $\lVert(1, 1, 0)\rVert = \sqrt{2} = 1.414$ 。辺までの距離です。
$\mathbf{p} = (0.5, 0, 0)$ ： $\mathbf{q} = (-0.5, -1, -1)$ なので $0 + \min(-0.5, 0) = -0.5$ 。内側で、+Xの面から0.5単位です。

どんなSDFでも定数を引くと、形が其の分だけ大きくなり角が丸くなるので（第6章）、sdBox(q, vec3(0.45)) - 0.08は角の丸い箱です。

トーラス

トーラスは、半径 $r$ の円をxz平面上の半径 $R$ の輪に沿って回した物です。
先ず、y軸と $\mathbf{p}$ を含む平面の中で、 $\mathbf{p}$ から輪までの距離を求めます。
横方向は $\sqrt{p_x^2 + p_z^2} - R$ 、縦方向は $p_y$ です。
するとトーラスは其の2D平面での円のSDFになります：

d_{\text{torus}}(\mathbf{p}) = \left\lVert \left(\sqrt{p_x^2 + p_z^2} - R,\; p_y\right) \right\rVert - r

float sdTorus(vec3 p, float R, float r) {
  vec2 q = vec2(length(p.xz) - R, p.y);
  return length(q) - r;
}

例： $R = 0.5$ 、 $r = 0.18$ 、 $\mathbf{p} = (1, 0.3, 0)$ なら $\mathbf{q} = (0.5, 0.3)$ で、 $\sqrt{0.34} - 0.18 = 0.583 - 0.18 = 0.403$ です。

組み合わせと変換

第6章と第7章の道具は3Dでも其のまま使えます：

// 滑らかな和集合（第6章のsminと同じ）
float smin(float a, float b, float k) {
  float h = clamp(0.5 + 0.5 * (b - a) / k, 0.0, 1.0);
  return mix(b, a, h) - k * h * (1.0 - h);
}

mat2 rot(float a) {
  float c = cos(a), s = sin(a);
  return mat2(c, s, -s, c);
}

mat2(c, s, -s, c)は行列を列毎に埋めるので、列は $(\cos a, \sin a)$ と $(-\sin a, \cos a)$ です。
線型代数の第8章の反時計回りの回転です。
3Dの物体をy軸周りに回すには、点のxz成分を回します。
q.xz = rot(a) * q.xzです。
第7章と同じく、点を変換すると形は逆向きに動きますが、アニメーションなら全く関係ありません。

新しいmap：

// シーン全体：x＝距離、y＝マテリアル番号
vec2 map(vec3 p) {
  // 床（y = 0の平面）
  vec2 res = vec2(p.y, 0.0);

  // 中央：二つの球を滑らかにくっつけた物
  float bob = 0.35 * sin(uTime * 1.5);
  float a = sdSphere(p - vec3(0.0, 0.7, 0.0), 0.6);
  float b = sdSphere(p - vec3(0.45, 1.3 + bob, 0.1), 0.35);
  float blob = smin(a, b, 0.3);
  if (blob < res.x) res = vec2(blob, 1.0);

  // 左：回転する角の丸い箱
  vec3 q = p - vec3(-1.9, 0.6, 0.3);
  q.xz = rot(uTime * 0.6) * q.xz;
  float box = sdBox(q, vec3(0.45)) - 0.08;
  if (box < res.x) res = vec2(box, 2.0);

  // 右：傾いたトーラス
  q = p - vec3(1.9, 0.75, 0.0);
  q.yz = rot(0.9) * q.yz;
  float torus = sdTorus(q, 0.5, 0.18);
  if (torus < res.x) res = vec2(torus, 3.0);

  return res;
}

materialColorには色を2つ追加します：

  if (id < 1.5) return vec3(0.9, 0.25, 0.35);
  if (id < 2.5) return vec3(0.95, 0.55, 0.1);
  return vec3(0.1, 0.7, 0.6);

3つの形：角の丸い箱、滑らかにくっついた2つの球、トーラス

小さい球は上下に揺れて大きい球に溶け込み、箱は回転します。
sminのkを変えて、球の間の首が太くなったり細くなったりするのを見てみて下さい。

正しさについて：スフィアトレーシングには、本当の距離より大きくならない距離が必要です。さもないと、1歩で薄い物体を通り抜けてしまう可能性があります。sminと回転は大丈夫です。
点を不均一に拡大縮小したり捻じったりすると、大抵は大丈夫ではありません。其の場合は距離に0.5の様な安全係数を掛けて下さい（小さく安全な歩幅になります）。

14.7 影

硬い影

此処でレイマーチングが輝きます。
影は唯のもう1本のレイです。
当たった点から光の方へマーチし、其のレイが何かに当たれば、点は影の中です。

// 硬い影：光の方へレイを飛ばし、何かに当たったら0
float hardShadow(vec3 ro, vec3 rd) {
  float t = 0.02;
  for (int i = 0; i < 64; i++) {
    float h = map(ro + rd * t).x;
    if (h < 0.0005) return 0.0;
    t += h;
    if (t > 12.0) break;
  }
  return 1.0;
}

2つの細かい点が大事です。
影のレイは表面から少し離れたp + n * 0.002から、t = 0.02で始めます。
当たった点其の物では距離が約0なので、レイが直ぐに自分の表面に「当たって」しまうからです（此のアーティファクトはどのレンダラーでもシャドウアクネと呼ばれます）。
そして影は直接光だけを暗くし、アンビエント（空）の部分は暗くしません：

    float shadow = hardShadow(p + n * 0.002, l);
    float specular = pow(max(dot(n, h), 0.0), 32.0) * diffuse;
    float sky = 0.5 + 0.5 * n.y;

    vec3 lit = albedo * vec3(1.3, 1.2, 1.0) * diffuse * shadow;
    lit += albedo * vec3(0.12, 0.15, 0.22) * sky;
    lit += vec3(0.8) * specular * shadow;

殆どタダの柔らかい影

本物の光には大きさがあるので、影には柔らかい縁、半影があります。
SDFは安くて有名な近似を与えてくれます（Inigo Quilezが広めた物です）。
光の方へマーチする間、各ステップで進んだ距離 $t$ と一番近い物体までの距離 $h$ が分かります。
比 $h/t$ は、点から見て隙間がどれ程広いか、大体レイが遮る物を外れる角度を表します。

影のレイはtだけ進んだ所で、遮る物から距離hの所を通る

比が小さいと光は殆ど遮られているので、レイ全体で一番小さい比を取ります：

s = \operatorname{clamp}\left(\min_i \frac{k\,h_i}{t_i},\; 0,\; 1\right)

// 柔らかい影：光へ向かうレイが物体の近くを通る程暗くする
float softShadow(vec3 ro, vec3 rd, float k) {
  float res = 1.0;
  float t = 0.02;
  for (int i = 0; i < 64; i++) {
    float h = map(ro + rd * t).x;
    if (h < 0.0005) return 0.0;
    res = min(res, k * h / t);
    t += clamp(h, 0.02, 0.4);
    if (t > 12.0) break;
  }
  return clamp(res, 0.0, 1.0);
}

$k$ は硬さを決めます。
大きい $k$ （32等）は鋭い影、小さい $k$ （4等）はとても柔らかい影になり、此処では8を使います。
歩幅は0.02〜0.4に制限します。
最小値はループが這う様に遅くなるのを防ぎ、最大値はレイが物体に一番近付く点を飛び越すのを防ぎます。
飛び越すと半影の一番暗い部分を見逃してしまいます。

左：硬い影。右：k = 8の柔らかい影。

左では全ての影の縁が同じ様に鋭いです。
右では影が箱に接する所の縁は鋭く、箱から離れる程柔らかくなり、本物の影と同じです。

14.8 アンビエントオクルージョン、霧、繰り返し

アンビエントオクルージョン

角や物体の真下の地面は、物体が空の一部を遮るので空からの光が少なくなります。
其れがアンビエントオクルージョン（AO）で、SDFなら法線に沿った数回のサンプルで近似出来ます。
開けた平らな表面の高さ $h$ の所では距離は丁度 $h$ のはずです。
其れより小さければ、近くに何かがあります：

\text{occ} = \sum_{i=1}^{5} w_i\left(h_i - d(\mathbf{p} + h_i\mathbf{n})\right), \qquad \text{ao} = \operatorname{clamp}(1 - 2\,\text{occ},\; 0,\; 1)

$h_i = 0.03 + 0.12i$ 、重みは $w_i = 0.6^{i-1}$ なので、近いサンプル程大きく数えます。

// アンビエントオクルージョン：法線方向に何回か距離場を調べる
float ambientOcclusion(vec3 p, vec3 n) {
  float occ = 0.0;
  float weight = 1.0;
  for (int i = 1; i <= 5; i++) {
    float h = 0.03 + 0.12 * float(i); // 表面からの距離
    float d = map(p + n * h).x;       // 実際に空いている距離
    occ += (h - d) * weight;
    weight *= 0.6;
  }
  return clamp(1.0 - 2.0 * occ, 0.0, 1.0);
}

AOはアンビエントの項だけに掛けます（lit += albedo * vec3(0.12, 0.15, 0.22) * sky * ao;）。
AOだけを色として表示すると：

アンビエントオクルージョンだけ：白＝開けている、暗い＝囲まれている。物体の下の接触部分が一番暗い。

照らした画像での効果は控えめですが重要です。
物体を床に「くっつけて」くれます。

modで無限に繰り返す

第7章のmodは空間を繰り返します。
x座標を幅2.5で0が中心の1つのセルに折り畳むと、其処で評価するSDFは全てのセルに1回ずつ、永遠に、1回分の値段で現れます：

  // 奥：modで無限に並ぶ柱
  q = p - vec3(0.0, 0.0, -8.0);
  q.x = mod(q.x + 1.25, 2.5) - 1.25;
  float pillar = sdBox(q, vec3(0.25, 2.5, 0.25));
  if (pillar < res.x) res = vec2(pillar, 4.0);

GLSLのmod(x, y)は $x - y\lfloor x/y \rfloor$ で、正の $y$ なら負の $x$ でも常に0〜 $y$ になるので（Cのfmodとは違います）、セルは左側にも綺麗に続きます。
+ 1.25と- 1.25はセルを0〜2.5ではなく−1.25〜1.25にずらし、柱をセルの中心に置きます。
形がセルの中に収まっている間だけ正確です。
幅2.5のセルに幅0.5の柱なら安全です。
柱はマテリアル4の灰色です。
if (id < 3.5) return vec3(0.1, 0.7, 0.6);の後にreturn vec3(0.55, 0.55, 0.6);です。

霧

無限の柱には、無限を隠す物が必要です。
指数関数の霧は、距離と共に空の色へ混ぜます：

f = 1 - e^{-\rho t^2}

$\rho = 0.004$ で $t = 5$ なら $1 - e^{-0.1} = 0.095$ で殆ど見えず、 $t = 20$ なら $1 - e^{-1.6} = 0.80$ で殆ど空です。
$t$ ではなく $t^2$ を使うと、近くの物体ははっきりしたまま、霧が急に迫って来ます。
exp(x)は $e^x$ です。

    // 霧：遠い程空の色に近づける
    float fog = 1.0 - exp(-0.004 * hit.x * hit.x);
    color = mix(lit, skyColor, fog);

周りを回るカメラ

最後に、カメラはシーンの周りを回り、マウスを左右に動かすと更に回ります：

  // カメラ：原点の周りを回る（マウスを左右に動かすと角度がずれる）
  float angle = 0.25 * uTime + 6.2831853 * (uMouse.x / uResolution.x - 0.5);
  vec3 ro = vec3(5.0 * sin(angle), 2.2, 5.0 * cos(angle)); // レイの原点＝カメラの位置

$(5\sin\alpha, 2.2, 5\cos\alpha)$ は高さ2.2の半径5の円で、基底は毎フレームroとtargetから作り直すので、カメラは常に中心を見ます。
マウスがウィンドウの真ん中なら追加の角度は0で、左端から右端までで丁度1回転（ $2\pi$ ）です。

完成したplayground.glsl全体は章の最後のDIFFにあります。
其のmainは全てを組み合わせます：

    float diffuse = max(dot(n, l), 0.0);
    float shadow = softShadow(p + n * 0.002, l, 8.0);
    float specular = pow(max(dot(n, h), 0.0), 32.0) * diffuse;
    float ao = ambientOcclusion(p, n);
    float sky = 0.5 + 0.5 * n.y;

    vec3 lit = albedo * vec3(1.3, 1.2, 1.0) * diffuse * shadow;
    lit += albedo * vec3(0.12, 0.15, 0.22) * sky * ao;
    lit += vec3(0.8) * specular * shadow;

14.9 パフォーマンス

レイマーチングは簡単ですが安くはありません。
1フレームのコストは大体：

\text{pixels} \times \text{steps} \times \text{cost of } \operatorname{map}

あたし達のピクセルは、マーチ1回（最大128ステップ）、法線の為のmap6回、影の最大64回、AOの5回を行います。
800×600では、1フレームで何千万回ものmapの呼び出しになり得ます。
効果がある物を、大体重要な順に：

何も起きない所でステップを減らす：適切なMAX_DIST、相対的な当たりの許容誤差、そして14.3のステップ数の表示で高いピクセルを探す事。
mapを安くする：複雑な物体の周りの安い外接形状（球や箱）を先に調べ、点が近い時だけ高いSDFを評価する。
追加のレイを減らす：ステップ数と最大距離を減らした影のレイ、6回ではなく4回の法線。
解像度を下げる：多くのデモは半分の解像度でレンダリングして拡大します。
分岐（ダイバージェンス）：GPUはピクセルをグループ（32又は64個ずつ）で実行し、グループは同じ命令を一緒に実行します。グループ内の1つのピクセルが120ステップ必要で他が10ステップなら、全員が遅いピクセルを待ちます。隣のレイ同士を似た物に保つと（地平線にノイズの多いシーンを置かない）、GPUを忙しく保てます。

macOS（OpenGL 4.1）でも此の章は全て其のまま動きます。
プレイグラウンドの機能しか使っていません。
統合GPUでは、大きなウィンドウだと完成したシーンが遅いかもしれません。
ウィンドウを小さくするか、MAX_STEPSを64に下げて下さい。

14.10 此の章の新しいGLSL

此の章で新しいGLSLは殆どありません。
既に知っている物の新しい使い方です。

機能	内容
`break`と`return`付きの`for`ループ	マーチングのループ。必ず最大回数を付ける
`exp(x)`	$e^x$ 、此処では指数関数の霧に
座標に`mod(x, y)`	$x - y\lfloor x/y \rfloor$ ：空間の無限の繰り返し
`mat2(c, s, -s, c)`	列毎のコンストラクター：2Dの回転
`q.xz = rot(a) * q.xz`	スウィズルの書き込みマスク：`vec3`の2成分を回す
`e.xyy`、`e.yxy`、`e.yyx`	3つの軸のオフセットを作るスウィズル
`vec2`の戻り値	距離とマテリアル番号を一緒に返す

14.11 実際のゲームやエンジンでの使い方

Shadertoyとデモシーン：数キロバイトのシェーダーコードで世界全体を作り、殆どがSDFのレイマーチングです。
エンジンの距離場の影とAO：Unreal Engineはメッシュの符号付き距離場を作り、遠くまでの柔らかい影とアンビエントオクルージョンに使います。此処と同じmin(k h / t)と法線方向のサンプルの考え方です。
SDFベースのゲーム：Media MoleculeのDreamsではプレイヤーがSDFの形で造形し、Claybookは変形する粘土の世界を距離場のトレースで描きました。
ボリュームの効果：雲、煙、霧、光芒は、密度場の中を一定の歩幅でレイマーチして描きます（スフィアトレーシングの近い親戚です）。
スクリーンスペースの反射とアンビエントオクルージョンは、SDFではなく深度バッファの中をマーチします。
衝突と物理：SDFは「何かからどれ位離れているか？」に1回の呼び出しで答えるので、パーティクルの衝突やキャラクターコントローラーに使われます。

14.12 まとめ

レイマーチングはフラグメントシェーダーで3Dのシーンを描きます。ピクセル毎に1本のレイ、 $\mathbf{p}(t) = \mathbf{r}_o + t\,\mathbf{r}_d$ 。
カメラのレイはlookAtの基底から作ります。 $\mathbf{r}_d = \operatorname{normalize}(u_x\mathbf{r} + u_y\mathbf{u} + F\mathbf{f})$ で、焦点距離 $F$ が視野角を決めます。
スフィアトレーシングはSDFの値だけ進み、此れは常に安全です。当たり（d < SURF_DIST * t）、MAX_DIST、MAX_STEPSで止めます。
map(p)がシーン全体を表します。min（又は比較）で物体を組み合わせ、sminで混ぜ、点を変換して動かし、modで繰り返します。
法線はSDFの勾配を正規化した物で、中心差分で測ります。
影は光の方へのレイで、途中のmin(k h / t)で柔らかい影になります。AOは法線に沿ってSDFをサンプルします。霧は $1 - e^{-\rho t^2}$ で空の色へ混ぜます。
コストはピクセル × ステップ × mapです。ステップ数の表示で時間が何処に掛かっているか分かり、掠めるレイが一番高く付きます。

14.13 課題

★ sdBoxの式で、 $\mathbf{b} = (1, 0.5, 2)$ として $\mathbf{p} = (3, 0, 0)$ 、 $(2, 1.5, 0)$ 、 $(0, 0, 1)$ の距離を計算して下さい（2、 $\sqrt{2}$ 、−0.5になるはずです）。
★ 14.3の計算例を、球 $B$ を $(0, 0, -6)$ に動かしてやり直して下さい。何ステップ掛かり、どの $t$ で当たりますか（丁度 $t = 5$ になるはずです）？
★★ カプセル（半径付きの線分、第6章の線分のSDFの3D版： $\mathbf{p} - \mathbf{a}$ を $\mathbf{b} - \mathbf{a}$ に射影し、0〜1に制限し、其の点までの距離から半径を引く）を追加し、5番目のマテリアルとして床に寝かせて下さい。
★★ calcNormalを4回の「四面体」版に置き換えて下さい。4つのベクトルk = $(1,-1,-1)$ 、 $(-1,-1,1)$ 、 $(-1,1,-1)$ 、 $(1,1,1)$ についてk * map(p + k * 0.001).xを足し、正規化します。mapの呼び出しが3分の1減っても、画像は同じに見えるはずです。
★★ カメラがマウスの両方向に付いて来る様にして下さい。uMouse.yでカメラの高さを0.5〜5の間で動かします。床がカメラより上に消えない様に、高さを制限して下さい。
★★ q.xzにmodを使って柱を $x$ と $z$ の両方に繰り返し（柱の森）、セル番号floor(...)で柱毎に高さを変えて下さい。柱がセルの境界に近付く場合は、距離を小さくするのを忘れずに。
★★★ 反射を追加して下さい。レイが床に当たったら、当たった点からreflect(rd, n)の方向（reflect(i, n) $= \mathbf{i} - 2(\mathbf{n} \cdot \mathbf{i})\mathbf{n}$ ）に2本目のレイをマーチし、其の照らした色の30%を床に混ぜます。チェッカー模様に物体が映って見えるはずです。

DIFFをダウンロード

ファイル： a/shader/playground.glsl | b/shader/playground.glsl

 uniform vec2 uResolution;  // 画面の大きさ（ピクセル）
 uniform vec2 uMouse;       // マウスの位置（ピクセル、左下が原点）
+const int MAX_STEPS = 128;   // レイ1本あたりの最大ステップ数
+const float MAX_DIST = 40.0; // 此れより遠くは「何もない」
+const float SURF_DIST = 0.001;
+//------------------------------------------
+// 距離関数（SDF）
+//------------------------------------------
+float sdSphere(vec3 p, float r) {
+  return length(p) - r;
+}
+float sdBox(vec3 p, vec3 b) {
+  vec3 q = abs(p) - b;
+  return length(max(q, 0.0)) + min(max(q.x, max(q.y, q.z)), 0.0);
+}
+float sdTorus(vec3 p, float R, float r) {
+  vec2 q = vec2(length(p.xz) - R, p.y);
+  return length(q) - r;
+}
+// 滑らかな和集合（第6章のsminと同じ）
+float smin(float a, float b, float k) {
+  float h = clamp(0.5 + 0.5 * (b - a) / k, 0.0, 1.0);
+  return mix(b, a, h) - k * h * (1.0 - h);
+}
+mat2 rot(float a) {
+  float c = cos(a), s = sin(a);
+  return mat2(c, s, -s, c);
+}
+// シーン全体：x＝距離、y＝マテリアル番号
+vec2 map(vec3 p) {
+  // 床（y = 0の平面）
+  vec2 res = vec2(p.y, 0.0);
+  // 中央：二つの球を滑らかにくっつけた物
+  float bob = 0.35 * sin(uTime * 1.5);
+  float a = sdSphere(p - vec3(0.0, 0.7, 0.0), 0.6);
+  float b = sdSphere(p - vec3(0.45, 1.3 + bob, 0.1), 0.35);
+  float blob = smin(a, b, 0.3);
+  if (blob < res.x) res = vec2(blob, 1.0);
+  // 左：回転する角の丸い箱
+  vec3 q = p - vec3(-1.9, 0.6, 0.3);
+  q.xz = rot(uTime * 0.6) * q.xz;
+  float box = sdBox(q, vec3(0.45)) - 0.08;
+  if (box < res.x) res = vec2(box, 2.0);
+  // 右：傾いたトーラス
+  q = p - vec3(1.9, 0.75, 0.0);
+  q.yz = rot(0.9) * q.yz;
+  float torus = sdTorus(q, 0.5, 0.18);
+  if (torus < res.x) res = vec2(torus, 3.0);
+  // 奥：modで無限に並ぶ柱
+  q = p - vec3(0.0, 0.0, -8.0);
+  q.x = mod(q.x + 1.25, 2.5) - 1.25;
+  float pillar = sdBox(q, vec3(0.25, 2.5, 0.25));
+  if (pillar < res.x) res = vec2(pillar, 4.0);
+  return res;
+}
+//------------------------------------------
+// レイマーチング
+//------------------------------------------
+// 戻り値：x＝当たった距離t、y＝マテリアル番号（何にも当たらなかったら-1）
+vec2 raymarch(vec3 ro, vec3 rd) {
+  float t = 0.0;
+  for (int i = 0; i < MAX_STEPS; i++) {
+    vec2 h = map(ro + rd * t);
+    if (h.x < SURF_DIST * t) return vec2(t, h.y); // 十分に近い：当たり
+    t += h.x;                                     // 安全な距離だけ進む
+    if (t > MAX_DIST) break;                      // 遠過ぎ：外れ
+  }
+  return vec2(t, -1.0);
+}
+// 法線＝距離場の勾配（中心差分）
+vec3 calcNormal(vec3 p) {
+  const vec2 e = vec2(0.001, 0.0);
+  return normalize(vec3(
+    map(p + e.xyy).x - map(p - e.xyy).x,
+    map(p + e.yxy).x - map(p - e.yxy).x,
+    map(p + e.yyx).x - map(p - e.yyx).x));
+}
+// 柔らかい影：光へ向かうレイが物体の近くを通る程暗くする
+float softShadow(vec3 ro, vec3 rd, float k) {
+  float res = 1.0;
+  float t = 0.02;
+  for (int i = 0; i < 64; i++) {
+    float h = map(ro + rd * t).x;
+    if (h < 0.0005) return 0.0;
+    res = min(res, k * h / t);
+    t += clamp(h, 0.02, 0.4);
+    if (t > 12.0) break;
+  }
+  return clamp(res, 0.0, 1.0);
+}
+// アンビエントオクルージョン：法線方向に何回か距離場を調べる
+float ambientOcclusion(vec3 p, vec3 n) {
+  float occ = 0.0;
+  float weight = 1.0;
+  for (int i = 1; i <= 5; i++) {
+    float h = 0.03 + 0.12 * float(i); // 表面からの距離
+    float d = map(p + n * h).x;       // 実際に空いている距離
+    occ += (h - d) * weight;
+    weight *= 0.6;
+  }
+  return clamp(1.0 - 2.0 * occ, 0.0, 1.0);
+}
+vec3 materialColor(float id, vec3 p) {
+  if (id < 0.5) {
+    float checker = mod(floor(p.x) + floor(p.z), 2.0);
+    return mix(vec3(0.08), vec3(0.35), checker);
+  }
+  if (id < 1.5) return vec3(0.9, 0.25, 0.35);
+  if (id < 2.5) return vec3(0.95, 0.55, 0.1);
+  if (id < 3.5) return vec3(0.1, 0.7, 0.6);
+  return vec3(0.55, 0.55, 0.6);
+}
 void main() {
-  // ピクセル座標を0〜1に正規化
+  // 中心が原点、縦が-1〜1の座標（第5章）
-  vec2 uv = gl_FragCoord.xy / uResolution;
+  vec2 uv = (2.0 * gl_FragCoord.xy - uResolution) / uResolution.y;
+  // カメラ：原点の周りを回る（マウスを左右に動かすと角度がずれる）
+  float angle = 0.25 * uTime + 6.2831853 * (uMouse.x / uResolution.x - 0.5);
+  vec3 ro = vec3(5.0 * sin(angle), 2.2, 5.0 * cos(angle)); // レイの原点＝カメラの位置
+  vec3 target = vec3(0.0, 0.6, 0.0);
+  // lookAtと同じ正規直交基底（線型代数の第11章）
+  vec3 forward = normalize(target - ro);
+  vec3 right = normalize(cross(forward, vec3(0.0, 1.0, 0.0)));
+  vec3 up = cross(right, forward);
+  const float FOCAL = 1.8; // 焦点距離：大きい程視野が狭い
+  vec3 rd = normalize(uv.x * right + uv.y * up + FOCAL * forward);
+  vec3 skyColor = vec3(0.45, 0.6, 0.85) - 0.4 * rd.y;
+  vec3 color = skyColor;
+  vec2 hit = raymarch(ro, rd);
+  if (hit.y >= 0.0) {
+    vec3 p = ro + rd * hit.x;
+    vec3 n = calcNormal(p);
+    vec3 albedo = materialColor(hit.y, p);
+    vec3 l = normalize(vec3(0.8, 0.6, 0.1)); // 光の方へ向かう方向
+    vec3 h = normalize(l - rd);              // ハーフベクトル（視線は-rd）
+    float diffuse = max(dot(n, l), 0.0);
+    float shadow = softShadow(p + n * 0.002, l, 8.0);
+    float specular = pow(max(dot(n, h), 0.0), 32.0) * diffuse;
+    float ao = ambientOcclusion(p, n);
+    float sky = 0.5 + 0.5 * n.y;
-  // 時間で動くグラデーション
+    vec3 lit = albedo * vec3(1.3, 1.2, 1.0) * diffuse * shadow;
-  vec3 color = 0.5 + 0.5 * cos(uTime + uv.xyx + vec3(0.0, 2.0, 4.0));
+    lit += albedo * vec3(0.12, 0.15, 0.22) * sky * ao;
+    lit += vec3(0.8) * specular * shadow;
-  // マウスの周りに白い円
+    // 霧：遠い程空の色に近づける
-  float d = length(gl_FragCoord.xy - uMouse);
+    float fog = 1.0 - exp(-0.004 * hit.x * hit.x);
-  color = mix(color, vec3(1.0), 1.0 - smoothstep(18.0, 20.0, d));
+    color = mix(lit, skyColor, fog);
+  }
-  out_color = vec4(color, 1.0);
+  out_color = vec4(pow(color, vec3(1.0 / 2.2)), 1.0); // ガンマ補正
 }

コメント（0）

まだコメントはありません。