寫在前面

這兩天微信打算換個頭像，于是打算弄個光追的小渲染器來渲染新頭像，，，下面是結果：

在這里插入圖片描述

多貼幾張不同渲染滴圖片：

在這里插入圖片描述

要實作光線追蹤可不容易，最近我也在網上看了一些他人的博客和資料，基本上要么是籠統的理論介紹，貼幾個高大上的數學公式，或者是 ppt 的截圖，然后就是你復制我我復制你的內容，最后貼幾張渲染的結果，我看完他們寫的之后基本都是一臉懵 β 不知所措，我太菜了

…

于是我打算自己寫一篇博客，以【把大象裝進冰箱】進行比喻：如果說網路上的資料是【打開門，送大象進去，關門】的話，那本篇博客就是討論【怎么裝】的，比如冰箱門把手在哪里，朝哪個方向拉才能打開冰箱，要用香蕉還是鞭子把大象趕進去，怎么關上冰箱的門…

批話不多說，直接開始正文的內容!

光線追蹤簡介

在傳統的計算機圖形學中，通常使用光柵化的方法來生成像素，因為光柵化之后我們只選取我們看到的圖形，視野之外的幾何資訊被裁剪，導致場景的全域資訊丟失，那么就會造成一些渲染上的 “不正確”，比如水面只能反射我們螢屏看到的像素：

此外，一些基于全域資訊的渲染效果不能很好的運行，比如全域光照，在光柵管線下我們往往需要通過各種近似演算法去模擬真實世界的光影現象，非常麻煩并且效果 “還湊合”

光線追蹤超越了傳統的圖形流水線，是一種現代的渲染，光線追蹤是一種基于物理的渲染方法，通過模擬光在介質中的真實表現來輸出足以逼近現實的影像，光線從光源出發，經過各種反射折射，進入攝像機：

因為從光源出發有 無數條 到達相機的路線，正向地從燈泡開始模擬光的傳播變得不是很現實，于是曲線救國，轉而從相機方向開始向場景投射射線，試圖找出到達光源的可能路徑：

本篇博客代碼主要實作的是路徑追蹤技術，路徑追蹤技術是光線追蹤的一個分支，通過向場景中投射光線，并且模擬光線的行徑，在物體和光源之間找出一條條可行的路徑，最后回傳積累的顏色，

和光柵化不同，光線追蹤的第一步是向世界空間中投射光線，這一步叫做 RayCasting，投射的光線碰到哪個物體，我們對應位置的像素就代表該物體，同時要給對應位置的像素賦以一定的顏色：

在這里插入圖片描述

假設我們發射的光線命中了一個點，我們就回傳這個點的顏色作為最終像素的顏色，

那么怎樣描述一點的顏色呢？我們引入渲染方程 ↓

渲染方程

故事要從 ~~元和二年~~ 1986 年開始說起，科學家首次提出渲染方程，以描述一個點受到的光照的影響，渲染一個像素的顏色就是求解渲染方程的程序，大致意思是這樣的：一個點的光照由 2 部分組成，分別是：

該點自己發出的光
來自任意方向上的光打到該點積累形成的光

注意這個任意方向，允許該點接收來自其 法向半球 范圍內方向的入射光：

那么想要計算來自法向半球內任意方向的入射光的積累，我們必須計算一個 半球積分！ 這里給出渲染方程的簡化形式：

L ( p ) = E ( p ) + ∫ L ( q ) ? c o s ( θ ) d ω i L\left(p\right)=E\left(p\right)+\int_{} L\left(q\right) * cos(\theta) \ \mathrm{d} \omega_{i} L(p)=E(p)+∫?L(q)?cos(θ) dωi?

其中符號的解釋如下：

L ( x ) → x 點的光強度 E ( x ) → x 點發出的光 q → 從 p 點出發 , 方向為 ω i 的光命中了 q 點的物體 θ → ω i 與 q 點法向量的夾角 \begin{array}{l} L(x) \rightarrow x \ 點的光強度 \\ E(x) \rightarrow x \ 點發出的光 \\ q \rightarrow 從\ p \ 點出發,方向為 \ \omega_{i} \ 的光命中了 \ q \ 點的物體 \\ \theta \rightarrow \omega_{i} \ 與 \ q \ 點法向量的夾角 \end{array} L(x)→x 點的光強度E(x)→x 點發出的光q→從 p 點出發,方向為 ωi? 的光命中了 q 點的物體θ→ωi? 與 q 點法向量的夾角?

此外， ∫ \int_{} ∫? 就是對法向半球入（出）射光線方向的積分，注意到渲染方程的遞回形式，要想求解 p 點的光照我們必須先遞回計算 q 點的光照值，如圖：

遞回式告訴我們：q 點也反射了光到 p 點，反射光的強度等于 q 點的光強度乘以反射到 p 點的百分比 cos(θ）

了解了渲染方程的成分，就可以進行渲染方程的求解了，在高等數學中，積分的計算需要找到被積函式的原函式，和積分變數，可是渲染方程是一個困難積分，無法精確地計算其原函式，于是需要找尋其他方法對困難積分進行計算，

那么怎樣計算一個困難積分的？我們引入蒙特卡洛方法 ↓

蒙特卡洛方法

我更愿意稱之為丟豆子，考慮最直觀的情況，欲求區間上一函式的積分，我們往 x 區間上面丟豆子，并且計算豆子命中的位置的 y 的值，最后把他們加起來作為積分的估計：

丟豆子的程序稱之為【采樣】，如果我們使用 均勻分布 的 x 進行丟豆子，就能得到上圖等寬度的柱狀圖近似，

事實上在實際問題中，豆子的位置不會總是服從均勻分布，那么每一個豆子的貢獻，除了豆子命中位置的 y 值，還取決于豆子 命中該位置的概率，蒙特卡洛方法允許我們使用 x 的任意的概率密度函式，來對積分進行估計，

假設采樣 x 的概率密度函式為 P D F ( x ) PDF(x) PDF(x)，被積函式為 f ( x ) f(x) f(x)，那么 x 點采樣的貢獻就是：

f ( x ) P D F ( x ) \frac{f(x)}{PDF(x)} PDF(x)f(x)?

所以積分計算就很簡單了，首先按照產生一堆符合概率密度函式 PDF 分布的隨機變數 x，然后對每一個 x 都計算 f ( x ) P D F ( x ) \frac{f(x)}{PDF(x)} PDF(x)f(x)? 最后求他們的均值即可，現在回過頭來看剛剛的均勻分布的丟豆子，其中 P D F ( x ) = 1 b ? a PDF(x) = \frac{1}{ b-a} PDF(x)=b?a1?，那么我們估計 x 2 x^2 x2 的積分就可以這么計算：

在這里插入圖片描述

可以看到僅 5 次采樣就可以獲得還不錯的結果，我們和真實的積分值十分逼近了！采樣的次數越多，差異就越少，所以蒙特卡洛方法可以做到對積分結果的無偏估計，這是好特性，

知曉了困難積分的近似求解方式，我們開始正式求解渲染方程 ↓

渲染方程求解偽代碼

渲染方程是對于被求解點 p 的法向半球的積分，那么我們的隨機變數就是在法向半球內的射線的方向，假設我們取法向半球內 均勻分布的射線 方向，那么就有 P D F ( x ) = 1 2 π PDF(x) = \frac{1}{2\pi} PDF(x)=2π1?，因為半球面積就是 2 π 2\pi 2π，如圖：

在這里插入圖片描述

于是有如下的偽代碼：

pathTracing(p)
{
	L = 0.0
	for(i in SAMPLE)
	{
		wi = random()	// 隨機選擇一個方向
		if(wi hit q)	// 射線 wi 擊中 q 點
		{
			cosine = dot(nq, wi)	// q 點法向量 nq 和 wi 的夾角余弦值
			L += cosine * pathTracing(q) / PDF(wi)
		}
	}
	return L / SAMPLE	// 回傳均值
}

其中 SAMPLE 就是我們的采樣次數，那么問題就來了，這個遞回的復雜度是指數，復雜度非常高帶來的就是極大的計算資源的消耗，因為光線會有爆炸般的增長，以 SAMPLE=100 為例：

那么我們只有限制 SAMPLE=1 才能防止指數增長，而一次采樣的估計肯定是不準確的，于是我們對每個像素，發射多條光線，然后平均他們的結果，每個像素的光線數目叫做 SPP，即（Sample Pre Pixel），下圖演示了 SPP=3 的情況，我們找尋了 3 條到光源的路徑，并且平均他們的貢獻：

注意只有在第一次采樣時發射若干條光線，其余的時候我們只隨機選取一個方向發射光線并且遞回計算，那么偽代碼就改成：

// 追蹤一條光線
pathTracing(p)
{
	L = 0.0
	wi = random()	// 隨機選擇一個方向
	if(wi hit q)	// 射線 wi 擊中 q 點
	{
		cosine = dot(nq, wi)	// q 點法向量 nq 和 wi 的夾角余弦值
		L += cosine * pathTracing(q) / PDF(wi)
	}

	return L
}

// 對一個像素投射 SPP 條光線
L = 0.0
for(i in SPP)
{
	wi = random()	// 隨機選擇一個方向
	if(wi hit q)	// 射線 wi 擊中 q 點
	{
		L += pathTracing(q) / PDF(wi)
	}
}
L /= SPP

這一步也沒啥特別的，就是向每一個像素投射光線，然后求解渲染方程，沒了，，，

渲染方程的偽代碼有了，我們通過 c++ 實作它 ↓

編程前的準備

著手撰寫一個在 Windows 10 下運行的 x64 程式，程式以圖片的形式輸出場景的渲染結果，我們以 Vusial Studio 2019 作為 IDE，此外我們還需要額外的幫助庫，

數學運算庫

首先是數學運算庫，我們需要一個能夠表示 三維向量，并且對向量進行加減乘除點積叉乘等操作的幫助庫，你可以自己寫一個簡易的 class，也可以使用現成的第三方庫，這里我使用的是 glm，它的網站在這里，你也可以從它的 GitHub 上面獲取，此外，也可以通過 vcpkg 包管理工具來下載，只需要運行命令：

vcpkg install glm

如果在安裝時遇到任何困難，可以參考我以前的博客：傳送門①，傳送門②

影像輸出

你可以使用任何流行的影像處理的庫來進行影像輸出，他們可以是 Opencv，Qt，甚至是 OpenGL，但是這里我使用非常輕量級的 svpng，svpng 不是一個 c++ 的第三方庫，它僅是一個 inc檔案：

你只需要把它放在你的工程目錄下，然后再 #include "svpng.inc" 即可呼叫它，svpng 就一個非常簡單的功能，就可以幫我們保存 png 影像，呼叫 svpng 函式即可，函式的原型長這樣：

void svpng(FILE* fp, unsigned w, unsigned h, const unsigned char* img, int alpha)

其中 FILE 是檔案指標，w 和 h 是圖片的寬度和高度，img 是影像的像素值陣列，alpha 是透明度，我們填 0 即可，通過如下的代碼就可以將一個范圍在 [0, 1] 之間的 double 浮點數 buffer 輸出到圖片上：

// 輸出 SRC 陣列中的資料到影像
void imshow(double* SRC)
{
    
    unsigned char* image = new unsigned char[WIDTH * HEIGHT * 3];// 影像buffer
    unsigned char* p = image;
    double* S = SRC;    // 源資料

    FILE* fp;
    fopen_s(&fp, "image.png", "wb");

    for (int i = 0; i < HEIGHT; i++)
    {
        for (int j = 0; j < WIDTH; j++)
        {
            *p++ = (unsigned char)clamp((*S++) * 255, 0.0, 255.0);  // R 通道
            *p++ = (unsigned char)clamp((*S++) * 255, 0.0, 255.0);  // G 通道
            *p++ = (unsigned char)clamp((*S++) * 255, 0.0, 255.0);  // B 通道
        }
    }

    svpng(fp, WIDTH, HEIGHT, image, 0);
}

clamp 是截斷函式，glm 庫帶的，如果報錯那么您可以刪掉它并且換成您自己的 clamp，只是時刻注意 SRC 是 [0, 1] 范圍的 double，我們習慣這么表示顏色，同時方便計算，不容易被截斷，此外，svpng 默認影像的 RGB 通道是相鄰的，我們直接利用指標進行遍歷即可，

隨便在 SRC 中寫點什么，比如輸出 xy 的值作為 rg 通道，如果你看到如下的圖片被生成，那么很成功！

在這里插入圖片描述

如果找不到 svpng.inc 那么檢查你的 vs 工程是否配置正確，將包含 svpng 的目錄添加到 vs 的 include 目錄：

在這里插入圖片描述

多執行緒加速

光線追蹤運算量巨大，單靠簡單的單執行緒程式無法高效執行，但是因為 每個光線的采樣是相互獨立的，于是我們可以利用多執行緒加速，Visual Studio 有自帶多執行緒加速的 openmp 庫，無需自己手動下載，只需要引入：

#include <omp.h>    // openmp多執行緒加速

同時在專案設定中，允許 vs 使用多執行緒：

在這里插入圖片描述

然后在需要并行執行的 for 回圈之前加上：

omp_set_num_threads(50); // 執行緒個數
#pragma omp parallel for

for()
{
	...
}

就可以享受多執行緒加速的福利了，此外，我墻裂建議你打開 O2 優化同時將運行模式調整到 Release，以獲取最大運行速度：

在這里插入圖片描述

一切就緒？我們準備進入光與粒的世界 ↓

相機配置與光線生成

光線追蹤的第一步是投射光線，我們模擬相機投影與成像的規則，指定一個 [-1, 1] 范圍內的投影平面和一個視點，然后根據輸出圖片的像素位置，計算其對應投影平面上的坐標，最后用坐標減去視點坐標，得到 視線的方向向量，如圖：

值得注意的是，圖片的 xy 軸原點是在圖片左上方，而實際投影我們需要一個在左下方的原點（即平面幾何坐標系），所以 y 要做一次 flip，此外，在世界坐標系下，我們確定相機的位置和投影平面的位置，讓相機看向 z 軸負方向：

相機配置就緒，我們嘗試 輸出相機的射線投射方向，其中 imshow 是上面撰寫的顯示圖片的函式：

double* image = new double[WIDTH * HEIGHT * 3];
memset(image, 0.0, sizeof(double) * WIDTH * HEIGHT * 3);
double* p = image;
for (int i = 0; i < HEIGHT; i++)
{
    for (int j = 0; j < WIDTH; j++)
    {
        // 像素坐標轉投影平面坐標
        double x = 2.0 * double(j) / double(WIDTH) - 1.0;
        double y = 2.0 * double(HEIGHT - i) / double(HEIGHT) - 1.0;

        vec3 coord = vec3(x, y, SCREEN_Z);          // 計算投影平面坐標
        vec3 direction = normalize(coord - EYE);    // 計算光線投射方向
        
        vec3 color = direction;

        *p = color.x; p++;  // R 通道
        *p = color.y; p++;  // G 通道
        *p = color.z; p++;  // B 通道
    }
}

imshow(image);

輸出如下結果說明相機有在作業：
在這里插入圖片描述

相機投射了光線，光線和場景物體相交，我們希望描述這一程序 ↓

三角形與光線求交

在計算機圖形學中通常使用三角形來描述任意形狀的物體，因為三角形具有很好的幾何特征，并且易于進行求交，裁剪等操作，在開始之前，我們先做一些規范化的定義：

結構定義

假設我們用起點（start）和方向（direction）來描述一個射線：

// 光線
typedef struct Ray
{
    vec3 startPoint = vec3(0, 0, 0);    // 起點
    vec3 direction = vec3(0, 0, 0);     // 方向
}Ray;

在開始撰寫三角形類之前，我們先確定求交操作到底要回傳那些資訊：

是否相交
交點位置，用于作為我們下一次彈射的起點
相交位置的表面屬性：比如法向量，表面顏色，材質屬性，發光度，粗糙度等

那么首先我們有表面屬性的定義，使用結構體能很好的幫我們組織資料，同時易于拓展：

// 物體表面材質定義
typedef struct Material
{
    bool isEmissive = false;        // 是否發光
    vec3 normal = vec3(0, 0, 0);    // 法向量
    vec3 color = vec3(0, 0, 0);     // 顏色
}Material;

然后是光線求交結果的定義：

// 光線求交結果
typedef struct HitResult
{
    bool isHit = false;             // 是否命中
    double distance = 0.0f;         // 與交點的距離
    vec3 hitPoint = vec3(0, 0, 0);  // 光線命中點
    Material material;              // 命中點的表面材質
}HitResult;

然后是三角形的 class 的定義：

class Shape
{
public:
    Shape(){}
    virtual HitResult intersect(Ray ray) { return HitResult(); }
};

// 三角形
class Triangle : public Shape
{
public:
    Triangle(){}
    Triangle(vec3 P1, vec3 P2, vec3 P3, vec3 C) 
    { 
        p1 = P1, p2 = P2, p3 = P3; 
        material.normal = normalize(cross(p2 - p1, p3 - p1)); material.color = C;
    }
    vec3 p1, p2, p3;    // 三頂點
    Material material;  // 材質

    // 與光線求交
    HitResult intersect(Ray ray) 
    { 
        HitResult res;
        // ...
        return res; 
    };
};

這里我墻裂建議使用 虛函式+指標+繼承 的編程習慣，因為我們光線和任意圖形求交，都有一致的回傳結果，即 HitResult 結構體，我們使 c++ 的指標特性，可以通過一套代碼，完成多種復雜圖元的求交，此外，在添加一種新圖元的時候，主代碼不需要任何的改動！（雖然現在我們只有三角形，，，

求交計算

和三角形的求交分為兩個步驟：首先判斷光線和三角形所在平面是否有交點，然后再判斷交點是否在三角形內部，思路很清晰，我們先來看光線與三角形所在平面相交的判斷：

在這里插入圖片描述

其中 t 表示了射線起點到交點的距離，如果 t 小于 0 那么表示三角形在攝像機之后！ ，然后開始判斷點是否在三角形中，我們連接頂點與 P 點，然后判斷連線與邊的叉乘方向是否與 法向量 一致，如果三個頂點的判斷都通過，說明 P 在三角形中，否則不在：

在這里插入圖片描述

注意此處法向量 N 是通過 p 1 p 2 ? × p 1 p 3 ? \vec{p_1p_2} \times \vec{p_1p_3} p1?p2? ?×p1?p3? ? 得到的，上圖的 N 垂直螢屏向外，于是有 Triangle 類的求交代碼：

// 與光線求交
HitResult intersect(Ray ray) 
{ 
    HitResult res;

    vec3 S = ray.startPoint;        // 射線起點
    vec3 d = ray.direction;         // 射線方向
    vec3 N = material.normal;       // 法向量
    if (dot(N, d) > 0.0f) N = -N;   // 獲取正確的法向量

    // 如果視線和三角形平行
    if (fabs(dot(N, d)) < 0.00001f) return res;

    // 距離
    float t = (dot(N, p1) - dot(S, N)) / dot(d, N);
    if (t < 0.0005f) return res;    // 如果三角形在相機背面

    // 交點計算
    vec3 P = S + d * t;

    // 判斷交點是否在三角形中
    vec3 c1 = cross(p2 - p1, P - p1);
    vec3 c2 = cross(p3 - p2, P - p2);
    vec3 c3 = cross(p1 - p3, P - p3);
    vec3 n = material.normal;   // 需要 "原生法向量" 來判斷
    if (dot(c1, n) < 0 || dot(c2, n) < 0 || dot(c3, n) < 0) return res;

    // 裝填回傳結果
    res.isHit = true;
    res.distance = t;
    res.hitPoint = P;
    res.material = material;
    res.material.normal = N;    // 要回傳正確的法向
    return res; 
};

代碼比想象中的要復雜，因為我們考慮了更多的因數，比如視線和三角形平行，此外，法向量可能與我們視線的方向相同，那么三角形是背向我們的，.于是要改寫法向 N 的方向：

此外，判斷點是否在三角形中，我們要用 改寫之前 的法向量來計算，因為原生法向量取決于頂點定義的順序（p1, p2, p3），我們也是按照 p1, p2, p3 的順序來進行叉乘的

注：
你可能注意到上面代碼的那個 if (t < 0.0005f) 了
為啥是 0.0005 呢？是為了防止在三角形 T 上彈射的光線再次命中三角形 T （我打我自己
因為浮點精度不足，交點可能出現在原三角形的里側，那么彈射時就會自己打到自己，如圖：

光線投射與求交都準備就緒，我們開始第一個光線追蹤程式 ↓

Triangle, Again

我們在場景中添加一個三角形：

const vec3 RED(1, 0.5, 0.5);

...

vector<Shape*> shapes;  // 幾何物體的集合
shapes.push_back(new Triangle(vec3(-0.5, -0.5, 0), vec3(0.0, 0.5, 0), vec3(0.5, -0.5, 0), RED));

隨后遍歷場景，逐個求交并且回傳 hit 的結果，我們撰寫一幫助函式來回傳最近距離的交點及其屬性：

// 回傳距離最近 hit 的結果
HitResult shoot(vector<Shape*>& shapes, Ray ray)
{
    HitResult res, r;
    res.distance = 1145141919.810f; // inf

    // 遍歷所有圖形，求最近交點
    for (auto& shape : shapes)
    {
        r = shape->intersect(ray);
        if (r.isHit && r.distance < res.distance) res = r;  // 記錄距離最近的求交結果
    }

    return res;
}

然后我們逐像素的投射光線，并且輸出光線第一個碰到的交點的顏色：

double* image = new double[WIDTH * HEIGHT * 3];
memset(image, 0.0, sizeof(double) * WIDTH * HEIGHT * 3);
double* p = image;
for (int i = 0; i < HEIGHT; i++)
{
    for (int j = 0; j < WIDTH; j++)
    {
        // 像素坐標轉投影平面坐標
        double x = 2.0 * double(j) / double(WIDTH) - 1.0;
        double y = 2.0 * double(HEIGHT - i) / double(HEIGHT) - 1.0;

        vec3 coord = vec3(x, y, SCREEN_Z);          // 計算投影平面坐標
        vec3 direction = normalize(coord - EYE);    // 計算光線投射方向

        // 生成光線
        Ray ray;
        ray.startPoint = coord;
        ray.direction = direction;

        // 找交點并輸出交點的顏色
        HitResult res = shoot(shapes, ray);
        vec3 color = res.material.color;

        *p = color.x; p++;  // R 通道
        *p = color.y; p++;  // G 通道
        *p = color.z; p++;  // B 通道
    }
}

imshow(image);

程式輸出了如下的圖片：
在這里插入圖片描述
這和 Learn OpenGL 中的 hello world 一致，只是我們以另一種方式來實作它！此外，最好嘗試多種三角形的組合，比如遮擋關系或輸出法向量，以便我們 debug，要 確保求交程式的正確性

第一個 hello world 程式運行良好，現在來準備進行光線追蹤 ↓

球面隨機向量

現在還差最后一塊拼圖：亂數，在渲染方程的求解中，我們在法向半球上隨機選取一個方向作為光線的彈射方向，首先獲取一個 [0-1] 范圍的隨機浮點數：

// 0-1 亂數生成
std::uniform_real_distribution<> dis(0.0, 1.0);
random_device rd;
mt19937 gen(rd());
double randf()
{
    return dis(gen);
}

然后我們隨機生成 3 個坐標 xyz，如果坐標在單位球內，我們拒絕并且重新選取 xyz，從而產生 均勻分布 的球面隨機向量：

// 單位球內的隨機向量
vec3 randomVec3()
{
    vec3 d;
    do
    {
        d = 2.0f * vec3(randf(), randf(), randf()) - vec3(1, 1, 1);
    } while (dot(d, d) > 1.0);
    return normalize(d);
}

在此之后我們還要根據碰撞點的表面，生成分布在法向半球的隨機向量，一種可行的策略是使用仍然拒絕法，一旦隨機向量不在法向半球內，我們就拒絕它，同時再產生一個新的隨機向量，代碼如下：

// 法向半球隨機向量
vec3 randomDirection(vec3 n)
{
    // 法向半球
    vec3 d;
    do
    {
        d = randomVec3();
    } while (dot(d, n) < 0.0f);
    return d;
}

值得注意的是，在 Peter Shirley 寫的《Ray Tracing in One Weekend Book Series》系列中，有一種更加簡潔的求法向半球隨機向量的方法，就是以法向量的終點為球心，產生單位球面上的隨機向量，然后連接法向量起點和隨機向量的終點就是最終的隨機方向 d，如圖：

在這里插入圖片描述
寫成代碼也很簡單：

// 法向半球隨機向量
vec3 randomDirection(vec3 n)
{
    return normalize(randomVec3() + n);
}

這種方法很奇怪，但是不得不承認它的效果更好，并且更快，下圖是兩種隨機方法的比較：

在這里插入圖片描述

個人猜測是因為方法 2 生成的向量在 法向量 的方向上具有更高的出現概率：

在這里插入圖片描述
最后一塊拼圖已經被補齊了，正式開始光線追蹤 ↓

路徑追蹤，干就完了！

還記得剛剛我們的 hello world 是如何輸出顏色的嗎？我們直接輸出了碰到的物體的顏色，接下來我們改變策略，對碰到的物體，我們要求其渲染方程下的顏色，我們定義一個函式，它接收整個場景的資訊，和一條光線，然后根據路徑追蹤，回傳該光線最侄訓累的顏色：

// 路徑追蹤
vec3 pathTracing(vector<Shape*>& shapes, Ray ray)
{
    ...
    return xxx;
}

直接光照

我們考慮直接光照，如果射線碰到光源，我們回傳光源的顏色，否則我們回傳純黑：

于是路徑追蹤的函式變得非常簡單，僅 4 行代碼：

// 路徑追蹤
vec3 pathTracing(vector<Shape*>& shapes, Ray ray)
{
    HitResult res = shoot(shapes, ray);

    if (!res.isHit) return vec3(0); // 未命中

    // 如果發光則回傳顏色
    if (res.material.isEmissive) return res.material.color;

    // 否則直接回傳
    return vec3(0);
}

然后我們修改主函式，我們添加一些三角形：

// 采樣次數
const int SAMPLE = 128;
// 每次采樣的亮度
const double BRIGHTNESS = (2.0f * 3.1415926f) * (1.0f / double(SAMPLE));

...

vector<Shape*> shapes;  // 幾何物體的集合
// 三角形
shapes.push_back(new Triangle(vec3(-0.5, -0.5, -0.5), vec3(0.5, -0.5, -0.5), vec3(0, -0.5, 0.5), CYAN));
// 底部平面
shapes.push_back(new Triangle(vec3(10, -1, 10), vec3(-10, -1, -10), vec3(-10, -1, 10), WHITE));
shapes.push_back(new Triangle(vec3(10, -1, 10), vec3(10, -1, -10), vec3(-10, -1, -10), WHITE));
// 光源
Triangle l1 = Triangle(vec3(0.6, 0.99, 0.4), vec3(-0.2, 0.99, -0.4), vec3(-0.2, 0.99, 0.4), WHITE);
Triangle l2 = Triangle(vec3(0.6, 0.99, 0.4), vec3(0.6, 0.99, -0.4), vec3(-0.2, 0.99, -0.4), WHITE);
l1.material.isEmissive = true;
l2.material.isEmissive = true;
shapes.push_back(&l1);
shapes.push_back(&l2);

場景長這樣：上方的是光源三角形組成的四邊形，中間的淡藍色三角形是不發光的物體，而底部則是一個很大的平面，

在這里插入圖片描述

然后使用多執行緒進行采樣，最后使其輸出直接光照下的場景：

double* image = new double[WIDTH * HEIGHT * 3];
memset(image, 0.0, sizeof(double) * WIDTH * HEIGHT * 3);

omp_set_num_threads(50); // 執行緒個數
#pragma omp parallel for
for (int k = 0; k < SAMPLE; k++)
{
    double* p = image;
    for (int i = 0; i < HEIGHT; i++)
    {
        for (int j = 0; j < WIDTH; j++)
        {
            // 像素坐標轉投影平面坐標
            double x = 2.0 * double(j) / double(WIDTH) - 1.0;
            double y = 2.0 * double(HEIGHT - i) / double(HEIGHT) - 1.0;

            vec3 coord = vec3(x, y, SCREEN_Z);          // 計算投影平面坐標
            vec3 direction = normalize(coord - EYE);    // 計算光線投射方向

            // 生成光線
            Ray ray;
            ray.startPoint = coord;
            ray.direction = direction;

            // 與場景的交點
            HitResult res = shoot(shapes, ray);
            vec3 color = vec3(0, 0, 0);

            if (res.isHit)
            {
                // 命中光源直接回傳光源顏色
                if (res.material.isEmissive)
                {
                    color = res.material.color;
                }
                // 命中物體則選擇一個隨機方向重新發射光線并且進行路徑追蹤
                else
                {
                    // 根據交點處法向量生成交點處反射的隨機半球向量
                    Ray randomRay;
                    randomRay.startPoint = res.hitPoint;
                    randomRay.direction = randomDirection(res.material.normal);
                    
                    // 顏色積累
                    vec3 srcColor = res.material.color;
                    vec3 ptColor = pathTracing(shapes, reflectRay);
                    color = ptColor * srcColor;    // 和原顏色混合
                    color *= BRIGHTNESS;
                }
            }

            *p += color.x; p++;  // R 通道
            *p += color.y; p++;  // G 通道
            *p += color.z; p++;  // B 通道
        }
    }
}

imshow(image);

可以看到輸出的圖片具有非常真實的特性（盡管有很多噪點）：

在這里插入圖片描述
可以看到在僅有直接光照的情況下，就能夠實作非常多的特效，比如光照，軟陰影，并且是基于物理的！而這些特效在傳統光柵管線中都是代價及其昂貴的特效，此外，增大每個像素的采樣次數（SPP，代碼中的 SAMPLE 引數）能提升品質：

在這里插入圖片描述

間接光照

直接光照僅考慮了渲染方程的自發光項，事實上除了來自光源的直接光照，還有來自其他物體反射的光：

因為渲染方程已經給了間接光照的計算公式，我們直接遞回計算，但是值得注意的是遞回的出口，我們可以簡單的使用一個遞回深度 depth 來控制，更加巧妙的方法是每次搖一個亂數 P，如果 P 小于某個閾值就結束，這個方法有一個很霸氣的名字，叫做毛子輪盤：

但是我們仍然要通過深度保證不會出現死遞回，此外每次回傳時應該將顏色除以 P 以保證顏色的 期望值 始終不變，于是有：

// 路徑追蹤
vec3 pathTracing(vector<Shape*>& shapes, Ray ray, int depth)
{
    if (depth > 8) return vec3(0);
    HitResult res = shoot(shapes, ray);

    if (!res.isHit) return vec3(0); // 未命中

    // 如果發光則回傳顏色
    if (res.material.isEmissive) return res.material.color;

    // 有 P 的概率終止
    double r = randf();
    float P = 0.8;
    if (r > P) return vec3(0);

    // 否則繼續
    Ray randomRay;
    randomRay.startPoint = res.hitPoint;
    randomRay.direction = randomDirection(res.material.normal);

    float cosine = dot(-ray.direction, res.material.normal);
    vec3 srcColor = res.material.color;
    vec3 ptColor = pathTracing(shapes, randomRay, depth+1) * cosine;
    vec3 color = ptColor * srcColor;    // 和原顏色混合

    return color / P;
}

我們重新布置一下場景，讓三角形的背光面朝向鏡頭，現在間接光照能夠被加入計算了：

在這里插入圖片描述

注意和直接光照最大的區別就是背光面能夠被間接地照亮，注意上圖 左右兩個三角形的背面分別反射了紅藍兩種顏色，而僅有直接光照的渲染是無法照亮背光面的：

在這里插入圖片描述

這里為了看的更清楚，順手加了個伽馬矯正，

到這里路徑追蹤的基本框架就已經完成了，下面是一些錦上添花的效果 ↓

繪制球體

我們試圖往場景中添加球體，球體的求交也比較簡單，不多說了，上車：

在這里插入圖片描述
這里給出 Sphere 類代碼：

// 球
class Sphere : public Shape
{
public:
    Sphere(){}
    Sphere(vec3 o, double r, vec3 c) { O = o; R = r; material.color = c; }
    vec3 O;             // 圓心
    double R;           // 半徑
    Material material;  // 材質

    // 與光線求交
    HitResult intersect(Ray ray)
    {
        HitResult res;

        vec3 S = ray.startPoint;        // 射線起點
        vec3 d = ray.direction;         // 射線方向

        float OS = length(O - S);
        float SH = dot(O - S, d);
        float OH = sqrt(pow(OS, 2) - pow(SH, 2));

        if (OH > R) return res; // OH大于半徑則不相交

        float PH = sqrt(pow(R, 2) - pow(OH, 2));

        float t1 = length(SH) - PH;
        float t2 = length(SH) + PH;
        float t = (t1 < 0) ? (t2) : (t1);   // 最近距離
        vec3 P = S + t * d;     // 交點

        // 防止自己交自己
        if (fabs(t1) < 0.0005f || fabs(t2) < 0.0005f) return res;

        // 裝填回傳結果
        res.isHit = true;
        res.distance = t;
        res.hitPoint = P;
        res.material = material;
        res.material.normal = normalize(P - O); // 要回傳正確的法向
        return res;
    }
};

一樣的上面的 0.0005f 也是為了防止自己交自己的情況出現，撒，這時候你就會發現為啥一開始我要不厭其煩地定義一個 Shape 類并且使用基類指標與虛函式，因為這允許你 添加任何的圖形而不用改動主代碼，比如我們加三個球：

vector<Shape*> shapes;  // 幾何物體的集合

// 球
shapes.push_back(new Sphere(vec3(-0.6, -0.8, 0.6), 0.2, WHITE));
shapes.push_back(new Sphere(vec3(-0.1, -0.7, 0.2), 0.3, WHITE));
shapes.push_back(new Sphere(vec3(0.5, -0.6, -0.5), 0.4, WHITE));

重新運行程式并且得到一個驚艷的結果：

在這里插入圖片描述
好吧，，，我終于知道為啥幾乎所有的光線追蹤的博客都會拿球來做例子，因為確實求交簡單，而且表現力非常強！

鏡面反射

光打到材質上，有一部分發生漫反射，有一部分發生鏡面反射，這取決于材質的屬性，我們在材質中新定義一個屬性叫做反射率 s，

// 物體表面材質定義
typedef struct Material
{
    ...
    double specularRate = 0.0f;      // 反射光占比
}Material;

入射光有 s 的概率被反射，否則繼續漫反射，那么我們的光追要用如下的流程：先搖一個亂數，如果其小于反射率，那么我們光線被反射，于是通過入射光和法向量的夾角計算反射光線，并且繼續遞回，否則我們正常地隨機取一個方向投射光線，

于是修改路徑追蹤的函式，注意主代碼中 ray casting 投射光線的時候也要做同樣的修改：

// 路徑追蹤
vec3 pathTracing(vector<Shape*>& shapes, Ray ray, int depth)
{
	前半部分和之前一樣
	
    ...
    
    Ray randomRay;
    randomRay.startPoint = res.hitPoint;
    randomRay.direction = randomDirection(res.material.normal);
    
    vec3 color = vec3(0);
    float cosine = fabs(dot(-ray.direction, res.material.normal));

    // 根據反射率決定光線最終的方向
    r = randf();
    if (r < res.material.specularRate)  // 鏡面反射
    {
        randomRay.direction = normalize(reflect(ray.direction, res.material.normal));
        color = pathTracing(shapes, randomRay, depth + 1) * cosine;
    }
    else    // 漫反射
    {
        vec3 srcColor = res.material.color;
        vec3 ptColor = pathTracing(shapes, randomRay, depth+1) * cosine;
        color = ptColor * srcColor;    // 和原顏色混合
    }

    return color / P;
}

在主代碼中添加對應不同反射率的球，然后重新運行程式，可以看到反射非常漂亮，而且可以看到螢屏空間內不存在的內容：

在這里插入圖片描述
我們也可修改反射率，下圖從左到右分別是 0.3，0.6，0.9 的反射率：

在這里插入圖片描述
注意到我們的球非常光滑，如果想模擬粗糙的反射，那么在生成反射光線方向的時候，加入隨機向量的擾動即可，而擾動的程度取決于材質的粗糙度，這也是一個材質屬性，我們加入它：

// 物體表面材質定義
typedef struct Material
{
    ...
    double roughness = 1.0f;        // 粗糙程度
}Material;

然后我們反射的時候不再按照反射光線的方向，而是根據粗糙度，在隨機向量和反射光線的方向做一個 線性插值 以決定最終反射的方向：

if (r < res.material.specularRate)  // 鏡面反射
{
    vec3 ref = normalize(reflect(ray.direction, res.material.normal));
    randomRay.direction = mix(ref, randomRay.direction, res.material.roughness);
    color = pathTracing(shapes, randomRay, depth + 1) * cosine;
}

注意主代碼中 ray casting 投射光線的部分也要做同樣的改動，下面分別是 0.9 反射率下，0.2，0.4，0.6 粗糙度下的效果，當然你也可以嘗試自己自定義這些引數：

在這里插入圖片描述

折射

光線通過介質發生折射，折射角取決于入射方向和物體表面法線，和反射類似，我們直接計算折射角即可，折射也有發生的概率，我們在材質結構體中添加一些欄位：

// 物體表面材質定義
typedef struct Material
{
    ...
    double refractRate = 0.0f;      // 折射光占比
    double refractAngle = 1.0f;     // 折射率
    double refractRoughness = 0.0f; // 折射粗糙度
}Material;

值得注意的是我們的概率計算：當亂數小于 reflectRate 的時候發生反射，亂數在 reflectRate 和 refractRate 之間發生折射，亂數大于 refractRate 的時候才是漫反射：

然后 path tracing 中對應的代碼也要添加多一個 if elsl 選項，同時主函式中投射光線時第一次彈射也要做同樣的修改，此外，和鏡面反射一樣，我們也加入了粗糙度來表示折射的擾動，但是擾動的方向是法向半球的 負方向，看代碼：

if (r < res.material.specularRate)  // 鏡面反射
{
    ...
}
else if (res.material.specularRate <= r && r <= res.material.refractRate)    // 折射
{
    vec3 ref = normalize(refract(ray.direction, res.material.normal, float(res.material.refractAngle)));
    randomRay.direction = mix(ref, -randomRay.direction, res.material.refractRoughness);
    color = pathTracing(shapes, randomRay, depth + 1) * cosine;
}
else    // 漫反射
{
    ...
}

結果非常漂亮，我們可以獲得很多免費的效果，比如地上的焦散光斑：

在這里插入圖片描述

注：
這段代碼嚴格意義上是錯的，因為沒有考慮射入球和射出球的兩種不同的情況
上圖使用 0.1 的折射角（refract 函式的 eta 引數）
因為我沒有學過光學相關的課程，我也不知道這個引數該怎么取，事實上我是隨便取的

通過改變粗糙度來調整折射光，這里使用 1.0 的折射角，表示光不發生折射，而是直接穿透，從左到右分別是 0.1，0.2，0.3 的粗糙度，我們可以得到毛玻璃的效果：

在這里插入圖片描述

抗鋸齒

你可能會發現場景中有很多鋸齒：
在這里插入圖片描述
在光線追蹤渲染器中使用抗鋸齒非常簡單，我們可以在發射光線的時候，在光線的方向上加一個小的偏移量，以實作一個像素多個方向的采樣，就好比光柵管線里面的 MSAA 一樣：

// MSAA
x += (randf() - 0.5f) / double(WIDTH);
y += (randf() - 0.5f) / double(HEIGHT);

效果：

在這里插入圖片描述

完整代碼

#include <iostream>
#include <vector>
#include <random>
#include <stdlib.h>
#include <glm/glm.hpp>  // 數學庫支持
#include "svpng.inc"    // png輸出 ref: https://github.com/miloyip/svpng
#include <omp.h>    // openmp多執行緒加速

using namespace glm;
using namespace std;

// --------------------- end of include --------------------- //

// 采樣次數
const int SAMPLE = 4096;

// 每次采樣的亮度
const double BRIGHTNESS = (2.0f * 3.1415926f) * (1.0f / double(SAMPLE));

// 輸出影像解析度
const int WIDTH = 256;
const int HEIGHT = 256;

// 相機引數
const double SCREEN_Z = 1.1;        // 視平面 z 坐標
const vec3 EYE = vec3(0, 0, 4.0);   // 相機位置

// 顏色
const vec3 RED(1, 0.5, 0.5);
const vec3 GREEN(0.5, 1, 0.5);
const vec3 BLUE(0.5, 0.5, 1);
const vec3 YELLOW(1.0, 1.0, 0.1);
const vec3 CYAN(0.1, 1.0, 1.0);
const vec3 MAGENTA(1.0, 0.1, 1.0);
const vec3 GRAY(0.5, 0.5, 0.5);
const vec3 WHITE(1, 1, 1);

// --------------- end of global variable definition --------------- //

// 光線
typedef struct Ray
{
    vec3 startPoint = vec3(0, 0, 0);    // 起點
    vec3 direction = vec3(0, 0, 0);     // 方向
}Ray;

// 物體表面材質定義
typedef struct Material
{
    bool isEmissive = false;        // 是否發光
    vec3 normal = vec3(0, 0, 0);    // 法向量
    vec3 color = vec3(0, 0, 0);     // 顏色
    double specularRate = 0.0f;     // 反射光占比
    double roughness = 1.0f;        // 粗糙程度
    double refractRate = 0.0f;      // 折射光占比
    double refractAngle = 1.0f;     // 折射率
    double refractRoughness = 0.0f; // 折射粗糙度
}Material;

// 光線求交結果
typedef struct HitResult
{
    bool isHit = false;             // 是否命中
    double distance = 0.0f;         // 與交點的距離
    vec3 hitPoint = vec3(0, 0, 0);  // 光線命中點
    Material material;              // 命中點的表面材質
}HitResult;

class Shape
{
public:
    Shape(){}
    virtual HitResult intersect(Ray ray) { return HitResult(); }
};

// 三角形
class Triangle : public Shape
{
public:
    Triangle(){}
    Triangle(vec3 P1, vec3 P2, vec3 P3, vec3 C) 
    { 
        p1 = P1, p2 = P2, p3 = P3; 
        material.normal = normalize(cross(p2 - p1, p3 - p1)); material.color = C;
    }
    vec3 p1, p2, p3;    // 三頂點
    Material material;  // 材質

    // 與光線求交
    HitResult intersect(Ray ray) 
    { 
        HitResult res;

        vec3 S = ray.startPoint;        // 射線起點
        vec3 d = ray.direction;         // 射線方向
        vec3 N = material.normal;       // 法向量
        if (dot(N, d) > 0.0f) N = -N;   // 獲取正確的法向量

        // 如果視線和三角形平行
        if (fabs(dot(N, d)) < 0.00001f) return res;

        // 距離
        float t = (dot(N, p1) - dot(S, N)) / dot(d, N);
        if (t < 0.0005f) return res;    // 如果三角形在相機背面

        // 交點計算
        vec3 P = S + d * t;

        // 判斷交點是否在三角形中
        vec3 c1 = cross(p2 - p1, P - p1);
        vec3 c2 = cross(p3 - p2, P - p2);
        vec3 c3 = cross(p1 - p3, P - p3);
        vec3 n = material.normal;   // 需要 "原生法向量" 來判斷
        if (dot(c1, n) < 0 || dot(c2, n) < 0 || dot(c3, n) < 0) return res;

        // 裝填回傳結果
        res.isHit = true;
        res.distance = t;
        res.hitPoint = P;
        res.material = material;
        res.material.normal = N;    // 要回傳正確的法向
        return res; 
    };
};

// 球
class Sphere : public Shape
{
public:
    Sphere(){}
    Sphere(vec3 o, double r, vec3 c) { O = o; R = r; material.color = c; }
    vec3 O;             // 圓心
    double R;           // 半徑
    Material material;  // 材質

    // 與光線求交
    HitResult intersect(Ray ray)
    {
        HitResult res;

        vec3 S = ray.startPoint;        // 射線起點
        vec3 d = ray.direction;         // 射線方向

        float OS = length(O - S);
        float SH = dot(O - S, d);
        float OH = sqrt(pow(OS, 2) - pow(SH, 2));

        if (OH > R) return res; // OH大于半徑則不相交

        float PH = sqrt(pow(R, 2) - pow(OH, 2));

        float t1 = length(SH) - PH;
        float t2 = length(SH) + PH;
        float t = (t1 < 0) ? (t2) : (t1);   // 最近距離
        vec3 P = S + t * d;     // 交點

        // 防止自己交自己
        if (fabs(t1) < 0.0005f || fabs(t2) < 0.0005f) return res;

        // 裝填回傳結果
        res.isHit = true;
        res.distance = t;
        res.hitPoint = P;
        res.material = material;
        res.material.normal = normalize(P - O); // 要回傳正確的法向
        return res;
    }
};

// ---------------------------- end of class definition ---------------------------- //

// 輸出 SRC 陣列中的資料到影像
void imshow(double* SRC)
{
    
    unsigned char* image = new unsigned char[WIDTH * HEIGHT * 3];// 影像buffer
    unsigned char* p = image;
    double* S = SRC;    // 源資料

    FILE* fp;
    fopen_s(&fp, "image.png", "wb");

    for (int i = 0; i < HEIGHT; i++)
    {
        for (int j = 0; j < WIDTH; j++)
        {
            *p++ = (unsigned char)clamp(pow(*S++, 1.0f / 2.2f) * 255, 0.0, 255.0);  // R 通道
            *p++ = (unsigned char)clamp(pow(*S++, 1.0f / 2.2f) * 255, 0.0, 255.0);  // G 通道
            *p++ = (unsigned char)clamp(pow(*S++, 1.0f / 2.2f) * 255, 0.0, 255.0);  // B 通道
        }
    }

    svpng(fp, WIDTH, HEIGHT, image, 0);
}

// 回傳距離最近 hit 的結果
HitResult shoot(vector<Shape*>& shapes, Ray ray)
{
    HitResult res, r;
    res.distance = 1145141919.810f; // inf

    // 遍歷所有圖形，求最近交點
    for (auto& shape : shapes)
    {
        r = shape->intersect(ray);
        if (r.isHit && r.distance < res.distance) res = r;  // 記錄距離最近的求交結果
    }

    return res;
}

// 0-1 亂數生成
std::uniform_real_distribution<> dis(0.0, 1.0);
random_device rd;
mt19937 gen(rd());
double randf()
{
    return dis(gen);
}

// 單位球內的隨機向量
vec3 randomVec3()
{
    
    vec3 d;
    do
    {
        d = 2.0f * vec3(randf(), randf(), randf()) - vec3(1, 1, 1);
    } while (dot(d, d) > 1.0);
    return normalize(d);
    /*
    double r1 = randf(), r2 = randf();
    double z = sqrt(1.0f - r2);
    double phi = 2 * 3.1415926 * r1;
    float x = cos(phi) * sqrt(r2);
    float y = sin(phi) * sqrt(r2);
    return normalize(vec3(x, y, z));
    */
}

// 法向半球隨機向量
vec3 randomDirection(vec3 n)
{
    /*
    // 法向半球
    vec3 d;
    do
    {
        d = randomVec3();
    } while (dot(d, n) < 0.0f);
    return d;
    */
    // 法向球
    return normalize(randomVec3() + n);
}

// 路徑追蹤
vec3 pathTracing(vector<Shape*>& shapes, Ray ray, int depth)
{
    if (depth > 8) return vec3(0);
    HitResult res = shoot(shapes, ray);

    if (!res.isHit) return vec3(0); // 未命中

    // 如果發光則回傳顏色
    if (res.material.isEmissive) return res.material.color;
   
    // 有 P 的概率終止
    double r = randf();
    float P = 0.8;
    if (r > P) return vec3(0);
    
    // 否則繼續
    Ray randomRay;
    randomRay.startPoint = res.hitPoint;
    randomRay.direction = randomDirection(res.material.normal);
    
    vec3 color = vec3(0);
    float cosine = fabs(dot(-ray.direction, res.material.normal));

    // 根據反射率決定光線最終的方向
    r = randf();
    if (r < res.material.specularRate)  // 鏡面反射
    {
        vec3 ref = normalize(reflect(ray.direction, res.material.normal));
        randomRay.direction = mix(ref, randomRay.direction, res.material.roughness);
        color = pathTracing(shapes, randomRay, depth + 1) * cosine;
    }
    else if (res.material.specularRate <= r && r <= res.material.refractRate)    // 折射
    {
        vec3 ref = normalize(refract(ray.direction, res.material.normal, float(res.material.refractAngle)));
        randomRay.direction = mix(ref, -randomRay.direction, res.material.refractRoughness);
        color = pathTracing(shapes, randomRay, depth + 1) * cosine;
    }
    else    // 漫反射
    {
        vec3 srcColor = res.material.color;
        vec3 ptColor = pathTracing(shapes, randomRay, depth+1) * cosine;
        color = ptColor * srcColor;    // 和原顏色混合
    }

    return color / P;
}

// ---------------------------- end of functions ---------------------------- //

int main()
{
    vector<Shape*> shapes;  // 幾何物體的集合

    Sphere s1 = Sphere(vec3(-0.65, -0.7, 0.0), 0.3, GREEN);
    Sphere s2 = Sphere(vec3(0.0, -0.3, 0.0), 0.4, WHITE);
    Sphere s3 = Sphere(vec3(0.65, 0.1, 0.0), 0.3, BLUE);
    s1.material.specularRate = 0.3;
    s1.material.roughness = 0.1;

    s2.material.specularRate = 0.3;
    s2.material.refractRate = 0.95;
    s2.material.refractAngle = 0.1;
    //s2.material.refractRoughness = 0.05;

    s3.material.specularRate = 0.3;

    shapes.push_back(&s1);
    shapes.push_back(&s2);
    shapes.push_back(&s3);

    shapes.push_back(new Triangle(vec3(-0.15, 0.4, -0.6), vec3(-0.15, -0.95, -0.6), vec3(0.15, 0.4, -0.6), YELLOW));
    shapes.push_back(new Triangle(vec3(0.15, 0.4, -0.6), vec3(-0.15, -0.95, -0.6), vec3(0.15, -0.95, -0.6), YELLOW));

    Triangle tt = Triangle(vec3(-0.2, -0.2, -0.95), vec3(0.2, -0.2, -0.95), vec3(-0.0, -0.9, 0.4), YELLOW);
    //tt.material.specularRate = 0.1;
    //tt.material.refractRate = 0.85;
    //tt.material.refractRoughness = 0.3;
    //shapes.push_back(&tt);
    
    // 發光物
    Triangle l1 = Triangle(vec3(0.4, 0.99, 0.4), vec3(-0.4, 0.99, -0.4), vec3(-0.4, 0.99, 0.4), WHITE);
    Triangle l2 = Triangle(vec3(0.4, 0.99, 0.4), vec3(0.4, 0.99, -0.4), vec3(-0.4, 0.99, -0.4), WHITE);
    l1.material.isEmissive = true;
    l2.material.isEmissive = true;
    shapes.push_back(&l1);
    shapes.push_back(&l2);

    // 背景盒子
    // bottom
    shapes.push_back(new Triangle(vec3(1, -1, 1), vec3(-1, -1, -1), vec3(-1, -1, 1), WHITE));
    shapes.push_back(new Triangle(vec3(1, -1, 1), vec3(1, -1, -1), vec3(-1, -1, -1), WHITE));
    // top
    shapes.push_back(new Triangle(vec3(1, 1, 1), vec3(-1, 1, 1), vec3(-1, 1, -1), WHITE));
    shapes.push_back(new Triangle(vec3(1, 1, 1), vec3(-1, 1, -1), vec3(1, 1, -1), WHITE));
    // back
    shapes.push_back(new Triangle(vec3(1, -1, -1), vec3(-1, 1, -1), vec3(-1, -1, -1), CYAN));
    shapes.push_back(new Triangle(vec3(1, -1, -1), vec3(1, 1, -1), vec3(-1, 1, -1), CYAN));
    // left
    shapes.push_back(new Triangle(vec3(-1, -1, -1), vec3(-1, 1, 1), vec3(-1, -1, 1), BLUE));
    shapes.push_back(new Triangle(vec3(-1, -1, -1), vec3(-1, 1, -1), vec3(-1, 1, 1), BLUE));
    // right
    shapes.push_back(new Triangle(vec3(1, 1, 1), vec3(1, -1, -1), vec3(1, -1, 1), RED));
    shapes.push_back(new Triangle(vec3(1, -1, -1), vec3(1, 1, 1), vec3(1, 1, -1), RED));
    
    
    double* image = new double[WIDTH * HEIGHT * 3];
    memset(image, 0.0, sizeof(double) * WIDTH * HEIGHT * 3);

    omp_set_num_threads(50); // 執行緒個數
    #pragma omp parallel for
    for (int k = 0; k < SAMPLE; k++)
    {
        double* p = image;
        for (int i = 0; i < HEIGHT; i++)
        {
            for (int j = 0; j < WIDTH; j++)
            {
                // 像素坐標轉投影平面坐標
                double x = 2.0 * double(j) / double(WIDTH) - 1.0;
                double y = 2.0 * double(HEIGHT - i) / double(HEIGHT) - 1.0;

                // MSAA
                x += (randf() - 0.5f) / double(WIDTH);
                y += (randf() - 0.5f) / double(HEIGHT);

                vec3 coord = vec3(x, y, SCREEN_Z);          // 計算投影平面坐標
                vec3 direction = normalize(coord - EYE);    // 計算光線投射方向

                // 生成光線
                Ray ray;
                ray.startPoint = coord;
                ray.direction = direction;

                // 與場景的交點
                HitResult res = shoot(shapes, ray);
                vec3 color = vec3(0, 0, 0);

                if (res.isHit)
                {
                    // 命中光源直接回傳光源顏色
                    if (res.material.isEmissive)
                    {
                        color = res.material.color;
                    }
                    // 命中物體則選擇一個隨機方向重新發射光線并且進行路徑追蹤
                    else
                    {
                        // 根據交點處法向量生成交點處反射的隨機半球向量
                        Ray randomRay;
                        randomRay.startPoint = res.hitPoint;
                        randomRay.direction = randomDirection(res.material.normal);

                        // 根據反射率決定光線最終的方向
                        double r = randf();
                        if (r < res.material.specularRate)  // 鏡面反射
                        {
                            vec3 ref = normalize(reflect(ray.direction, res.material.normal));
                            randomRay.direction = mix(ref, randomRay.direction, res.material.roughness);
                            color = pathTracing(shapes, randomRay, 0);
                        }
                        else if (res.material.specularRate <= r && r <= res.material.refractRate)    // 折射
                        {
                            vec3 ref = normalize(refract(ray.direction, res.material.normal, float(res.material.refractAngle)));
                            randomRay.direction = mix(ref, -randomRay.direction, res.material.refractRoughness);
                            color = pathTracing(shapes, randomRay, 0);
                        }
                        else    // 漫反射
                        {
                            vec3 srcColor = res.material.color;
                            vec3 ptColor = pathTracing(shapes, randomRay, 0);
                            color = ptColor * srcColor;    // 和原顏色混合
                        }
                        color *= BRIGHTNESS;
                    }
                }

                *p += color.x; p++;  // R 通道
                *p += color.y; p++;  // G 通道
                *p += color.z; p++;  // B 通道
            }
        }
    }
    
    imshow(image);

    return 0;
}

后記與總結

hhh 這是我第一次寫這么長的文章，可能寫的太多有點啰嗦了，主要是原理的部分過于拖沓，，，不過這樣也有好處，畢竟能讓自己更清楚的了解每個像素背后的一切，是通過精妙的演算法和物理公式得到的，而不是魔理沙的月球魔法，，，

立個 flag 這是我今年寫的最認真的文章了（事實上我一一共寫了 3 版代碼，博客放的是最終版），通過撰寫這篇文章我們得到了一個簡易的光線追蹤離線渲染器，誠然它有很多不足，但是起碼通過它我們能對光線追蹤這一新的世界有一個管中窺豹的認識（翻譯：擺爛了我不寫了）

光追所面臨的種種問題，老黃在 18 年 RTX 新顯卡發布的時候就給出了非常棒的解決方案，涉及到非常多的優化方法及復雜的數學原理，能夠讓光追邁入實時渲染，這不，mc 的光追版已經出來了嘛（可惜我買不起顯卡），新的方法不斷出現，規范也在被逐漸完善，所以說，學習的路還很長啊，，，

在這里插入圖片描述

轉載請註明出處，本文鏈接：https://www.uj5u.com/qita/256317.html

標籤：AI

上一篇：人工智能本身也是一種新型基礎設施

下一篇：作業經驗｜lambada處理集合的常用10種實戰騷操作，我都記錄下來了

光線追蹤渲染實戰：蒙特卡洛路徑追蹤及其c++實作

目錄

寫在前面

光線追蹤簡介

渲染方程

蒙特卡洛方法

渲染方程求解偽代碼

編程前的準備

數學運算庫

影像輸出

多執行緒加速

相機配置與光線生成

三角形與光線求交

結構定義

求交計算

Triangle, Again

球面隨機向量

路徑追蹤，干就完了！

直接光照

間接光照

繪制球體

鏡面反射

折射

抗鋸齒

完整代碼

后記與總結