Трасування променів з нуля — перетини, затінення та BVH

1. Огляд конвеєра трасування променів

Модель трасування променів Уіттеда (1980) визначає рекурсивний конвеєр. Кожен піксель генерує один первинний промінь; кожне влучання породжує щонайбільше два вторинні промені (тіньовий + дзеркальний/ заломлений). Справжні трасувальники шляхів (path tracers) випускають багато променів на піксель і усереднюють результати — ця стаття охоплює модель Уіттеда як концептуальну основу.

Результатом є 2D-масив значень кольору — відрендерений на <canvas> через ImageData. Продуктивність становить O(W × H × N × log N) на кадр, де N — кількість геометрії; BVH у розділі 7 зменшує це до O(W × H × log N).

2. Генерація первинних променів

Промінь — це параметричний напівпрямий: P(t) = origin + t × direction, t > 0. Налаштування камери перетворює кожен піксель (px, py) на промінь у світовому просторі за допомогою піраміди видимості (frustum):

function generateRay(px, py, cam) {
  const s_x = (2 * (px + 0.5) / cam.width  - 1) * cam.tanFov * cam.aspect;
  const s_y = (1 - 2 * (py + 0.5) / cam.height) * cam.tanFov;
  const dir = normalise(
    add(
      add(scale(cam.right,   s_x),
          scale(cam.up,      s_y)),
      cam.forward
    )
  );
  return { origin: cam.position, dir };
}

3. Перетин променя зі сферою

Перетин променя P(t) = o + t·d зі сферою (центр c, радіус r) знаходять підстановкою у рівняння сфери |P − c|² = r²:

function intersectSphere(ray, sphere) {
  const oc = sub(ray.origin, sphere.center);
  const a  = dot(ray.dir, ray.dir);
  const b  = 2 * dot(oc, ray.dir);
  const c  = dot(oc, oc) - sphere.radius * sphere.radius;
  const D  = b * b - 4 * a * c;
  if (D < 0) return null;
  const sqrtD = Math.sqrt(D);
  const t1 = (-b - sqrtD) / (2 * a);
  const t2 = (-b + sqrtD) / (2 * a);
  const t  = t1 > 0.001 ? t1 : t2 > 0.001 ? t2 : null;
  if (t === null) return null;
  const point  = add(ray.origin, scale(ray.dir, t));
  const normal = normalise(sub(point, sphere.center));
  return { t, point, normal, mat: sphere.mat };
}

4. Перетин променя з трикутником (Меллера-Трумбора)

Алгоритм Меллера-Трумбора (1997) безпосередньо обчислює барицентричні координати, не знаходячи спершу площину трикутника. Це стандартний алгоритм, який використовують у кожному промисловому трасувальнику променів.

function intersectTriangle(ray, v0, v1, v2) {
  const e1  = sub(v1, v0), e2 = sub(v2, v0);
  const h   = cross(ray.dir, e2);
  const a   = dot(e1, h);
  if (Math.abs(a) < 1e-7) return null;  // паралельно
  const f   = 1 / a;
  const s   = sub(ray.origin, v0);
  const u   = f * dot(s, h);
  if (u < 0 || u > 1) return null;
  const q   = cross(s, e1);
  const v   = f * dot(ray.dir, q);
  if (v < 0 || u + v > 1) return null;
  const t   = f * dot(e2, q);
  if (t < 0.001) return null;
  const point  = add(ray.origin, scale(ray.dir, t));
  const normal = normalise(cross(e1, e2));  // плоске затінення
  return { t, point, normal, u, v };
}

Для гладкого затінення інтерполюйте нормалі вершин за допомогою барицентричних координат: N = u·n1 + v·n2 + (1−u−v)·n0.

5. Затінення за Фонгом і тіні

У кожній точці влучання модель Фонга поєднує фонову (ambient), дифузну та дзеркальну (specular) складові від кожного джерела світла. Тіньовий промінь перевіряє затінення перед додаванням будь-якого прямого освітлення:

function shade(hit, ray, scene, depth) {
  const { point, normal, mat } = hit;
  let colour = scale(mat.albedo, mat.ka);  // фонове

  for (const light of scene.lights) {
    const toLight  = sub(light.pos, point);
    const lightDist = length(toLight);
    const L        = scale(toLight, 1 / lightDist);

    // Тіньовий промінь
    const shadowRay = { origin: add(point, scale(normal, 0.001)), dir: L };
    const shadowHit = closestHit(shadowRay, scene, lightDist);
    if (shadowHit) continue;

    const diff = Math.max(dot(normal, L), 0);
    const R    = sub(scale(normal, 2 * dot(normal, L)), L);
    const V    = scale(ray.dir, -1);
    const spec = Math.pow(Math.max(dot(R, V), 0), mat.shininess);

    colour = add(colour,
      add(scale(mulV(mat.albedo, light.colour), mat.kd * diff),
          scale(light.colour, mat.ks * spec)));
  }
  return colour;
}

6. Відбиття та заломлення

Ідеальне дзеркальне відбиття

function reflect(d, n) {
  return sub(d, scale(n, 2 * dot(d, n)));
}
// У shade(), якщо mat.reflectivity > 0 та depth < MAX_DEPTH:
const rDir = reflect(ray.dir, normal);
const rCol = traceRay({ origin: add(point, scale(normal, 0.001)), dir: rDir }, scene, depth + 1);
colour = add(colour, scale(rCol, mat.reflectivity));

Заломлення за законом Снеліуса

function refract(d, n, ni, nt) {
  const nr   = ni / nt;
  const cosI = -dot(d, n);
  const sinT2 = nr * nr * (1 - cosI * cosI);
  if (sinT2 > 1) return null;  // повне внутрішнє відбиття
  return add(scale(d, nr), scale(n, nr * cosI - Math.sqrt(1 - sinT2)));
}
// Френель (наближення Шліка):
function schlick(cosI, n1, n2) {
  let r0 = (n1 - n2) / (n1 + n2); r0 = r0 * r0;
  return r0 + (1 - r0) * Math.pow(1 - cosI, 5);
}

Для скла (n₂ ≈ 1.5) коефіцієнт Френеля визначає співвідношення між відбитим та заломленим променями. Наближення Шліка швидке й точне в межах 1% для більшості кутів.

7. Ієрархія обмежувальних об'ємів (BVH)

Наївна перевірка кожного променя проти кожного примітива має складність O(N) на промінь. BVH організовує примітиви у двійкове дерево вирівняних за осями обмежувальних паралелепіпедів (AABB). Кожен вузол зберігає AABB, що охоплює всіх його нащадків; під час обходу промінь, який не влучає в AABB вузла, може пропустити все піддерево.

Евристика площі поверхні (SAH) — це стандартна функція вартості для вибору способу розділення примітивів вузла:

class BVHNode {
  constructor(prims) {
    this.aabb = computeAABB(prims);
    if (prims.length <= 4) { this.prims = prims; return; }

    // Розділяємо вздовж найширшої осі за медіаною центроїдів
    const axis = widestAxis(this.aabb);
    prims.sort((a, b) => centroid(a)[axis] - centroid(b)[axis]);
    const mid  = prims.length >> 1;
    this.left  = new BVHNode(prims.slice(0, mid));
    this.right = new BVHNode(prims.slice(mid));
  }

  intersect(ray) {
    if (!intersectAABB(ray, this.aabb)) return null;
    if (this.prims) return closestPrim(ray, this.prims);
    const hL = this.left.intersect(ray);
    const hR = this.right.intersect(ray);
    if (!hL) return hR;
    if (!hR) return hL;
    return hL.t < hR.t ? hL : hR;
  }
}

Для сцени з 10 000 трикутників добре побудоване BVH зменшує середню вартість перетину променя з 10 000 тестів до приблизно 25–40 тестів (близько 8–9 рівнів обходу дерева).

8. Повний трасувальник променів на JavaScript

Повний цикл рендерингу зв'язує все докупи. Усю важку роботу (попіксельний цикл) можна винести у Web Worker, щоб не блокувати головний потік.

// ── Мінімальні допоміжні функції для 3D-векторів ───────────────
const add  = (a, b) => ({ x: a.x+b.x, y: a.y+b.y, z: a.z+b.z });
const sub  = (a, b) => ({ x: a.x-b.x, y: a.y-b.y, z: a.z-b.z });
const scale = (v, s) => ({ x: v.x*s,   y: v.y*s,   z: v.z*s });
const dot  = (a, b) => a.x*b.x + a.y*b.y + a.z*b.z;
const cross = (a,b) => ({ x:a.y*b.z-a.z*b.y, y:a.z*b.x-a.x*b.z, z:a.x*b.y-a.y*b.x });
const length = v => Math.sqrt(dot(v, v));
const normalise = v => scale(v, 1 / (length(v) || 1e-9));
const mulV = (a,b) => ({ x:a.x*b.x, y:a.y*b.y, z:a.z*b.z });

// ── Сцена ──────────────────────────────────────────────────────
const scene = {
  spheres: [
    { center:{x:0,y:0,z:-5},   radius:1,   mat:{albedo:{x:0.9,y:0.2,z:0.2}, ka:0.05, kd:0.7, ks:0.4, shininess:40, reflectivity:0.3} },
    { center:{x:-2.5,y:0,z:-6}, radius:1,   mat:{albedo:{x:0.2,y:0.6,z:0.9}, ka:0.05, kd:0.7, ks:0.6, shininess:80, reflectivity:0.0} },
    { center:{x: 2.5,y:0,z:-6}, radius:1,   mat:{albedo:{x:0.9,y:0.9,z:0.9}, ka:0.05, kd:0.1, ks:0.9, shininess:120, reflectivity:0.9} },
    { center:{x:0,y:-101,z:-5}, radius:100, mat:{albedo:{x:0.6,y:0.6,z:0.6}, ka:0.05, kd:0.8, ks:0.1, shininess:10,  reflectivity:0.0} },
  ],
  lights: [
    { pos:{x:5,y:8,z:-2},  colour:{x:1,y:0.95,z:0.85} },
    { pos:{x:-4,y:3,z:-4}, colour:{x:0.3,y:0.4,z:0.8} },
  ],
};

// ── Тест сцени (BVH замінив би це) ─────────────────────────────
function closestHit(ray, maxDist = Infinity) {
  let best = null;
  for (const s of scene.spheres) {
    const h = intersectSphere(ray, s);
    if (h && h.t < (best?.t ?? maxDist)) best = h;
  }
  return best;
}

// ── Рекурсивне трасування ──────────────────────────────────────
const MAX_DEPTH = 4;
const SKY = {x:0.05,y:0.08,z:0.13};

function traceRay(ray, depth = 0) {
  if (depth >= MAX_DEPTH) return SKY;
  const hit = closestHit(ray);
  if (!hit) return SKY;

  const col = shade(hit, ray, scene, depth);

  if (hit.mat.reflectivity > 0) {
    const rDir = reflect(ray.dir, hit.normal);
    const rCol = traceRay({ origin: add(hit.point, scale(hit.normal, 0.001)), dir: rDir }, depth + 1);
    return add(scale(col, 1 - hit.mat.reflectivity), scale(rCol, hit.mat.reflectivity));
  }
  return col;
}

// ── Рендеринг на canvas ────────────────────────────────────────
function render(canvas) {
  const ctx  = canvas.getContext('2d');
  const img  = ctx.createImageData(canvas.width, canvas.height);
  const d    = img.data;
  const cam  = {
    position: {x:0,y:1,z:1},
    forward:  normalise({x:0,y:-0.1,z:-1}),
    right:    {x:1,y:0,z:0},
    up:       {x:0,y:1,z:0},
    width: canvas.width, height: canvas.height,
    tanFov: Math.tan(Math.PI / 4),
    aspect: canvas.width / canvas.height,
  };

  for (let py = 0; py < canvas.height; py++) {
    for (let px = 0; px < canvas.width; px++) {
      const ray = generateRay(px, py, cam);
      const col = traceRay(ray);
      const i   = (py * canvas.width + px) * 4;
      d[i]   = Math.min(col.x * 255, 255);
      d[i+1] = Math.min(col.y * 255, 255);
      d[i+2] = Math.min(col.z * 255, 255);
      d[i+3] = 255;
    }
  }
  ctx.putImageData(img, 0, 0);
}