FPGA開発日記

カテゴリ別記事インデックス https://msyksphinz.github.io/github_pages , English Version https://fpgadevdiary.hatenadiary.com/

レイトレーシングアルゴリズムについて基本を学ぶ(3. 実ソースコードの読み解き)

前回オリジナルLLVMバックエンドを使ってレイトレーシングのプログラムをコンパイルして動かしたのだが、実は「レイトレーシングとは何なのか」良く知らない。コンピュータグラフィックスの授業などは殆ど取ったことが無く、これまでの経験でも使用することが無かったため良く分からない。そこで資料を読んで勉強してみることにした。

  • Introduction to Ray Tracing: a Simple Method for Creating 3D Images

https://www.scratchapixel.com/lessons/3d-basic-rendering/introduction-to-ray-tracing

参考にしたのはまさにレイトレーシングのプログラムを拝借してきたサイトから。レイトレーシングアルゴリズムを丁寧に解説してある。これを読んでみよう。 結構長いのでしっかりGoogle Translateに活躍してもらった。Google Translateの結果を読みながら、要点をつまみ出していく。

また、画像も分かりやすいので本文中から引用させて頂いた。

今回は実ソースコードをCode Readingしてみる。一度独自LLVMコンパイルさせたコードではあるが、中身を読まずにコンパイルしたので一応概要を理解しておきたい。


ソースコードの読み解き

今回はコンパイルとシミュレーションに使用したレイトレースプログラムの読み解きを行ってみる。これまでの解説を読んだ知識で読み解くことができるか。

  • Vec3クラス:頂点を示すクラス。例えば球体の中心を表す場合などに使用する。
template<typename T>
class Vec3
{
public:
    T x, y, z;
    Vec3() : x(T(0)), y(T(0)), z(T(0)) {}
    Vec3(T xx) : x(xx), y(xx), z(xx) {}
    Vec3(T xx, T yy, T zz) : x(xx), y(yy), z(zz) {}
...
  • Sphereクラス。球体を表現する。centerVec3fクラスで表現され、あとは半径、半径の2乗、surfacecolorは表面の色(RGBのはずだがこれもVec3fで表現するのか)、放射色、透明度と反射率を含んでいる。
class Sphere
{
public:
    Vec3f center;                           /// position of the sphere
    float radius, radius2;                  /// sphere radius and radius^2
    Vec3f surfaceColor, emissionColor;      /// surface color and emission (light)
    float transparency, reflection;         /// surface transparency and reflectivity
...
  • main関数。やっていることは大きく2つに分けられる。球体の登録と光源の登録、最後にレンダリングである。ecに0以外を設定すると光源となり、そこから発光される。
//[comment]
// In the main function, we will create the scene which is composed of 5 spheres
// and 1 light (which is also a sphere). Then, once the scene description is complete
// we render that scene, by calling the render() function.
//[/comment]
int main(int argc, char **argv)
{
    srand48(13);
    std::vector<Sphere> spheres;
    // position, radius, surface color, reflectivity, transparency, emission color
    spheres.push_back(Sphere(Vec3f( 0.0f, -10004, -20), 10000, Vec3f(0.20f, 0.20f, 0.20f), 0, 0.0f));
    spheres.push_back(Sphere(Vec3f( 0.0f,      0, -20),     4, Vec3f(1.00f, 0.32f, 0.36f), 1, 0.5f));
    spheres.push_back(Sphere(Vec3f( 5.0f,     -1, -15),     2, Vec3f(0.90f, 0.76f, 0.46f), 1, 0.0f));
    spheres.push_back(Sphere(Vec3f( 5.0f,      0, -25),     3, Vec3f(0.65f, 0.77f, 0.97f), 1, 0.0f));
    spheres.push_back(Sphere(Vec3f(-5.5f,      0, -15),     3, Vec3f(0.90f, 0.90f, 0.90f), 1, 0.0f));
    // light
    spheres.push_back(Sphere(Vec3f( 0.0f,     20, -30),     3, Vec3f(0.00f, 0.00f, 0.00f), 0, 0.0f, Vec3f(3)));
    render(spheres);

    return 0;
}

ベクトルspheresに球体を5つ登録し、最後に光源を登録する。renderを呼び出しレンダリングを行う訳だ。

renderは大きく2つに分かれているらしい。光源のトレースをするためにすべての画素に対してtrace()を呼び出す部分と実際に画像を出力する部分だ。

    unsigned width = 640, height = 480;
    Vec3f *image = new Vec3f[width * height], *pixel = image;
    float invWidth = 1.0f / float(width), invHeight = 1.0f / float(height);
    float fov = 30, aspectratio = width / float(height);
    float angle = tanf(M_PI * 0.5f * fov / 180.0f);
    // Trace rays
    for (unsigned y = 0; y < height; ++y) {
        for (unsigned x = 0; x < width; ++x, ++pixel) {
            float xx = (2 * ((x + 0.5f) * invWidth) - 1) * angle * aspectratio;
            float yy = (1 - 2 * ((y + 0.5f) * invHeight)) * angle;
            Vec3f raydir(xx, yy, -1);
            raydir.normalize();
            trace(Vec3f(0), raydir, spheres, 0, pixel);
        }
    }

んんん?raydirというは何だ?視点からの方向を示しているのかなあ。。。とりあえずtrace()の中に入ってみる。

void trace(
    const Vec3f &rayorig,
    const Vec3f &raydir,
    const std::vector<Sphere> &spheres,
    const int &depth,
    Vec3f *ret_vec3f)
{

まず視点からスタートして、全ての球体に対して交点を探す。確か見つかった交点の中で一番視点に対して近いものを選択するのだったかな。

   // find intersection of this ray with the sphere in the scene
    for (unsigned i = 0; i < spheres.size(); ++i) {
        float t0 = INFINITY, t1 = INFINITY;
        if (spheres[i].intersect(rayorig, raydir, t0, t1)) {
            if (t0 < 0) t0 = t1;
            if (t0 < tnear) {
                tnear = t0;
                sphere = &spheres[i];
            }
        }
    }

もしどの交点にもぶつからなかったら、それは背景なのですぐに戻る。

    // if there's no intersection return black or background color
    if (!sphere) {
      *ret_vec3f = Vec3f(2);
      return;
    }

交差した物体が透過する物質または反射する物体ならば、さらに交点から新たな光の計算をしてtraceを再開しなければならない。

    if ((sphere->transparency > 0 || sphere->reflection > 0) && depth < MAX_RAY_DEPTH) {
        float facingratio = -raydir.dot(nhit);
        // change the mix value to tweak the effect
        float fresneleffect = mix(powf(1 - facingratio, 3), 1, 0.1f);
        // compute reflection direction (not need to normalize because all vectors
        // are already normalized)
        Vec3f refldir = raydir - nhit * 2 * raydir.dot(nhit);
        refldir.normalize();
        Vec3f reflection;
        trace(phit + nhit * bias, refldir, spheres, depth + 1, &reflection);
...

trace再帰的に呼び出しているが、これだといつまでたっても終わらない可能性があるのである程度まで再帰的に計算したら強制的に打ち切りにするためにdepth < MAX_RAY_DEPTHという条件が追加されている。

もうひとつは透過する場合の屈折の計算。これも同様にtrace()を呼び出して再帰的に計算する。

        // if the sphere is also transparent compute refraction ray (transmission)
        if (sphere->transparency) {
            float ior = 1.1, eta = (inside) ? ior : 1.0f / ior; // are we inside or outside the surface?
            float cosi = -nhit.dot(raydir);
            float k = 1.0f - eta * eta * (1.0f - cosi * cosi);
            Vec3f refrdir = raydir * eta + nhit * (eta *  cosi - sqrtf(k));
            refrdir.normalize();
            trace(phit - nhit * bias, refrdir, spheres, depth + 1, &refraction);
        }

最後に反射によって得た表面の色と屈折によって得た表面の色を混ぜ合わせて最終的なピクセルの色とする。

        // the result is a mix of reflection and refraction (if the sphere is transparent)
        surfaceColor = (
            reflection * fresneleffect +
            refraction * (1 - fresneleffect) * sphere->transparency) * sphere->surfaceColor;

屈折率と透過率が定義されている場合について計算したが、そうでない場合は単純。新たに光源からの光線を生成して対象となる物体に対して光を発射する。もし発射した途中で別の物体にぶつかった場合には、対象となる物体に対して光は当たらないので影となる。

        // it's a diffuse object, no need to raytrace any further
        for (unsigned i = 0; i < spheres.size(); ++i) {
            if (spheres[i].emissionColor.x > 0) {
                // this is a light
                Vec3f transmission = 1;
                Vec3f lightDirection = spheres[i].center - phit;
                lightDirection.normalize();
                for (unsigned j = 0; j < spheres.size(); ++j) {
                    if (i != j) {
                        float t0, t1;
                        if (spheres[j].intersect(phit + nhit * bias, lightDirection, t0, t1)) {
                            transmission = 0;
                            break;
                        }
                    }
                }
                surfaceColor += sphere->surfaceColor * transmission *
                std::max(float(0), nhit.dot(lightDirection)) * spheres[i].emissionColor;
            }
        }

なるほど、アルゴリズムの外観は分かった。レイトレースの基本的なプログラムの流れは理解できた気がする。

このアルゴリズムの高速化についてだが、全てピクセル単位での計算なのでまずはピクセル毎に並列に実行できそうな気がする。レイトレースの中身の演算でも、好転の計算などはすべての物体に対して独立並列に行われているので、並列化の要素はたくさんありそうだ。